Showing papers in "SIAM Journal on Matrix Analysis and Applications in 1991"

Rational iterative methods for the matrix sign function

Abstract: In this paper an analysis of rational iterations for the matrix sign function is presented. This analysis is based on Pade approximations of a certain hypergeometric function and it is shown that l...

A black box generalized conjugate gradient solver with inner iterations and variable-step preconditioning

Abstract: The generalized conjugate gradient method proposed by Axelsson is studied in the case when a variable-step preconditioning is used. This can be the case when the preconditioned system is solved app...

The restricted singular value decomposition: properties and applications

Abstract: The restricted singular value decomposition (RSVD) is the factorization of a given matrix, relative to two other given matrices. It can be interpreted as the ordinary singular value decomposition with different inner products in row and column spaces. Its properties and structure are investigated in detail as well as its connection to generalized eigenvalue problems, canonical correlation analysis and other generalizations of the singular value decomposition. Applications that are discussed include the analysis of the extended shorted operator, unitarily invariant norm minimization with rank constraints, rank minimization in matrix balls, the analysis and solution of linear matrix equations, rank minimization of a partitioned matrix and the connection with generalized Schur complements, constrained linear and total linear least squares problems, with mixed exact and noisy data, including a generalized Gauss-Markov estimation scheme. Two constructive proofs of the RSVD in terms of other generalizations of the ordinary singular value decomposition are provided as well.

The restricted singular value decomposition of matrix triplets

Abstract: In this paper the concept of restricted singular values of matrix triplets is introduced. A decomposition theorem concerning the general matrix triplet $( A,B,C )$, where $A \in \mathcal{C}^{m \times n} ,B \in \mathcal{C}^{m \times p} $, and $C \in \mathcal{C}^{q \times n} $, which is called the restricted singular value decomposition (RSVD), is proposed. This result generalizes the well-known singular value decomposition, the generalized singular value decomposition, and the recently proposed product-induced singular value decomposition. Connection of restricted singular values with the problem of determination of matrix rank under restricted perturbation is also discussed.

Exponential nonnegativity on the ice cream cone

Abstract: Let $K_n $ denote the n-dimensional ice cream cone. This paper investigates the structure of those matrices A such that $e^{tA} K_n \subset K_n $ for all $t\geqq 0$. The characterizations extend to general ellipsoidal cones.

The restricted total least squares problem: formulation, algorithm, and properties

Abstract: The Total Least Squares (TLS) method has been devised as a more global fitting technique than the ordinary least squares technique for solving overdetermined sets of linear equations AX ≈ B when errors occur in all data. If the errors on the measurements A and B are uncorrelated with zero mean and equal variance, TLS is able to compute a strongly consistent estimate of the true solution of the corresponding unperturbed set A 0 X = B 0. In this paper the TLS computations are generalized in order to maintain consistency of the solution in the following cases: first of all, some columns of A may be error-free and secondly, the errors on the remaining data may be correlated and not equally sized. Hereto, a numerically reliable Generalized TLS algorithm GTLS, based on the Generalized Singular Value Decomposition (GSVD), is developed. Additionally, the equivalence between the GTLS solution and alternative expressions of consistent estimators, described in literature, is proven. These relations allow to deduce the main statistical properties of the GTLS solution.

Updating and downdating of orthogonal polynomials with data fitting applications

Abstract: We derive and assess new methods for updating and downdating least squares polynomial fits to discrete data using polynomials orthogonal on all the data points being used. Rather than fixing on one basis throughout, the methods adaptively update and downdate both the least squares fit and the polynomial basis. This is achieved by performing similarity transformations on the tridiagonal Jacobi matrices representing the basis. Although downdating is potentially unstable, experimental results show that the methods give satisfactory results for low degree fits. The most economical of the algorithms implementing the methods needs 14n+ O(1) flops and 2n square roots to update a fit of order n.

Nested epsilon decompositions of linear systems: weakly coupled and overlapping blocks

