Back to Xu's home page.

PDF Format

We present weighted Golub-Kahan-Lanczos algorithms for computing a bidiagonal form of a product of two invertible matrices. We give convergence analysis for the algorithms when they are used for computing extreme eigenvalues of a product of two symmetric positive definite matrices. We apply the algorithms for computing extreme eigenvalues of a double-sized matrix arising from linear response problem and provide convergence results. Numerical testing results are also reported. As another application we show a connection between the proposed algorithms and the conjugate gradient (CG) algorithm for positive definite linear systems.

Preprint 32-2015 (TU Berlin)

The sign characteristics of Hermitian matrix polynomials are discussed, and in particular an appropriate definition of the sign characteristics associated with the eigenvalue infinity. The concept of sign characteristic arises in different forms in many scientific fields, and is essential for the stability analysis in Hamiltonian systems or the perturbation behavior of eigenvalues under structured perturbations. We extend classical results by Gohberg, Lancaster, and Rodman to the case of infinite eigenvalues. We derive a systematic approach, studying how sign characteristics behave after an analytic change of variables, including the important special case of Mobius transformations, and we prove a signature constraint theorem. We also show that the sign characteristic at infinity stays invariant in a neighborhood under perturbations for even degree Hermitian matrix polynomials, while it may change for odd degree matrix polynomials. We argue that the non-uniformity can be resolved by introducing an extra zero leading matrix coefficient.

PDF Format

We show that a Schur form of a real orthogonal matrix can be obtained from a full CS decomposition. Based on this fact a CS decomposition-based orthogonal eigenvalue method is developed. We also describe an algorithm for orthogonal similarity transformation of an orthogonal matrix to a condensed product form, and an algorithm for full CS decomposition. The latter uses mixed shifted and zero-shift iterations for high accuracy. Numerical examples are presented.

PDF Format

In this work, we are interested in effects of a simple profile of permanent charges on ionic flows. We determine when a permanent charge produces current reversal. We adopt the classical PoissonÐNernstÐPlanck (PNP) models of ionic flows for this study. The starting point of our analysis is the recently developed geometric singular perturbation approach for PNP models. Under the setting in the paper for case studies, we are able to identify a single governing equation for the existence and the value of the permanent charge for a current reversal. A number of interesting features are established. The related topic on reversal potential can be viewed as a dual problem and is briefly examined in this work too.

PDF Format

In this work, we examine the stationary one-dimensional classical PoissonÐNernstÐPlanck(cPNP) model for ionic flow Ð a singularly perturbed boundary value problem (BVP). For the case of zero permanent charge, we provide a complete answer concerning the existence and uniqueness of the BVP. The analysis relies on a number of ingredients: a geometric singular perturbation framework for a reduction to a singular BVP, a reduction of the singular BVP to a matrix eigenvalue problem, a relation between the matrix eigenvalues and zeros of a meromorphic function, and an application of the Cauchy Argument Principle for identifying zeros of the meromorphic function. Once the zeros of the meromorphic function in a stripe are determined, an explicit solution of the singular BVP is available. It is expected that this work would be useful for studies of other PNP systems.

PDF Format

The long standing problem is discussed of how to deflate the part associated with the eigenvalue infinity in a structured matrix pencil using structure preserving unitary transformations. We derive such a deflation procedure and apply this new technique to symmetric, Hermitian or alternating pencils and in a modified form to (anti)-palindromic pencils. We present a detailed error and perturbation analysis of this and other deflation procedures and demonstrate the properties of the new algorithm with several numerical examples.

PDF Format

Bounds are developed for the condition number of the linear finite element equations of an anisotropic diffusion problem with arbitrary meshes. They depend on three factors. The first factor is proportional to a power of the number of mesh elements and represents the condition number of the linear finite element equations for the Laplacian operator on a uniform mesh. The other two factors arise from the mesh nonuniformity viewed in the Euclidean metric and in the metric defined by the diffusion matrix. The new bounds reveal that the conditioning of the finite element equations with adaptive anisotropic meshes is much better than what is commonly assumed. Diagonal scaling for the linear system and its effects on the conditioning are also studied. It is shown that the Jacobi preconditioning, which is an optimal diagonal scaling for a symmetric positive definite sparse matrix, can eliminate the effects of mesh nonuniformity viewed in the Euclidean metric and reduce those effects of the mesh viewed in the metric defined by the diffusion matrix. Tight bounds on the extreme eigenvalues of the stiffness and mass matrices are obtained. Numerical examples are given.

