A finite element solution for the fractional advection–dispersion equation

The fractional advection–dispersion equation (FADE) known as its non-local dispersion, has been proven to be a promising tool to simulate anomalous solute transport in groundwater. We present an unconditionally stable finite element (FEM) approach to solve the one-dimensional FADE based on the Caputo definition of the fractional derivative with considering its singularity at the boundaries. The stability and accuracy of the FEM solution is verified against the analytical solution, and the sensitivity of the FEM solution to the fractional order α and the skewness parameter β is analyzed. We find that the proposed numerical approach converge to the numerical solution of the advection–dispersion equation (ADE) as the fractional order α equals 2. The problem caused by using the first- or third-kind boundary with an integral-order derivative at the inlet is remedied by using the third-kind boundary with a fractional-order derivative there. The problems for concentration estimation at boundaries caused by the singularity of the fractional derivative can be solved by using the concept of transition probability conservation. The FEM solution of this study has smaller numerical dispersion than that of the FD solution by Meerschaert and Tadjeran (J Comput Appl Math 2004). For a given α, the spatial distribution of concentration exhibits a symmetric non-Fickian behavior when β = 0. The spatial distribution of concentration shows a Fickian behavior on the left-hand side of the spatial domain and a notable non-Fickian behavior on the right-hand side of the spatial domain when β = 1, whereas when β = −1 the spatial distribution of concentration is the opposite of that of β = 1. Finally, the numerical approach is applied to simulate the atrazine transport in a saturated soil column and the results indicat that the FEM solution of the FADE could better simulate the atrazine transport process than that of the ADE, especially at the tail of the breakthrough curves.     

A central difference method with low numerical dispersion for solving the scalar wave equation

In this paper, we propose a nearly‐analytic central difference method, which is an improved version of the central difference method. The new method is fourth‐order accurate with respect to both space and time but uses only three grid points in spatial directions. The stability criteria and numerical dispersion for the new scheme are analysed in detail. We also apply the nearly‐analytic central difference method to 1D and 2D cases to compute synthetic seismograms. For comparison, the fourth‐order Lax‐Wendroff correction scheme and the fourth‐order staggered‐grid finite‐difference method are used to model acoustic wavefields. Numerical results indicate that the nearly‐analytic central difference method can be used to solve large‐scale problems because it effectively suppresses numerical dispersion caused by discretizing the scalar wave equation when too coarse grids are used. Meanwhile, numerical results show that the minimum sampling rate of the nearly‐analytic central difference method is about 2.5 points per minimal wavelength for eliminating numerical dispersion, resulting that the nearly‐analytic central difference method can save greatly both computational costs and storage space as contrasted to other high‐order finite‐difference methods such as the fourth‐order Lax‐Wendroff correction scheme and the fourth‐order staggered‐grid finite‐difference method.     

用于声波方程数值模拟的时间-空间域有限差分系数确定新方法

声波方程数值模拟已广泛应用于理论地震计算,同时构成了地震逆时偏移成像技术的基础.对于有限差分法而言,在满足一定的稳定性条件时,普遍存在着因网格化而形成的数值频散效应.如何有效地缓解或压制数值频散是有限差分方法研究的关键所在.为精确求解空间偏导数,相继发展了高阶差分格式优化方法和伪谱方法.近期,为更好地缓解数值频散,提出了时间-空间域有限差分方法,该方法采用了泰勒展开近似方法来确定有限差分格式系数,因而只能保证在一定的小范围内很好的拟合波场传播规律.为进一步压制数值频散效应,本文引入了时间-空间域特定波数点满足频散关系的方法,根据震源、波速和网格间距确定波数范围,同时考虑了多个传播角度,然后建立方程确定了相应的有限差分格式系数,使得差分系数能在更大范围符合波场传播规律.通过频散分析和正演模拟,验证了本文方法的有效性.     

有限差分的时间和空间频散及其对逆时偏移和全波形反演的影响

