Trigonometric transform splitting methods for real symmetric Toeplitz systems

In this paper we study efficient iterative methods for real symmetric Toeplitz systems based on the trigonometric transformation splitting (TTS) of the real symmetric Toeplitz matrix $A$. Theoretical analyses show that if the generating function $f$ of the $n\times n$ Toeplitz matrix $A$ is a real positive even function, then the TTS iterative methods converge to the unique solution of the linear system of equations for sufficient large $n$. Moreover, we derive an upper bound of the contraction factor of the TTS iteration which is dependent solely on the spectra of the two TTS matrices involved.Different from the CSCS iterative method in Ng (2003) in which all operations counts concern complex operations when the DFTs are employed, even if the Toeplitz matrix $A$ is real and symmetric, our method only involves real arithmetics when the DCTs and DSTs are used. The numerical experiments show that our method works better than CSCS iterative method and much better than the positive definite and skew-symmetric splitting (PSS) iterative method in Bai et al. (2005) and the symmetric Gauss–Seidel (SGS) iterative method.     

并行J-变量块Cholesky分解算法的仿真研究

该文提出一个针对大型实对称正定稠密方程组或复对称非Hermitian稠密方程组线性求解器的并行分布式算法。它使用了不同于ScaLAPACK的J-变量块Cholesky分解算法和一维块循环列数据分配。该算法以MPI作为消息传递库,在最多可达16个处理器的集群上针对实对称正定稠密方程组可提供与ScaLAPACK近似的浮点操作性能,并可解决一些涉及复对称非Hermitian稠密方程组的电磁场散射问题。该算法的优点是执行Cholesky分解所需的存储量只是标准并行库ScaLAPACK的一半。仿真的数值结果表明该算法是正确、有效的。     

Iterative Toeplitz solvers with local quadratic convergence

We study an iterative, locally quadratically convergent algorithm for solving Toeplitz systems of equations from [R. P. Brent, F. G. Gustavson and D. Y. Y. Yun. “Fast solution of Toeplitz systems of equations and computation of Padé approximations”,J. Algorithms, 1:259–295, 1980]. We introduce a new iterative algorithm that is locally quadratically convergent when used to solve symmetric positive definite Toeplitz systems. We present a set of numerical experiments on randomly generated symmetric positive definite Toeplitz matrices. In these experiments, our algorithm performed significantly better than the previously proposed algorithm.     

Modified HSS iteration methods for a class of complex symmetric linear systems

In this paper, we introduce and analyze a modification of the Hermitian and skew-Hermitian splitting iteration method for solving a broad class of complex symmetric linear systems. We show that the modified Hermitian and skew-Hermitian splitting (MHSS) iteration method is unconditionally convergent. Each iteration of this method requires the solution of two linear systems with real symmetric positive definite coefficient matrices. These two systems can be solved inexactly. We consider acceleration of the MHSS iteration by Krylov subspace methods. Numerical experiments on a few model problems are used to illustrate the performance of the new method.     

Existence and Stability of Choleski Q.I.F. for Symmetric Linear Systems

Recently, Evans [1] has given a variant of the classical Choleski factorization, called Choleski Q.I.F., for the solution of symmetric linear systems. In the present paper we address the questions of existence and stability of the factorization. We show that Evans' Choleski Q.I.F. exists if the coefficient matrix is symmetric positive definite. Our proof of existence of the factorization is via obtaining the diagonal elements of the factor matrices in a closed-form in terms of certain principal minors of the coefficient matrix. By interpreting the Choleski Q.I.F. as an elimination procedure in which eliminations are carried out, alternately, in columns from the left and the right, we show that the Choleski Q.I.F. is also (numerically) strongly stable if the coefficient matrix is symmetric positive definite. We include arithmetical operations counts for the Choleski Q.I.F.; these agree with the counts for the classical Choleski factorization.     

分块带边结构线性规划并行算法

