Conjugate gradient methods and ILU preconditioning of non-symmetric matrix systems with arbitrary sparsity patterns

Preconditioning techniques based on incomplete Gaussian elimination for large, sparse, non-symmetric matrix systems are described. A certain level of fill-in may be specified in the incomplete factorizations. All methods considered may be applied to matrices with arbitrary sparsity patterns, for instance those associated with the general preprocessor algorithms or adaptive mesh techniques. The preconditioners have been combined with five conjugate gradient-like methods and tested on finite element discretized scalar convection-diffusion equations in 2D and 3D. It is found from numerical experiments that an amount of fill-in corresponding to about 50% of the number of original non-zero matrix entries is the optimal choice for this class of preconditioners. The preconditioners show almost no sensitivity to grid distortion. In problems with significantly variable coefficients or anisotropy the preconditioners stabilize the basic iterative schemes in addition to reducing the computational work substantially, mostly by more than 90%. The modified preconditioning technique, where fill-in is added on the main diagonal, performs in general better than the standard incomplete LU factorization, but is inferior to the latter in 3D problems and for matrix systems with complicated sparsity patterns.     

An ILU preconditioner with coupled node fill-in for iterative solution of the mixed finite element formulation of the 2D and 3D Navier-Stokes equations

In the present paper, preconditioning of iterative equation solvers for the Navier-Stokes equations is investigated. The Navier-Stokes equations are solved for the mixed finite element formulation. The linear equation solvers used are the orthomin and the Bi-CGSTAB algorithms. The storage structure of the equation matrix is given special attention in order to avoid swapping and thereby increase the speed of the preconditioner. The preconditioners considered are Jacobian, SSOR and incomplete LU preconditioning of the matrix associated with the velocities. A new incomplete LU preconditioning with fill-in for the pressure matrix at locations in the matrix where the corner nodes are coupled is designed. For all preconditioners, inner iterations are investigated for possible improvement of the preconditioning. Numerical experiments are executed both in two and three dimensions.     

ROBUST NUMERICAL METHODS FOR TRANSONIC FLOWS

In this paper, numerical methods for solving the transonic full potential equation are developed. The governing equation is discretized by a flux-biasing finite volume method. The resulting non-linear algebraic system is solved by using a continuation method with full Newton iteration. The continuation method is based on solving a highly ‘upstream-weighted’ discretization and then gradually reducing the upstream weighting. A general PCG-like sparse matrix iterative solver is used to solve the Jacobians at each non-linear step. Various types of incomplete LU (ILU) preconditioners and ordering techniques are compared. Numerical results are presented to demonstrate that these methods are efficient and robust for solving the transonic potential equation in the workstation computing environment. © 1997 by John Wiley & Sons, Ltd.     

COMPUTATION OF TRANSIENT NONLINEAR SHIP WAVES USING AN ADAPTIVE ALGORITHM

An indirect boundary integral method is used to solve transient nonlinear ship wave problems. A resulting mixed boundary value problem is solved at each time-step using a mixed Eulerian– Lagrangian time integration technique. Two dynamic node allocation techniques, which basically distribute nodes on an ever changing body surface, are presented. Both two-sided hyperbolic tangent and variational grid generation algorithms are developed and compared on station curves. A ship hull form is generated in parametric space using a B-spline surface representation. Two-sided hyperbolic tangent and variational adaptive curve grid-generation methods are then applied on the hull station curves to generate effective node placement. The numerical algorithm, in the first method, used two stretching parameters. In the second method, a conservative form of the parametric variational Euler–Lagrange equations is used the perform an adaptive gridding on each station. The resulting unsymmetrical influence coefficient matrix is solved using both a restarted version of GMRES based on the modified Gram–Schmidt procedure and a line Jacobi method based on LU decomposition. The convergence rates of both matrix iteration techniques are improved with specially devised preconditioners. Numerical examples of node placements on typical hull cross-sections using both techniques are discussed and fully nonlinear ship wave patterns and wave resistance computations are presented.     

A comparison of preconditioners for incompressible Navier–Stokes solvers

