Explicit preconditioned iterative methods for solving large unsymmetric finite element systems

A class of Generalized Approximate Inverse Matrix (GAIM) techniques, based on the concept of LU-sparse factorization procedures, is introduced for computing explicitly approximate inverses of large sparse unsymmetric matrices of irregular structure, without inverting the decomposition factors. Explicit preconditioned iterative methods, in conjunction with modified forms of the GAIM techniques, are presented for solving numerically initial/boundary value problems on multiprocessor systems. Application of the new methods on linear boundary-value problems is discussed and numerical results are given.     

Vectorizable preconditioners for mixed finite element solution of second-order elliptic problems

We compare the performance of the CG (conjugate gradient) method, the point-wise incomplete factorization preconditioned CG method, and the block-incomplete factorization preconditioned CG method for solving problems arising in mixed finite element discretization of second order elliptic differential equations. The robustness and vectorizable properties of these methods are illustrated by a large set of numerical tests.     

GPU-accelerated preconditioned iterative linear solvers

This work is an overview of our preliminary experience in developing a high-performance iterative linear solver accelerated by GPU coprocessors. Our goal is to illustrate the advantages and difficulties encountered when deploying GPU technology to perform sparse linear algebra computations. Techniques for speeding up sparse matrix-vector product (SpMV) kernels and finding suitable preconditioning methods are discussed. Our experiments with an NVIDIA TESLA M2070 show that for unstructured matrices SpMV kernels can be up to 8 times faster on the GPU than the Intel MKL on the host Intel Xeon X5675 Processor. Overall performance of the GPU-accelerated Incomplete Cholesky (IC) factorization preconditioned CG method can outperform its CPU counterpart by a smaller factor, up to 3, and GPU-accelerated The incomplete LU (ILU) factorization preconditioned GMRES method can achieve a speed-up nearing 4. However, with better suited preconditioning techniques for GPUs, this performance can be further improved.     

Parallel methods for optimality criteria-based topology optimization

Topology optimization problems require the repeated solution of finite element problems that are often extremely ill-conditioned due to highly heterogeneous material distributions. This makes the use of iterative linear solvers inefficient unless appropriate preconditioning is used. Even then, the solution time for topology optimization problems is typically very high. These problems are addressed by considering the use of non-overlapping domain decomposition-based parallel methods for the solution of topology optimization problems. The parallel algorithms presented here are based on the solid isotropic material with penalization (SIMP) formulation of the topology optimization problem and use the optimality criteria method for iterative optimization. We consider three parallel linear solvers to solve the equilibrium problem at each step of the iterative optimization procedure. These include two preconditioned conjugate gradient (PCG) methods: one using a diagonal preconditioner and one using an incomplete LU factorization preconditioner with a drop tolerance. A third substructuring solver that employs a hybrid of direct and iterative (PCG) techniques is also studied. This solver is found to be the most effective of the three solvers studied, both in terms of parallel efficiency and in terms of its ability to mitigate the effects of ill-conditioning. In addition to examining parallel linear solvers, we consider the parallelization of the iterative optimality criteria method. To tackle checkerboarding and mesh dependence, we propose a multi-pass filtering technique that limits the number of “ghost” elements that need to be exchanged across interprocessor boundaries.     

The preconditioned simultaneous displacement method (PSD method) for elliptic difference equations

This paper introduces the Preconditioned Simultaneous Displacement iterative method (PSD method) in a new “computable” form for the numerical solution of linear systems of the form Au=b, where the matrix A is large and sparse. The convergence properties of the method are analysed under certain assumptions on the matrix A. Moreover, “good” values (near the optimum) for the involved parameters are determined in terms of bounds on the eigenvalues of certain matrices. Bounds on the reciprocal rate of convergence of the PSD method are also given. The method is shown to be superior over the well known Symmetric Successive Overrelaxation method (SSOR method) (at the optimum stage PSD is shown to converge approximately two times faster than SSOR) and in certain cases over the Successive Overrelaxation method (SOR method).     

