对流占优扩散问题的并行计算

１．引言 在刻画流体运动的某些物理现象,以及研究热的传导、粒子的扩散等问题时,都会归结到求解对流扩散方程．用有限差分方法求解该方程,若采用显式方法,计算格式简单,但它们都是条件稳定的,时间步长必须取得非常小;若采用隐式方法,方法是无条件稳定的,但要解代数方程组,求解比较困难．Ｄ．Ｊ．ＥＶＡＮＳ和Ａ．Ｒ．ＡＨＭＡＤ在文［２］中提出了用显式交替方向法求解定态椭圆型方程,对Ｌａｐｌａｃｅ方程做了数值实验．本文将这个方法推广到了时间依赖的问题,而且适用于对流占优扩散问题的求解．基于二阶迎风格式［１］;本…     

Hydrodynamics in Porous Media: A Finite Volume Lattice Boltzmann Study

Fluid flow through porous media is of great importance for many natural systems, such as transport of groundwater flow, pollution transport and mineral processing. In this paper, we propose and validate a novel finite volume formulation of the lattice Boltzmann method for porous flows, based on the Brinkman–Forchheimer equation. The porous media effect is incorporated as a force term in the lattice Boltzmann equation, which is numerically solved through a cell-centered finite volume scheme. Correction factors are introduced to improve the numerical stability. The method is tested against fully porous Poiseuille, Couette and lid-driven cavity flows. Upon comparing the results with well-documented data available in literature, a satisfactory agreement is observed. The method is then applied to simulate the flow in partially porous channels, in order to verify its potential application to fractured porous conduits, and assess the influence of the main porous media parameters, such as Darcy number, porosity and porous media thickness.     

一维非定常对流扩散方程的高阶组合紧致迎风格式

通过将对流项采用四五阶组合迎风紧致格式离散,扩散项采用四阶对称紧致格式离散之后,对得到的半离散格式在时间方向采用四阶龙格库塔方法求解,从而得到了一种求解非定常对流扩散方程问题的高精度组合紧致有限差分格式,其收敛阶为O(h~4+τ~4).经Fourier精度分析和数值验证,证实了格式的良好性能.三个数值算例包括线性常系数问题,矩形波问题和非线性问题,数值结果表明:该格式具有很高的分辨率,且适用于对高雷诺数问题的数值模拟.     

A posteriori error analysis for a fractional differential equation

Numerical treatment for a fractional differential equation (FDE) is proposed and analysed. The solution of the FDE may be singular near certain domain boundaries, which leads to numerical difficulty. We apply the upwind finite difference method to the FDE. The stability properties and a posteriori error analysis for the discrete scheme are given. Then, a posteriori adapted mesh based on a posteriori error analysis is established by equidistributing arc-length monitor function. Numerical experiments illustrate that the upwind finite difference method on a posteriori adapted mesh is more accurate than the method on uniform mesh.     

一维对流扩散反应方程的高精度数值方法

通过将原方程变换为对流扩散方程,将所得方程的对流项采用四阶组合紧致迎风格式离散,扩散项采用四阶对称紧致格式离散之后,对得到的空间半离散格式采用四阶龙格库塔方法进行时间推进,得到了一种求解非定常对流扩散反应问题的高精度方法,其收敛阶为O(h4+τ4).经数值实验并与文献结果进行对比,表明该格式适用于对流占优问题的数值模拟,验证了格式的良好性能.     

Approximate boundary-flux calculations

A technique for determining derivatives (fluxes or stresses) from finite element solutions is developed. The approach is a generalization to higher dimensions of a procedure known to give highly accurate results in one dimension. Numerical experiments demonstrate that certain difficulties are associated with corners in the higher-dimensional extensions and two variants of the method are examined. We consider both triangular and quadrilateral elements and observe some interesting differences in the numerical rates of convergence. Finally, this post-processing scheme is tested for nonlinear problems.     

Finite Volume Approximation of One-Dimensional Stiff Convection-Diffusion Equations

In this work, we present a novel method to approximate stiff problems using a finite volume (FV) discretization. The stiffness is caused by the existence of a small parameter in the equation which introduces a boundary layer. The proposed semi-analytic method consists in adding in the finite volume space the boundary layer corrector which encompasses the singularities of the problem. We verify the stability and convergence of our finite volume schemes which take into account the boundary layer structures. A major feature of the proposed scheme is that it produces an efficient stable second order scheme to be compared with the usual stable upwind schemes of order one or the usual costly second order schemes demanding fine meshes.     