We consider solution methods for large systems of linear equations that arise from the finite element discretization of the incompressible Navier–Stokes equations. These systems are of the so‐called saddle point type, which means that there is a large block of zeros on the main diagonal. To solve these types of systems efficiently, several block preconditioners have been published. These types of preconditioners require adaptation of standard finite element packages. The alternative is to apply a standard ILU preconditioner in combination with a suitable renumbering of unknowns. We introduce a reordering technique for the degrees of freedom that makes the application of ILU relatively fast. We compare the performance of this technique with some block preconditioners. The performance appears to depend on grid size, Reynolds number and quality of the mesh. For medium‐sized problems, which are of practical interest, we show that the reordering technique is competitive with the block preconditioners. Its simple implementation makes it worthwhile to implement it in the standard finite element method software. Copyright © 2007 John Wiley & Sons, Ltd.     

Numerical perturbation method for approximate solution of poisson's equation on a moderately deforming grid

In problems such as the computation of incompressible flows with moving boundaries, it may be necessary to solve Poisson's equation on a large sequence of related grids. In this paper the LU decomposition of the matrix A ₀ representing Poisson's equation discretized on one grid is used to efficiently obtain an approximate solution on a perturbation of that grid. Instead of doing an LU decomposition of the new matrix A , the RHS is perturbed by a Taylor expansion of A ^?1 about A ₀. Each term in the resulting series requires one ‘backsolve’ using the original LU . Tests using Laplace's equation on a square/rectangle deformation look promising; three and seven correction terms for deformations of 20% and 40% respectively yielded better than 1% accuracy. As another test, Poisson's equation was solved in an ellipse (fully developed flow in a duct) of aspect ratio 2/3 by perturbing about a circle; one correction term yielded better than 1% accuracy. Envisioned applications other than the computation of unsteady incompressible flow include: three-dimensional parabolic problems in tubes of varying cross-section, use of ‘elimination’ techniques other than LU decomposition, and the solution of PDEs other than Poisson's equation.     

Evaluation of the deflated preconditioned CG method to solve bubbly and porous media flow problems on GPU and CPU

In both bubbly and porous media flow, the jumps in coefficients may yield an ill‐conditioned linear system. The solution of this system using an iterative technique like the conjugate gradient (CG) is delayed because of the presence of small eigenvalues in the spectrum of the coefficient matrix. To accelerate the convergence, we use two levels of preconditioning. For the first level, we choose between out‐of‐the‐box incomplete LU decomposition, sparse approximate inverse, and truncated Neumann series‐based preconditioner. For the second level, we use deflation. Through our experiments, we show that it is possible to achieve a computationally fast solver on a graphics processing unit. The preconditioners discussed in this work exhibit fine‐grained parallelism. We show that the graphics processing unit version of the two‐level preconditioned CG can be up to two times faster than a dual quad core CPU implementation. John Wiley & Sons, Ltd.     

ILUS: An incomplete LU preconditioner in sparse skyline format

Incomplete LU factorizations are among the most effective preconditioners for solving general large, sparse linear systems arising from practical engineering problems. This paper shows how an ILU factorization may be easily computed in sparse skyline storage format, as opposed to traditional row-by-row schemes. This organization of the factorization has many advantages, including its amenability when the original matrix is in skyline format, the ability to dynamically monitor the stability of the factorization and the fact that factorizations may be produced with symmetric structure. Numerical results are presented for Galerkin finite element matrices arising from the standard square lid-driven cavity problem. © 1997 John Wiley & Sons, Ltd.     

Numerical study on the vortex motion patterns around a rotating circular cylinder and their critical characters

A hybrid finite difference and vortex method (HFDV), based on the domain decomposition method (DDM), is used for calculating the flow around a rotating circular cylinder at Reynolds number Re=1000, 200 and the angular‐to‐rectilinear speed ratio α∈(0.5, 3.25) respectively. A fully implicit third‐order eccentric finite difference scheme is adopted in the finite difference method, and the deduced large broad band sparse matrix equations are solved by a highly efficient modified incomplete LU decomposition conjugate gradient method (MILU‐CG). The long‐time, fully developed features about the variations of the vortex patterns in the wake, as well as the drag and lift forces on the cylinder, are given. The calculated streamline contours are in good agreement with the experimentally visualized flow pictures. The existence of the critical state is confirmed again, and the single side shed vortex pattern at the critical state is shown for the first time. Also, the optimized lift‐to‐drag force ratio is obtained near the critical state. Copyright © 1999 John Wiley & Sons, Ltd.     

ℋ︁‐LU factorization in preconditioners for augmented Lagrangian and grad‐div stabilized saddle point systems

