On the impact of triangle shapes for boundary layer problems using high-order finite element discretization

The impact of triangle shapes, including angle sizes and aspect ratios, on accuracy and stiffness is investigated for simulations of highly anisotropic problems. The results indicate that for high-order discretizations, large angles do not have an adverse impact on solution accuracy. However, a correct aspect ratio is critical for accuracy for both linear and high-order discretizations. Large angles are also found to be not problematic for the conditioning of the linear systems arising from the discretizations. Further, when choosing preconditioning strategies, coupling strengths among elements rather than element angle sizes should be taken into account. With an appropriate preconditioner, solutions on meshes with and without large angles can be achieved within a comparable time.     

一种可扩展的三维半导体器件并行数值模拟算法

在基于漂移-扩散模型的三维半导体器件数值模拟中,通过有限体积法进行数值离散,采用完全耦合的牛顿迭代求解非线性代数方程组,并使用基于代数多重网格预条件子的GMRES方法求解牛顿迭代中的线性方程组,构造一种稳健且高度可扩展的非结构四面体网格上求解半导体方程的并行算法.基于PHG平台实现该算法的并行计算程序,并对PN结和MOS场效应晶体管等问题进行了最大网格规模达到5亿单元、最大并行规模达到1 024进程的大规模数值模拟实验,结果表明,该算法计算效率高,可扩展性好.     

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.     

Preconditioning methods for discontinuous Galerkin solutions of the Navier–Stokes equations

A Newton–Krylov method is developed for the solution of the steady compressible Navier–Stokes equations using a discontinuous Galerkin (DG) discretization on unstructured meshes. Steady-state solutions are obtained using a Newton–Krylov approach where the linear system at each iteration is solved using a restarted GMRES algorithm. Several different preconditioners are examined to achieve fast convergence of the GMRES algorithm. An element Line-Jacobi preconditioner is presented which solves a block-tridiagonal system along lines of maximum coupling in the flow. An incomplete block-LU factorization (Block-ILU(0)) is also presented as a preconditioner, where the factorization is performed using a reordering of elements based upon the lines of maximum coupling. This reordering is shown to be superior to standard reordering techniques (Nested Dissection, One-way Dissection, Quotient Minimum Degree, Reverse Cuthill–Mckee) especially for viscous test cases. The Block-ILU(0) factorization is performed in-place and an algorithm is presented for the application of the linearization which reduces both the memory and CPU time over the traditional dual matrix storage format. Additionally, a linear p-multigrid preconditioner is also considered, where Block-Jacobi, Line-Jacobi and Block-ILU(0) are used as smoothers. The linear multigrid preconditioner is shown to significantly improve convergence in term of number of iterations and CPU time compared to a single-level Block-Jacobi or Line-Jacobi preconditioner. Similarly the linear multigrid preconditioner with Block-ILU smoothing is shown to reduce the number of linear iterations to achieve convergence over a single-level Block-ILU(0) preconditioner, though no appreciable improvement in CPU time is shown.     

Extension of the spectral volume method to high-order boundary representation

In this paper, the spectral volume method is extended to the two-dimensional Euler equations with curved boundaries. It is well-known that high-order methods can achieve higher accuracy on coarser meshes than low-order methods. In order to realize the advantage of the high-order spectral volume method over the low order finite volume method, it is critical that solid wall boundaries be represented with high-order polynomials compatible with the order of the interpolation for the state variables. Otherwise, numerical errors generated by the low-order boundary representation may overwhelm any potential accuracy gains offered by high-order methods. Therefore, more general types of spectral volumes (or elements) with curved edges are used near solid walls to approximate the boundaries with high fidelity. The importance of this high-order boundary representation is demonstrated with several well-know inviscid flow test cases, and through comparisons with a second-order finite volume method.     

磁流体方腔槽道流整体线性稳定性数值方法研究

将Chebyshev谱配置法和基于非均匀网格的高阶FD-q差分格式运用于磁流体方腔槽道流整体线性稳定性研究,比较两类数值方法的优缺点.Chebyshev谱配置法收敛快且精度高,但需要构造非常庞大的满矩阵,存储量和计算开销巨大;高阶FD-q差分格式采用了基于Kosloff-Tal-Ezer变换的Chebyshev谱配置点作为离散网格,在保持较高网格收敛精度的同时,差分格式可以采用稀疏矩阵进行存储,显著降低了存储量和计算开销.相比传统的谱配置法,基于非均匀网格的高阶FD-q差分格式计算效率得到显著的提升,将高阶FD-q差分格式运用于非正则模线性最优瞬态增长的计算,计算效果良好.     

