First-order Spatial Derivatives Research Articles

A comparison of linear consistent correction methods for first-order SPH derivatives

A well-known issue with the widely used Smoothed Particle Hydrodynamics (SPH) method is the neighborhood deficiency. Near the surface, the SPH interpolant fails to accurately capture the underlying fields due to a lack of neighboring particles. These errors may introduce ghost forces or other visual artifacts into the simulation. In this work we investigate three different popular methods to correct the first-order spatial derivative SPH operators up to linear accuracy, namely the Kernel Gradient Correction (KGC), Moving Least Squares (MLS) and Reproducing Kernel Particle Method (RKPM). We provide a thorough, theoretical comparison in which we remark strong resemblance between the aforementioned methods. We support this by an analysis using synthetic test scenarios. Additionally, we apply the correction methods in simulations with boundary handling, viscosity, surface tension, vorticity and elastic solids to showcase the reduction or elimination of common numerical artifacts like ghost forces. Lastly, we show that incorporating the correction algorithms in a state-of-the-art SPH solver only incurs a negligible reduction in computational performance.

Large-time behavior of global solutions to the Cauchy problem of one-dimensional viscous and heat-conducting ionized gas

This paper is concerned with the Cauchy problem of one-dimensional viscous and heat-conducting ionized gas when the viscosity is a positive constant or depends on the density. According to Saha's ionization equation, the equations of state for the basic thermodynamic variables of an ionized gas depend also on the degree of ionization, which leads to the loss of concavity of the physical entropy in some small bounded domain. Such a property of the ionized gas together with the unboundedness of the domain under our consideration makes it hard to derive the desired dissipative estimates on the first-order spatial derivative of both the bulk velocity and the absolute temperature, thus making the problem more challenging.In this paper, for a class of constant non-vacuum equilibrium states with positive temperature, we succeed in deducing the desired uniform-in-time bounds on the above mentioned dissipative estimates on both the bulk velocity and the absolute temperature. Based on such an estimate, we can then establish the global existence of solutions to the Cauchy problem of one-dimensional viscous and heat-conducting ionized gas with positive constant or density-dependent viscosity and show that such a class of constant non-vacuum equilibrium states is time asymptotically nonlinear stable.

Global existence for a class of non-equilibrium reaction–diffusion systems with flux limitation

For a class of general reaction–diffusion systems with flux limitation in one dimension, when homogeneous Neumann boundary conditions are considered, we will prove global upper and lower boundness estimates for the positive classical solutions and boundness for their first-order spatial derivatives. In particular, the non-equilibrium radiation flux-limited diffusion coupled with material heat conduction model is included as a special case of the general reaction–diffusion systems. The key point of our method is to at first prove that the second-order spatial derivative of the solution of the second diffusion equation can be linearly controlled by the first-order spatial derivative of the solution of the flux-limited diffusion equation. Then, with this control, we can apply the Bernstein method to the flux-limited diffusion equation, while a stratifying time technique is introduced to get a complete first order spatial derivative estimate. With these global prior estimates, the global existence and uniqueness of positive classical solutions for the general flux-limited diffusion systems are proved.

Large Time Behavior of Solutions to the 3D Rotating Navier–Stokes Equations

We consider the large time behavior of the solutions for the initial value problem of the Navier–Stokes equations with the Coriolis force in the three-dimensional whole space. We show the L^p temporal decay estimates with the dispersion effect of the Coriolis force for the global solutions. Moreover, we prove the large time asymptotic expansion of the solutions behaving like the first-order spatial derivatives of the integral kernel of the corresponding linear solution.

Novel hybrid compact schemes for stream function–velocity formulation of the incompressible Navier–Stokes equations