Abstract: A graph-theoretic algorithm is proposed for decomposition of a linear system of equations into subsystems having a prescribed size of mutual interactions The algorithm generates a whole range of nested decompositions with an apparent trade-off between levels of coupling and sizes of the subsystems Both disjoint and overlapping subsystems are considered Having a linear time complexity for a selected strength of coupling, the algorithm is suitable for conditioning large systems and achieving fast convergence rates in block-iterative computations via parallel multiprocessor schemes

Stable solvers and block elimination for bordered systems

Abstract: Linear systems with a fairly well conditioned matrix M of the form \[ \begin{gathered} \begin{pmatrix} A & b \\ c & d \end{pmatrix} \begin{matrix} n \\ 1 \end{matrix}, \\ \begin{matrix} n & 1 \end{matrix} \end{gathered} \] for which a “black-box” solver for A is available, are considered. To solve systems with M, a mixed block elimination algorithm, called BEM, is proposed. It has the following advantages: (1) It is easier to understand and to program than the widely accepted deflated block elimination (DBE) proposed by Chan, yet allows the same broad class of solvers and has comparable accuracy. (2) It requires one less solve with A. (3) It allows a rigorous error analysis that shows why it may fail in exceptional cases (all other black-box methods known to us also fail in these cases).BEM is also compared to iterative refinement of Crout block elimination (BEC) introduced by Pryce and Govaerts. BEC allows a more restricted class of solvers than BEM but is faster in cases where a solver is given not for ...

Chasing algorithms for the eigenvalues problem

Abstract: A generic chasing algorithm for the matrix eigenvalue problem is introduced and studied. This algorithm includes, as special cases, the implicit, multiple-step $QR$ and $LR$ algorithms and similar bulge-chasing algorithms for the standard eigenvalue problem. The scope of the generic chasing algorithm is quite broad; it encompasses a number of chasing algorithms that cannot be analyzed by the traditional (e.g., implicit Q theorem) approach. These include the $LR$ algorithm with partial pivoting and other chasing algorithms that employ pivoting for stability, as well as hybrid algorithms that combine elements of the $LR$ and $QR$ algorithms. The main result is that each step of the generic chasing algorithm amounts to one step of the generic $GR$ algorithm. Therefore the convergence theorems for $GR$ algorithms that were proven in a previous work [D. S. Watkins and L. Elsner, Linear Algebra Appl., 143 (1991), pp. 19–47] also apply to the generic chasing algorithm.

A variant of the Gohberg-Semencul formula involving circulant matrices

Abstract: The Gohberg–Semencul formula expresses the inverse of a Toeplitz matrix as the difference of products of lower triangular and upper triangular Toeplitz matrices. In this paper the idea of cyclic displacement structure is used to show that the upper triangular matrices in this formula can be replaced by circulant matrices. The resulting computational savings afforded by this modified formula is discussed.

On random correlation matrices

Abstract: This report contains a detailed study of random correlation matrices, including algebraic, statistical, and historical background. Such matrices are of particular interest because they serve to model “average signals” for simulation testing of signal processing algorithms. The statistical behavior of spectral functions of the two major types of random correlation matrices is extensively discussed in the latter half, from both theoretical and empirical aspects. The emphasis is on eigenvalue distribution and condition number behavior. Actual application to algorithm testing will be described in a subsequent report.

Observability of linear time-varying descriptor systems

Abstract: A characterization of observability for linear time-varying descriptor systems $E( t )x' ( t ) + F ( t ) x ( t ) = B (t) u (t), y ( t ) = C ( t ) x ( t )$, is given. E is not required to have constant rank. The characterization is designed to reduce symbolic computation and has potential advantages even when E is nonsingular. It is also shown that all observable analytic descriptor systems are smoothly observable even if they are not uniformly observable. Finally, the external behavior of time-varying descriptor systems is characterized.

On growth in Gaussian elimination with complete pivoting

Abstract: It has been conjectured that when Gaussian elimination with complete pivoting is applied to a real n-by-n matrix, the maximum possible growth is n. In this note, a 13-by-13 matrix is given, for which the growth is 13.0205. The matrix was constructed by solving a large nonlinear programming problem. Growth larger than n has also been observed for matrices of orders 14, 15, and 16.