有限差分是最常用的地震波方程数值模拟方法,但时间和空间离散会产生数值频散.正演是逆时偏移和全波形反演的基本单元,成像和反演的精度很大程度依赖于所采用的数值模拟算法.本文研究了有限差分的时间和空间频散特性及其对逆时偏移和全波形反演的影响.通过时间有限差分+伪谱法、时间频散校正+空间有限差分、时间频散校正+伪谱法获取时间频散、空间频散和无频散数据;发展了抗时间频散、抗空间频散、抗时间+空间频散的逆时偏移和全波形反演方法;采用理论模型和实际资料对提出的方法进行了测试.数值结果表明：逆时偏移同时受时间和空间频散影响,时间频散导致同相轴不聚焦、成像位置偏离,空间频散会产生高频噪声和虚假反射界面;全波形反演在低频大尺度反演中几乎不受时间和空间频散影响,高频精细反演中时间频散引起波形相移、降低反演精度,空间频散增加多解性、导致反演不收敛;抗频散方法可以有效缓解时间和空间频散影响,获得高质量的偏移剖面和反演结果.

     

用于弹性波方程数值模拟的有限差分系数确定方法

有限差分方法是波场数值模拟的一个重要方法,交错网格差分格式比规则网格差分格式稳定性更好,但方法本身都存在因网格化而形成的数值频散效应,这会降低波场模拟的精度与分辨率.为了缓解有限差分算子的数值频散效应,精确求解空间偏导数,本文把求解波动方程的线性化方法推广到用于求解弹性波方程交错网格有限差分系数;同时应用最大最小准则作为模拟退火(SA)优化算法求解差分系数的数值频散误差判定标准来求解有限差分系数.通过上述两种方法,分别利用均匀各向同性介质和复杂构造模型进行了数值正演模拟和数值频散分析,并与传统泰勒展开算法、最小二乘算法进行比较,验证了线性化方法和模拟退火方法都能有效压制数值频散,并比较了各个算法的特点.     

Extreme value theory based on the r largest annual events: a robust approach

In Smith (1986, J. Hydrol. 86, 27–43), a family of statistical distributions and estimators for extreme values based on a fixed number r > = 1 of the largest annual events are presented. The method of estimation was numerical maximum likelihood. In this paper, we consider the robust estimation of parameters in such families of distributions. The estimation technique, which is based on optimal B-robust estimates, will assign weights to each observation and give estimates of the parameters based on the data which are well modeled by the distribution. Thus, observations which are not consistent with the proposed distribution can be identified and the validity of the model can be assessed. The method is illustrated on Venice sea level data.     

黏弹介质波动方程有限差分解的稳定性研究

稳定性问题是地震波数值模拟的一个重要问题.基于地震波传播理论,从黏弹介质本构方程出发,对矩形网格下不同黏弹模型波动方程有限差分解的稳定性进行了理论分析,导出了Kelvin-Voigt黏弹模型和Maxwell黏弹模型在任意空间差分精度下稳定性条件的表达式;给出了品质因子Q≥5时的简化式,并通过数值算例验证了理论研究结论的正确性;总结了地震波速度、频率、空间网格大小、差分系数以及品质因子与稳定性条件的关系;通过误差分析给出了近似公式的使用条件.     

Stability condition of finite difference solution for viscoelastic wave equations

The stability problem is a very important aspect in seismic wave numerical modeling. Based on the theory of seismic waves and constitutive equations of viscoelastic models, the stability problems of finite difference scheme for KelvinVoigt and Maxwell models with rectangular grids are analyzed. Expressions of stability conditions with arbitrary spatial accuracies for two viscoelastic models are derived. With approximation of quality factor Q≥?5, simplified expressions are developed and some numerical models are given to verify the validity of the corresponding theoretical results. Then this paper summarizes the influences of seismic wave velocity, frequency, size of grid and difference coefficients, as well as quality factor on stability condition. Finally the prerequisite conditions of the simplified stability equations are given with error analysis.     

二维地震波场的组合型紧致有限差分数值模拟

地震波场数值模拟在地球物理勘探和地震学中具有重要的支撑作用.本文将组合型紧致差分格式用于声波和弹性波方程的数值模拟中.根据泰勒级数展开和声波方程,建立了位移场时间四阶离散格式,并将组合型紧致差分格式用于位移场空间导数的求取,然后对该差分格式进行了精度分析、误差分析、频散分析和稳定性分析.理论研究结果表明：①该差分格式为时间四阶、空间六阶精度,与常规七点六阶中心差分和五点六阶紧致差分相比,具有更小的截断误差和更高的模拟精度;②每个波长仅需要5.6个采样点,且满足稳定性条件的库郎数为0.792,可以使用粗网格和较大时间步长进行计算.所以该方法具有占用内存少、计算效率高和低数值频散等优势.最后,本文进行了二维各向同性完全弹性介质的声波和弹性波方程的数值模拟,实验结果表明本文提出的方法具有更高的计算精度,能够大幅度的节约计算量和内存需求,对于三维大尺度模型问题具有更好的适应性.     