In this paper, we have developed new hybrid compact schemes for the simulation of stream function–velocity formulation of a two-dimensional incompressible Navier–Stokes equation. The first-order spatial derivatives are approximated by an optimized upwind compact scheme, and the Laplacian and biharmonic operators are discretized using fourth-order hybrid compact schemes. Moreover, we have also performed Fourier analysis to assess the resolution and added numerical diffusion properties of numerical schemes for stream function–velocity formulation of the linear Navier–Stokes equation. For time discretization, we have used an explicit fourth-stage fourth-order Runge–Kutta method and hybrid filters. Furthermore, to validate the accuracy and efficiency of the schemes, several fluid flow problems, including a test problem with a non-homogeneous source term and a lid-driven cavity problem, are considered. Numerical results exhibit a great match to the results reported in the literature at lower computational cost with hybrid filters.

On Solving Initial Boundary Value Problems for a Parabolic Equation in Dini Classes

The method of potentials is one of the classical methods for solving initial boundary value problems for parabolic equations and systems. The method is based on studying the smoothness of parabolic potentials in various functional spaces and studying the solvability of the corresponding integral equations and systems to which the problems are reduced. This method allows one to find a constructive form of solutions to initial-boundary value problems and directly investigate the smoothness of solutions. This article continues a series of works on constructing the solvability theory of parabolic initial boundary value problems in domains with nonsmooth lateral boundaries in the time variable. The article considers the first and second initial boundary value problems with zero initial condition for a uniformly parabolic equation that is one-dimensional in the spatial variable. It is assumed that the equation coefficients are bounded and uniformly continuous with a modulus of continuity satisfying twice the Dini condition. The specific feature of the study is that the right side of the equation can grow to infinity in a certain way when approaching the initial data assignment line. The solvability (in the classical sense) of these problems is established in a semi-bounded curvilinear domain in a plane with a nonsmooth lateral boundary from the Dini-Hölder class, admitting, in particular, so-called "beaks". The solution is constructed as a sum of the plane parabolic potential and the single layer potential, the kernel of which is the fundamental solution of the equation. The smoothness of the solution of these problems in the Dini-Hölder classes is investigated, and the corresponding correctness estimates are obtained. These estimates characterize the behavior of solutions and their first-order spatial derivative in the closure of the domain. The results obtained can be used to study heat and mass transfer processes. For the proof, the smoothness in the Dini-Hölder space of a plane parabolic potential with an unlimited density of the Dini class in a curvilinear domain with a nonsmooth lateral boundary is investigated. The obtained result can be used in solving other initial boundary value problems for an inhomogeneous parabolic equation.

Wave Equation Numerical Simulation and RTM With Mixed Staggered-Grid Finite-Difference Schemes

For the conventional staggered-grid finite-difference scheme (C-SFD), although the spatial finite-difference (FD) operator can reach 2Mth-order accuracy, the FD discrete wave equation is the only second-order accuracy, leading to low modeling accuracy and poor stability. We proposed a new mixed staggered-grid finite-difference scheme (M-SFD) by constructing the spatial FD operator using axial and off-axial grid points jointly to approximate the first-order spatial partial derivative. This scheme is suitable for modeling the stress–velocity acoustic and elastic wave equation. Then, based on the time–space domain dispersion relation and the Taylor series expansion, we derived the analytical expression of the FD coefficients. Theoretically, the FD discrete acoustic wave equation and P- or S-wave in the FD discrete elastic wave equation given by M-SFD can reach the arbitrary even-order accuracy. For acoustic wave modeling, with almost identical computational costs, M-SFD can achieve higher modeling accuracy than C-SFD. Moreover, with a larger time step used in M-SFD than that used in C-SFD, M-SFD can achieve higher computational efficiency and reach higher modeling accuracy. For elastic wave simulation, compared to C-SFD, M-SFD can obtain higher modeling accuracy with almost the same computational efficiency when the FD coefficients are calculated based on the S-wave time–space domain dispersion relation. Solving the split elastic wave equation with M-SFD can further improve the modeling accuracy but will decrease the efficiency and increase the memory usage as well. Stability analysis shows that M-SFD has better stability than C-SFD for both acoustic and elastic wave simulations. Applying M-SFD to reverse time migration (RTM), the imaging artifacts caused by the numerical dispersion are effectively eliminated, which improves the imaging accuracy and resolution of deep formation.