The equation AXB+CYD = E over a principal ideal domain

Abstract: This paper considers the equation $AXB + CYD = E$ in matrices over a principal ideal domain. Under the assumption that $[ A:C ]$ is left invertible and $[ B':D' ]^\prime $ is right invertible in the domain, it is shown that this equation is solvable if and only if both $AX + YD = E$ and $XB + CY = E$ are solvable. The set of all solutions to the equation are shown to be in bijective correspondence with the set of solutions to the latter two equations modulo an equivalence relation.

Fast QR decomposition of Vandermonde-like matrices and polynomial least squares approximation

Abstract: Let f and g be functions defined at the real and distinct nodes $x_k $, and consider the inner product $( f,g ): = \sum_{k = 1}^m f ( x_k ) g ( x_k ) w_k^2 $ with positive weights $w_k^2 $. The present paper discusses the computation of orthonormal polynomials $\pi _0 ,\pi _1 , \cdots ,\pi _{n - 1} ,n\leqq m$, with respect to this inner product, and the use of these polynomials in a fast scheme for computing a QR decomposition of the transpose of Vandermonde-like matrices. Two methods are compared for computing the recurrence coefficients for the polynomials $\pi _j $ and their values at the nodes $x_k $: the Stieltjes procedure and a method in which an inverse eigenvalue problem for a tridiagonal symmetric matrix is solved by an algorithm proposed by Rutishauser, Gragg, and Harrod. The latter method is found to generally yield higher accuracy than the Stieltjes procedure if n is close to m, and roughly the same accuracy otherwise. This method for solving an inverse eigenvalue problem is applied in an alg...

An algorithm for Ax = lBx with symmetric and positive-definite A and B

Abstract: An algorithm is given for computing the solution of the eigenvalues of $Ax = \lambda Bx$ with symmetric and positive-definite A and B. It reduces $Ax = \lambda Bx$ to the generalized singular value problem $LL^T x = \lambda ( L_B L_B^T )x$ by the Cholesky decompositions $A = LL^T $ and $B = L_B L_B^T $, and then reduces the generalized singular value decomposition of $L^T $ and $L_B^T $ to the CS decomposition of Q by the QR decomposition $( L,L_B )^T = QR$. Finally, it reduces A and B to diagonal forms by singular value decompositions. The algorithm provided is stable and, what is more, faster than the QZ algorithm. Numerical examples are also presented.

Reducing the computations of the singular value decomposition array given by Brent and Luk

Abstract: A new, efficient, two-plane rotation (TPR) method for computing two-sided rotations involved in singular value decomposition (SVD) is presented. It is shown that a two-sided rotation can be evaluated by only two plane rotations and a few additions. This leads to significantly reduced computations. Moreover, if coordinate rotation digital computer(CORDIC) processors are used for realizing the processing elements (PEs) of the SVD array given by Brent and Luk, the computational overhead of the diagonal PEs due to angle calculations can be avoided. The resultingSVD array has a homogeneous structure with identical diagonal and off-diagonal PEs. Similar results can also be obtained if the TPR method is applied to Luk's triangular SVD array and to Stewart's Schur decomposition array. In this paper, we develop a two-plane rotation (TPR) method for computing two- sided rotations. We show that the above computational complexity can be reduced sig- nificantly because each two-sided rotation can be evaluated by only two plane rotations and a few additions. Moreover, the SVD array given by Brent and Luk becomes ho- mogeneous with identical diagonal and off-diagonal PEs when CORDIC processors are

Tridiagonal approach to the algebraic environment of Toeplitz matrices, part I: basic results

Abstract: This paper contains a thorough investigation of a family of symmetric "predictor polynomials" associated with a nonnegative-definite Toeplitz matrix. These polynomials are constructed from the classical predictors and from the values assumed by some dual predictors in a fixed point of unit modulus; the appropriate duality is induced by changing the sequence of reflection coefficients into its conjugate mirror image, within a unit modulus factor. The central theme of the paper is a well-defined three-term recurrence relation satisfied by these symmetric polynomials; it motivates the "tridiagonal" terminology. The properties of the recurrence are studied in detail; special attention is paid to the important issue of computing the recurrence coefficients from the reflection coefficients. It is shown how this three-term recurrence formula produces an efficient solution method, called the split Levinson algorithm, for the linear prediction problem.