PDF Format

We present formulas for the construction of optimal H_{&infin}
controllers that can be implemented in a numerically robust way. We
base the formulas on the &gamma-iteration developed in
[Benner, Byers, Mehrmann, and Xu, Linear Algebra Appl. 425, 548 - 570, 2007].
The controller formulas proposed here avoid the solution of algebraic
Riccati equations with their problematic matrix inverses and matrix products. They are also applicable in the neighborhood
of the optimal &gamma , where the classical formulas may call for the inverse of
singular or ill-conditioned matrices. The advantages
of the new formulas are demonstrated by several numerical examples.

PDF Format

Formulas are given for functions of a matrix *A* in
terms of Krylov matrices of *A*, under the assumption
that *A* is nonderogatory.
Relations between the
coefficients of a polynomial of *A* and the generating
vector of a Krylov matrix of *A* are provided.
With the formulas,
linear transformations between
Krylov matrices and functions of *A* are introduced,
and associated algebraic properties are derived.
Hessenberg reduction forms are revisited equipped
with appropriate
inner products and related matrix factorizations
are given.
Properties about functions, Krylov matrices,
and nonderogatory conditions are also provided.

PostScript Format or PDF Format

The classical singular value decomposition for a matrix
A in C^{m,n} is a
canonical form for A that also displays the eigenvalues of the Hermitian
matrices AA^{*} and A^{*} A. In this paper, we develop a
corresponding
decomposition for A that provides the Jordan canonical forms for the
complex symmetric matrices AA^{T} and A^{T}A. More
generally, we consider the matrix triple (A,G_{1},G_{2}),
where G_{1} in C ^{m,m},
G_{2} in C^{n,n} are invertible and
either complex symmetric and complex skew-symmetric,
and we provide a canonical form under transformations of the form
(A,G_{1},G_{2}) -- (X^{T} A Y,
X^{T} G_{1}X, Y^{T} G_{2}Y),
where X,Y are nonsingular.

PDF Format

