SIAM Digital Library
 
 
 

SIAM J. on Matrix Analysis and Applications

Year Range: 

All Journal content published prior to 1997 is part of LOCUS.

Search Issue | RSS Feeds RSS
Next Issue

2012

Volume 33, Issue 1, pp. 1-289


Stable Computation of the CS Decomposition: Simultaneous Bidiagonalization

Brian D. Sutton

SIAM. J. Matrix Anal. & Appl. 33, pp. 1-21 (21 pages)

Online Publication Date: January 05, 2012

Full Text: | Download PDF

Show Abstract
Since its discovery in 1977, the CS decomposition (CSD) has resisted computation, even though it is a sibling of the well-understood eigenvalue and singular value decompositions. Several algorithms have been developed for the reduced 2-by-1 form of the decomposition, but none have been extended to the complete 2-by-2 form of the decomposition in Stewart's original paper. In this article, we present an algorithm for simultaneously bidiagonalizing the four blocks of a unitary matrix partitioned into a 2-by-2 block structure. This serves as the first, direct phase of a two-stage algorithm for the CSD, much as Golub–Kahan–Reinsch bidiagonalization serves as the first stage in computing the singular value decomposition. Backward stability is proved.

dqds with Aggressive Early Deflation

Yuji Nakatsukasa, Kensuke Aishima, and Ichitaro Yamazaki

SIAM. J. Matrix Anal. & Appl. 33, pp. 22-51 (30 pages)

Online Publication Date: January 05, 2012

Full Text: | Download PDF

Show Abstract
The dqds algorithm computes all the singular values of an $n \times n$ bidiagonal matrix to high relative accuracy in $O(n^2)$ cost. Its efficient implementation is now available as a LAPACK subroutine and is the preferred algorithm for this purpose. In this paper we incorporate into dqds a technique called aggressive early deflation, which has been applied successfully to the Hessenberg QR algorithm. Extensive numerical experiments show that aggressive early deflation often reduces the dqds runtime significantly. In addition, our theoretical analysis suggests that with aggressive early deflation, the performance of dqds is largely independent of the shift strategy. We confirm through experiments that the zero-shift version is often as fast as the shifted version. We give a detailed error analysis to prove that with our proposed deflation strategy, dqds computes all the singular values to high relative accuracy.

Difference Filter Preconditioning for Large Covariance Matrices

Michael L. Stein, Jie Chen, and Mihai Anitescu

SIAM. J. Matrix Anal. & Appl. 33, pp. 52-72 (21 pages)

Online Publication Date: January 05, 2012

Full Text: | Download PDF

Show Abstract
In many statistical applications one must solve linear systems involving large, dense, and possibly irregularly structured covariance matrices. These matrices are often ill-conditioned; for example, the condition number increases at least linearly with respect to the size of the matrix when observations of a random process are obtained from a fixed domain. This paper discusses a preconditioning technique based on a differencing approach such that the preconditioned covariance matrix has a bounded condition number independent of the size of the matrix for some important process classes. When used in large scale simulations of random processes, significant improvement is observed for solving these linear systems with an iterative method.

Compact Fourier Analysis for Multigrid Methods based on Block Symbols

Thomas K. Huckle and Christos Kravvaritis

SIAM. J. Matrix Anal. & Appl. 33, pp. 73-96 (24 pages)

Online Publication Date: January 05, 2012

Full Text: | Download PDF

Show Abstract
The notion of compact Fourier analysis (CFA) is discussed. CFA allows description of multigrid (MG) in a nutshell and offers a clear overview on all MG components. The principal idea of CFA is to model the MG mechanisms by means of scalar generating functions and matrix functions (block symbols). The formalism of the CFA approach is presented by describing the symbols of the fine and coarse grid problems, the prolongation and restriction, the smoother, and the coarse grid correction, resp., smoothing corrections. CFA uses matrix functions and their features (e.g., product, inverse, adjugate, norm, spectral radius, eigenvectors, eigenvalues of multilevel $\omega$-circulant matrices), and scalar functions and their roots. This leads to an elementary description and allows for an easy analysis of MG algorithms. A first application is to utilize CFA for deriving MG as a direct solver, i.e., an MG cycle that will converge in just one iteration step. Necessary and sufficient conditions that have to be fulfilled by the MG components are given for obtaining MG functioning as a direct solver. Furthermore, new general and practical smoothers and transfer operators that lead to efficient MG methods are introduced. In addition, we study sparse approximations of the Galerkin coarse grid operator yielding efficient and practicable MG algorithms (approximately direct solvers). Numerical experiments demonstrate the theoretical results.

