Application of a block modified Chebyshev algorithm to the iterative solution of symmetric linear systems with multiple right hand side vectors

Summary. An adaptive Richardson iteration method is described for the solution of large sparse symmetric positive definite linear systems of equations with multiple right-hand side vectors. This scheme ``learns' about the linear system to be solved by computing inner products of residual matrices during the iterations. These inner products are interpreted as block modified moments. A block version of the modified Chebyshev algorithm is presented which yields a block tridiagonal matrix from the block modified moments and the recursion coefficients of the residual polynomials. The eigenvalues of this block tridiagonal matrix define an interval, which determines the choice of relaxation parameters for Richardson iteration. Only minor modifications are necessary in order to obtain a scheme for the solution of symmetric indefinite linear systems with multiple right-hand side vectors. We outline the changes required. Received April 22, 1993     

Iterative Methods for a Linearly Perturbed Algebraic Matrix Riccati Equation Arising in Stochastic Control

We start with a discussion of coupled algebraic Riccati equations arising in the study of linear-quadratic optimal control problems for Markov jump linear systems. Under suitable assumptions, this system of equations has a unique positive semidefinite solution, which is the solution of practical interest. The coupled equations can be rewritten as a single linearly perturbed matrix Riccati equation with special structures. We study the linearly perturbed Riccati equation in a more general setting and obtain a class of iterative methods from different splittings of a positive operator involved in the Riccati equation. We prove some special properties of the sequences generated by these methods and determine and compare the convergence rates of these methods. Our results are then applied to the coupled Riccati equations of jump linear systems. We obtain linear convergence of the Lyapunov iteration and the modified Lyapunov iteration, and confirm that the modified Lyapunov iteration indeed has faster convergence than the original Lyapunov iteration.     

The use of approximate factorization in stiff ODE solvers

We consider implicit integration methods for the numerical solution of stiff initial-value problems. In applying such methods, the implicit relations are usually solved by Newton iteration. However, it often happens that in subintervals of the integration interval the problem is nonstiff or mildly stiff with respect to the stepsize. In these nonstiff subintervals, we do not need the (expensive) Newton iteration process. This motivated us to look for an iteration process that converges in mildly stiff situations and is less costly than Newton iteration. The process we have in mind uses modified Newton iteration as the outer iteration process and a linear solver for solving the linear Newton systems as an inner iteration process. This linear solver is based on an approximate factorization of the Newton system matrix by splitting this matrix into its lower and upper triangular part. The purpose of this paper is to combine fixed point iteration, approximate factorization iteration and Newton iteration into one iteration process for use in initial-value problems where the degree of stiffness is changing during the integration.     

Numerical solution of positive control problem via linear programming

The reachability problem for linear time-invariant discrete-time control systems with sign-restricted input is considered. The time-optimal control is constructed by an iterative procedure. Each step of the iteration is defined as a linear programming problem. This problem is solved by the simplex algorithm. The initial feasible solution for the simplex algorithm is provided by the preceding step of the iteration. The inversion of the basis matrix is reduced to a bordering procedure. The structural stability of the solution is investigated.     

A hybrid iterative method for symmetric positive definite linear systems

A hybrid iterative scheme that combines the Conjugate Gradient (CG) method with Richardson iteration is presented. This scheme is designed for the solution of linear systems of equations with a large sparse symmetric positive definite matrix. The purpose of the CG iterations is to improve an available approximate solution, as well as to determine an interval that contains all, or at least most, of the eigenvalues of the matrix. This interval is used to compute iteration parameters for Richardson iteration. The attraction of the hybrid scheme is that most of the iterations are carried out by the Richardson method, the simplicity of which makes efficient implementation on modern computers possible. Moreover, the hybrid scheme yields, at no additional computational cost, accurate estimates of the extreme eigenvalues of the matrix. Knowledge of these eigenvalues is essential in some applications.Research supported in part by NSF grant DMS-9409422.Research supported in part by NSF grant DMS-9205531.     

Coloring update methods

For linear update methods (such as SOR), a coloring method is introduced for which the multicolor iteration matrix has the same spectrum as the original iteration matrix. It applies to general linear systems, not necessarily arising from PDEs. When the iteration matrices are nonsingular, it is shown that they are similar to each other.     

Parallel iterative solvers for boundary value methods

A parallel variant of the block Gauss-Seidel iteration for the solution of block-banded linear systems is presented. The coefficient matrix is partitioned among the processors as in the domain decomposition methods and then it is split so that the resulting iterative method has the same spectral properties of the block Gauss-Seidel iteration.The parallel algorithm is applied to the solution of block-banded linear systems arising from the numerical discretization of initial value problems by means of Boundary Value Methods (BVMs). BVMs define a new approach for the solution of ordinary differential equations and seem to be attractive for their interesting stability properties and a possible parallel implementation. In this paper, we refer to BVMs based on the extended trapezoidal rules.     

A fast LU update for linear programming