Minimax polynomial preconditioning for Hermitian linear systems

Abstract: This paper explores the use of polynomial preconditioning for Hermitian positive definite and indefinite linear systems $Ax = b$. Unlike preconditioners based on incomplete factorizations, polynomial preconditioners are easy to employ and well suited to vector and/or parallel architectures. It is shown that any polynomial iterative method may be used to define a preconditioning polynomial, and that the optimum polynomial preconditioner is obtained from a minimax approximation problem. A variety of preconditioning polynomials are then surveyed, including the Chebyshev, de Boor and Rice, Grcar, and bilevel polynomials. Adaptive procedures for each of these polynomials are also discussed, and it is shown that the new bilevel polynomial is particularly well suited for use in adaptive CG algorithms.

Tridiagonal approach to the algebraic environment of Toeplitz matrices, part II.: zero and eigenvalue problems

Abstract: The central subject of this paper is the three-term recurrence formula satisfied by the symmetric (first-kind) and antisymmetric (second-kind) polynomials relative to a given sequence of reflection coefficients, the last element of which has unit modulus. The theory of these polynomials is shown to have interesting analogies with the classical theory of orthogonal polynomials on the real line. In the case of real data, the former is equivalent to a special case of the latter (by a change of variable). The main application considered here is the problem of computing the zeros of the highest-degree symmetric polynomial, which is identified as a predictor polynomial. This problem occurs not only in some modelling techniques for digital signal processing, but can also be interpreted as the eigenvalue problem for a unitary Hessenberg matrix. Attractive solution methods are derived from the "tridiagonal approach," based on the three-term recurrence relation.

Divide-and-conquer solutions of least-squares problems for matrices with displacement structure

Abstract: A divide-and-conquer implementation of a generalized Schur algorithm enables (exact and) least-squares solutions of various block-Toeplitz or Toeplitz-block systems of equations with $O ( \alpha ^3 n\log ^2 n )$ operations to be obtained, where the displacement rank $\alpha $ is a small constant (typically between two to four for scalar near-Toeplitz matrices) independent of the size of the matrices.

Reduction of a general matrix to tridiagonal form

Abstract: An algorithm for reducing a nonsymmetric matrix to tridiagonal form as a first step toward finding its eigenvalues is described. The algorithm uses a variation of threshold pivoting, where at each step, the pivot is chosen to minimize the maximum entry in the transformation matrix that reduces the next column and row of the matrix. Situations are given where the tridiagonalization process breaks down, and two recovery methods are presented for these situations. Although no existing tridiagonalization algorithm is guaranteed to succeed, this algorithm is found to be very robust and fast in practice. A gradual loss of similarity is also observed as the order of the matrix increases.

Least squares approximation by real normal matrices with specified spectrum

Abstract: The problem of best approximating a given real matrix in the Frobenius norm by real, normal matrices subject to a prescribed spectrum is considered. The approach is based on using the projected gradient method. The projected gradient of the objective function on the manifold of constraints can be formulated explicitly. This gives rise to a descent flow that can be followed numerically. The explicit form also facilitates the computation of the second-order optimality condition from which some interesting properties of the stationary points are related to the well-known Wielandt–Hoffman theorem.

Linear operators preserving certain equivalence relations on matrices

Abstract: Using a uniform approach, characterizations are obtained of linear operators on matrix spaces that preserve certain equivalence relations such as consimilarity, $ * $-congruence, nonsingular equivalence, and unitary equivalence.

A class of arbitrarily ill conditioned floating-point matrices