On Iterative Solution for Linear Complementarity Problem with an $H_{+}$-Matrix

A. Hadjidimos, M. Lapidakis, and M. Tzoumas

SIAM. J. Matrix Anal. & Appl. 33, pp. 97-110 (14 pages)

Online Publication Date: January 10, 2012

Full Text: | Download PDF

Show Abstract
The numerous applications of the linear complementarity problem (LCP) in, e.g., the solution of linear and convex quadratic programming, free boundary value problems of fluid mechanics, and moving boundary value problems of economics make its efficient numerical solution a very imperative and interesting area of research. For the solution of the LCP, many iterative methods have been proposed, especially, when the matrix of the problem is a real positive definite or an $H_{+}$-matrix. In this work we assume that the real matrix of the LCP is an $H_{+}$-matrix and solve it by using a new method, the scaled extrapolated block modulus algorithm, as well as an improved version of the very recently introduced modulus-based matrix splitting modified AOR iteration method. As is shown by numerical examples, the two new methods are very effective and competitive with each other. (A corrected PDF is attached to this article.)

Uni-mode and Partial Uniqueness Conditions for CANDECOMP/PARAFAC of Three-Way Arrays with Linearly Dependent Loadings

Xijing Guo, Sebastian Miron, David Brie, and Alwin Stegeman

SIAM. J. Matrix Anal. & Appl. 33, pp. 111-129 (19 pages)

Online Publication Date: January 13, 2012

Full Text: | Download PDF

Show Abstract
In this paper, three sufficient conditions are derived for the three-way CANDECOMP/PARAFAC (CP) model, which ensure uniqueness in one of the three modes (“uni-mode-uniqueness”). Based on these conditions, a partial uniqueness condition is proposed which allows collinear loadings in only one mode. We prove that if there is uniqueness in one mode, then the initial CP model can be uniquely decomposed in a sum of lower-rank tensors for which identifiability can be independently assessed. This condition is simpler and easier to check than other similar conditions existing in the specialized literature. These theoretical results are illustrated by numerical examples.

Verified Bounds for Least Squares Problems and Underdetermined Linear Systems

Siegfried M. Rump

SIAM. J. Matrix Anal. & Appl. 33, pp. 130-148 (19 pages)

Online Publication Date: January 13, 2012

Full Text: | Download PDF

Show Abstract
New algorithms are presented for computing verified error bounds for least squares problems and underdetermined linear systems. In contrast to previous approaches the new methods do not rely on normal equations and are applicable to sparse matrices. Computational results demonstrate that the new methods are faster than existing ones.

Block Tensor Unfoldings

Stefan Ragnarsson and Charles F. Van Loan

SIAM. J. Matrix Anal. & Appl. 33, pp. 149-169 (21 pages)

Online Publication Date: January 13, 2012

Full Text: | Download PDF

Show Abstract
Within the field of numerical multilinear algebra, block tensors are increasingly important. Accordingly, it is appropriate to develop an infrastructure that supports reasoning about block tensor computation. In this paper we establish concise notation that is suitable for the analysis and development of block tensor algorithms, prove several useful block tensor identities, and make precise the notion of a block tensor unfolding.

Alternating-directional Doubling Algorithm for $M$-Matrix Algebraic Riccati Equations

Wei-guo Wang, Wei-chao Wang, and Ren-Cang Li

SIAM. J. Matrix Anal. & Appl. 33, pp. 170-194 (25 pages)

Online Publication Date: March 08, 2012

Full Text: | Download PDF

Show Abstract
A new doubling algorithm—the alternating-directional doubling algorithm (ADDA)—is developed for computing the unique minimal nonnegative solution of an $M$-matrix algebraic Riccati equation (MARE). It is argued by both theoretical analysis and numerical experiments that ADDA is always faster than two existing doubling algorithms: SDA of Guo, Lin, and Xu (Numer. Math., 103 (2006), pp. 393–412) and SDA-ss of Bini, Meini, and Poloni (Numer. Math., 116 (2010), pp. 553–578) for the same purpose. Also demonstrated is that all three methods are capable of delivering minimal nonnegative solutions with entrywise relative accuracies as warranted by the defining coefficient matrices of a MARE. The three doubling algorithms, differing only in their initial setups, correspond to three special cases of the general bilinear (also called Möbius) transformation. It is explained that ADDA is the best among all possible doubling algorithms resulted from all bilinear transformations.

A Posteriori Error Analysis of Parameterized Linear Systems Using Spectral Methods