Canonical forms for matrix triples (A,G,{\hat G}), where A is arbitrary rectangular and G, {\hat G} are either real symmetric or skew symmetric, or complex Hermitian or skew Hermitian, are derived. These forms generalize classical product Schur forms as well as singular value decompositions. An new proof for the complex case is given, where there is no need to distinguish whether G and {\hat G} are Hermitian or skew Hermitian. This proof is independent from the results in [Bolschakov and Reichstein'95], where a similar canonical form has been obtained for the complex case, and it allows generalization to the real case. Here, the three cases, i.e., that G and {\hat G} are both symmetric, both skew symmetric or one each, are treated separately.

PostScript Format or PDF Format

We give several different formulations for the discrete-time linear-quadratic control problem in terms of structured eigenvalue problems, and discuss the relationships among the associated structured objects: symplectic matrices and pencils, BVD-pencils and polynomials, and the recently introduced classes of palindromic pencils and matrix polynomials. We show how these structured objects can be transformed into each other, and also how their eigenvalues, eigenvectors and invariant/deflating subspaces are related.

PostScript Format or PDF Format

We discuss the eigenvalue problem for general and structured matrix polynomials which may be singular and may have eigenvalues at infinity. We derive condensed forms that allow (partial) deflation of the infinite eigenvalue and singular structure of the matrix polynomial. The remaining reduced order staircase form leads to new types of linearizations which determine the finite eigenvalues and corresponding eigenvectors. The new linearizations also simplify the construction of structure preserving linearizations.

PostScript Format or PDF Format

We derive formulas for the minimal positive solution of a particular non-symmetric Riccati equation arising in transport theory. The formulas are based on the eigenvalues of an associated matrix. We use the formulas to explore some new properties of the minimal positive solution and to derive fast and highly accurate numerical methods. Some numerical tests demonstrate the properties of the new methods.

PostScript Format or PDF Format

We propose a scaling scheme for Newton's iteration
for calculating the polar decomposition. The scaling
factors are generated by a simple scalar iteration in
which the initial value depends only on estimates of
the extreme singular values of the original matrix, which
can for example be the Frobenius norms of the matrix and its
inverse. In exact arithmetic, for matrices with condition number
no greater than 10^{16}, with this scaling scheme,
no more than 9 iterations are needed for convergence to the unitary
polar factor with a convergence tolerance roughly equal to
10^{-16}. It is proved that if matrix inverses computed in
finite precision arithmetic satisfy a backward-forward error model
then the numerical method is backward stable. It is also proved
that Newton's method with Higham's scaling or with Frobenius norm
scaling is backward stable.

PostScript Format or PDF Format

We discuss the perturbation theory for purely imaginary eigenvalues of Hamiltonian matrices under Hamiltonian and non-Hamiltonian perturbations. We show that there is a substantial difference in the behavior under these perturbations. We also discuss the perturbation of real eigenvalues of real skew-Hamiltonian matrices under structured perturbations and use these results to analyze the properties of the URV method of computing the eigenvalues of Hamiltonian matrices.

PostScript Format or PDF Format

We present a numerical method for the
solution of the optimal H_{&infin}
control problem based on the gamma-iteration and
a novel extended matrix pencil formulation
of the state-space solution to the (sub)optimal H_{&infin}
control problem. In particular, instead of algebraic Riccati
equations or unstructured matrix pencils, our approach is solely
based on solving even generalized eigenproblems. The enhanced numerical
robustness of the method is
derived from the fact that using the structure of the problem,
spectral symmetries are preserved. Moreover,
these methods are also applicable even if the pencil has
eigenvalues on the imaginary axis.
We compare the new method with conventional methods and present
several examples.

PostScript Format or PDF Format

We introduce a transformation that connects the discrete-time and continuous-time algebraic Riccati equations. We show that under mild conditions one algebraic Riccati equation can be derived from another by the transformation, and both algebraic Riccati equations share common Hermitian solutions. Moreover, the properties that are parallelly imposed, commonly in system and control setting, to the coefficient matrices and Hermitian solutions of the discrete-time and continuous-times algebraic Riccati equations are equivalently related. The transformation is simple and all the relations can be easily derived. We also introduce a generalized transformation that needs weaker conditions. The proposed transformations may provide a unified tool to develop the theories and numerical methods for the algebraic Riccati equations and the associated system and control problems. They also give a different way to understand and interpret the algebraic Riccati equations.

PostScript Format or PDF Format

We present structure preserving algorithms for the numerical computation of structured staircase forms of skew-symmetric/symmetric matrix pencils along with the Kronecker indices of the associated skew-symmetric/symmetric Kronecker-like canonical form. These methods allow deflation of the singular structure and deflation of infinite eigenvalues with index greater than one. Two algorithms are proposed: one for general skew-symmetric/symmetric pencils and one for pencils in which the skew-symmetric matrix is a direct sum of 0 and ${\cal J}=[0, I;-I, 0]$. We show how to use the structured staircase form to solve boundary value problems arising in control applications and present numerical examples.

PostScript Format or PDF Format

Preprint No 175 (old version), MATHEON, DFG Research Center Mathematics for key technologies in Berlin, TU Berlin, Str. des 17. Juni 136, D-10623 Berlin, Germany, 2004.

We introduce a transformation between the generalized symplectic pencils and the skew-Hermitian/Hermitian pencils. Under the transformation the regularity of the matrix pencils is preserved, and the equivalence relations about their eigenvalues and deflating subspaces are established. The eigenvalue problems of the generalized symplectic pencils and skew-Hermitian/Hermitian pencils are strongly related to the discrete-time and continuous-time robust control problems, respectively. With the transformation a simple connection between these two types of robust control problems is made. The connection may help to develop unified methods for solving the robust control problems.

PostScript Format or PDF Format

We present a numerical method to compute the SVD-like decomposition B=QDS^{-1}, where Q is orthogonal, S is symplectic and D is a permuted diagonal matrix. The method can be applied directly to compute the canonical form of the Hamiltonian matrices of the form JB^TB, where J=[0,I;-I,0]. It can also be applied to solve the related application problems such as the gyroscopic systems and linear Hamiltonian systems. Error analysis and numerical examples show that the eigenvalues of JB^TB computed by this method are more accurate than that computed by the methods working on the explicit product JB^TB or BJB^T.

PostScript Format or PDF Format

We present a numerical method for solving the indefinite least squares problem. We first normalize the coefficient matrix. Then we compute the hyperbolic QR factorization of the normalized matrix. Finally we compute the solution by solving several triangular systems. We give the first order error analysis to show that the method is backward stable. The method is more efficient than the backward stable method proposed by Chandrasekaran, Gu and Sayed.

PostScript Format or PDF Format

We discuss matrix pencils with a double symmetry structure that arise in the Hartree-Fock model in quantum chemistry. We derive anti-triangular condensed forms from which the eigenvalues can be read off. Ideally these would be condensed forms under unitary equivalence transformations that can be turned into stable (structure preserving) numerical methods. For the pencils under consideration this is a difficult task and not always possible. We present necessary and sufficient conditions when this is possible. If this is not possible then we show how we can include other transformations that allow to reduce the pencil to an almost anti-triangular form.

PostScript Format or PDF Format

A matrix S \in {\mathbb C}^{2m \times 2m} is symplectic if S J S^\ast = J, where J=[0, I; -I, 0]. Symplectic matrices play an important role in the analysis and numerical solution of matrix problems involving the indefinite inner product x^\ast (iJ) y. In this paper we provide several matrix factorizations related to symplectic matrices. We introduce a singular value-like decomposition B = Q D S^{-1} for any real matrix B \in {\mathbb R}^{n \times 2m}, where Q is real orthogonal, S is real symplectic, and D is permuted diagonal. We show the relation between this decomposition and the canonical form of real skew-symmetric matrices and a class of Hamiltonian matrices. We also show that if S is symplectic it has the structured singular value decomposition S=U D V^\ast, where U, V are unitary and symplectic, D = diag(\Omega, \Omega^{-1}) and \Omega is positive diagonal. We study the B J B^T factorization of real skew-symmetric matrices. The B J B^T factorization has the applications in solving the skew-symmetric systems of linear equations, and the eigenvalue problem for skew-symmetric/symmetric pencils. The B J B^T factorization is not unique, and in numerical application one requires the factor B with small norm and condition number to improve the numerical stability. By employing the singular value-like decomposition and the singular value decomposition of symplectic matrices we give the general formula for B with minimal norm and condition number.

PostScript Format or PDF Format

We discuss the numerical solution of structured generalized eigenvalue problems that arise from linear-quadratic optimal control problems, H_infinity optimization, multibody systems, and many other areas of applied mathematics, physics, and chemistry. The classical approach for these problems requires computing invariant and deflating subspaces of matrices and matrix pencils with Hamiltonian and/or skew-Hamiltonian structure. We extend the recently developed methods for Hamiltonian matrices to the general case of skew-Hamiltonian/Hamiltonian pencils. The algorithms circumvent problems with skew-Hamiltonian/Hamiltonian matrix pencils that lack structured Schur forms by embedding them into matrix pencils that always admit a structured Schur form. The rounding error analysis of the resulting algorithms is favorable. For the embedded matrix pencils, the algorithms use structure preserving unitary matrix computations and are strongly backwards stable, i.e., they compute the exact structured Schur form of a nearby matrix pencil with the same structure.

PostScript Format or PDF Format

We study the perturbation theory for the eigenvalue problem of a formal matrix product A1^(s1) ... Ap^(sp) , where all Ak are square and sk = 1 or -1. We generalize the classical perturbation results for matrices and matrix pencils to perturbation results for generalized deflating subspaces and eigenvalues of such formal matrix products. As an application we then extend the structured perturbation theory for the eigenvalue problem of Hamiltonian matrices to Hamiltonian/skew-Hamiltonian pencils.

PostScript Format or PDF Format

The existence, uniqueness and parametrization of Lagrangian invariant subspaces for Hamiltonian matrices is studied. Necessary and sufficient conditions and a complete parametrization are given. Some necessary and sufficient conditions for the existence of Hermitian solutions of algebraic Riccati equations follow as simple corollaries.

PostScript Format or PDF Format

In this paper we derive canonical forms under structure preserving equivalence transformations for matrices and matrix pencils that have a multiple structure, which is either an H-selfadjoint or H-skew-adjoint structure, where the matrix H is a complex nonsingular Hermitian or skew-Hermitian matrix. Matrices and pencils of such multiple structures arise for example in quantum chemistry in Hartree-Fock models or random phase approximation.

PostScript Format or PDF Format

We study classical control problems like pole assignment, stabilization, linear quadratic control and H-infty control from a numerical analysis point of view. We present several examples that show the difficulties with classical approaches and suggest re-formulations of the problems in a more general framework. We also discuss some new algorithmic approaches.

PostScript Format or PDF Format

For inner products defined by a symmetric indefinite matrix Sigma_{p,q} = [Ip, 0; 0 -Iq], we study canonical forms for real or complex Sigma_{p,q}-Hermitian matrices, Sigma_{p,q}-skew Hermitian matrices and Sigma_{p,q}-unitary matrices under equivalence transformations which keep the class invariant.

PostScript Format or PDF Format

We study canonical forms for Hamiltonian and symplectic matrices or pencils under equivalence transformations which keep the class invariant. In contrast to other canonical forms our forms are as close as possible to a triangular structure in the same class. We give necessary and sufficient conditions for the existence of Hamiltonian and symplectic triangular Jordan, Kronecker and Schur forms. The presented results generalize results of Lin and Ho and simplify the proofs presented there.

PostScript Format or PDF Format

In this paper we describe a simple observation that can be used to extend two recently proposed structure preserving methods for the eigenvalue problem for real Hamiltonian matrices to the case of complex Hamiltonian and skew-Hamiltonian matrices.

PostScript Format or PDF Format

We present a constructive existence proof that every real skew-Hamiltonian matrix W has a real Hamiltonian square root. The key step in this construction shows how one may bring any such W into a real quasi-Jordan canonical form via symplectic similarity. We show further that every W has infinitely many real Hamiltonian square roots, and give a lower bound on the dimension of the set of all such square roots. Some extensions to complex matrices are also presented.

PostScript Format or PDF Format

A new method is presented for the numerical computation of the generalized eigenvalues 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 square root of the machine precision, the new method computes the eigenvalues to full possible accuracy.

PostScript Format or PDF Format

We discuss the single-input pole placement problem (SIPP) and analyze how the conditioning of the problem can be estimated and improved if the poles are allowed to vary in specific regions in the complex plane. Under certain assumptions we give formulas as well as bounds for the norm of the feedback gain and the condition number of the closed loop matrix. Via several numerical examples we demonstrate how these results can be used to estimate the condition number of a given SIPP problem and also how to select the poles to improve the conditioning.

PostScript Format or PDF Format

For an Hermitian matrix the QR transform is diagonally similar to two steps of the LR transform. Even for non Hermitian matrices the QR transform may be written in rational form.

PostScript Format or PDF Format

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; H, 0] and makes use of the symplectic URV-like decomposition that was recently introduced by the authors.