A hybrid explicit implicit staggered grid finite-difference scheme for the first-order acoustic wave equation modeling

Implicit staggered-grid finite-difference (SGFD) methods are widely used for the first-order acoustic wave-equation modeling. The identical implicit SGFD operator is commonly used for all of the first-order spatial derivatives in the first-order acoustic wave-equation. In this paper, we propose a hybrid explicit implicit SGFD (HEI-SGFD) scheme which could simultaneously preserve the wave-equation simulation accuracy and increase the wave-equation simulation speed. We use a second-order explicit SGFD operator for half of the first-order spatial derivatives in the first-order acoustic wave-equation. At the same time, we use the implicit SGFD operator with added points in the diagonal direction for the other first-order spatial derivatives in the first-order acoustic wave-equation. The proposed HEI-SGFD scheme nearly doubles the wave-equation simulation speed compared to the implicit SGFD schemes. In essence, the proposed HEI-SGFD scheme is equivalent to the second-order FD scheme with ordinary grid format. We then determine the HEI-SGFD coefficients in the time–space domain by minimizing the phase velocity error using the least-squares method. Finally, the effectiveness of the proposed method is demonstrated by dispersion analysis and numerical simulations.

Dual-pairing summation by parts finite difference methods for large scale elastic wave simulations in 3D complex geometries

High-order accurate summation-by-parts (SBP) finite difference (FD) methods constitute efficient numerical methods for simulating large-scale hyperbolic wave propagation problems. Traditional SBP FD operators that approximate first-order spatial derivatives with central-difference stencils often have spurious unresolved numerical wave-modes in their computed solutions. Recently derived high order accurate dual pair SBP operators based non-central (upwind) FD stencils have the potential to suppress these poisonous spurious wave-modes on marginally resolved computational grids. In this paper, we demonstrate that not all high order dual pair SBP FD operators are applicable. Numerical dispersion relation analysis shows that odd-order dual pair SBP FD operators also support spurious unresolved high-frequency wave-modes on marginally resolved meshes. Meanwhile, even-order dual pair SBP FD operators (of order 2,4,6) do not support spurious unresolved high frequency wave-modes and also have better numerical dispersion properties.For all the dual pair SBP FD operators we discretise the three space dimensional (3D) elastic wave equation on boundary-conforming curvilinear meshes. Using the energy method we prove that the semi-discrete approximation is stable and energy-conserving. We derive a priori error estimate and prove the convergence of the numerical error for smooth solutions. Numerical experiments for the 3D elastic wave equation in complex geometries corroborate the theoretical analysis. Numerical simulations of the 3D elastic wave equation in heterogeneous media with complex non-planar free surface topography are given, including numerical simulations of community developed seismological benchmark problems. Computational results show that even-order dual pair SBP FD operators are more efficient, robust and less prone to numerical dispersion errors on marginally resolved meshes when compared to the odd-order dual pair and traditional SBP FD operators. Finally, scaling tests demonstrate nearly perfect strong scaling and verify the efficiency of our parallel implementation of the high order upwind SBP methods.

Hexagonal Grid Computation of the Derivatives of the Solution to the Heat Equation by Using Fourth-Order Accurate Two-Stage Implicit Methods