High-order adaptive methods for parabolic systems

We consider the adaptive solution of parabolic partial differential systems in one and two space dimensions by finite element procedures that automatically refine and coarsen computational meshes, vary the degree of the piecewise polynomial basis and, in one dimension, move the computational mesh. Two-dimensional meshes of triangular, quadrilateral, or a mixture of triangular and quadrilateral elements are generated using a finite quadtree procedure that is also used for data management. A posteriori estimates, used to control adaptive enrichment, are generated from the hierarchical polynomial basis. Temporal integration, within a method-of-lines framework, uses either backward difference methods or a variant of the singly implicit Runge-Kutta (SIRK) methods. A high-level user interface facilitates use of the adaptive software.     

An anisotropic mesh adaptation method for the finite element solution of heterogeneous anisotropic diffusion problems

Heterogeneous anisotropic diffusion problems arise in the various areas of science and engineering including plasma physics, petroleum engineering, and image processing. Standard numerical methods can produce spurious oscillations when they are used to solve those problems. A common approach to avoid this difficulty is to design a proper numerical scheme and/or a proper mesh so that the numerical solution validates the discrete counterpart (DMP) of the maximum principle satisfied by the continuous solution. A well known mesh condition for the DMP satisfaction by the linear finite element solution of isotropic diffusion problems is the non-obtuse angle condition that requires the dihedral angles of mesh elements to be non-obtuse. In this paper, a generalization of the condition, the so-called anisotropic non-obtuse angle condition, is developed for the finite element solution of heterogeneous anisotropic diffusion problems. The new condition is essentially the same as the existing one except that the dihedral angles are now measured in a metric depending on the diffusion matrix of the underlying problem. Several variants of the new condition are obtained. Based on one of them, two metric tensors for use in anisotropic mesh generation are developed to account for DMP satisfaction and the combination of DMP satisfaction and mesh adaptivity. Numerical examples are given to demonstrate the features of the linear finite element method for anisotropic meshes generated with the metric tensors.     

A field expansions method for scattering by periodic multilayered media

The interaction of acoustic and electromagnetic waves with periodic structures plays an important role in a wide range of problems of scientific and technological interest. This contribution focuses upon the robust and high-order numerical simulation of a model for the interaction of pressure waves generated within the earth incident upon layers of sediment near the surface. Herein described is a boundary perturbation method for the numerical simulation of scattering returns from irregularly shaped periodic layered media. The method requires only the discretization of the layer interfaces (so that the number of unknowns is an order of magnitude smaller than finite difference and finite element simulations), while it avoids not only the need for specialized quadrature rules but also the dense linear systems characteristic of boundary integral/element methods. The approach is a generalization to multiple layers of Bruno and Reitich's "Method of Field Expansions" for dielectric structures with two layers. By simply considering the entire structure simultaneously, rather than solving in individual layers separately, the full field can be recovered in time proportional to the number of interfaces. As with the original field expansions method, this approach is extremely efficient and spectrally accurate.     

Smoothed aggregation multigrid solvers for high-order discontinuous Galerkin methods for elliptic problems

We develop a smoothed aggregation-based algebraic multigrid solver for high-order discontinuous Galerkin discretizations of the Poisson problem. Algebraic multigrid is a popular and effective method for solving the sparse linear systems that arise from discretizing partial differential equations. However, high-order discontinuous Galerkin discretizations have proved challenging for algebraic multigrid. The increasing condition number of the matrix and loss of locality in the matrix stencil as p increases, in addition to the effect of weakly enforced Dirichlet boundary conditions all contribute to the challenging algebraic setting.     

Multigrid algorithms for high-order discontinuous Galerkin discretizations of the compressible Navier–Stokes equations

Multigrid algorithms are developed for systems arising from high-order discontinuous Galerkin discretizations of the compressible Navier–Stokes equations on unstructured meshes. The algorithms are based on coupling both p- and h-multigrid (ph-multigrid) methods which are used in nonlinear or linear forms, and either directly as solvers or as preconditioners to a Newton–Krylov method.The performance of the algorithms are examined in solving the laminar flow over an airfoil configuration. It is shown that the choice of the cycling strategy is crucial in achieving efficient and scalable solvers. For the multigrid solvers, while the order-independent convergence rate is obtained with a proper cycle type, the mesh-independent performance is achieved only if the coarsest problem is solved to a sufficient accuracy. On the other hand, the multigrid preconditioned Newton–GMRES solver appears to be insensitive to this condition and mesh-independent convergence is achieved under the desirable condition that the coarsest problem is solved using a fixed number of multigrid cycles regardless of the size of the problem.It is concluded that the Newton–GMRES solver with the multigrid preconditioning yields the most efficient and robust algorithm among those studied.     