T. Butler, P. Constantine, and T. Wildey

SIAM. J. Matrix Anal. & Appl. 33, pp. 195-209 (15 pages)

Online Publication Date: March 13, 2012

Full Text: | Download PDF

Show Abstract
We develop computable a posteriori error estimates for the pointwise evaluation of linear functionals of a solution to a parameterized linear system of equations. These error estimates are based on a variational analysis applied to polynomial spectral methods for forward and adjoint problems. We also use this error estimate to define an improved linear functional and we prove that this improved functional converges at a much faster rate than the original linear functional given a pointwise convergence assumption on the forward and adjoint solutions. The advantage of this method is that we are able to use low order spectral representations for the forward and adjoint systems to cheaply produce linear functionals with the accuracy of a higher order spectral representation. The method presented in this paper also applies to the case where only the convergence of the spectral approximation to the adjoint solution is guaranteed. We present numerical examples showing that the error in this improved functional is often orders of magnitude smaller. We also demonstrate that in higher dimensions, the computational cost required to achieve a given accuracy is much lower using the improved linear functional.

Riemannian Optimization on Tensor Products of Grassmann Manifolds: Applications to Generalized Rayleigh-Quotients

O. Curtef, G. Dirr, and U. Helmke

SIAM. J. Matrix Anal. & Appl. 33, pp. 210-234 (25 pages)

Online Publication Date: March 13, 2012

Full Text: | Download PDF

Show Abstract
We introduce a generalized Rayleigh-quotient $\rho_A$ on the direct product of Grassmannians $\mathrm{Gr}({\bf m},{\bf n})$ enabling a unified approach to well-known optimization tasks from different areas of numerical linear algebra, such as best low-rank approximations of tensors (data compression), geometric measures of entanglement (quantum computing), and subspace clustering (image processing). We compute the Riemannian gradient of $\rho_A$, characterize its critical points, and prove that they are generically nondegenerated. Moreover, we derive an explicit necessary condition for the nondegeneracy of the Hessian. Finally, we present two intrinsic methods for optimizing $\rho_A$—a Newton-like and a conjugated gradient—and compare our algorithms tailored to the above-mentioned applications with established ones from the literature.

On a Nonlinear Matrix Equation Arising in Nano Research

Chun-Hua Guo, Yueh-Cheng Kuo, and Wen-Wei Lin

SIAM. J. Matrix Anal. & Appl. 33, pp. 235-262 (28 pages)

Online Publication Date: March 15, 2012

Full Text: | Download PDF

Show Abstract
The matrix equation $X+A^{\top}X^{-1}A=Q$ arises in Green's function calculations in nano research, where $A$ is a real square matrix and $Q$ is a real symmetric matrix dependent on a parameter and is usually indefinite. In practice one is mainly interested in those values of the parameter for which the matrix equation has no stabilizing solutions. The solution of interest in this case is a special weakly stabilizing complex symmetric solution $X_*$, which is the limit of the unique stabilizing solution $X_{\eta}$ of the perturbed equation $X+A^{\top}X^{-1}A=Q+i\eta I$, as $\eta\to 0^+$. It has been shown that a doubling algorithm can be used to compute $X_{\eta}$ efficiently even for very small values of $\eta$, thus providing good approximations to $X_*$. It has been observed by nano scientists that a modified fixed-point method can sometimes be quite useful, particularly for computing $X_{\eta}$ for many different values of the parameter. We provide a rigorous analysis of this modified fixed-point method and its variant and of their generalizations. We also show that the imaginary part $X_I$ of the matrix $X_*$ is positive semidefinite and we determine the rank of $X_I$ in terms of the number of unimodular eigenvalues of the quadratic pencil $\lambda^2 A^{\top}-\lambda Q+A$. Finally we present a new structure-preserving algorithm that is applied directly on the equation $X+A^{\top}X^{-1}A=Q$. In doing so, we work with real arithmetic most of the time.

On the Design of Deterministic Matrices for Fast Recovery of Fourier Compressible Functions

J. Bailey, M. A. Iwen, and C. V. Spencer

SIAM. J. Matrix Anal. & Appl. 33, pp. 263-289 (27 pages)

Online Publication Date: March 27, 2012

Full Text: | Download PDF

Show Abstract
We present a general class of compressed sensing matrices which are then demonstrated to have associated sublinear-time sparse approximation algorithms. We then develop methods for constructing specialized matrices from this class which are sparse when multiplied with a discrete Fourier transform matrix. Ultimately, these considerations improve previous sampling requirements for deterministic sparse Fourier transform methods.
Close

close