On the indefinite Helmholtz equation: Complex stretched absorbing boundary layers,iterative analysis,and preconditioning

This paper studies and analyzes a preconditioned Krylov solver for Helmholtz problems that are formulated with absorbing boundary layers based on complex coordinate stretching. The preconditioner problem is a Helmholtz problem where not only the coordinates in the absorbing layer have an imaginary part, but also the coordinates in the interior region. This results into a preconditioner problem that is invertible with a multigrid cycle. We give a numerical analysis based on the eigenvalues and evaluate the performance with several numerical experiments. The method is an alternative to the complex shifted Laplacian and it gives a comparable performance for the studied model problems.     

A sweeping preconditioner for time-harmonic Maxwell’s equations with finite elements

This paper is concerned with preconditioning the stiffness matrix resulting from finite element discretizations of Maxwell’s equations in the high frequency regime. The moving PML sweeping preconditioner, first introduced for the Helmholtz equation on a Cartesian finite difference grid, is generalized to an unstructured mesh with finite elements. The method dramatically reduces the number of GMRES iterations necessary for convergence, resulting in an almost linear complexity solver. Numerical examples including electromagnetic cloaking simulations are presented to demonstrate the efficiency of the proposed method.     

二维声学多层快速多极子边界元及其应用

与模型自由度的平方成正比的存储量和计算量,使传统边界元无法应用到大型模型的计算。为此,发展了一种二维声学多层快速多极子边界元算法。通过二维Helmholtz核函数展开理论的简要介绍,推导了源点矩计算、源点矩转移、源点矩至本地展开转移、本地展开转移公式,并详细描述了二维声学快速多极子边界元算法的具体实现步骤。使用快速傅里叶插值进行源点矩和本地展开系数的多层传递。采用对角左预处理方法,改善边界方程的条件数,减少迭代求解次数。最后通过数值算例,验证了所发展的二维声学快速多极子算法的正确性和高效性。      

An Iterative Two-Grid Method of a Finite Element PML Approximation for the Two Dimensional Maxwell Problem

In this paper, we propose an iterative two-grid method for the edge finite element discretizations (a saddle-point system) of Perfectly Matched Layer (PML) equations to the Maxwell scattering problem in two dimensions. Firstly, we use a fine space to solve a discrete saddle-point system of $H(grad)$ variational problems, denoted by auxiliary system 1. Secondly, we use a coarse space to solve the original saddle-point system. Then, we use a fine space again to solve a discrete$\boldsymbol{H}(curl)$-elliptic variational problems, denoted by auxiliary system 2. Furthermore, we develop a regularization diagonal block preconditioner for auxiliary system 1 and use $H$-$X$ preconditioner for auxiliary system 2. Hence we essentially transform the original problem in a fine space to a corresponding (but much smaller) problem on a coarse space, due to the fact that the above two preconditioners are efficient and stable. Compared with some existing iterative methods for solving saddle-point systems, such as PMinres, numerical experiments show the competitive performance of our iterative two-grid method.     

Performance of a parallel algebraic multilevel preconditioner for stabilized finite element semiconductor device modeling

In this study results are presented for the large-scale parallel performance of an algebraic multilevel preconditioner for solution of the drift-diffusion model for semiconductor devices. The preconditioner is the key numerical procedure determining the robustness, efficiency and scalability of the fully-coupled Newton–Krylov based, nonlinear solution method that is employed for this system of equations. The coupled system is comprised of a source term dominated Poisson equation for the electric potential, and two convection–diffusion-reaction type equations for the electron and hole concentration. The governing PDEs are discretized in space by a stabilized finite element method. Solution of the discrete system is obtained through a fully-implicit time integrator, a fully-coupled Newton-based nonlinear solver, and a restarted GMRES Krylov linear system solver. The algebraic multilevel preconditioner is based on an aggressive coarsening graph partitioning of the nonzero block structure of the Jacobian matrix. Representative performance results are presented for various choices of multigrid V-cycles and W-cycles and parameter variations for smoothers based on incomplete factorizations. Parallel scalability results are presented for solution of up to 10⁸${10}^8$ unknowns on 4096 processors of a Cray XT3/4 and an IBM POWER eServer system.     

A dispersion minimizing finite difference scheme and preconditioned solver for the 3D Helmholtz equation