The (mixed finite element) discretization of the linearized Navier–Stokes equations leads to a linear system of equations of saddle point type. The iterative solution of this linear system requires the construction of suitable preconditioners, especially in the case of high Reynolds numbers. In the past, a stabilizing approach has been suggested which does not change the exact solution but influences the accuracy of the discrete solution as well as the effectiveness of iterative solvers. This stabilization technique can be performed on the continuous side before the discretization, where it is known as ‘grad‐div’ (GD) stabilization, as well as on the discrete side where it is known as an ‘augmented Lagrangian’ (AL) technique (and does not change the discrete solution). In this paper, we study the applicability of ??‐LU factorizations to solve the arising subproblems in the different variants of stabilized saddle point systems. We consider both the saddle point systems that arise from the stabilization in the continuous as well as on the discrete setting. Recently, a modified AL preconditioner has been proposed for the system resulting from the discrete stabilization. We provide a straightforward generalization of this approach to the GD stabilization. We conclude the paper with numerical tests for a variety of problems to illustrate the behavior of the considered preconditioners as well as the suitability of ??‐LU factorization in the preconditioners. Copyright © 2010 John Wiley & Sons, Ltd.     

ELEMENT－BY－ELEMENT MATRIX DECOMPOSITION ANDSTEP－BY－STEP INTEGRATION METHOD FOR TRANSIENTDYNAMIC PROBLEMS

ＥＬＥＭＥＮＴ－ＢＹ－ＥＬＥＭＥＮＴＭＡＴＲＩＸＤＥＣＯＭＰＯＳＩＴＩＯＮＡＮＤＳＴＥＰ－ＢＹ－ＳＴＥＰＩＮＴＥＧＲＡＴＩＯＮＭＥＴＨＯＤＦＯＲＴＲＡＮＳＩＥＮＴＤＹＮＡＭＩＣＰＲＯＢＬＥＭＳＷａｎｇＨｕａｉｚｈｏｎｇ（王怀忠）（ＲｅｃｅｉｖｅｄＪｕ...     

Fast boundary‐domain integral algorithm for the computation of incompressible fluid flow problems

The paper deals with the numerical solution of fluid dynamics using the boundary‐domain integral method (BDIM). A velocity–vorticity formulation of the Navier–Stokes equations is adopted, where the kinematic equation is written in its parabolic form. Computational aspects of the numerical simulation of two‐dimensional flows is described in detail. In order to lower the computational cost, the subdomain technique is applied. A preconditioned Krylov subspace method (PKSM) is used for the solution of systems of linear equations. Level‐based fill‐in incomplete lower upper decomposition (ILU) preconditioners are developed and their performance is examined. Scaling of stopping criteria is applied to minimize the number of iterations for the PKSM. The effectiveness of the proposed method is tested on several benchmark test problems. Copyright © 1999 John Wiley & Sons, Ltd.     

Flame computations by a Chebyshev multidomain method

A Chebyshev collocation method is proposed for the computation of laminar flame propagation in a two-dimensional gaseous medium. The method is based on a domain decomposition technique associated with co-ordinate transforms to map the infinite physical subdomains into finite computational ones. The influence matrix method is used to handle the patching conditions at the interfaces. This technique is particularly efficient since at each time step only matrix products have to be performed. The method is tested first on an elliptic model problem; it is then applied to laminar flame computations, including calculations of cellular instabilities of flame fronts.     

Preconditioned conjugate gradient methods for the incompressible Navier-Stokes equations

A robust technique for solving primitive variable formulations of the incompressible Navier-Stokes equations is to use Newton iteration for the fully implicit non-linear equations. A direct sparse matrix method can be used to solve the Jacobian but is costly for large problems; an alternative is to use an iterative matrix method. This paper investigates effective ways of using a conjugate-gradient-type method with an incomplete LU factorization preconditioner for two-dimensional incompressible viscous flow problems. Special attention is paid to the ordering of unknowns, with emphasis on a minimum updating matrix (MUM) ordering. Numerical results are given for several test problems.     

带源参数的二维热传导反问题的无网格方法

利用无网格有限点法求解带源参数的二维热传导反问题,推导了相应的离散方程. 与 其它基于网格的方法相比,有限点法采用移动最小二乘法构造形函数,只需要节点信息,不 需要划分网格,用配点法离散控制方程,可以直接施加边界条件,不需要在区域内部求积分. 用有限点法求解二维热传导反问题具有数值实现简单、计算量小、可以任意布置节点等优点. 最后通过算例验证了该方法的有效性.     