A full discrete two-grid finite-volume method for a nonlinear parabolic problem

A fully discrete two-grid finite-volume method (FVM) for a nonlinear parabolic problem is studied in this paper. This method involves solving a nonlinear parabolic equation on coarse mesh space and a linearized parabolic equation on fine grid. Both L ² and H ¹ norm error estimates of the standard FVM for the nonlinear parabolic problem are derived. Compared with the standard FVM, the two-level method is of the same order as the one-level method in the H ¹-norm as long as the mesh sizes satisfy h=𝒪(H ^3/2). However, the two-level method involves much less work than the standard method. Numerical results are provided to demonstrate the effectiveness of our algorithm.     

磁轴承变饱和柔性变结构控制

采用磁轴承取代传统机械轴承,可消除飞轮储能装置机械摩擦,提高转子临界转速.针对磁轴承控制系统响应快、运行稳定的要求,设计了一种变饱和柔性变结构控制器.根据磁轴承的基本结构与悬浮力模型,建立了系统的状态方程;通过有限元仿真,分析了磁轴承耦合性,确定了系统基本控制方案;将柔性变结构和变饱和项结合,当系统存在较大偏差时只通过线性控制器调节,在中等偏差时通过线性控制器和变饱和状态控制器一起调节,在小偏差时由线性控制器和变饱和状态控制器转变的线性控制器一起调节.仿真与实验结果表明,变饱和柔性变结构控制满足实时控制的要求,磁轴承控制系统具有较快的响应速度和较高的运行稳定性.     

A Second-Order Difference Scheme for the Penalized Black–Scholes Equation Governing American Put Option Pricing

In this paper we present a stable finite difference scheme on a piecewise uniform mesh along with a power penalty method for solving the American put option problem. By adding a power penalty term the linear complementarity problem arising from pricing American put options is transformed into a nonlinear parabolic partial differential equation. Then a finite difference scheme is proposed to solve the penalized nonlinear PDE, which combines a central difference scheme on a piecewise uniform mesh with respect to the spatial variable with an implicit time stepping technique. It is proved that the scheme is stable for arbitrary volatility and arbitrary interest rate without any extra conditions and is second-order convergent with respect to the spatial variable. Furthermore, our method can efficiently treats the singularities of the non-smooth payoff function. Numerical results support the theoretical results.     

Multi-level adaptive simulation of transient two-phase flow in heterogeneous porous media

An implicit pressure and explicit saturation (IMPES) finite element method (FEM) incorporating a multi-level shock-type adaptive refinement technique is presented and applied to investigate transient two-phase flow in porous media. Local adaptive mesh refinement is implemented seamlessly with state-of-the-art artificial diffusion stabilization allowing simulations that achieve both high resolution and high accuracy. Two benchmark problems, modelling a single crack and a random porous medium, are used to demonstrate the robustness of the method and illustrate the capabilities of the adaptive refinement technique in resolving the saturation field and the complex interaction (transport phenomena) between two fluids in heterogeneous media.     

Numerical simulation of contaminant biodegradation by higher order methods and adaptive time stepping

In this work we present and analyze a reliable and robust approximation scheme for biochemically reacting transport in the subsurface following Monod type kinetics. Water flow is modeled by the Richards equation. The proposed scheme is based on higher order finite element methods for the spatial discretization and the two step backward differentiation formula for the temporal one. The resulting nonlinear algebraic systems of equations are solved by a damped version of Newtons method. For the linear problems of the Newton iteration Krylov space methods are used. In computational experiments conducted for realistic subsurface (groundwater) contamination scenarios we show that the higher order approximation scheme significantly reduces the amount of inherent numerical diffusion compared to lower order ones. Thereby an artificial transverse mixing of the species leading to a strong overestimation of the biodegradation process is avoided. Finally, we present a robust adaptive time stepping technique for the coupled flow and transport problem which allows efficient long-term predictions of biodegradation processes.     

Implementation of variably saturated flow into PHREEQC for the simulation of biogeochemical reactions in the vadose zone