Abstract: Let $\mathbb{F}$ be a floating-point number system with basis $\beta \geqq 2$ and an exponent range consisting of at least the exponents 1 and 2. A class of arbitrarily ill conditioned matrices is described, the coefficients of which are elements of $\mathbb{F}$. Due to the very rapidly increasing sensitivity of those matrices, they might be regarded as “almost” ill posed problems.The condition of those matrices and their sensitivity with respect to inversion is given by means of a closed formula. The condition is rapidly increasing with the dimension. For example, in the double precision of the IEEE 754 floating-point standard (base 2, 53 bits in the mantissa including implicit 1), matrices with $2n$ rows and columns are given with a condition number of approximately $4\cdot10^{32n} $.

Monotonicity properties of the Toda flow, the QR-Flow, and subspace iteration

Abstract: Let ${\bf X}( t )$ denote the Toda flow on the space of $n \times n$ matrices, with ${\bf X}( 0 )$ a symmetric matrix, and let ${\bf X}_r ( t )$ denote the $r \times r$ upper left corner principal submatrix of ${\bf X}( t )$, i.e., ${\bf X}_r ( t ) = {\bf E}_r^T {\bf X}( t ){\bf E}_r $ where ${\bf E}_r = \begin{bmatrix} {I_r } \\ 0 \end{bmatrix}$. Then the r ordered eigenvalues $\lambda _1 ( {\bf X}_r ( t ) ) \geqq \lambda _2 ( {\bf X}_r ( t ) )\geqq \cdots \geqq \lambda_r ( {\bf X}_r ( t ) )$ of ${\bf X}_r ( t )$ are each a nondecreasing function of t, for $1\leqq r\leqq n$. A similar result is proved for the QR-flow ${\bf Y} ( t ) = \exp ( {\bf X} ( t ) )$, for the eigenvalues of ${\bf Y}_r ( t ) = {\bf E}_r^T {\bf Y} ( t ){\bf E}_r $. For any generalized Toda flow $f( {\bf X} ( t ) )$ with $f ( \cdot )$ a nondecreasing function, it is shown that $\operatorname{Tr} ( {\bf E}_r^T f ( {\bf X} ( t ) ){\bf E}_r )$ is a nondecreasing function of t. The QR-flow inequalities are used to show that the Ritz valu...

Computation of the Euler angles of a symmetric 3X3 matrix

Abstract: Closed form formulas for computing the eigenvectors of a symmetric $3 \times 3$ matrix are presented. The matrix of the eigenvectors is computed as a product of three rotations through Euler angles. The formulas require approximately 90 arithmetic operations, six trigonometric evaluations, and two root evaluations. These formulas may be applied as a subroutine in a parallel one-sided Jacobi-type method in which three rather than two columns, as is the case in the standard Jacobi method, are operated on in each step.

On the Smith normal form of structured polynomial matrices

Abstract: The Smith normal form of a polynomial matrix $D( s ) = Q( s ) + T( s )$ is investigated, where $D( s )$ is “structured” in the sense that (i) the coefficients of the entries of $Q( s )$ belong to a field ${\bf K}$, (ii) the nonzero coefficients of the entries of $T( s )$ are algebraically independent over ${\bf K}$, and (iii) every minor of $Q( s )$ is a monomial in s. Such matrices have been useful in the structural approach in control theory. It is shown that all the invariant polynomials except for the last are monomials in s and the last invariant polynomial is expressed in terms of the combinatorial canonical form (CCF) of a layered mixed matrix associated with $D( s )$. On the basis of this, the Smith form of $D( s )$ can be computed by means of an efficient (polynomial-time) matroid-theoretic algorithm that involves arithmetic operations in the base field ${\bf K}$ only.

On the identification of local minimizers in inertia-controlling methods for quadratic programming

Abstract: The verification of a local minimizes of a general (i.e., nonconvex) quadratic program is in general an NP-hard problem. The difficulty concerns the optimality of certain points (which we call dead points) at which the first-order necessary conditions for optimality are satisfied, but strict complementarity does not hold. Inertia-controlling quadratic programming (ICQP) methods form an important class of methods for solving general quadratic programs. We derive a computational scheme for proceeding at a dead point that is appropriate for a general ICQP method.