Efficient iterative solvers for structural dynamics problems

Direct solvers are commonly used in implicit finite element codes for structural dynamics problems. This study explores an alternative approach to solving the resulting linear systems by using preconditioned iterative schemes such as the preconditioned conjugate gradient algorithm and its variants. Preconditioners used in this study include approximate Cholesky factorization, block Jacobi, and the symmetric Gauss-Seidel over-relaxation scheme. The effects of various preconditioners and ordering schemes on the solution time and required storage are investigated. Performance results of these iterative solver are compared against the MA57 direct solver routine of the commercial Harwell Software Library.     

On Accelerated Iterative Methods for the Solution of Systems of Linear Equations

In literature, it has been reported that the convergence of some preconditioned stationary iterative methods using certain type upper triangular matrices as preconditioners are faster than the basic iterative methods. In this paper, a new preconditioned iterative method for the numerical solution of linear systems has been introduced, and the convergence analysis of the proposed method and an existing one have been done. Some numerical examples have also been given, which show the effectiveness of both of the methods.     

Explicit approximate inverse preconditioning techniques

Summary  The numerical treatment and the production of related software for solving large sparse linear systems of algebraic equations, derived mainly from the discretization of partial differential equation, by preconditioning techniques has attracted the attention of many researchers. In this paper we give an overview of explicit approximate inverse matrix techniques for computing explicitly various families of approximate inverses based on Choleski and LU—type approximate factorization procedures for solving sparse linear systems, which are derived from the finite difference, finite element and the domain decomposition discretization of elliptic and parabolic partial differential equations. Composite iterative schemes, using inner-outer schemes in conjunction with Picard and Newton method, based on approximate inverse matrix techniques for solving non-linear boundary value problems, are presented. Additionally, isomorphic iterative methods are introduced for the efficient solution of non-linear systems. Explicit preconditioned conjugate gradient—type schemes in conjunction with approximate inverse matrix techniques are presented for the efficient solution of linear and non-linear system of algebraic equations. Theoretical estimates on the rate of convergence and computational complexity of the explicit preconditioned conjugate gradient method are also presented. Applications of the proposed methods on characteristic linear and non-linear problems are discussed and numerical results are given.     

一种针对大波数Helmholtz方程的高性能并行预条件迭代求解算法

针对传统串行迭代法求解大波数Helmholtz方程存在效率低下且受限于单机内存的问题,提出了一种基于消息传递接口(Message Passing Interface,MPI) 的并行预条件迭代法。该算法利用复移位拉普拉斯算子对Helmholtz方程进行预条件处理,联合稳定双共轭梯度法和基于矩阵的多重网格法来求解预条件方程离散后的大规模线性系统,在Linux集群系统上基于 MPI环境实现了求解算法的并行计算,重点解决了多重网格的并行划分、信息传递和多重网格组件的构建问题。数值实验表明,对于大波数问题,提出的算法具有良好的并行加速比,相较于串行算法极大地提高了计算效率。     

On sparse lu factorization procedures for the solution of parabolic differential equations in three space dimensions

The numerical solution of partial differential equations in 3 dimensions by finite difference methods leads to the problem of solving large order sparse structured linear systems.

In this paper, a factorization procedure in algorithmic form is derived yielding direct and iterative methods of solution of some interesting boundary value problems in physics and engineering.     

基于曙光并行机的超大规模非线性方程组并行算法研究

该文讨论了一类求解大规模非线性方程组算法的并行性能及其在曙光并行机上的实现过程，与传统的算法不同之处是用一个块对角矩阵作为迭代矩阵，且该矩阵可由一个仅包含向量内积和矩阵与向量乘积的递推关系简便计算得到，在对算法进行描述之后，分析了算法的并行加速比和存储需求，讨论了算法在基于消息传递的MPI并行环境下的实现流程，数值计算表明理论分析与数值结果相比，算法在分布式并行环境下具有有较好的并行主攻较低的存储要求，可适用于大规模科学与工程的高性能计算。     