We give fourth-order accurate implicit methods for the computation of the first-order spatial derivatives and second-order mixed derivatives involving the time derivative of the solution of first type boundary value problem of two dimensional heat equation. The methods are constructed based on two stages: At the first stage of the methods, the solution and its derivative with respect to time variable are approximated by using the implicit scheme in Buranay and Arshad in 2020. Therefore, Oh4+τ of convergence on constructed hexagonal grids is obtained that the step sizes in the space variables x1, x2 and in time variable are indicated by h, 32h and τ, respectively. Special difference boundary value problems on hexagonal grids are constructed at the second stages to approximate the first order spatial derivatives and the second order mixed derivatives of the solution. Further, Oh4+τ order of uniform convergence of these schemes are shown for r=ωτh2≥116, ω>0. Additionally, the methods are applied on two sample problems.

High-order weighted compact nonlinear scheme for one- and two-dimensional Hamilton-Jacobi equations

This paper designs a fifth-order weighted compact nonlinear scheme (WCNS) on a five-point stencil to solve one- and two-dimensional Hamilton-Jacobi equations. The five-point WCNS is used to compute the left and right limits of first-order spatial derivatives of the HJ equations in the Lax-Friedrichs monotone numerical Hamiltonian. The WENO-Z type interpolation for cell-edge values of the solutions is used to suppress numerical oscillations which may appear near discontinuities. Five- and seven-point WENO-Z schemes for Hamilton-Jacobi equations are also designed for comparisons. The performance of the WCNS and the two WENO-Z schemes is demonstrated by several numerical examples in one-dimensional and two-dimensional cases.

Classical solution for a nonlinear hybrid system modeling combustion in a multilayer porous medium

Combustion processes in porous media have been used by the petroleum engineering industry to extract heavy oil from reservoirs. This study focuses on a one-dimensional nonlinear hybrid system consisting of n reaction–diffusion–convection equations coupled with n ordinary differential equations, which models a combustion front moving through a porous medium with n parallel layers. The state variables are the temperature and fuel concentration in each layer. Coupling occurs in both the reaction function and differential operator coefficients. We prove the existence of a classical solution, first locally and then globally over time, to an initial and boundary value problem for the corresponding system. The proof uses a new approach for combustion problems in porous media. The local solution is obtained by defining an operator in a set of Hölder continuous functions and using Schauder’s fixed-point theorem to find a fixed point as the desired solution. Using Zorn’s lemma, we extend the local solution to a global solution, proving that the first-order spatial derivative of the temperature in each layer is a bounded function.

Combining First- and Second-Order Continuity Constraints in Ultrasound Elastography.

Ultrasound elastography is a prominent noninvasive medical imaging technique that estimates tissue elastic properties to detect abnormalities in an organ. A common approximation to tissue elastic modulus is tissue strain induced after mechanical stimulation. To compute tissue strain, ultrasound radio frequency (RF) data can be processed using energy-based algorithms. These algorithms suffer from ill-posedness to tackle. A continuity constraint along with the data amplitude similarity is imposed to obtain a unique solution to the time-delay estimation (TDE) problem. Existing energy-based methods exploit the first-order spatial derivative of the displacement field to construct a regularizer. This first-order regularization scheme alone is not fully consistent with the mechanics of tissue deformation while perturbed with an external force. As a consequence, state-of-the-art techniques suffer from two crucial drawbacks. First, the strain map is not sufficiently smooth in uniform tissue regions. Second, the edges of the hard or soft inclusions are not well-defined in the image. Herein, we address these issues by formulating a novel regularizer taking both first- and second-order derivatives of the displacement field into account. The second-order constraint, which is the principal novelty of this work, contributes both to background continuity and edge sharpness by suppressing spurious noisy edges and enhancing strong boundaries. We name the proposed technique: Second-Order Ultrasound eLastography (SOUL). Comparative assessment of qualitative and quantitative results shows that SOUL substantially outperforms three recently developed TDE algorithms called Hybrid, GLUE, and MPWC-Net++. SOUL yields 27.72%, 62.56%, and 81.37% improvements of the signal-to-noise ratio (SNR) and 72.35%, 54.03%, and 65.17% improvements of the contrast-to-noise ratio (CNR) over GLUE with data pertaining to simulation, phantom, and in vivo tissue, respectively. The SOUL code can be downloaded from code.sonography.ai.