In this paper, a new 27-point finite difference method is presented for solving the 3D Helmholtz equation with perfectly matched layer (PML), which is a second order scheme and pointwise consistent with the equation. An error analysis is made between the numerical wavenumber and the exact wavenumber, and a refined choice strategy based on minimizing the numerical dispersion is proposed for choosing weight parameters. A full-coarsening multigrid-based preconditioned Bi-CGSTAB method is developed for solving the linear system stemming from the Helmholtz equation with PML by the finite difference scheme. The shifted-Laplacian is extended to precondition the 3D Helmholtz equation, and a spectral analysis is given. The discrete preconditioned system is solved by the Bi-CGSTAB method, with a multigrid method used to invert the preconditioner approximately. Full-coarsening multigrid is employed, and a new matrix-based prolongation operator is constructed accordingly. Numerical results are presented to demonstrate the efficiency of both the new 27-point finite difference scheme with refined parameters, and the preconditioned Bi-CGSTAB method with the 3D full-coarsening multigrid.     

An efficient preconditioner for monolithically-coupled large-displacement fluid–structure interaction problems with pseudo-solid mesh updates

We present a block preconditioner for the efficient solution of the linear systems that arise when employing Newton’s method to solve monolithically-coupled large-displacement fluid–structure interaction problems in which the update of the moving fluid mesh is performed by the equations of large-displacement elasticity. Following a theoretical analysis of the preconditioner, we propose an efficient implementation that yields a solver with near-optimal computational cost, in the sense that the time for the solution of the linear systems scales approximately linearly with the number of unknowns. We evaluate the performance of the preconditioner in selected two- and three-dimensional test problems.     

Efficient algebraic multigrid for migration–diffusion–convection–reaction systems arising in electrochemical simulations

The article discusses components and performance of an algebraic multigrid (AMG) preconditioner for the fully coupled multi-ion transport and reaction model (MITReM) with nonlinear boundary conditions, important for electrochemical modeling. The governing partial differential equations (PDEs) are discretized in space by a combined finite element and residual distribution method. Solution of the discrete system is obtained by means of a Newton-based nonlinear solver, and an AMG-preconditioned BICGSTAB Krylov linear solver. The presented AMG preconditioner is based on so-called point-based classical AMG. The linear solver is compared to a standard direct and several one-level iterative solvers for a range of geometries and chemical systems with scientific and industrial relevance. The results indicate that point-based AMG methods, carefully designed, are an attractive alternative to more commonly employed numerical methods for the simulation of complex electrochemical processes.     

An eigen-based high-order expansion basis for structured spectral elements

We present an eigen-based high-order expansion basis for the spectral element approach with structured elements. The new basis exhibits a numerical efficiency significantly superior, in terms of the conditioning of coefficient matrices and the number of iterations to convergence for the conjugate gradient solver, to the commonly-used Jacobi polynomial-based expansion basis. This basis results in extremely sparse mass matrices, and it is very amenable to the diagonal preconditioning. Ample numerical experiments demonstrate that with the new basis and a simple diagonal preconditioner the number of conjugate gradient iterations to convergence has essentially no dependence or only a very weak dependence on the element order. The expansion bases are constructed by a tensor product of a set of special one-dimensional (1D) basis functions. The 1D interior modes are constructed such that the interior mass and stiffness matrices are simultaneously diagonal and have identical condition numbers. The 1D vertex modes are constructed to be orthogonal to all the interior modes. The performance of the new basis has been investigated and compared with other expansion bases.     

Helmholtz and parabolic equation solutions to a benchmark problem in ocean acoustics

The Helmholtz equation (HE) describes wave propagation in applications such as acoustics and electromagnetics. For realistic problems, solving the HE is often too expensive. Instead, approximations like the parabolic wave equation (PE) are used. For low-frequency shallow-water environments, one persistent problem is to assess the accuracy of the PE model. In this work, a recently developed HE solver that can handle a smoothly varying bathymetry, variable material properties, and layered materials, is used for an investigation of the errors in PE solutions. In the HE solver, a preconditioned Krylov subspace method is applied to the discretized equations. The preconditioner combines domain decomposition and fast transform techniques. A benchmark problem with upslope-downslope propagation over a penetrable lossy seamount is solved. The numerical experiments show that, for the same bathymetry, a soft and slow bottom gives very similar HE and PE solutions, whereas the PE model is far from accurate for a hard and fast bottom. A first attempt to estimate the error is made by computing the relative deviation from the energy balance for the PE solution. This measure gives an indication of the magnitude of the error, but cannot be used as a strict error bound.     