This paper discusses sparse matrix kernels of simplex-based linear programming software. State-of-the-art implementations of the simplex method maintain an LU factorization of the basis matrix which is updated at each iteration. The LU factorization is used to solve two sparse sets of linear equations at each iteration. We present new implementation techniques for a modified Forrest-Tomlin LU update which reduce the time complexity of the update and the solution of the associated sparse linear systems. We present numerical results on Netlib and other real-life LP models.     

Parallel linear system solvers for Runge-Kutta methods

If the nonlinear systems arising in implicit Runge-Kutta methods like the Radau IIA methods are iterated by (modified) Newton, then we have to solve linear systems whose matrix of coefficients is of the form I-A hJ with A the Runge-Kutta matrix and J an approximation to the Jacobian of the righthand side function of the system of differential equations. For larger systems of differential equations, the solution of these linear systems by a direct linear solver is very costly, mainly because of the LU-decompositions. We try to reduce these costs by solving the linear systems by a second (inner) iteration process. This inner iteration process is such that each inner iteration again requires the solution of a linear system. However, the matrix of coefficients in these new linear systems is of the form I - B hJ where B is similar to a diagonal matrix with positive diagonal entries. Hence, after performing a similarity transformation, the linear systems are decoupled into s subsystems, so that the costs of the LU-decomposition are reduced to the costs of s LU-decompositions of dimension d. Since these LU-decompositions can be computed in parallel, the effective LU-costs on a parallel computer system are reduced by a factor s ³ . It will be shown that matrices B can be constructed such that the inner iterations converge whenever A and J have their eigenvalues in the positive and nonpositive halfplane, respectively. The theoretical results will be illustrated by a few numerical examples. A parallel implementation on the four-processor Cray-C98/4256 shows a speed-up ranging from at least 2.4 until at least 3.1 with respect to RADAU5 applied in one-processor mode.     

A new single-step iteration method for solving complex symmetric linear systems

For solving a class of complex symmetric linear systems, we introduce a new single-step iteration method, which can be taken as a fixed-point iteration adding the asymptotical error (FPAE). In order to accelerate the convergence, we further develop the parameterized variant of the FPAE (PFPAE) iteration method. Each iteration of the FPAE and the PFPAE methods requires the solution of only one linear system with a real symmetric positive definite coefficient matrix. Under suitable conditions, we derive the spectral radius of the FPAE and the PFPAE iteration matrices, and discuss the quasi-optimal parameters which minimize the above spectral radius. Numerical tests support the contention that the PFPAE iteration method has comparable advantage over some other commonly used iteration methods, particularly when the experimental optimal parameters are not used.     

SEQUENTIAL SYSTEMS OF LINEAR EQUATIONS ALGORITHM FOR NONLINEAR OPTIMIZATION PROBLEMS-INEQUALITY CONSTRAINED PROBLEMS

AbstractIn this paper, a new superlinearly convergent algorithm of sequential systems of linear equations (SSLE) for nonlinear optimization problems with inequality constraints is proposed. Since the new algorithm only needs to solve several systems of linear equations having a same coefficient matrix per iteration, the computation amount of the algorithm is much less than that of the existing SQP algorithms per iteration. Moreover, for the SQP type algorithms, there exist so-called inconsistent problems, i.e., quadratic programming subproblems of the SQP algorithms may not have a solution at some iterations, but this phenomenon will not occur with the SSLE algorithms because the related systems of linear equations always have solutions. Some numerical results are reported.     

The Cayley transform in the numerical solution of unitary differential systems

In recent years some numerical methods have been developed to integrate matrix differential systems whose solutions are unitary matrices. In this paper we propose a new approach that transforms the original problem into a skew-Hermitian differential system by means of the Cayley transform. The new methods are semi-explicit, that is, no iteration is required but the solution of a certain number of linear matrix systems at each step is needed. Several numerical comparisons with known unitary integrators are reported. This revised version was published online in June 2006 with corrections to the Cover Date.     

Generalized skew‐Hermitian triangular splitting iteration methods for saddle‐point linear systems

A generalized skew‐Hermitian triangular splitting iteration method is presented for solving non‐Hermitian linear systems with strong skew‐Hermitian parts. We study the convergence of the generalized skew‐Hermitian triangular splitting iteration methods for non‐Hermitian positive definite linear systems, as well as spectrum distribution of the preconditioned matrix with respect to the preconditioner induced from the generalized skew‐Hermitian triangular splitting. Then the generalized skew‐Hermitian triangular splitting iteration method is applied to non‐Hermitian positive semidefinite saddle‐point linear systems, and we prove its convergence under suitable restrictions on the iteration parameters. By specially choosing the values of the iteration parameters, we obtain a few of the existing iteration methods in the literature. Numerical results show that the generalized skew‐Hermitian triangular splitting iteration methods are effective for solving non‐Hermitian saddle‐point linear systems with strong skew‐Hermitian parts. Copyright © 2013 John Wiley & Sons, Ltd.     

A hybrid iterative method for symmetric indefinite linear systems