A software tool for the simulation of one-dimensional unsaturated flow and solute transport together with biogeochemical reactions in the vadose zone was developed by integrating a numerical solution of Richards' equation into the geochemical modelling framework PHREEQC. Unlike existing software, the new simulation tool is fully based on existing capabilities of PHREEQC without source code modifications or coupling to external software packages. Because of the direct integration, the full set of PHREEQC's geochemical reactions with connected databases for reaction constants is immediately accessible. Thus, unsaturated flow and solute transport can be simulated together with complex solution speciation, equilibrium and kinetic mineral reactions, redox reactions, ion exchange reactions and surface adsorption including diffuse double layer calculations. Liquid phase flow is solved as a result of element advection, where the Darcy flux is calculated according to Richards' equation. For accurate representation of advection-dominated transport, a total-variation-diminishing scheme was implemented. Dispersion was simulated with a standard approach; however, PHREEQC's multicomponent transport capabilities can be used to simulate species-dependent diffusion. Since liquid phase saturation is recalculated after every reaction step, biogeochemical processes that modify liquid phase saturation, such as dissolution or precipitation of hydrated minerals are considered a priori. Implementation of the numerical schemes for flow and solute transport have been described, along with examples of the extensive code verification, before illustrating more advanced application examples. These include (i) the simulation of infiltration with saturation-dependent cation exchange capacity, (ii) changes in hydraulic properties due to mineral reactions and (iii) transport of mobile organic substances with complexation of heavy metals in varying geochemical conditions.     

Time-splitting pseudo-spectral domain decomposition method for the soliton solutions of the one- and multi-dimensional nonlinear Schrödinger equations

In this paper, we study the simulation of nonlinear Schrödinger equation in one, two and three dimensions. The proposed method is based on a time-splitting method that decomposes the original problem into two parts, a linear equation and a nonlinear equation. The linear equation in one dimension is approximated with the Chebyshev pseudo-spectral collocation method in space variable and the Crank–Nicolson method in time; while the nonlinear equation with constant coefficients can be solved exactly. As the goal of the present paper is to study the nonlinear Schrödinger equation in the large finite domain, we propose a domain decomposition method. In comparison with the single-domain, the multi-domain methods can produce a sparse differentiation matrix with fewer memory space and less computations. In this study, we choose an overlapping multi-domain scheme. By applying the alternating direction implicit technique, we extend this efficient method to solve the nonlinear Schrödinger equation both in two and three dimensions, while for the solution at each time step, it only needs to solve a sequence of linear partial differential equations in one dimension, respectively. Several examples for one- and multi-dimensional nonlinear Schrödinger equations are presented to demonstrate high accuracy and capability of the proposed method. Some numerical experiments are reported which show that this scheme preserves the conservation laws of charge and energy.     

The dual reciprocity boundary element formulation for nonlinear diffusion problems

This paper presents an extension of the dual reciprocity boundary element method (DRBEM) to deal with nonlinear diffusion problems in which thermal conductivity, specific heat, and density coefficients are all functions of temperature. The DRBEM, recently applied to the solution of problems governed by parabolic and hyperbolic equations, consists in the transformation of the differential equation into an integral equation involving boundary integrals only, the solution of which is achieved by employing a standard boundary element discretization coupled with a two-level finite difference time integration scheme. Contrary to previous formulations for the diffusion equation, the dual reciprocity BEM utilizes the well-known fundamental solution to Laplace's equation, which is space-dependent only. This avoids complex time integrations that normally appear in formulations employing time-dependent fundamental solutions, and permits accurate numerical solutions to be obtained in an efficient way. For nonlinear problems, the integral of conductivity is introduced as a new variable to obtain a linear diffusion equation in the Kirchhoff transform space. This equation involves a modified time variable which is itself a function of position. The problem is solved in an iterative way by using an efficient Newton-Raphson technique which is shown to be rapidly convergent.     

An energy-preserving algorithm for nonlinear Hamiltonian wave equations with Neumann boundary conditions

In this paper, a novel energy-preserving numerical scheme for nonlinear Hamiltonian wave equations with Neumann boundary conditions is proposed and analyzed based on the blend of spatial discretization by finite element method (FEM) and time discretization by Average Vector Field (AVF) approach. We first use the finite element discretization in space, which leads to a system of Hamiltonian ODEs whose Hamiltonian can be thought of as the semi-discrete energy of the original continuous system. The stability of the semi-discrete finite element scheme is analyzed. We then apply the AVF approach to the Hamiltonian ODEs to yield a new and efficient fully discrete scheme, which can preserve exactly (machine precision) the semi-discrete energy. The blend of FEM and AVF approach derives a new and efficient numerical scheme for nonlinear Hamiltonian wave equations. The numerical results on a single-soliton problem and a sine-Gordon equation are presented to demonstrate the remarkable energy-preserving property of the proposed numerical scheme.     

An analysis and comparison of conservative upwind difference schemes for ideal and non-ideal gases