Accelerating a FFT-based solver for numerical homogenization of periodic media by conjugate gradients

In this short note, we present a new technique to accelerate the convergence of a FFT-based solver for numerical homogenization of complex periodic media proposed by Moulinec and Suquet [1]. The approach proceeds from discretization of the governing integral equation by the trigonometric collocation method due to Vainikko [2], to give a linear system which can be efficiently solved by conjugate gradient methods. Computational experiments confirm robustness of the algorithm with respect to its internal parameters and demonstrate significant increase of the convergence rate for problems with high-contrast coefficients at a low overhead per iteration.     

A Fast High Order Iterative Solver for the Electromagnetic Scattering by Open Cavities Filled with the Inhomogeneous Media

The scattering of the open cavity filled with the inhomogeneous media is studied. The problem is discretized with a fourth order finite difference scheme and the immersed interface method, resulting in a linear system of equations with the high order accurate solutions in the whole computational domain. To solve the system of equations, we design an efficient iterative solver, which is based on the fast Fourier transformation, and provides an ideal preconditioner for Krylov subspace method. Numerical experiments demonstrate the capability of the proposed fast high order iterative solver.     

A Boundary Meshless Method for Solving Heat Transfer Problems Using the Fourier Transform

Fourier transform is applied to remove the time-dependent variable in the diffusion equation. Under non-harmonic initial conditions this gives rise to a non-homogeneous Helmholtz equation, which is solved by the method of fundamental solutions and the method of particular solutions. The particular solution of Helmholtz equation is available as shown in [4, 15]. The approximate solution in frequency domain is then inverted numerically using the inverse Fourier transform algorithm. Complex frequencies are used in order to avoid aliasing phenomena and to allow the computation of the static response. Two numerical examples are given to illustrate the effectiveness of the proposed approach for solving 2-D diffusion equations.     

地球外核热对流运动的区域分解多重网格并行数值模拟

旋转球层中热对流运动的数值模拟是地球发电机模型的重要组成部分,对研究地球发电机作用机理具有重要意义.本文设计一个基于国产超级计算平台并行性能良好的地球外核热对流运动并行数值模型.时间积分方案采用与Crank-Nicolson格式和二阶Adams-Bashford公式相结合的近似分解分步法,空间离散基于立方球网格的二阶精度有限体积格式.所得到的两个大规模稀疏线性代数方程组采用带预处理的Krylov子空间迭代法进行求解.为加速迭代求解过程及提高并行性能,迭代过程采用区域分解多重网格的多层限制型加法Schwarz预处理子,减少了求解程序的计算时间,提高了数值模型的并行性能,模型被很好地扩展到上万处理器核数.数值模拟结果与基准模型算例0的参考值吻合得很好.     

The equilibrium limit of a constitutive model for two-phase granular mixtures and its numerical approximation

In this paper we analyze the equilibrium limit of the constitutive model for two-phase granular mixtures introduced in Papalexandris (2004) [13], and develop an algorithm for its numerical approximation. At, equilibrium, the constitutive model reduces to a strongly coupled, overdetermined system of quasilinear elliptic partial differential equations with respect to the pressure and the volume fraction of the solid granular phase. First we carry a perturbation analysis based on standard hydrostatic-type scaling arguments which reduces the complexity of the coupling of the equations. The perturbed system is then supplemented by an appropriate compatibility condition which arises from the properties of the gradient operator. Further, based on the Helmholtz decomposition and Ladyzhenskaya’s decomposition theorem, we develop a projection-type, Successive-Over-Relaxation numerical method. This method is general enough and can be applied to a variety of continuum models of complex mixtures and mixtures with micro-structure. We also prove that this method is both stable and consistent hence, under standard assumptions, convergent. The paper concludes with the presentation of representative numerical results.     

Preconditioning bandgap eigenvalue problems in three-dimensional photonic crystals simulations

