High-order Finite Difference Method Research Articles

Global SH-wavefield calculation for a two-dimensional whole-Earth model with the parallel hybrid PSM/FDM algorithm

We present a parallel hybrid algorithm based on pseudospectral method (PSM) and finite difference method (FDM) for two-dimensional (2-D) global SH-wavefield simulation. The whole-Earth model is taken as a cross section of spherical Earth, and corresponding wave equations are defined in 2-D cylindrical coordinates. Spatial derivatives in the wave equations are approximated with efficient and high accuracy PSM in the lateral and high-order FDM in the radial direction on staggered grids. This algorithm allows us to divide the whole-Earth into sub-domains in radial direction and implement efficient parallel computing on PC cluster, while retains high accuracy and efficiency of PSM in lateral direction. A transformation of moment tensor between 3-D spherical Earth and our 2-D model was proposed to give corresponding moment tensor components used in 2-D modeling. Comparison of modeling results with those obtained by direct solution method shows very good accuracy of our algorithm. We also demonstrate its feasibility with a lateral heterogeneous whole-Earth model with localized velocity perturbation.

Smoothed Particle Method for Real-Space Electronic Structure Calculations

We present a meshfree particle method for real-space electronic structure calculations based on symmetric smoothed particle hydrodynamics (SSPH). As a simple example, we applied SSPH to the Schrödinger equation for a one-dimensional harmonic oscillator. The results calculated using SSPH are in good agreement with analytical solutions. We compared results for different cases of kernel function, smoothing length, number of particles, and the order of the Taylor series expansion. The results using SSPH and the higher-order finite-difference method are also compared. Our results indicate that SSPH can be applied to real-space electronic structure calculations that require high accuracy.

Viscoacoustic waveform inversion of velocity structures in the time domain

Because of the conversion of elastic energy into heat, seismic waves are attenuated and dispersed as they propagate. The attenuation effects can reduce the resolution of velocity models obtained from waveform inversion or even cause the inversion to produce incorrect results. Using a viscoacoustic model consisting of a single standard linear solid, we discovered a theoretical framework of viscoacoustic waveform inversion in the time domain for velocity estimation. We derived and found the viscoacoustic wave equations for forward modeling and their adjoint to compensate for the attenuation effects in viscoacoustic waveform inversion. The wave equations were numerically solved by high-order finite-difference methods on centered grids to extrapolate seismic wavefields. The finite-difference methods were implemented satisfying stability conditions, which are also presented. Numerical examples proved that the forward viscoacoustic wave equation can simulate attenuative behaviors very well in amplitude attenuation and phase dispersion. We tested acoustic and viscoacoustic waveform inversions with a modified Marmousi model and a 3D field data set from the deep-water Gulf of Mexico for comparison. The tests with the modified Marmousi model illustrated that the seismic attenuation can have large effects on waveform inversion and that choosing the most suitable inversion method was important to obtain the best inversion results for a specific seismic data volume. The tests with the field data set indicated that the inverted velocity models determined from the acoustic and viscoacoustic inversions were helpful to improve images and offset gathers obtained from migration. Compared to the acoustic inversion, viscoacoustic inversion is a realistic approach for real earth materials because the attenuation effects are compensated.

Stable, high order accurate adaptive schemes for long time, highly intermittent geophysics problems

Many geophysical phenomena are characterized by properties that evolve over a wide range of scales which introduce difficulties when attempting to model these features in one computational method. We have developed a high-order finite difference method for the elastic wave equation that is able to efficiently handle varying temporal and spatial scales in a single, stand-alone framework. We apply this method to earthquake cycle models characterized by extremely long interseismic periods interspersed with abrupt, short periods of dynamic rupture. Through the use of summation-by-parts operators and weak enforcement of boundary conditions we derive a provably stable discretization. Time stepping is achieved through the implicit θ-method which allows us to take large time steps during the intermittent period between earthquakes and adapts appropriately to fully resolve rupture.

An unconditionally stable compact ADI method for three-dimensional time-fractional convection–diffusion equation