适用于声波方程数值模拟的时间-空间域隐式有限差分算子优化方法

声波方程数值模拟已广泛应用于理论地震计算,同时构成了地震逆时偏移成像技术的基础.对于有限差分法而言,在满足一定的稳定性条件时,普遍存在着因网格化而形成的数值频散效应.如何有效地缓解或压制数值频散是有限差分方法研究的关键所在.有限差分格式分为显式有限差分和隐式有限差分.隐式有限差分能够进一步压制数值频散效应.因此本文提出了给定频率范围满足时间-空间域隐式有限差分频散关系的方法,并根据震源频率、波速和网格间距确定波数范围,在此基础上建立方程确定了相应的隐式有限差分系数,使得差分系数能在更大频率范围符合波场传播规律.通过频散分析和正演模拟,验证了本文方法的有效性.     

三角网格有限元法声波与弹性波模拟频散分析

本文对声波与弹性波方程进行有限元法离散,构造有限元法频散关系的一般特征值问题,分析了时间离散格式为中心差分的三角网格有限元法声波与弹性波模拟的频散特性. 比较了三种质量矩阵即分布式质量矩阵、集中质量矩阵和混合质量矩阵对有限元法频散的影响;选取四种典型三角网格,分析了混合质量矩阵有限元（MFEM）频散的方向各向异性;数值频散、方向各向异性随插值阶数的增加逐渐减弱,当空间为三阶插值时,频散主要表现为随采样率的变化而几乎无明显方向各向异性, 其频散幅值也较小. 控制其他影响因素不变的情况下,研究了不同波速比介质中弹性波的数值频散. 最后给出了三角网格MFEM的数值耗散性.     

Truncation errors of selected finite difference methods for two-dimensional advection-diffusion equation with mixed derivatives

The spread of a passive contaminant in an open-channel reach is considered with use of a two-dimensional advection-diffusion equation with the included off-diagonal dispersion coefficients. This paper presents the calculation of truncation errors, namely numerical diffusion and numerical dispersion for various finite difference schemes. The accuracy of the considered finite-difference approximations is analysed by deriving and studying the relevant modified partial differential equation.     

Evaluation of the Ross fast solution of Richards’ equation in unfavourable conditions for standard finite element methods

Ross [Ross PJ. Modeling soil water and solute transport – fast, simplified numerical solutions. Agron J 2003;95:1352–61] developed a fast, simplified method for solving Richards’ equation. This non-iterative 1D approach, using Brooks and Corey [Brooks RH, Corey AT. Hydraulic properties of porous media. Hydrol. papers, Colorado St. Univ., Fort Collins; 1964] hydraulic functions, allows a significant reduction in computing time while maintaining the accuracy of the results. The first aim of this work is to confirm these results in a more extensive set of problems, including those that would lead to serious numerical difficulties for the standard numerical method. The second aim is to validate a generalisation of the Ross method to other mathematical representations of hydraulic functions.     

Optimal implicit staggered‐grid finite‐difference schemes based on the sampling approximation method for seismic modelling

We propose new implicit staggered‐grid finite‐difference schemes with optimal coefficients based on the sampling approximation method to improve the numerical solution accuracy for seismic modelling. We first derive the optimized implicit staggered‐grid finite‐difference coefficients of arbitrary even‐order accuracy for the first‐order spatial derivatives using the plane‐wave theory and the direct sampling approximation method. Then, the implicit staggered‐grid finite‐difference coefficients based on sampling approximation, which can widen the range of wavenumber with great accuracy, are used to solve the first‐order spatial derivatives. By comparing the numerical dispersion of the implicit staggered‐grid finite‐difference schemes based on sampling approximation, Taylor series expansion, and least squares, we find that the optimal implicit staggered‐grid finite‐difference scheme based on sampling approximation achieves greater precision than that based on Taylor series expansion over a wider range of wavenumbers, although it has similar accuracy to that based on least squares. Finally, we apply the implicit staggered‐grid finite difference based on sampling approximation to numerical modelling. The modelling results demonstrate that the new optimal method can efficiently suppress numerical dispersion and lead to greater accuracy compared with the implicit staggered‐grid finite difference based on Taylor series expansion. In addition, the results also indicate the computational cost of the implicit staggered‐grid finite difference based on sampling approximation is almost the same as the implicit staggered‐grid finite difference based on Taylor series expansion.     