To explore band structures of three-dimensional photonic crystals numerically, we need to solve the eigenvalue problems derived from the governing Maxwell equations. The solutions of these eigenvalue problems cannot be computed effectively unless a suitable combination of eigenvalue solver and preconditioner is chosen. Taking eigenvalue problems due to Yee’s scheme as examples, we propose using Krylov–Schur method and Jacobi–Davidson method to solve the resulting eigenvalue problems. For preconditioning, we derive several novel preconditioning schemes based on various preconditioners, including a preconditioner that can be solved by Fast Fourier Transform efficiently. We then conduct intensive numerical experiments for various combinations of eigenvalue solvers and preconditioning schemes. We find that the Krylov–Schur method associated with the Fast Fourier Transform based preconditioner is very efficient. It remarkably outperforms all other eigenvalue solvers with common preconditioners like Jacobi, Symmetric Successive Over Relaxation, and incomplete factorizations. This promising solver can benefit applications like photonic crystal structure optimization.     

On Convergence of a Least-Squares Kansa's Method for the Modified Helmholtz Equations

We analyze a least-squares asymmetric radial basis function collocation method for solving the modified Helmholtz equations. In the theoretical part, we proved the convergence of the proposed method providing that the collocation points are sufficiently dense. For numerical verification, direct solver and a subspace selection process for the trial space (the so-called adaptive greedy algorithm) is employed, respectively, for small and large scale problems.     

A Riemann solver for single-phase and two-phase shallow flow models based on relaxation. Relations with Roe and VFRoe solvers

We present a Riemann solver derived by a relaxation technique for classical single-phase shallow flow equations and for a two-phase shallow flow model describing a mixture of solid granular material and fluid. Our primary interest is the numerical approximation of this two-phase solid/fluid model, whose complexity poses numerical difficulties that cannot be efficiently addressed by existing solvers. In particular, we are concerned with ensuring a robust treatment of dry bed states. The relaxation system used by the proposed solver is formulated by introducing auxiliary variables that replace the momenta in the spatial gradients of the original model systems. The resulting relaxation solver is related to Roe solver in that its Riemann solution for the flow height and relaxation variables is formally computed as Roe’s Riemann solution. The relaxation solver has the advantage of a certain degree of freedom in the specification of the wave structure through the choice of the relaxation parameters. This flexibility can be exploited to handle robustly vacuum states, which is a well known difficulty of standard Roe’s method, while maintaining Roe’s low diffusivity. For the single-phase model positivity of flow height is rigorously preserved. For the two-phase model positivity of volume fractions in general is not ensured, and a suitable restriction on the CFL number might be needed. Nonetheless, numerical experiments suggest that the proposed two-phase flow solver efficiently models wet/dry fronts and vacuum formation for a large range of flow conditions.As a corollary of our study, we show that for single-phase shallow flow equations the relaxation solver is formally equivalent to the VFRoe solver with conservative variables of Gallouët and Masella [T. Gallouët, J.-M. Masella, Un schéma de Godunov approché C.R. Acad. Sci. Paris, Série I, 323 (1996) 77–84]. The relaxation interpretation allows establishing positivity conditions for this VFRoe method.     

An h-adaptive finite element solver for the calculations of the electronic structures

In this paper, a framework of using h-adaptive finite element method for the Kohn–Sham equation on the tetrahedron mesh is presented. The Kohn–Sham equation is discretized by the finite element method, and the h-adaptive technique is adopted to optimize the accuracy and the efficiency of the algorithm. The locally optimal block preconditioned conjugate gradient method is employed for solving the generalized eigenvalue problem, and an algebraic multigrid preconditioner is used to accelerate the solver. A variety of numerical experiments demonstrate the effectiveness of our algorithm for both the all-electron and the pseudo-potential calculations.     

An efficient block variant of robust structured multifrontal factorization method

Based on the two-dimensional three-temperature （2D3T） radiation diffusion equations and its discrete system, using the block diagonal structure of the three-temperature matrix, the reordering and symbolic decomposition parts of the RSMF method are replaced with corresponding block operation in order to improve the solution efficiency. We call this block form method block RSMF （in brief, BRSMF） method. The new BRSMF method not only makes the reordering and symbolic decomposition become more effective, but also keeps the cost of numerical factorization from increasing and ensures the precision of solution very well. The theoretical analysis of the computation complexity about the new BRSMF method shows that the solution efficiency about the BRSMF method is higher than the original RSMF method. The numerical experiments also show that the new BRSMF method is more effective than the original RSMF method.