A partial preconditioned conjugate gradient method for large eigenproblems

In this paper a preconditioned iterative method suitable for the solution of the generalized eigenvalue problem Ax = λBx is presented. The proposed method is suitable for the determination of extreme eigenvalues and their corresponding eigenvectors of large sparse eigenproblems derived from the finite element discretization method. The solution is obtained through the minimization of the Rayleigh quotient by a preconditioned conjugate gradient (CG) method. The proposed triangular splitting preconditioner combines Evans' SSOR preconditioner with a drop-off tolerance criterion and is called partial preconditioner (PPR). The PPR is attractive in a large FE framework because it is simple and can be implemented at the element level, as opposed to incomplete Choleski preconditioners, which require a sparse assembly. Because of the renewed interest in CG techniques for FE work on microprocessors and parallel computers it is believed that this improved approach to the generalized eigenvalue problem, through the minimization of the Rayleigh quotient, is likely to be very promising.     

Systolic Preconditioned Algorithms for the Iterative Solution of Linear Systems

Systolic algorithms for preconditioned iterative procedures are described and discussed in relation to existing systolic arrays for iterative methods applied to the systems of linear equations arising from the solution of partial differential equations. In boundary value problems it is shown that if the cost of the preconditioned systolic arrays in terms of hardware is related to the (standard) iterative methods, then savings in the number of array cells can be made when the system involved is large and sparse (narrow bandwidth) with a significant improvement in convergence rate.     

一种LU分解与迭代法的结合策略及算法实现

在矩阵求解算法中，直接法或迭代法都不能有效地求解大规模稀疏或病态矩阵,因此提出一种LU分解与迭代法结合的策略。采用LU分解对矩阵进行预处理，以提高迭代法的收敛性,并采用一种判断策略使矩阵的LU分解结果可最大限度地重复利用。此结合策略应用于两种共轭梯度（CG）法，得到CLUCG和CLUTCG两种算法。它们已应用于模拟和混合信号电路模拟器ZeniVDE中。大量实验结果表明此结合策略是很有效的，得到的两种算法具有较快的速度和较好的收敛性。     

Efficient enforcement of hard articulation constraints in the presence of closed loops and contacts

In rigid body simulation, one must distinguish between contacts (so‐called unilateral constraints) and articulations (bilateral constraints). For contacts and friction, iterative solution methods have proven most useful for interactive applications, often in combination with Shock‐Propagation in cases with strong interactions between contacts (such as stacks), prioritizing performance and plausibility over accuracy. For articulation constraints, direct solution methods are preferred, because one can rely on a factorization with linear time complexity for tree‐like systems, even in ill‐conditioned cases caused by large mass‐ratios or high complexity. Despite recent advances, combining the advantages of direct and iterative solution methods wrt. performance has proven difficult and the intricacy of articulations in interactive applications is often limited by the convergence speed of the iterative solution method in the presence of closed kinematic loops (i.e. auxiliary constraints) and contacts. We identify common performance bottlenecks in the dynamic simulation of unilateral and bilateral constraints and are able to present a simulation method, that scales well in the number of constraints even in ill‐conditioned cases with frictional contacts, collisions and closed loops in the kinematic graph. For cases where many joints are connected to a single body, we propose a technique to increase the sparsity of the positive definite linear system. A solution to these bottlenecks is presented in this paper to make the simulation of a wider range of mechanisms possible in real‐time without extensive parameter tuning.     

Numerical solution of nonlinear elliptic partial differential equations by a generalized conjugate gradient method