Estimating finite difference block equivalent hydraulic conductivity for numerically solving the Richards' equation

A method is proposed for calculating the equivalent hydraulic conductivity (EHC) within a finite difference block (FDB). Application of the constant‐flux assumption of Darcy's Law, the EHC equals to the integration of effective hydraulic conductivity (K_w) as a function of pressure head (h_w) divided by the head difference at the ends of the FDB. Error analysis show that the constant‐flux (CF) EHC estimates are better than those computed by the commonly used arithmetic‐mean (AM), geometric‐mean (GM), and harmonic‐mean (HM) techniques. CF EHC results are even more superior at larger interblock head difference situations. Simulations of water infiltration experiments show that simulations using the CF EHC or AM or GM weighting technique have only slight difference while applying the Neumann type boundary condition at the ground surface. In case of the Dirichlet type boundary condition, however, the CF EHC is superior to the other two in correctly estimating the depth of infiltration while enlarging the grid size. Therefore, it is recommended to adopt the CF EHC with a larger grid size to the more stable and more efficient results. Copyright © 2007 John Wiley & Sons, Ltd.     

一种优化的频率域三维声波有限差分模拟方法

为提高频率域有限差分（FD,finite-difference）正演模拟技术的计算精度和效率,基于旋转坐标系统的优化差分格式被广泛应用,但是只应用于正方形网格的情况.基于平均导数法（ADM）的优化差分格式,应用于正方形和长方形网格模拟.这些频率域有限差分算子,各自具有不同的差分格式和对应的优化系数求解表达式.本文基于三维声波方程发展了一种新的优化方法,只要给定FD模板形式,可直接构造频散方程,求取FD模板上各节点的优化系数.此方法的优点在于频率域FD算子的优化系数对应各个节点,可扩展优化其他格式.运用此优化方法,计算得到了不同空间采样间距比情况下27点和7点格式的优化系数.数值实验表明,优化27点格式与ADM 27点格式具有相同的精度,优化7点格式比经典的7点格式具有更小的数值频散.

     

A finite volume upwind scheme for the solution of the linear advection–diffusion equation with sharp gradients in multiple dimensions

A finite volume upwind numerical scheme for the solution of the linear advection equation in multiple dimensions on Cartesian grids is presented. The small-stencil, Modified Discontinuous Profile Method (MDPM) uses a sub-cell piecewise constant reconstruction and additional information at the cell interfaces, rather than a spatial extension of the stencil as in usual methods. This paper presents the MDPM profile reconstruction method in one dimension and its generalization and algorithm to two- and three-dimensional problems. The method is extended to the advection–diffusion equation in multiple dimensions. The MDPM is tested against the MUSCL scheme on two- and three-dimensional test cases. It is shown to give high-quality results for sharp gradients problems, although some scattering appears. For smooth gradients, extreme values are best preserved with the MDPM than with the MUSCL scheme, while the MDPM does not maintain the smoothness of the original shape as well as the MUSCL scheme. However the MDPM is proved to be more efficient on coarse grids in terms of error and CPU time, while on fine grids the MUSCL scheme provides a better accuracy at a lower CPU.     

昆仑山断裂带围陷波的有限差分数值模拟解释

利用交错网格有限差分方法对昆仑山断裂带人工爆破产生的围陷波进行了三维数值模拟解释.为提高断裂带最终模型的可信度,在围陷波模拟的同时考虑了人工爆破记录的三个分量.对昆仑山断裂带围陷波的模拟结果表明,影响围陷波特性的断裂带深度主要在1.0 km以内.S波速度和断裂带宽度对围陷波的到时、频率、振幅和相位影响较大.数值模拟解释获得的昆仑山断裂带的细结构参数是：浅部断裂带宽度为300 m,深部为250 m;深度在400 m以上断裂带内S波速度为0.98 km/s,外部围岩S波速度为1.70 km/s,Q值为13.8;S波速度和Q值随着深度的增加而增加;1000 m以下断裂带内S波速度为2.80 km/s,围岩S波速度为3.3 km/s.