Making sparse matrices sparser: Computational results

Many optimization algorithms involve repeated processing of a fixed set of linear constraints. If we pre-process the constraint matrixA to be sparser, then algebraic operations onA will become faster. We consider the problem of making a given matrix as sparse as possible, theSparsity Problem (SP). In a companion paper with S. Frank Chang, we developed some theoretical algorithms for SP under a non-degeneracy assumption (McCormick and Chang, 1988). Here we investigate what must be done to make those algorithms applicable in practice. We report encouraging computational results in making linear programming constraint matrices sparser. We also find that the Simplex Algorithm can solve the reduced LPs faster. Comparisons are made to a heuristic algorithm for SP of Adler et al. (1989).This work was partially supported by NSF Grants ECS-84-04350 and CDR-84-21402, and by ONR Contract N0014-87-K0214.     

A spectral algorithm for envelope reduction of sparse matrices

The problem of reordering a sparse symmetric matrix to reduce its envelope size is considered. A new spectral algorithm for computing an envelope-reducing reordering is obtained by associating a Laplacian matrix with the given matrix and then sorting the components of a specified eigenvector of the Laplacian. This Laplacian eigenvector solves a continuous relaxation of a discrete problem related to envelope minimization called the minimum 2-sum problem. The permutation vector computed by the spectral algorithm is a closest permutation vector to the specified Laplacian eigenvector. Numerical results show that the new reording algorithm usually computes smaller envelope sizes than those obtained from the current standards such as the Gibbs—Poole—Stockmeyer (GPS) algorithm or the reverse Cuthill—McKee (RCM) algorithm in SPARSPAK, in some cases reducing the envelope by more than a factor of two.     

A bilinear algorithm for sparse representations

We consider the following sparse representation problem: represent a given matrix X∈ℝ ^m×N as a multiplication X=AS of two matrices A∈ℝ ^m×n (m≤n<N) and S∈ℝ ^n×N , under requirements that all m×m submatrices of A are nonsingular, and S is sparse in sense that each column of S has at least n−m+1 zero elements. It is known that under some mild additional assumptions, such representation is unique, up to scaling and permutation of the rows of S. We show that finding A (which is the most difficult part of such representation) can be reduced to a hyperplane clustering problem. We present a bilinear algorithm for such clustering, which is robust to outliers. A computer simulation example is presented showing the robustness of our algorithm.     

A transfer algorithm for hierarchical clustering

Transfer algorithms are usually used to optimize an objective function that is defined on the set of partitions of a finite set X. In this paper we define an equivalence relation ? on the set of fuzzy equivalence relations on X and establish a bijection from the set of hierarchies on X to the set of equivalence classes with respect to ?. Thus, hierarchies can be identified with fuzzy equivalence relations and the transfer algorithm can be modified in order to optimize an objective function that is defined on the set of hierarchies on X.     

A note on quasi-newton formulae for sparse second derivative matrices

In order to apply quasi-Newton methods to solve unconstrained minimization calculations when the number of variables is very large, it is usually necessary to make use of any sparsity in the second derivative matrix of the objective function. Therefore, it is important to extend to the sparse case the updating formulae that occur in variable metric algorithms to revise the estimate of the second derivative matrix. Suitable extensions suggest themselves when the updating formulae are derived by variational methods [1, 3]. The purpose of the present paper is to give a new proof of a theorem of Dennis and Schnabel [1], that shows the effect of sparsity on updating formulae for second derivative estimates.     

Solving the eigenvalue problem for sparse matrices

A modification of the Danilewski method is presented, permitting the solution of the eigenvalue problem for a constant sparse matrix of large order to be reduced to the solution of the same problem for a polynomial matrix of lower order. Certain solution algorithms are proposed for a partial eigenvalue problem for the polynomial matrix. Questions of the realization of the algorithms on a model PRORAB computer are examined.Translated from Zapiski Nauchnykh Seminarov Leningradskogo Otdeleniya Matematicheskogo Instituta im. V. A. Steklova AN SSSR, Vol. 58, pp. 92–110, 1976.     

Determinantal formulae for matrices with sparse inverses

The determinant of a matrix is expressed in terms of certain of its principal minors by a formula which can be “read off” from the graph of the inverse of the matrix. The only information used is the zero pattern of the inverse, and each zero pattern yields one or more corresponding formulae for the determinant.     

A sparse nonlinear optimization algorithm

One of the most effective numerical techniques for solving nonlinear programming problems is the sequential quadratic programming approach. Many large nonlinear programming problems arise naturally in data fitting and when discretization techniques are applied to systems described by ordinary or partial differential equations. Problems of this type are characterized by matrices which are large and sparse. This paper describes a nonlinear programming algorithm which exploits the matrix sparsity produced by these applications. Numerical experience is reported for a collection of trajectory optimization problems with nonlinear equality and inequality constraints.The authors wish to acknowledge the insightful contributions of Dr. William Huffman.     

A linearly convergent algorithm for sparse signal reconstruction

For the sparse signal reconstruction problem in compressive sensing, we propose a projection-type algorithm without any backtracking line search based on a new formulation of the problem. Under suitable conditions, global convergence and its linear convergence of the designed algorithm are established. The efficiency of the algorithm is illustrated through some numerical experiments on some sparse signal reconstruction problem.     

A conjugate gradient algorithm for sparse linear inequalities

