A new backward stable, structure preserving method of complexity O(n^3) is presented for computing the stable invariant subspace of a real Hamiltonian matrix and the stabilizing solution of the continuous-time algebraic Riccati equation. The new method is based on the relationship between the invariant subspaces of the Hamiltonian matrix H and the extended matrix /0 H\ and makes use \H 0/ of the symplectic URV-like decomposition that was recently introduced by the authors.

Rank-revealing ULV and URV factorizations are useful tools to determine the rank and to compute bases for null-spaces of a matrix. However, in the practical ULV (resp. URV ) factorization each left (resp. right) null vector is recomputed from its corresponding right (resp. left) null vector via triangular solves. Triangular solves are required at initial factorization, refinement and updating. As a result, algorithms based on these factorizations may be expensive, especially on parallel computers where triangular solves are expensive. In this paper we propose an alternative approach. Our new rank-revealing ULV factorization, which we call ¨top-down¨ ULV factorization ( TDULV -factorization) is based on right null vectors of lower triangular matrices and therefore no triangular solves are required. Right null vectors are easy to estimate accurately using condition estimators such as incremental condition estimator (ICE). The TDULV factorization is shown to be equivalent to the URV factorization with the advantage of circumventing triangular solves.

A new backward stable, structure preserving method of complexity O(n^3) is presented for computing the stable invariant subspace of a real Hamiltonian matrix and the stabilizing solution of the continuous-time algebraic Riccati equation. The new method is based on the relationship between the invariant subspaces of the Hamiltonian matrix H and the extended matrix /0 H\ and makes use \H 0/ of the symplectic URV-like decomposition that was recently introduced by the authors.

Rank-revealing ULV and URV factorizations are useful tools to determine the rank and to compute bases for null-spaces of a matrix. However, in the practical ULV (resp. URV ) factorization each left (resp. right) null vector is recomputed from its corresponding right (resp. left) null vector via triangular solves. Triangular solves are required at initial factorization, refinement and updating. As a result, algorithms based on these factorizations may be expensive, especially on parallel computers where triangular solves are expensive. In this paper we propose an alternative approach. Our new rank-revealing ULV factorization, which we call ¨top-down¨ ULV factorization ( TDULV -factorization) is based on right null vectors of lower triangular matrices and therefore no triangular solves are required. Right null vectors are easy to estimate accurately using condition estimators such as incremental condition estimator (ICE). The TDULV factorization is shown to be equivalent to the URV factorization with the advantage of circumventing triangular solves.

A perturbation analysis shows that if a numerically stable procedure is used to compute the matrix sign function, then it is competitive with conventional methods for computing invariant subspaces. Stability analysis of the Newton iteration improves an earlier result of Byers and confirms that ill-conditioned iterates may cause numerical instability. Numerical examples demonstrate the theoretical results.

A new method is presented for the numerical computation of the generalized eigen- values of real Hamiltonian or symplectic pencils and matrices. The method is strongly backward stable, i.e., it is numerically backward stable and preserves the structure (i.e., Hamiltonian or symplectic). In the case of a Hamiltonian matrix the method is closely related to the square reduced method of Van Loan, but in contrast to that method which may suffer from a loss of accuracy of order sqrt(epsilon), where epsilon is the machine precision, the new method computes the eigenvalues to full possible accuracy.

Two efficient methods for solving generalized Lyapunov equations and their implementations in FORTRAN 77 are presented. The first one is a generalization of the Bartels--Stewart method and the second is an extension of Hammarling's method to generalized Lyapunov equations. Our LAPACK based subroutines are implemented in a quite flexible way. They can handle the transposed equations and provide scaling to avoid overflow in the solution. Moreover, the Bartels--Stewart subroutine offers the optional estimation of the separation and the reciprocal condition number. A brief description of both algorithms is given. The performance of the software is demonstrated by numerical experiments.

This paper describes LAPACK-based Fortran 77 subroutines for the reduction of a Hamiltonian matrix to square-reduced form and the approximation of all its eigenvalues using the implicit version of Van Loan's method. The transformation of the Hamilto- nian matrix to a square-reduced Hamiltonian uses only orthogonal symplectic similarity transformations. The eigenvalues can then be determined by applying the Hessenberg QR iteration to a matrix of half the order of the Hamiltonian matrix and taking the square roots of the computed values. Using scaling strategies similar to those suggested for algebraic Riccati equations can in some cases improve the accuracy of the computed eigenvalues. We demonstrate the performance of the subroutines for several examples and show how they can be used to solve some control-theoretic problems.

A new method is presented for the numerical computation of the generalized eigen- values of real Hamiltonian or symplectic pencils and matrices. The method is strongly backward stable, i.e., it is numerically backward stable and preserves the structure (i.e., Hamiltonian or symplectic). In the case of a Hamiltonian matrix the method is closely related to the square reduced method of Van Loan, but in contrast to that method which may suffer from a loss of accuracy of order sqrt(epsilon), where epsilon is the machine precision, the new method computes the eigenvalues to full possible accuracy.

A restarted symplectic Lanczos method for the Hamiltonian eigenvalue problem is presented. The Lanczos vectors are constructed to form a symplectic basis. Breakdowns and near-breakdowns are overcome by inexpensive implicit restarts. The method is used to compute eigenvalues, eigenvectors and invariant subspaces of large and sparse Hamiltonian matrices and low rank approximations to the solution of continuous-time algebraic Riccati equations with large and sparse coefficient matrices.