Hybrid iterative methods that combine a conjugate direction method with a simpler iteration scheme, such as Chebyshev or Richardson iteration, were first proposed in the 1950s. The ease with which Chebyshev and Richardson iteration can be implemented efficiently on a large variety of computer architectures has in recent years lead to renewed interest in iterative methods that use Chebyshev or Richardson iteration. This paper presents a new hybrid iterative method for the solution of linear systems of equations with a symmetric indefinite matrix. Our method combines the conjugate residual method with Richardson iteration. Special attention is paid to the determination of two real intervals, one on each side of the origin, that contain most of the eigenvalues of the matrix. These intervals are used to compute suitable iteration parameters for Richardson iteration. We also discuss when to switch between the methods. The hybrid scheme typically uses the Richardson method for most iterations, and this reduces the number of arithmetic vector operations significantly compared with the number of arithmetic vector operations required when only the conjugate residual method is used. Computed examples illustrate the competitiveness of the hybrid scheme.     

Inner and outer bounds for the solution set of parametric linear systems

Consider linear systems involving affine-linear dependencies on interval parameters. Presented is a free C-XSC software implementing a generalized parametric fixed-point iteration method for verified enclosure of the parametric solution set. Some specific features of the corresponding algorithm concerning sharp enclosure of the contracting matrix and inner approximation of the solution enclosure are discussed.     

A two-step matrix splitting iteration paradigm based on one single splitting for solving systems of linear equations

For solving large sparse systems of linear equations, we construct a paradigm of two-step matrix splitting iteration methods and analyze its convergence property for the nonsingular and the positive-definite matrix class. This two-step matrix splitting iteration paradigm adopts only one single splitting of the coefficient matrix, together with several arbitrary iteration parameters. Hence, it can be constructed easily in actual applications, and can also recover a number of representatives of the existing two-step matrix splitting iteration methods. This result provides systematic treatment for the two-step matrix splitting iteration methods, establishes rigorous theory for their asymptotic convergence, and enriches algorithmic family of the linear iteration solvers, for the iterative solutions of large sparse linear systems.     

Block Two-stage Methods for Singular Systems and Markov Chains

The use of block two-stage methods for the iterative solution of consistent singular linear systems is studied. In these methods, suitable for parallel computations, different blocks, i.e., smaller linear systems, can be solved concurrently by different processors. Each of these smaller systems are solved by an (inner) iterative method. Hypotheses are provided for the convergence of non-stationary methods, i.e., when the number of inner iterations may vary from block to block and from one outer iteration to another. It is shown that the iteration matrix corresponding to one step of the block method is convergent, i.e., that its powers converge to a limit matrix. A theorem on the convergence of the infinite product of matrices with the same eigenspace corresponding to the eigenvalue 1 is proved, and later used as a tool in the convergence analysis of the block method. The methods studied can be used to solve any consistent singular system, including discretizations of certain differential equations. They can also be used to find stationary probability distribution of Markov chains. This last application is considered in detail.     

A single-step HSS method for non-Hermitian positive definite linear systems

In this paper, based on the Hermitian and skew-Hermitian splitting (HSS) iteration method, a single-step HSS (SHSS) iteration method is introduced to solve the non-Hermitian positive definite linear systems. Theoretical analysis shows that, under a loose restriction on the iteration parameter, the SHSS method is convergent to the unique solution of the linear system. Furthermore, we derive an upper bound for the spectral radius of the SHSS iteration matrix, and the quasi-optimal parameter is obtained to minimize the above upper bound. Numerical experiments are reported to the efficiency of the SHSS method; numerical comparisons show that the proposed SHSS method is superior to the HSS method under certain conditions.     

A splitting preconditioner for saddle point problems

For large sparse systems of linear equations iterative techniques are attractive. In this paper, we study a splitting method for an important class of symmetric and indefinite system. Theoretical analyses show that this method converges to the unique solution of the system of linear equations for all t>0 (t is the parameter). Moreover, all the eigenvalues of the iteration matrix are real and nonnegative and the spectral radius of the iteration matrix is decreasing with respect to the parameter t. Besides, a preconditioning strategy based on the splitting of the symmetric and indefinite coefficient matrices is proposed. The eigensolution of the preconditioned matrix is described and an upper bound of the degree of the minimal polynomials for the preconditioned matrix is obtained. Numerical experiments of a model Stokes problem and a least‐squares problem with linear constraints presented to illustrate the effectiveness of the method. Copyright © 2011 John Wiley & Sons, Ltd.     

Randomized Preprocessing of Homogeneous Linear Systems of Equations

Our randomized preprocessing enables pivoting-free and orthogonalization-free solution of homogeneous linear systems of equations. In the case of Toeplitz inputs, we decrease the estimated solution time from quadratic to nearly linear, and our tests show dramatic decrease of the CPU time as well. We prove numerical stability of our approach and extend it to solving nonsingular linear systems, inversion and generalized (Moore-Penrose) inversion of general and structured matrices by means of Newton’s iteration, approximation of a matrix by a nearby matrix that has a smaller rank or a smaller displacement rank, matrix eigen-solving, and root-finding for polynomial and secular equations and for polynomial systems of equations. Some by-products and extensions of our study can be of independent technical intersest, e.g., our extensions of the Sherman-Morrison-Woodbury formula for matrix inversion, our estimates for the condition number of randomized matrix products, and preprocessing via augmentation.