This paper presents a conjugate gradient method for solving systems of linear inequalities. The method is of dual optimization type and consists of two phases which can be implemented in a common framework. Phase 1 either finds the minimum-norm solution of the system or detects the inconsistency of the system. In the latter event, the method proceeds to Phase 2 in which an approximate least-squares solution to the system is obtained. The method is particularly suitable to large scale problems because it preserves the sparsity structure of the problem. Its efficiency is shown by computational comparisons with an SOR type method.     

A perfect matching algorithm for sparse bipartite graphs

This paper considers an algorithm for finding a perfect matching, if there is one, in a bipartite graph G. It is shown that the search for a perfect matching in G may be carried on separately in the strongly connected components of appropriate directed graphs. The algorithm may be particularly useful for block triangularization of very large, sparse, nonsingular matrices.     

Spectral methods with sparse matrices

Summary Spectral methods employ global polynomials for approximation. Hence they give very accurate approximations for smooth solutions. Unfortunately, for Dirichlet problems the matrices involved are dense and have condition numbers growing asO(N ⁴) for polynomials of degree N in each variable. We propose a new spectral method for the Helmholtz equation with a symmetric and sparse matrix whose condition number grows only asO(N ²). Certain algebraic spectral multigrid methods can be efficiently used for solving the resulting system. Numerical results are presented which show that we have probably found the most effective solver for spectral systems.     

The preconditioned inverse iteration for hierarchical matrices

The preconditioned inverse iteration is an efficient method to compute the smallest eigenpair of a symmetric positive definite matrix M. Here we use this method to find the smallest eigenvalues of a hierarchical matrix. The storage complexity of the data‐sparse ‐matrices is almost linear. We use ‐arithmetic to precondition with an approximate inverse of M or an approximate Cholesky decomposition of M. In general, ‐arithmetic is of linear‐polylogarithmic complexity, so the computation of one eigenvalue is cheap. We extend the ideas to the computation of inner eigenvalues by computing an invariant subspace S of (M ? μI)² by subspace preconditioned inverse iteration. The eigenvalues of the generalized matrix Rayleigh quotient μ_M(S) are the desired inner eigenvalues of M. The idea of using (M ? μI)² instead of M is known as the folded spectrum method. As we rely on the positive definiteness of the shifted matrix, we cannot simply apply shifted inverse iteration therefor. Numerical results substantiate the convergence properties and show that the computation of the eigenvalues is superior to existing algorithms for non‐sparse matrices.Copyright © 2012 John Wiley & Sons, Ltd.     

A sensitivity analysis algorithm for hierarchical decision models

In this paper, a comprehensive algorithm is developed to analyze the sensitivity of hierarchical decision models (HDM), including the analytic hierarchy process and its variants, to single and multiple changes in the local contribution matrices at any level of the decision hierarchy. The algorithm is applicable to all HDM that use an additive function to derive the overall contribution vector. It is independent of pairwise comparison scales, judgment quantification techniques and group opinion combining methods. The allowable range/region of perturbations, contribution tolerance, operating point sensitivity coefficient, total sensitivity coefficient and the most critical decision element at a certain level are identified in the HDM SA algorithm. An example is given to demonstrate the application of the algorithm and show that HDM SA can reveal information more significant and useful than simply knowing the rank order of the decision alternatives.     

A sparse sequential quadratic programming algorithm

Described here is the structure and theory for a sequential quadratic programming algorithm for solving sparse nonlinear optimization problems. Also provided are the details of a computer implementation of the algorithm along with test results. The algorithm maintains a sparse approximation to the Cholesky factor of the Hessian of the Lagrangian. The solution to the quadratic program generated at each step is obtained by solving a dual quadratic program using a projected conjugate gradient algorithm. An updating procedure is employed that does not destroy sparsity.     

A facial reduction algorithm for finding sparse SOS representations

The facial reduction algorithm reduces the size of the positive semidefinite cone in SDP. The elimination method for a sparse SOS polynomial [M. Kojima, S. Kim, H. Waki, Sparsity in sums of squares of polynomials, Math. Program. 103 (2005) 45-62] removes monomials which do not appear in any SOS representations. In this paper, we establish a relationship between a facial reduction algorithm and the elimination method for a sparse SOS polynomial.     

A recursive algorithm for constructing complicated Dixon matrices

Whether the determinant of the Dixon matrix equals zero or not is used for determining if a system of n + 1 polynomial equations in n variables has a common root, and is a very efficient quantifier elimination approach too. But for a complicated polynomial system, it is not easy to construct the Dixon matrix. In this paper, a recursive algorithm to construct the Dixon matrix is proposed by which some problems that cannot be tackled by other methods can be solved on the same computer platform. A dynamic programming algorithm based on the recursive formula is developed and compared for speed and efficiency to the recursive algorithm.     

A look-ahead Bareiss algorithm for general Toeplitz matrices

Summary. The Bareiss algorithm is one of the classical fast solvers for systems of linear equations with Toeplitz coefficient matrices. The method takes advantage of the special structure, and it computes the solution of a Toeplitz system of order~ with only~ arithmetic operations, instead of~ operations. However, the original Bareiss algorithm requires that all leading principal submatrices be nonsingular, and the algorithm is numerically unstable if singular or ill-conditioned submatrices occur. In this paper, an extension of the Bareiss algorithm to general Toeplitz systems is presented. Using look-ahead techniques, the proposed algorithm can skip over arbitrary blocks of singular or ill-conditioned submatrices, and at the same time, it still fully exploits the Toeplitz structure. Implementation details and operations counts are given, and numerical experiments are reported. We also discuss special versions of the proposed look-ahead Bareiss algorithm for Hermitian indefinite Toeplitz systems and banded Toeplitz systems. Received August 27, 1993 / Revised version received March 1994