PostScript Format or PDF Format

Two results about the matrix exponential are given. One is to characterize the matrices A which satisfy e^A e^(A*) = e^(A*)e^A, another is about the upper bounds of trace(e^Ae^(A*)). When A is stable, the bounds preserve the asymptotic stability.

PostScript Format or PDF Format

For the solution of the multi-input pole placement problem we derive explicit formulas for the subspace from which the feedback gain matrix can be chosen and for the feedback gain as well as the eigenvector matrix of the closed-loop system. We discuss which Jordan structures can be assigned and also when diagonalizability can be achieved. Based on these formulas we study the conditioning of the pole-placement problem in terms of perturbations in the data and show how the conditioning depends on the condition number of the closed loop eigenvector matrix, the norm of the feedback matrix and the distance to uncontrollability.

PostScript Format or PDF Format

For the solution of the single-input pole placement problem we derive explicit expressions for the feedback gain matrix as well as the eigenvector matrix of the closed-loop system. Based on these formulas we study the conditioning of the pole-placement problem in terms of perturbations in the data and show how the conditioning depends on the condition number of the closed loop eigenvector matrix, which is a similar to a generalized Cauchy matrix, the norm of the feedback vector and the distance to uncontrollability.

We discuss some properties of a quadratic matrix equation with some restrictions, then use these results on the algebraic Riccati equation to get a new algorithm. The algorithm sufficiently takes account of the structure of the associated matrix; hence it is very effective.