An optimal method for a staggered-grid finite-difference solution of elastic wave equations including rotational deformation

The staggered-grid finite-difference (SGFD) method is one of the most popular means used in seismic numerical modelling because of its computational efficiency and it is easy to implement. However, it will lead to serious numerical errors while using coarse grids or wavelet with high-frequency. In this work, we proposed an optimal spatial SGFD method based on a new improved particle swarm optimization (PSO) algorithm. Then the dispersion of numerical solution of the first-order spatial derivatives is analyzed. Meanwhile, numerical modelling of elastic wave equations and elastic wave equations including rotational deformation are performed with the optimal spatial SGFD method based on the improved PSO algorithm. Numerical dispersion analysis and numerical modeling results indicate that the optimal spatial SGFD method based on the improved PSO algorithm has high modeling accuracy and can efficiently suppress the numerical dispersion.

Boundary treatment of implicit-explicit Runge-Kutta method for hyperbolic systems with source terms

In this paper, we develop a high order finite difference boundary treatment method for the implicit-explicit (IMEX) Runge-Kutta (RK) schemes solving hyperbolic systems with possibly stiff source terms on a Cartesian mesh. The main challenge is how to obtain the solutions at ghost points resulting from the wide stencil of the interior high order scheme. We address this problem by combining the idea of using the RK schemes at the boundary and an inverse Lax-Wendroff (ILW) procedure in [29,30]. Specifically, we only apply the ILW procedure in the starting stage of the RK method. In the intermediate stages, we solve out the intermediate solutions as well as their first-order spatial derivatives at the boundary by using the RK schemes, which are then used to compute solutions at ghost points by Taylor expansions. Our method is different from the widely used approach for the explicit RK schemes by imposing boundary conditions at intermediate stages [14,29,30] which does not apply to IMEX schemes. In addition, the intermediate boundary conditions are available for explicit RK schemes only up to third order while our method applies to IMEX and explicit RK schemes of arbitrary order. The stability and the desired accuracy of our boundary treatment with IMEX RK method up to fourth order are demonstrated numerically through 1D examples and 2D reactive Euler equations.

Numerical study on the performance of semicircular and rectangular submerged breakwaters

A systematic numerical comparative study of the performance of semicircular and rectangular submerged breakwaters interacting with solitary waves is the basis of this paper. To accomplish this task, Nwogu\'s extended Boussinesq model equations are employed to simulate the interaction of the wave with breakwaters. The finite difference technique has been used to discretize the spatial terms while a fourth-order predictor-corrector method is employed for time discretization in our numerical model. The proposed computational scheme uses a staggered-grid system where the first-order spatial derivatives have been discretized with fourth-order accuracy. For validation purposes, five test cases are considered and numerical results have been successfully compared with the existing analytical and experimental results. The performances of the rectangular and semicircular breakwaters have been examined in terms of the wave reflection, transmission, and dissipation coefficients (RTD coefficients) denoted by K_R, K_T, K_D. The latter coefficient K_D emerges due to the non-energy conserving K_R and K_T,. Our computational results and graphical illustrations show that the rectangular breakwater has higher reflection coefficients than semicircular breakwater for a fixed crest height, but as the wave height increases, the two reflection coefficients approach each other. On the other hand, the rectangular breakwater has larger dissipation coefficients compared to that of the semicircular breakwater and the difference between them increases as the height of the crest increases. However, the transmission coefficient for the semicircular breakwater is greater than that of the rectangular breakwater and the difference in their transmission coefficients increases with the crest height. Quantitatively, for rectangular breakwaters the reflection coefficients K_R are 5-15% higher while the diffusion coefficients K_D are 3-23% higher than that for the semicircular breakwaters, respectively. The transmission coefficients K_T for rectangular breakwater shows the better performance up to 2.47% than that for the semicircular breakwaters. Based on our computational results, one may conclude that the rectangular breakwater has a better overall performance than the semicircular breakwater. Although the model equations are non-dissipative, the non-energy conserving transmission and reflection coefficients due to wave-breakwater interactions lead to dissipation type contribution.