We have studied previously a generalized conjugate gradient method for solving sparse positive-definite systems of linear equations arising from the discretization of elliptic partial-differential boundary-value problems. Here, extensions to the nonlinear case are considered. We split the original discretized operator into the sum of two operators, one of which corresponds to a more easily solvable system of equations, and accelerate the associated iteration based on this splitting by (nonlinear) conjugate gradients. The behavior of the method is illustrated for the minimal surface equation with splittings corresponding to nonlinear SSOR, to approximate factorization of the Jacobian matrix, and to elliptic operators suitable for use with fast direct methods. The results of numerical experiments are given as well for a mildy nonlinear example, for which, in the corresponding linear case, the finite termination property of the conjugate gradient algorithm is crucial.     

New updates of incomplete LU factorizations and applications to large nonlinear systems

In this paper, we address the problem of preconditioning sequences of large sparse indefinite systems of linear equations and present two new strategies to construct approximate updates of factorized preconditioners. Both updates are based on the availability of an incomplete factorization for one matrix of the sequence and differ in the approximation of the so-called ideal update. For a general treatment, an incomplete LU (ILU) factorization is considered, but the proposed approaches apply to incomplete factorizations of symmetric matrices as well. The first strategy is an approximate diagonal update of the ILU factorization; the second strategy relies on banded approximations of the factors in the ideal update. The efficiency and reliability of the proposed preconditioners are shown in the solution of nonlinear systems of equations by preconditioned Newton–Krylov methods. Nearly matrix-free implementations of the updating strategy are provided, and numerical experiments are carried out on application problems.     

三维对流扩散方程的稀疏存储及预条件迭代

基于四阶紧致格式对三维对流扩散方程进行离散,并给出所得到的离散线性方程组的块三角稀疏矩阵形式。以带双阈值的不完全因子化LU分解[(ILUT(τ,s))]作为预条件子,分别用FGMRES、BICGSTAB和TFQMR作为迭代加速器,对离散线性方程组进行求解验证了格式精度并比较了不同迭代法的CPU时间和迭代步。此外,通过比较传统迭代法和预条件迭代法的计算效率,表明预条件迭代法不仅能够保证格式的四阶精度,还能极大地提高收敛效率。     

Shape reanalysis and sensitivities utilizing preconditioned iterative boundary solvers

The computational advantages associated with the utilization of preconditioned iterative equation solvers are quantified for the reanalysis of perturbed shapes using continuum structural boundary element analysis (BEA). Both single- and multizone three-dimensional problems are examined. Significant redutions in computer time are obtained by making use of previously computed solution vectors and preconditioners in subsequent analyses. The effectiveness of this technique is demonstrated for the computation of shape response sensitivities required in shape optimization. Computer times and accuracies achieved using the preconditioned iterative solvers are compared with those obtained via direct solvers and implicit differentiation of the boundary integral equations. It is concluded that this approach employing preconditioned iterative equation solvers in reanalysis and sensitivity analysis can be competitive with if not superior to those involving direct solvers.     

Optimal left and right additive Schwarz preconditioning for minimal residual methods with Euclidean and energy norms

For the solution of non-symmetric or indefinite linear systems arising from discretizations of elliptic problems, two-level additive Schwarz preconditioners are known to be optimal in the sense that convergence bounds for the preconditioned problem are independent of the mesh and the number of subdomains. These bounds are based on some kind of energy norm. However, in practice, iterative methods which minimize the Euclidean norm of the residual are used, despite the fact that the usual bounds are non-optimal, i.e., the quantities appearing in the bounds may depend on the mesh size; see [X.-C. Cai, J. Zou, Some observations on the l² convergence of the additive Schwarz preconditioned GMRES method, Numer. Linear Algebra Appl. 9 (2002) 379-397]. In this paper, iterative methods are presented which minimize the same energy norm in which the optimal Schwarz bounds are derived, thus maintaining the Schwarz optimality. As a consequence, bounds for the Euclidean norm minimization are also derived, thus providing a theoretical justification for the practical use of Euclidean norm minimization methods preconditioned with additive Schwarz. Both left and right preconditioners are considered, and relations between them are derived. Numerical experiments illustrate the theoretical developments.