A high-order compact finite difference method is presented for solving the three-dimensional (3D) time-fractional convection–diffusion equation (of order α∈(1,2)). The original equation is first transformed to a fractional diffusion-wave equation, then using fourth-order Padé approximation for spatial derivatives and the center difference method for time derivative respectively, a fully discrete implicit compact scheme is obtained. Furthermore, based on different splitting terms, three unconditionally stable ADI compact schemes with optimal convergence order are developed respectively. The resulting schemes in each ADI solution step corresponding to a strictly diagonally dominant matrix equation can be solved using the 1D tridiagonal Thomas algorithm with a considerable saving in computing time. Numerical experiments show that these schemes can significantly improve the time accuracy.

A Finite Difference Method of High Order Accuracy for the Solution of Two-Point Boundary Value Problems

We present a new high order finite difference method for second order differential equation y''=f(x,y) subject to boundary conditions y(a)=alpha and y(b)= beta. The method is based on rational function approximation and its development is based on power series expansions. Under appropriate conditions, local truncation error calculated and order of method estimated six. Our finite difference method leads to nonlinear system of equations. Numerical examples are given to illustrate the effectiveness, efficiency and high order accuracy of the method.

Review of summation-by-parts schemes for initial–boundary-value problems

High-order finite difference methods are efficient, easy to program, scale well in multiple dimensions and can be modified locally for various reasons (such as shock treatment for example). The main drawback has been the complicated and sometimes even mysterious stability treatment at boundaries and interfaces required for a stable scheme. The research on summation-by-parts operators and weak boundary conditions during the last 20 years has removed this drawback and now reached a mature state. It is now possible to construct stable and high order accurate multi-block finite difference schemes in a systematic building-block-like manner. In this paper we will review this development, point out the main contributions and speculate about the next lines of research in this area.

Viscoacoustic prestack reverse time migration based on the optimal time-space domain high-order finite-difference method

Prestack reverse time migration (RTM) is an accurate imaging method of subsurface media. The viscoacoustic prestack RTM is of practical significance because it considers the viscosity of the subsurface media. One of the steps of RTM is solving the wave equation and extrapolating the wave field forward and backward; therefore, solving accurately and efficiently the wave equation affects the imaging results and the efficiency of RTM. In this study, we use the optimal time-space domain dispersion high-order finite-difference (FD) method to solve the viscoacoustic wave equation. Dispersion analysis and numerical simulations show that the optimal time-space domain FD method is more accurate and suppresses the numerical dispersion. We use hybrid absorbing boundary conditions to handle the boundary reflection. We also use source-normalized cross-correlation imaging conditions for migration and apply Laplace filtering to remove the low-frequency noise. Numerical modeling suggests that the viscoacoustic wave equation RTM has higher imaging resolution than the acoustic wave equation RTM when the viscosity of the subsurface is considered. In addition, for the wave field extrapolation, we use the adaptive variable-length FD operator to calculate the spatial derivatives and improve the computational efficiency without compromising the accuracy of the numerical solution.

A Runge Kutta Discontinuous Galerkin approach to solve reactive flows on conforming hybrid grids: the parabolic and source operators

A Runge–Kutta Discontinuous Galerkin method (RKDG) to solve the parabolic and source parts of reactive Navier–Stokes equations written in conservation form is presented. The parabolic operator uses a recent recovery method set up by van Leer for structured grids and a new projection method proposed by Borrel–Ryan for unstructured grids. The physical model involves complex chemistry and detailed transport. Transport coefficients are evaluated using algorithms which provide empirical expressions. In 1-D test cases the RKDG method is compared with a high order finite difference method. 2-D test cases in structured, unstructured and hybrid meshes are presented.

Higher order finite difference method for the reaction and anomalous-diffusion equation

In this paper, our aim is to study the high order finite difference method for the reaction and anomalous-diffusion equation. According to the equivalence of the Riemann–Liouville and Grünwald–Letnikov derivatives under the suitable smooth condition, a second-order difference approximation for the Riemann–Liouville fractional derivative is derived. A fourth-order compact difference approximation for second-order derivative in spatial is used. We analyze the solvability, conditional stability and convergence of the proposed scheme by using the Fourier method. Then we obtain that the convergence order is O(τ2+h4), where τ is the temporal step length and h is the spatial step length. Finally, numerical experiments are presented to show that the numerical results are in good agreement with the theoretical analysis.

High-fidelity numerical solution of the time-dependent Dirac equation