A HWENO reconstruction based high-order compact gas-kinetic scheme on unstructured mesh

As an extension of a fourth-order compact gas kinetic scheme (GKS) on structured mesh [24], this work is about the development of a third-order compact GKS on unstructured mesh for the compressible Euler and Navier-Stokes solutions. Based on the time accurate high-order gas-kinetic evolution solution, the time-dependent gas distribution function at a cell interface in GKS provides not only the flux function and its time derivative, but also the time accurate flow variables there at the next time level. As a result, besides updating the conservative flow variables inside each control volume through the interface fluxes, the cell averaged first-order spatial derivatives of flow variables can be obtained as well using the updated flow variables at the closed cell interfaces around that cell through the divergence theorem. Therefore, with the cell-averaged flow variables and their first-order spatial derivatives inside each cell, the Hermite WENO (HWENO) techniques can be naturally implemented for the compact high-order reconstruction at the beginning of the next time step. Following the reconstruction technique in [64], a new HWENO reconstruction on triangular mesh is designed in the current scheme. Combined with the two-stage temporal discretization and second-order time accurate flux function, a third-order compact scheme on unstructured mesh has been constructed. Accurate solutions can be obtained for both inviscid and viscous flows without sensitive dependence on the quality of triangular mesh. The robustness and accuracy of the scheme have been validated through many cases, including strong shocks in the hypersonic viscous flow and smooth Navier-Stokes solution.

On Multi-Quadric Based RBF–FD Method for Second-Order Diffusion Filters

The usual three-node approximations for the first-order spatial derivatives in the radial basis functions-based finite difference (RBF-FD) scheme allows relatively smaller time step size than the standard scheme (forward difference in time and central difference in space) in obtaining stable results. In this paper, to better mimic the continuous model, we use nine nodal points in RBF–FD approximation of first-order spatial derivatives. We then employ these approximations to solve a two-dimensional nonlinear diffusion equation from image processing. The constructed numerical scheme is shown to be stable and allows up to obtain four times bigger time step size than the standard scheme.

Consistent nonlinear deterministic and stochastic wave evolution equations from deep water to the breaking region

Modelling the evolution of the wave field in coastal waters is a complicated task, partly due to triad nonlinear wave interactions, which are one of the dominant mechanisms in this area. Stochastic formulations already implemented into large-scale operational wave models, whilst very efficient, are one-dimensional in nature and fail to account for the majority of the physical properties of the wave field evolution. This paper presents new two-dimensional (2-D) formulations for the triad interactions source term. A quasi-two-dimensional deterministic mild slope equation is improved by including dissipation and first-order spatial derivatives in the nonlinear part of equation, significantly enhancing the accuracy in the breaking zone. The newly defined deterministic model is used to derive an updated stochastic model consistent from deep waters to the breaking region. It is localized following the approach derived in Vrecica & Toledo (J. Fluid Mech., vol. 794, 2016, pp. 310–342), to which several improvements are also presented. The model is compared to measurements of breaking and non-breaking spectral evolution, showing good agreement in both cases. Finally, the model is used to analyse several interesting 2-D properties of the shoaling wave field including the evolution of directionally spread seas.

Method and Algorithm for Obtaining Elements of the Tensor of Spatial Derivatives of the Magnetic Induction Vector in the Problem of Searching for Magnetic Anomalies

The values of all the components of the magnetic induction vector and its first-order spatial derivatives are determined from the spatial distribution of the magnetic field parameter values at each point of the observation plane. The inverse problem is solved using the analytic eigenvector method. The execution of the proposed algorithm was simulated using real data of magnetometric studies in the geomagnetic field.