In this paper we gives a Bauer-Fike like perturbation result about the separation of two matrices. We show that the perturbation bounds depend on the eigenvalues, the size of Jordan blocks, and the condition numbers of the matrices.

PostScript Format or PDF Format

We discuss the numerical solution of linear quadratic optimal control problems and H infty control problems. A standard approach for these problems is the solution of algebraic Riccati equations. Recently for this approach new structure preserving methods have been developed which are faster than the currently used methods and give results of full possible accuracy by making use of the underlying structure. These methods can be used also for Riccati equations with an associated Hamiltonian that has eigenvalues on the imaginary axis.

PostScript Format or PDF Format

A unified deflating subspace approach is presented for the solution of a large class of matrix equations, including Lyapunov, Sylvester, Riccati and also some higher order polynomial matrix equations including matrix m-th roots and matrix sector functions. A numerical method for the computation of the desired deflating subspace is presented that is based on adapted versions of the periodic QZ algorithm.

PostScript Format

We present a constructive existence proof that every real skew-Hamiltonian matrix W has a real Hamiltonian square root. The key step in this construction shows how one may bring any such W into a real quasi-Jordan canonical form via symplectic similarity. We show further that every W has infinitely many real Hamiltonian square roots, and give a lower bound on the dimension of the set of all such square roots. Some extensions to complex matrices are also presented.