A stable high-order accurate finite difference method for the time-dependent Dirac equation is derived. Grid-convergence studies in 1-D and 3-D corroborate the analysis. The method is applied to time-resolved quantum tunneling where a comparison with the solution to the time-dependent Schrödinger equation in 1-D illustrates the differences between the two equations. In contrast to the conventional tunneling probability decay predicted by the Schrödinger equation, the Dirac equation exhibits Klein tunneling. Solving the time-dependent Dirac equation with a step potential in 3-D reveals that particle spin affects the tunneling process. The observed spin-dependent reflection allows for a new type of spin-selective measurements.

Stability characteristic of hypersonic flow over a blunt wedge under freestream pulse wave

AbstractTo investigate the stability characteristic of hypersonic flow under the action of a freestream pulse wave, a high-order finite difference method was employed to do direction numerical simulation (DNS) of hypersonic unsteady flow over an 8° half-wedge-angle blunt wedge with freestream slow acoustic wave. The evolution of disturbance wave modes in the boundary layer under a pulse wave and a continuous wave are compared, and the wall temperature effect on the hypersonic boundary layer stability for a pulse wave disturbance is discussed. Results show that, both for a pulse wave and a continuous wave in freestream, the disturbance waves inside the nose boundary layer are mainly a fundamental mode; the Fourier amplitude of pressure disturbance mode in the boundary layer for a pulse wave is far less than that for a continuous wave, and the band frequency of the former is wider than that of the latter. All disturbance modes decay rapidly along the streamwise in the nose boundary layer. In the non-nose boundary layer, the dominant mode is transferred from fundamental mode into second harmonic. The transformation of dominant mode for a pulse wave appears much earlier than that for a continuous wave. Different frequency disturbance modes present different changes along streamline in the boundary layer, and the frequency band narrows around the second harmonic mode along the streamwise. Keen competition and the transformation of energy exist among different modes in the boundary layer. Wall temperature modifies the stability characteristic of the hypersonic boundary layer, which presents little effect on the development of fundamental modes and cooling wall could accelerates the growth of the high frequency mode as well as the dominant mode transformation.

Numerical Investigation of the Viscous Dissipation Term on 2D Heat Transfer

In this paper are analyzed, using high-order finite difference method, applications in which the viscous dissipation term can be neglected or not in the heat transfer equation. Some examples using various numerical values for the velocity field show that the viscous dissipation does not affect significantly the temperature field. Using the L2 norm, the numerical solution is compared with some examples that have an analytical solution.

High-order explicit-implicit numerical methods for nonlinear anomalous diffusion equations

In this paper, the high-order finite difference/element methods for the nonlinear anomalous diffusion equations of subdiffusion and superdiffusion are developed, where the high-order finite difference methods are used to approximate the time-fractional derivatives and the finite element methods are used in the spatial domain. The stability and error estimates are proved for both cases of superdiffusion and subdiffusion. Numerical examples are provided to confirm the theoretical analysis.

Summation-by-parts in time

We develop a new high order accurate time-integration technique for initial value problems. We focus on problems that originate from a space approximation using high order finite difference methods on summation-by-parts form with weak boundary conditions, and extend that technique to the time-domain. The new time-integration method is global, high order accurate, unconditionally stable and together with the approximation in space, it generates optimally sharp fully discrete energy estimates. In particular, it is shown how stable fully discrete high order accurate approximations of the Maxwells’ equations, the elastic wave equations and the linearized Euler and Navier–Stokes equations can obtained. Even though we focus on finite difference approximations, we stress that the methodology is completely general and suitable for all semi-discrete energy-stable approximations. Numerical experiments show that the new technique is very accurate and has limited order reduction for stiff problems.

A fully-coupled upwind discontinuous Galerkin method for incompressible porous media flows: High-order computations of viscous fingering instabilities in complex geometry

We present a new approach to the simulation of viscous fingering instabilities in incompressible, miscible displacement flows in porous media. In the past, high resolution computational simulations of viscous fingering instabilities have always been performed using high-order finite difference or Fourier-spectral methods which do not posses the flexibility to compute very complex subsurface geometries. Our approach, instead, by means of a fully-coupled nonlinear implementation of the discontinuous Galerkin method, possesses a fundamental differentiating feature, in that it maintains high-order accuracy on fully unstructured meshes. In addition, the proposed method shows very low sensitivity to mesh orientation, in contrast with classical finite volume approximation used in porous media flow simulations. The robustness and accuracy of the method are demonstrated in a number of challenging computational problems.

