Xu' Publications and Preprints

RECENTLY SUBMITTED PAPERS

The modeling framework of port-Hamiltonian systems is systematically extended to linear constrained dynamical systems (descriptor systems, differential-algebraic equations) of arbitrary index and with time-varying constraints. A new algebraically and geometrically defined system structure is derived. It is shown that this structure is invariant under equivalence transformations, and that it is adequate also for the modeling of high-index descriptor systems. The regularization procedure for descriptor systems to make them suitable for simulation and control is modified to preserve the port-Hamiltonian form. The relevance of the new structure is demonstrated with several examples.

PUBLISHED PAPERS

In this paper, we investigate the generalized eigenvalue problem Ax = @lambda Bx arising from economic models. Under certain conditions, there is a simple generalized eigenvalue @rho (A,B) in the interval (0, 1) with a positive eigenvector. Based on the Noda iteration, a modified Noda iteration (MNI) and a generalized Noda iteration (GNI) are proposed for finding the generalized eigenvalue @rho (A,B) and the associated unit positive eigenvector. It is proved that the GNI method always converges and has a quadratic asymptotic convergence rate. So GNI has a similar convergence behavior as MNI. The efficiency of these algorithms is illustrated by numerical examples.

We present weighted Golub-Kahan-Lanczos algorithms. We show their applications to solving the eigenvalue problem of a product of two symmetric positive definite matrices and an eigenvalue problem arising from the linear response problem. Convergence analysis is provided and numerical testing results are reported. As another application we make a connection between the proposed algorithms and the preconditioned conjugate gradient (PCG) algorithm.

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.

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.

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.

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.

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.

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.

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.

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.

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^TA. More generally, we consider the matrix triple (A,G₁,G₂), where G₁ in C ^m,m, G₂ 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₁,G₂) -- (X^T A Y, X^T G₁X, Y^T G₂Y), where X,Y are nonsingular.

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.

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.

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.

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.

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¹⁶, 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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

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.

PREPRINT and REPRINT

The question when a general linear time invariant control system is equivalent to a port-Hamiltonian systems is answered. Several equivalent characterizations are derived which extend the characterizations of [38] to the general non-minimal case. An explicit construction of the transformation matrices is presented. The methods are applied in the stability analysis of disc brake squeal.

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.

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.

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.

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.

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.

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.