In a recent paper [P. Glaister, Conservative upwind difference schemes for the Euler equations, Comput. Math. Appl. 45 (2003) 1673–1682] a number of numerical schemes were presented for the Euler equations governing compressible flows of an ideal gas, the principal one of which is based on a conservative linearisation approach. This scheme was subsequently extended to encompass compressible flows of real gases where the equation of state allows for non-ideal gases [P. Glaister, Conservative upwind difference schemes for compressible flows of a real gas, Comput. Math. Appl. 48 (2004) 469–480]. These schemes use different parameter vectors in their construction and, consequently, the scheme in [P. Glaister, Conservative upwind difference schemes for compressible flows of a real gas, Comput. Math. Appl. 48 (2004) 469–480] when applied to the special case of an ideal gas is not identical to the principal ideal gas scheme in [P. Glaister, Conservative upwind difference schemes for the Euler equations, Comput. Math. Appl. 45 (2003) 1673–1682]. In this paper it is shown how these schemes are related, followed by a numerical comparison when each is applied to two standard test problems.     

Fast and robust level set update for 3D non-planar X-FEM crack propagation modelling

In the X-FEM framework, the need to represent a discontinuity independently of the structural mesh relies on the level set technique. Hence crack propagation can be simulated by an update of two distinct level sets, the evolution of which is described by differential equations. The aim of this paper is to analyse the resolution of these equations in order to formulate a robust and fast numerical process allowing 3D crack propagation simulations even in presence of high kink angles occurring in mixed mode propagation. The numerical integration is accomplished by means of a robust finite difference upwind scheme applied to an auxiliary regular grid. An alternative level set update equation and a fast localisation of the integration domain, specifically developed for crack propagation problems, are formulated and proposed in the paper in order to gain in stability, robustness and performance.     

Supercloseness of Continuous Interior Penalty Methods on Shishkin Triangular Meshes and Hybrid Meshes for Singularly Perturbed Problems with Characteristic Layers

A fast preconditioned penalty method is developed for a system of parabolic linear complementarity problems (LCPs) involving tempered fractional order partial derivatives governing the price of American options whose underlying asset follows a geometry Lévy process with multi-state regime switching. By means of the penalty method, the system of LCPs is approximated with a penalty term by a system of nonlinear tempered fractional partial differential equations (TFPDEs) coupled by a finite-state Markov chain. The system of nonlinear TFPDEs is discretized with the shifted Grünwald approximation by an upwind finite difference scheme which is shown to be unconditionally stable. Semi-smooth Newton’s method is utilized to solve the finite difference scheme as an outer iterative method in which the Jacobi matrix is found to possess Toeplitz-plus-diagonal structure. Consequently, the resulting linear system can be fast solved by the Krylov subspace method as an inner iterative method via fast Fourier transform (FFT). Furthermore, a novel preconditioner is proposed to speed up the convergence rate of the inner Krylov subspace iteration with theoretical analysis. With the above-mentioned preconditioning technique via FFT, under some mild conditions, the operation cost in each Newton’s step can be expected to be \(\mathcal{O}(N\mathrm{log}N)\), where N is the size of the coefficient matrix. Numerical examples are given to demonstrate the accuracy and efficiency of our proposed fast preconditioned penalty method.     

Parallel hybrid particle/finite volume algorithm for transported PDF methods employing sub-time stepping

A previously presented hybrid finite volume/particle method for the solution of the joint-velocity-frequency-composition probability density function (JPDF) transport equation in complex 3D geometries is extended for parallel computing. The parallelization strategy is based on domain decomposition. The finite volume method (FVM) and the particle method (PM) are parallelized separately and the algorithm is fully synchronous. For the FVM a standard method based on transferring data in ghost cells is used. Moreover, a subdomain interior decomposition algorithm to efficiently solve the implicit time integration for hyperbolic systems is described. The parallelization of the PM is more complicated due to the use of a sub-time stepping algorithm for the particle trajectory integration. Hereby, each particle obeys its local CFL criterion, and the covered distances per global time step can vary significantly. Therefore, an efficient algorithm which deals with this issue and has minimum communication effort was devised and implemented. Numerical tests to validate the parallel vs. the serial algorithm are presented, where also the effectiveness of the subdomain interior decomposition for the implicit time integration was investigated. A 3D dump-combustor configuration test case with about 2.5 × 10⁵ cells was used to demonstrate the good performance of the parallel algorithm. The hybrid algorithm scales well and the maximum speedup on 60 processors for this configuration was 50 (≈80% parallel efficiency).