Prestack reverse-time migration with a time-space domain adaptive high-order staggered-grid finite-difference method

With advanced computational power, prestack reverse-time migration (RTM) is being used increasingly in seismic imaging. The accuracy and efficiency of RTM strongly depends on the algorithms used for numerical solutions of wave equations. Hence, how to solve the wave equation accurately and rapidly is very important in the process of RTM. In this paper, in order to improve the accuracy of the numerical solution, we use a time-space domain staggered-grid finite-difference (SFD) method to solve the acoustic wave equation, and develop a new acoustic prestack RTM scheme based on this time-space domain high-order SFD. Synthetic and real data tests demonstrate that the RTM scheme improves the imaging quality significantly compared with the conventional SFD RTM. Meanwhile, in the process of wavefield extrapolation, we apply adaptive variable-length spatial operators to compute spatial derivatives to decrease computational costs effectively with little reduction of the accuracy of the numerical solutions.

Numerical modeling for analyzing microfluidic acoustophoretic motion of cells and particles with application to identification of vibro-acoustic properties

Microfluidic, acoustophoretic cell/particle separation has gained significant interest recently. In order to analyze the motion of cells/particles in the acoustophoretic separation, a one-dimensional (1-D) analytical model in a “static” fluid medium has been widely used, while the effects of acoustic streaming, viscous boundary layers, and 2-D and 3-D geometries are usually not considered. Therefore, it is not sufficient to accurately predict the cell/particle motion. Thus, a numerical modeling procedure for analyzing the acoustophoretic cell/particle motion is presented to include the aforementioned effects. Here, the first-order acoustic pressure and the second-order acoustic streaming velocity are first calculated by using a high-order finite difference method. Then, acoustophoretic force is calculated based on the acoustophoretic force equation proposed by Gorkov and is applied to the Newton’s equation of motion to calculate the motion of cells/particles. Through various simulations, the effects of acoustic streaming on the motion of cells/particles are studied. Since the acoustophoretic motion depends on the vibro-acoustic properties (e.g., density, compressibility, and size) of particles/cells, the vibro-acoustic properties can be estimated by optimally fitting the experimental and simulated trajectories. The properties obtained from experimental results with polystyrene beads show good agreement with the data reported in literature.

Assessment of stretched vortex subgrid-scale models for LES of incompressible inhomogeneous turbulent flow

The physical space version of the stretched vortex subgrid scale model [Phys. Fluids 12, 1810 (2000)] is tested in large eddy simulations (LES) of the turbulent lid driven cubic cavity flow. LES is carried out using a higher order finite-difference method [J. Comput. Phys. 229, 8802 (2010)]. The effects of different vortex orientation models and subgrid turbulence spectrums are assessed through comparisons of the LES predictions against direct numerical simulations (DNS) [Phys. Fluids 12, 1363 (2000)]. Three Reynolds numbers 12000, 18000, and 22000 are studied. Good agreement with the DNS data for the mean and fluctuating quantities is observed.

Acoustic Prestack Reverse Time Migration Using the Adaptive High‐Order Finite‐Difference Method in the Time‐Space Domain

Abstract:At present, prestack reverse time migration (RTM) is the most accurate method for seismic imaging in theory. One of its key steps is wavefield extrapolation. The essence of wavefield extrapolation is to solve the wave equation. Thus solving the wave equation fast yet accurately is very important. In this paper, we use a time‐space domain dispersion‐relation‐based finite‐difference (FD) method to solve the acoustic wave equation, then analyze the dispersion and stability of the numerical solution, perform the numerical modeling, and compare the results with the conventional FD method. The results of analysis and numerical modeling demonstrate the time‐space domain FD method has better accuracy and stability than the conventional method under the same discretization. We adopt the time‐space domain high‐order FD method to solve the acoustic wave equation in the process of wavefield extrapolation, then use the normalized cross‐correlation imaging condition to obtain imaging in RTM. The results of synthetic data processing show that the imaging accuracy has enhanced. Meanwhile, in the process of wavefield extrapolation, we adopt adaptive variable‐length spatial operators to compute spatial derivatives and thus to decrease computing cost significantly without reducing the accuracy of numerical modeling and imaging.