BAND STRUCTURE AND TRANSMISSION CHARACTERISTICS CALCULATION OF PHONONIC CRYSTALS BASED ON MULTI-LEVEL SUBSTRUCTURE TECHNIQUE

基于多重多级子结构方法提出一种快速的声子晶体能带与传输特性的计算策略. 主要思想是将声子晶体划分成多层级子结构有限元模型,在能带计算中采用静凝聚和子结构周游树技术将子结构的内部刚度阵凝聚到Bloch 边界上. 由于内部刚度阵并不随着简约波矢变化,所以这种计算策略可以大大降低求解规模并提高计算效率,并不对整体有限元模型引入近似. 在传输特性计算中同样采用该策略,由于声子晶体单胞具有周期性,所以各个单胞的系数矩阵是相同的,从而减少计算量,并且可以灵活地选择是否回代求解单胞内部自由度. 数值算例以三维局域共振型声子晶体和二维Bragg 散射型声子晶体为例,计算结果验证了这种求解策略的正确性和高效性,并适用于复杂声子晶体分析.     

基于多重多级子结构声子能带与传输特性分析

基于多重多级子结构方法提出一种快速的声子晶体能带与传输特性的计算策略. 主要思想是将声子晶体划分成多层级子结构有限元模型,在能带计算中采用静凝聚和子结构周游树技术将子结构的内部刚度阵凝聚到Bloch 边界上. 由于内部刚度阵并不随着简约波矢变化,所以这种计算策略可以大大降低求解规模并提高计算效率,并不对整体有限元模型引入近似. 在传输特性计算中同样采用该策略,由于声子晶体单胞具有周期性,所以各个单胞的系数矩阵是相同的,从而减少计算量,并且可以灵活地选择是否回代求解单胞内部自由度. 数值算例以三维局域共振型声子晶体和二维Bragg 散射型声子晶体为例,计算结果验证了这种求解策略的正确性和高效性,并适用于复杂声子晶体分析.      

Improved preconditioned conjugate gradient method and its application in F. E. A. for engineering

In this paper two theorems with theoretical and practical significance are given in respect to the preconditioned conjugate gradient method(PCCG).The theorems discuss respectively the qualitative property of the iterative solution and the construction principle of the iterative matrix.The authors put forward a new incompletely LU factorizing technique for non-M-matrix and the method of constructing the iterative matrix.This improved PCCG is used to calculate the ill-conditioned problems and large-scale three-dimensional finite element problems,and simultaneously contrasted with other methods.The abnormal phenomenon is analyzed when PCCG is used to solve the system of ill-conditioned equations,It is shown that the method proposed in this paper is quite effective in solving the system of large-scale finite element equations and the system of ill-conditioned equations.     

直接施加本质边界条件的FE／EFG耦合算法及其数值实现

提出一种可以直接施加本质边界条件的有限元与无网格Galerkin（FE／EFG）耦合算法。将问题域分成FE和EFG两种类型的子域,采用转换矩阵耦舍两子域的交界面;通过另一转换矩阵将无网格区域本质边界上的名义位移转换成真实位移,从而可在其上直接施加本质边界条件;采用二次转换实现两种转换矩阵之间的协调。提出全域统一采用单元...     

Finite element model updating using variable separation

In the field of structural dynamics, reliable finite element response predictions are becoming increasingly important to industry and there is a genuine interest to improve these in the light of measured frequency response functions. Unlike modal-based model updating formulations, response-based methods have been applied only with limited success due to incomplete measurements and numerical ill-conditioning problems. The least squares approximation method is one of the methods used but often poses a problem of pseudo inverse due to the number of incomplete measurements. The proposed algorithm is a modification and extension of a previously-developed nonlinear least squares method for damage detection and finite element model updating. The paper derives explicit expressions for the first and second order partial derivatives with respect to the correction parameters and for the Jacobian matrix used in the Newton–Raphson solution of the nonlinear set of equations in order to avoid the pseudo inverse and to build a symmetrical system. The proposed method, assigned to a frequency parameterization which considers the minimum distance to be minimized, shows a good numerical stability. The performance of the method in localizing structural damage and updating model is examined using simulated measurements.