基于内点算法((Interior Point Method,IPM)框架,导出具有分块带边结构系数矩阵的线性规划(Linear Programming, I_P)问题的简化和最简修正方程,并证明最简修正方程的对角分块具有正定性。结合正定矩阵的Cholcsky分解和解藕技术设计了修正方程的并行求解方法,给出了LP的并行内点算法结构。集群环境下的数值实验表明,所提算法具有很好的加速比和可扩展性,适合求解大规模结构化工尹问题。     

Highly scalable parallel algorithms for sparse matrix factorization

In this paper, we describe scalable parallel algorithms for symmetric sparse matrix factorization, analyze their performance and scalability, and present experimental results for up to 1,024 processors on a Gray T3D parallel computer. Through our analysis and experimental results, we demonstrate that our algorithms substantially improve the state of the art in parallel direct solution of sparse linear systems-both in terms of scalability and overall performance. It is a well known fact that dense matrix factorization scales well and can be implemented efficiently on parallel computers. In this paper, we present the first algorithms to factor a wide class of sparse matrices (including those arising from two- and three-dimensional finite element problems) that are asymptotically as scalable as dense matrix factorization algorithms on a variety of parallel architectures. Our algorithms incur less communication overhead and are more scalable than any previously known parallel formulation of sparse matrix factorization. Although, in this paper, we discuss Cholesky factorization of symmetric positive definite matrices, the algorithms can be adapted for solving sparse linear least squares problems and for Gaussian elimination of diagonally dominant matrices that are almost symmetric in structure. An implementation of one of our sparse Cholesky factorization algorithms delivers up to 20 GFlops on a Gray T3D for medium-size structural engineering and linear programming problems. To the best of our knowledge, this is the highest performance ever obtained for sparse Cholesky factorization on any supercomputer     

Parallel computation of the eigenvalues of symmetric Toeplitz matrices through iterative methods

This paper presents a new procedure to compute many or all of the eigenvalues and eigenvectors of symmetric Toeplitz matrices. The key to this algorithm is the use of the “Shift–and–Invert” technique applied with iterative methods, which allows the computation of the eigenvalues close to a given real number (the “shift”). Given an interval containing all the desired eigenvalues, this large interval can be divided in small intervals. Then, the “Shift–and–Invert” version of an iterative method (Lanczos method, in this paper) can be applied to each subinterval. Since the extraction of the eigenvalues of each subinterval is independent from the other subintervals, this method is highly suitable for implementation in parallel computers. This technique has been adapted to symmetric Toeplitz problems, using the symmetry exploiting Lanczos process proposed by Voss [H. Voss, A symmetry exploiting Lanczos method for symmetric Toeplitz matrices, Numerical Algorithms 25 (2000) 377–385] and using fast solvers for the Toeplitz linear systems that must be solved in each Lanczos iteration. The method compares favourably with ScaLAPACK routines, specially when not all the spectrum must be computed.     

Systolic arrays for matrix over-relaxation methods

This paper presents the systolic designs for the Jacobi Overrelaxation (J.O.R.), Successive Overrelaxation (S.O.R.), Accelerated Overrelaxation (A.O.R.), stationary and non-stationary second order Richardson methods for the iterative solution of large linear systems. The linear systems are obtained from the discretisation of a two and three dimensional Laplace equation by the finite difference method and the coefficient matrices are sparse symmetric and positive definite. We investigate the hardware implementation of a linear systolic array to achieve a low cost and optimal area efficient VLSI solution.     

A Method for Outer Interval Solution of Systems of Linear Equations Depending Linearly on Interval Parameters

Consider the systems of linear interval equations whose coefficients are linear functions of interval parameters. Such systems, called parametrized systems of linear interval equations, are encountered in many practical problems, e.g in electrical engineering and structure mechanics. A direct method for computing a tight enclosure for the solution set is proposed in this paper. It is proved that for systems with real matrix and interval right-hand vector the method generates the hull of the solution set. For such systems an explicit formula for the hull is also given. Finally some numerical examples are provided to demonstrate the usefulness of the method in structure mechanics.     

Classical foundations of algorithms for solving positive definite Toeplitz equations

The ongoing development and analysis of efficient algorithms for solving positive definite Toeplitz equations is motivated to a large extent by the importance of these equations in signal processing applications. The role of positive definite Toeplitz matrices in this and other areas of mathematics and engineering stems from Schur's study of bounded analytic functions on the unit disk, and Szegő's theory of polynomials orthogonal on the unit circle. These ideas underlie several Toeplitz solvers, and provide a useful framework for understanding the relationships among these algorithms. In this paper we give an overview of several direct algorithms for solving positive definite Toeplitz systems of linear equations from this classical viewpoint.     

A Simple Derivation of the Hansen-Bliek-Rohn-Ning-Kearfott Enclosure for Linear Interval Equations

Recently, Ning & Kearfott derived a formula for the interval enclosure of the solution set of linear systems of equations with uncertain data ranging in intervals, in the case when the coefficient matrix is an H-matrix. The enclosure is optimal when the midpoint matrix is diagonal, and when the midpoint is the identity, it reduces to the optimal method for enclosing preconditioned systems found by Hansen and Bliek and simplified by Rohn.An elementary proof of this formula is given using only simple properties of H-matrices and Schur complements. The new proof gives additional insight into why the theorem is true. It is also shown how to preserve rigor in the enclosure when finite precision arithmetic is used.     

Block pivoting implementation of a symmetric Toeplitz solver

Toeplitz matrices are characterized by a special structure that can be exploited in order to obtain fast linear system solvers. These solvers are difficult to parallelize due to their low computational cost and their closely coupled data operations. We propose to transform the Toeplitz system matrix into a Cauchy-like matrix since the latter can be divided into two independent matrices of half the size of the system matrix and each one of these smaller arising matrices can be factorized efficiently in multicore computers. We use OpenMP and store data in memory by blocks in consecutive positions yielding a simple and efficient algorithm. In addition, by exploiting the fact that diagonal pivoting does not destroy the special structure of Cauchy-like matrices, we introduce a local diagonal pivoting technique which improves the accuracy of the solution and the stability of the algorithm.     

Second international symposium on the stability of steel structures

The convergence rates of the point-Jacobi method and of the successive overrelaxation (including the Gauss-Seidel) method for symmetric positive definite linear equations are determined as functions of the spectrum of the matrix which is scaled to ones on the principal diagonal. It is seen that the convergence rate is strongly dependent upon the magnitude of the smallest eigenvalue and that for even mildly ill-conditioned matrices the convergence is very slow. It is also shown that the finite difference representation of the three-dimensional Laplace operator leads to remarkably well-conditioned matrices even when the number of equations is very large. Numerical experiments were performed supporting the analysis. It is concluded that the methods studied are virtually limited to large finite difference systems, but that they are very well adapted to such systems.     

Forward and back substitution algorithms on GPU: a case study on modified incomplete Cholesky Preconditioner for three-dimensional finite difference method

Forward and back substitution algorithms are widely used for solving linear systems of equations after performing LU decomposition on the coefficient matrix. They are also essential in the implementation of high performance preconditioners which improve the convergence properties of the various iterative methods. In this paper, we describe an efficient approach to implementing forward and back substitution algorithms on a GPU and provide the implementation details of these algorithms on a Modified Incomplete Cholesky Preconditioner for the Conjugate Gradient (CG) algorithm. The resulting forward and back substitution algorithms are then used on a Modified Incomplete Cholesky Preconditioned Conjugate Gradient method to solve the sparse, symmetric, positive definite and linear systems of equations arising from the discretization of three dimensional finite difference ground-water flow models. By utilizing multiple threads, the proposed method yields speedups up to 60 times on GeForce GTX 280 compared to CPU implementation and up to 4.8 times speedup compared to cuSPARSE library function optimized for GPU by NVIDIA.     

Two Ways to Extend the Cholesky Decomposition to Block Matrices with Interval Entries

We define two block interval Cholesky decompositions in order to bound the solutions of linear systems of equations with symmetric block matrices and right-hand sides both of which are varying within given intervals. Feasibility results are shown for both of them and compared to each other.     

Global optimization for robust control synthesis based on the Matrix Product Eigenvalue Problem

In this paper, we propose a new formulation for a class of optimization problems which occur in general robust control synthesis, called the Matrix Product Eigenvalue Problem (MPEP): Minimize the maximum eigenvalue of the product of two block‐diagonal positive‐definite symmetric matrices under convex constraints. This optimization class falls between methods of guaranteed low complexity such as the linear matrix inequality (LMI) optimization and methods known to be NP‐hard such as the bilinear matrix inequality (BMI) formulation, while still addressing most robust control synthesis problems involving BMIs encountered in applications. The objective of this paper is to provide an algorithm to find a global solution within any specified tolerance ε for the MPEP. We show that a finite number of LMI problems suffice to find the global solution and analyse its computational complexity in terms of the iteration number. We prove that the worst‐case iteration number grows no faster than a polynomial of the inverse of the tolerance given a fixed size of the block‐diagonal matrices in the eigenvalue condition. Copyright 2001 © John Wiley & Sons, Ltd.     

A Method for Outer Interval Solution of Linear Parametric Systems

The paper deals with the problem of determining an outer interval solution (interval enclosure of the solution set) of linear systems whose elements are affine functions of interval parameters. An iterative method for finding such a solution is suggested. A numerical example illustrating the new method is solved.     

Mitigation of symmetry condition in positive realness for adaptive control

Feasibility of nonlinear and adaptive control methodologies in multivariable linear time-invariant systems with state-space realization {A,B,C} is apparently limited by the standard strictly positive realness conditions that imply that the product CB must be positive definite symmetric. This paper expands the applicability of the strictly positive realness conditions used for the proofs of stability of adaptive control or control with uncertainty by showing that the not necessarily symmetric CB is only required to have a diagonal Jordan form and positive eigenvalues. The paper also shows that under the new condition any minimum-phase systems can be made strictly positive real via constant output feedback. The paper illustrates the usefulness of these extended properties with an adaptive control example.