PostScript Format or PDF Format

Several MEX-files are developed based on SLICOT Fortran subroutines. The MEX-files provide new tools for the numerical solution of some classical control problems, such as the solution of linear or Riccati matrix equations computations in the MATLAB environment. Numerical tests show that the resulting MEX-files are equally accurate and much more efficient than the corresponding MATLAB functions in the control system toolbox and the robust control toolbox. In order to increase user-friendlyness the related m-files are also developed so that the MEX-file interface to the corresponding SLICOT routines can be implemented directly and easily.

Preprint SFB393/98-25, Sonderforschungsbereich 393, `Numerische Simulation auf massiv parallelen Rechnern'. Fakultat feur Mathematik, Tech University Chemnitz, FR Germany, 1998.

PostScript Format

The existence and uniqueness of Lagrangian invariant subspaces of Hamiltonian matrices is studied. Necessary and sufficient conditions are given in terms of the Jordan structure and certain sign characteristics that give uniqueness of these subspaces even in the presence of purely imaginary eigenvalues. These results are applied to obtain in special cases existence and uniqueness results for Hermitian solutions of continuous time algebraic Riccati equations.

PostScript Format

We study canonical forms for Hamiltonian and symplectic matrices or pencils under equivalence transformations which keep the class invariant. In contrast to other canonical forms our forms are as close as possible to a triangular structure in the same class. We give necessary and sufficient conditions for the existence of Hamiltonian and symplectic triangular Jordan, Kronecker and Schur forms. The presented results generalize results of Lin and Ho and simplify the proofs presented there.

Hongguo Xu < xu@math.ku.edu > Last modified Jan. 13, 2016.