二维多流体域耦合声场的光滑有限元解法

针对声学有限元分析中四节点等参单元计算精度低,对网格质量敏感的问题,将光滑有限元法引入到多流体域耦合声场的数值分析中,提出了二维多流体域耦合声场的光滑有限元解法。该方法在Helmholtz控制方程与多流体域耦合界面的声压/质点法向速度连续条件的基础上,得到二维多流体耦合声场的离散控制方程,并采用光滑有限元的分区光滑技术将声学梯度矩阵形函数导数的域内积分转换形函数的域边界积分,避免了雅克比矩阵的计算。以管道二维多流体域耦合内声场为数值分析算例,研究结果表明,与标准有限元相比,对单元尺寸较大或扭曲严重的四边形网格模型,光滑有限元的计算精度更高。因此光滑有限元能很好地应用于大尺寸单元或扭曲严重的网格模型下二维多流体域耦合声场的预测,具有良好的工程应用前景。      

Mimetic finite difference method for the Stokes problem on polygonal meshes

Various approaches to extend finite element methods to non-traditional elements (general polygons, pyramids, polyhedra, etc.) have been developed over the last decade. The construction of basis functions for such elements is a challenging task and may require extensive geometrical analysis. The mimetic finite difference (MFD) method works on general polygonal meshes and has many similarities with low-order finite element methods. Both schemes try to preserve the fundamental properties of the underlying physical and mathematical models. The essential difference between the two schemes is that the MFD method uses only the surface representation of discrete unknowns to build the stiffness and mass matrices. Since no extension of basis functions inside the mesh elements is required, practical implementation of the MFD method is simple for polygonal meshes that may include degenerate and non-convex elements. In this article, we present a new MFD method for the Stokes problem on arbitrary polygonal meshes and analyze its stability. The method is developed for the general case of tensor coefficients, which allows us to apply it to a linear elasticity problem, as well. Numerical experiments show, for the velocity variable, second-order convergence in a discrete L² norm and first-order convergence in a discrete H¹ norm. For the pressure variable, first-order convergence is shown in the L² norm.     

Analysis of the monotonicity conditions in the mimetic finite difference method for elliptic problems

The maximum principle is one of the most important properties of solutions of partial differential equations. Its numerical analog, the discrete maximum principle (DMP), is one of the most difficult properties to achieve in numerical methods, especially when the computational mesh is distorted to adapt and conform to the physical domain or the problem coefficients are highly heterogeneous and anisotropic. Violation of the DMP may lead to numerical instabilities such as oscillations and to unphysical solutions such as heat flow from a cold material to a hot one. In this work, we investigate sufficient conditions to ensure the monotonicity of the mimetic finite difference (MFD) method on two- and three-dimensional meshes. These conditions result in a set of general inequalities for the elements of the mass matrix of every mesh element. Efficient solutions are devised for meshes consisting of simplexes, parallelograms and parallelepipeds, and orthogonal locally refined elements as those used in the AMR methodology. On simplicial meshes, it turns out that the MFD method coincides with the mixed-hybrid finite element methods based on the low-order Raviart–Thomas vector space. Thus, in this case we recover the well-established conventional angle conditions of such approximations. Instead, in the other cases a suitable design of the MFD method allows us to formulate a monotone discretization for which the existence of a DMP can be theoretically proved. Moreover, on meshes of parallelograms we establish a connection with a similar monotonicity condition proposed for the Multi-Point Flux Approximation (MPFA) methods. Numerical experiments confirm the effectiveness of the considered monotonicity conditions.     

Interpolation-free monotone finite volume method for diffusion equations on polygonal meshes

We developed a new monotone finite volume method for diffusion equations. The second-order linear methods, such as the multipoint flux approximation, mixed finite element and mimetic finite difference methods, are not monotone on strongly anisotropic meshes or for diffusion problems with strongly anisotropic coefficients. The finite volume (FV) method with linear two-point flux approximation is monotone but not even first-order accurate in these cases. The developed monotone method is based on a nonlinear two-point flux approximation. It does not require any interpolation scheme and thus differs from other nonlinear finite volume methods based on a two-point flux approximation. The second-order convergence rate is verified with numerical experiments.     

An analysis of time discretization in the finite element solution of hyperbolic problems

The problem of the time discretization of hyperbolic equations when finite elements are used to represent the spatial dependence is critically examined. A modified equation analysis reveals that the classical, second-order accurate, time-stepping algorithms, i.e., the Lax-Wendroff, leap-frog, and Crank-Nicolson methods, properly combine with piecewise linear finite elements in advection problems only for small values of the time step. On the contrary, as the Courant number increases, the numerical phase error does not decrease uniformly at all wavelengths so that the optimal stability limit and the unit CFL property are not achieved. These fundamental numerical properties can, however, be recovered, while still remaining in the standard Galerkin finite element setting, by increasing the order of accuracy of the time discretization. This is accomplished by exploiting the Taylor series expansion in the time increment up to the third order before performing the Galerkin spatial discretization using piecewise linear interpolations. As a result, it appears that the proper finite element equivalents of second-order finite difference schemes are implicit methods of incremental type having third- and fourth-order global accuracy on uniform meshes (Taylor-Galerkin methods). Numerical results for several linear examples are presented to illustrate the properties of the Taylor-Galerkin schemes in one- and two-dimensional calculations.     

Bad behavior of Godunov mixed methods for strongly anisotropic advection–dispersion equations

We study the performance of Godunov mixed methods, which combine a mixed-hybrid finite element solver and a Godunov-like shock-capturing solver, for the numerical treatment of the advection–dispersion equation with strong anisotropic tensor coefficients. It turns out that a mesh locking phenomenon may cause ill-conditioning and reduce the accuracy of the numerical approximation especially on coarse meshes. This problem may be partially alleviated by substituting the mixed-hybrid finite element solver used in the discretization of the dispersive (diffusive) term with a linear Galerkin finite element solver, which does not display such a strong ill conditioning. To illustrate the different mechanisms that come into play, we investigate the spectral properties of such numerical discretizations when applied to a strongly anisotropic diffusive term on a small regular mesh. A thorough comparison of the stiffness matrix eigenvalues reveals that the accuracy loss of the Godunov mixed method is a structural feature of the mixed-hybrid method. In fact, the varied response of the two methods is due to the different way the smallest and largest eigenvalues of the dispersion (diffusion) tensor influence the diagonal and off-diagonal terms of the final stiffness matrix. One and two dimensional test cases support our findings.     

A taxonomy and comparison of parallel block multi-level preconditioners for the incompressible Navier–Stokes equations

In recent years, considerable effort has been placed on developing efficient and robust solution algorithms for the incompressible Navier–Stokes equations based on preconditioned Krylov methods. These include physics-based methods, such as SIMPLE, and purely algebraic preconditioners based on the approximation of the Schur complement. All these techniques can be represented as approximate block factorization (ABF) type preconditioners. The goal is to decompose the application of the preconditioner into simplified sub-systems in which scalable multi-level type solvers can be applied. In this paper we develop a taxonomy of these ideas based on an adaptation of a generalized approximate factorization of the Navier–Stokes system first presented in [A. Quarteroni, F. Saleri, A. Veneziani, Factorization methods for the numerical approximation of Navier–Stokes equations, Computational Methods in Applied Mechanical Engineering 188 (2000) 505–526]. This taxonomy illuminates the similarities and differences among these preconditioners and the central role played by efficient approximation of certain Schur complement operators. We then present a parallel computational study that examines the performance of these methods and compares them to an additive Schwarz domain decomposition (DD) algorithm. Results are presented for two and three-dimensional steady state problems for enclosed domains and inflow/outflow systems on both structured and unstructured meshes. The numerical experiments are performed using MPSalsa, a stabilized finite element code.     

High-Order Mesh Generation for Discontinuous Galerkin Methods Based on Elastic Deformation

In this paper, a high-order curved mesh generation method for Discontinuous Galerkin methods is introduced. First, a regular mesh is generated. Second, the solid surface is re-constructed using cubic polynomial. Third, the elastic governing equations are solved using high-order finite element method to provide a fully or partly curved grid. Numerical tests indicate that the intersection between element boundaries can be avoided by carefully defining the elasticity modulus.     

A Relaxed Dimensional Factorization preconditioner for the incompressible Navier–Stokes equations

In this paper we introduce a Relaxed Dimensional Factorization (RDF) preconditioner for saddle point problems. Properties of the preconditioned matrix are analyzed and compared with those of the closely related Dimensional Splitting (DS) preconditioner recently introduced by Benzi and Guo [7]. Numerical results for a variety of finite element discretizations of both steady and unsteady incompressible flow problems indicate very good behavior of the RDF preconditioner with respect to both mesh size and viscosity.