Inverse and Direct Solutions of the Discrete Time Lyapunov Equation with System Matrix in Companion Form

The discrete time Lyapunov equation is used in many applications and there is interest in its inverse and direct solutions. New methods are proposed to obtain solutions for cases where the system matrix is in controllable canonical form. The approach is based on the relationship between the discrete Lyapunov equation and the entries of one of the stability tables presented by Jury. It is shown that the inverse solution, which is based on this stability table, can be obtained using LDLt decomposition. Also the direct solution of the discrete Lyapunov equation can be obtained directly from the entries of this stability table. The proposed algorithms are illustrated by numerical examples.


I. INTRODUCTION
T HE Discrete Time Lyapunov Matrix Equation arises in many applications, for instance in the study of discrete time system stability [15], in covariance calculation [36], [37], in the iterative solution of the matrix Riccati equation [19] and in the optimal constant output feedback problem [11], [27]. Other relevant applications include, for instance statistics [20], [21], probability [2] and spectral analysis [12]. The equation is also a computational tool in the design of control systems [22] and in the coprime matrix fraction description of linear systems [49]. It is known that if the pair {A, b} is controllable and U is any positive definite matrix, then the existence of a unique positive definite symmetric solution to (1) is sufficient and necessary for all eigenvalues of the matrix A to lie inside the unit disc. Given the state-space realization x(k + 1) = Ax(k) + bu(k) where u(k) is a zero-mean stationary white random process with covariance U and A is stable, then the steady-state covariance defined by satisfies the Lyapunov equation in (1).
Inverse problems are integrated in many applications and there is high research interest in these problems in last decades. The inverse solution of the Discrete time Lyapunov matrix equation has been considered in applications, such state-covariance assignment [39], the generation of q-Markov covers for discrete systems [40] and the design of Fuzzy controllers [7]. Few results for the inverse solution of the Lyapunov equation have been presented in the literature [39]. Hence there is interest in the derivation of new methods for the inverse solution of the Lyapunov equation. In this paper our aim is to propose a new method to obtain the inverse solution of the Discrete time Lyapunov matrix equation for the case where the system matrix is in controllable canonical form.
In this paper, new methods for obtaining the inverse and the direct solutions of the Discrete time Lyapunov matrix equation for linear discrete systems with the system matrix in controllable canonical form are presented. They are based on the properties of the Mansour matrix [1] and lead to solutions which are closely related to the ∆ j 's, the entries of the stability table presented in [16], [1]. The proposed approach offers insight in the relationship between the discrete Lyapunov equation and the stability table test of the characteristic polynomial.
In section II some preliminary results on the Lyapunov equation and on the Mansour matrix form are presented. Also the stability table presented in [16], [1] is presented in detail. In section III the new methods are presented and illustrated with examples in section IV.

II. PRELIMINARIES
Consider a single-input linear discrete system described by where x(k) ∈ R n , u(k) ∈ R. The problem considered in this paper is to obtain the inverse and the direct solutions of the Lyapunov equation assuming the system matrices A and b are in controllable canonical form [29].
The problem for obtaining the inverse solution of the Lyapunov equation in (4) is stated as follows: Given the positive definite and symmetric matrix X find the system matrix A in the form of (5).The problem for obtaining the direct solution of the Lyapunov equation in (4) is stated as follows: Given the system matrices in (5), (6) find the positive definite and symmetric matrix X which satisfies (4). If the controllablesystem matrix is not given in this form, it is easy to find a similarity transformation matrix that transform the given matrix to the controllable canonical form.
Consider now the basic statements of the proposed method. Let f (z) = z n + n i=1 α i z n−i be a nth degree polynomial with real coefficients. A sequence of polynomials f (z) = F n (z), F n−1 (z), F n−2 (z), . . . of degree , n, n − 1, n − 2, . . . is defined via the stability table presented in [16], [1]. Let with the table The entries of any odd numbered row of the table are calculated from the entries of the previous two rows and given by the formula and equivalently to Define It should be noted that the polynomials F j (z) can be obtained recursively using the ∆ j 's. After manipulations it follows from (10) that (11) or equivalently the coefficients α ij can be expressed as Using the relation between the coefficients α i of f (z) and the quantities ∆ j 's the Mansour form {Σ, b Σ } of the system has been defined and its similarity to a companion matrix has been established. The Mansour form {Σ, b Σ } of the system is represented by the matrices in (13), (14) The following theorem gives the similarity relation between the Mansour matrix and the matrix in companion form.  (5); with the ∆ j as defined above, let the Mansour matrix associated with this system is given by (13) and with the α ij as defined above, let Then The necessary and sufficient conditions for the eigenvalues of the matrix A to be inside the unit disc (stability condition) is where ∆ j are the appropriate entries of the stability table.
The controllability Gramian P of a discrete system in Mansour form has some interesting properties which are summarized in the following lemma.
Lemma 1: The solution P of the Lyapunov equation where {Σ, b Σ } is in Mansour form can be given by the diagonal matrix with Proof: The proof of this lemma is consequence of the discussions in [1]. In this paper the LDL t decomposition is an important method. Let X be a positive-definite symmetric matrix. Then X has unique decomposition X = LDL t where: • L is unit lower triangular matrix, • D is diagonal matrix with positive diagonal entries. This decomposition is discussed in detail in [9].

III. THE INVERSE AND DIRECT SOLUTIONS OF THE DISCRETE TIME LYAPUNOV EQUATION
The inverse and the direct solutions of the Lyapunov equation proposed in this section are based on the properties of the Mansour form summarized in the previous section.

A. THE INVERSE SOLUTION OF THE LYAPUNOV EQUATION
Consider a single-input discrete system {A c , b c } in controllable canonical form and the positive definite and symmetric matrix X c which is the solution of the Lyapunov equation The problem to find the inverse solution of the Lyapunov equation in (23) was also solved in [39]. In this section it will be solved using the facts (presented in lemma 1) that the controllability Gramian in Mansour form is diagonal and that these diagonal elements are related to the ∆ j , j = 1, 2, · · · , n with simple expressions. Consider the LDL t decomposition [9] of X c i.e where L is a lower triangular matrix with 1's in the diagonal and P is a diagonal matrix with positive entries. Using (24) one obtains the following equation from (23) It can be easily seen that L −1 corresponds to the transformation matrix in lemma 1.
The matrix P is diagonal and represents the controllability Gramian in Mansour form given in lemma 1 and therefore, the diagonal elements of the matrix P are related to the quantities ∆ j , j = 1, 2, · · · , n. This implies that the ∆ j , j = 1, 2, · · · , n can be computed from the diagonal matrix P obtained from the the LDL t decomposition of the matrix X c . Further, the entries of the transformation matrix T are also functions of ∆ j , j = 1, 2, · · · , n − 1. The last column of the matrix T is given by [∆ n−1 ∆ n−2 · · · ∆ 2 ∆ 1 1] t . Therefore the ∆ j , j = 1, 2, · · · , n − 1 can be also obtained from this column of the matrix T . Using (20) one obtains the ∆ n and obviously there are two solutions for ∆ n . Consequently this leads to two solutions of this problem. Once the ∆ j , j = 1, 2, · · · , n are known, the system representation {Σ, b Σ } can be obtained directly using lemma 1. These observations lead to the following algorithm: The Algorithm for the Inverse Solution of the Lyapunov Equation Given the positive definite and symmetric matrix X c which is the solution of the Lyapunov equation in (23) and assuming the single-input discrete system {A c , b c } in controllable canonical form, the inverse solution of the Lyapunov equation can be obtained as follows: 1) Compute the lower triangular matrix L with 1's in the diagonal using LDL t decomposition [9] so that X c = LP L t (28) VOLUME X, XXXX P = diag{p 1 , p 2 , · · · , p n } 2) Find ∆ n from This will lead to two solutions for ∆ n .
3) Compute the ∆ j , j = 1, 2, . . . , (n−1) using the following These equations lead to a positive and negative solution for ∆ j . The sign of ∆ j is decided from the last column of matrix T as discussed before.

4) Compute the system matrix representation
Since there are two solutions for ∆ n , this will also lead to two solutions for {A c , b c }.
It should be mentioned that for the inverse solution of the discrete Lyapunov equation few results have been presented in the literature [39]. The method presented there involves the solution of linear systems of equations. The method for the inverse solution of the Lyapunov equation proposed here is novel and based on the relation between the discrete Lyapunov equation and the stability table and can be obtained using LDL t decomposition.
Remarks: 1. It can be easily seen that ∆ j , j = 1, 2, . . . , (n − 1) can be also obtained from the last column of the matrix T . However using the matrix P instead of matrix T is numerically more accurate.

B. THE DIRECT SOLUTION OF THE LYAPUNOV EQUATION
Consider a single-input discrete system {A c , b c } in controllable canonical form. Our aim is to obtain the matrix X c which is the solution of the Lyapunov equation The solution of the equation is a symmetric matrix withe real entries if and only if the system matrix A has no reciprocal eigenvalues,i.e. λ i λ j = 0, i, j ∈ 1, 2, · · · , n. Also the solution is unique and has the so-called Toeplitz symmetric structure. [52], [53] The system in controllable canonical form {A c , b c } can be transformed to the Mansour form {Σ, b Σ } in (13), (14) using the transformation matrix T in (15). Then, the matrix X c can be expressed as where T and P are given in (15), (19)- (22) respectively. The entries of the matrix T and the entries of the controllability Gramian matrix P of a system in Mansour form can be expressed as functions of the ∆ j 's of the stability table. This implies that the matrix X c can be computed directly from the ∆ j , j = 1, 2, · · · , n. These observations lead to the following algorithm: The Algorithm for the Direct Solution of the Lyapunov equation Given a single-input discrete system {A c , b c } in controllable canonical form, the symmetric matrix X c which is the solution of the Lyapunov equation in (4) can be obtained as follows: 1) Find the quantities ∆ j , j = 1, 2, · · · , n using the stability table in [15], [1].
3) Compute the solution of the Lyapunov equation for the discrete system in controllable canonical form as Clearly, the matrix P is positive definite if the matrix A is stable.
As mentioned before some methods for the direct solution of the Lyapunov equation with the system matrix in general form are based on the Kronecker product while others require the solution of a linear system of equations. The solution of the discrete Lyapunov equation with the system matrix in companion form has been studied in [4] and [5]. In [4] the solution involves the solution of a linear system of equations and the computational cost requires in general n 3 flops. In [5] the solution of the discrete Lyapunov equation is the inverse of the Schur-Cohn matrix associated with the characteristic polynomial of the system matrix. In the method proposed in this paper the solution of the Lyapunov equation is the product of three matrices, a diagonal P and two other matrices T and T −t which are triangular with 1's in the diagonal.
The solution involves the inversion of a triangular matrix and the computational cost for that requires n 2 flops. The proposed method requires obtaining the matrices T and P from the entries of the stability table and inverting the upper triangular matrix T with 1's in the diagonal which is simpler than inversion of general and/or positive definite matrices. The proposed algorithms will be illustrated with two examples in the next section.

Example 1
In this example the inverse solution of the Lyapunov equation is studied. Consider the fourth-order system in controllable canonical form described by the state equation (3) and the solution of the Lyapunov equation is given by The LDL t decomposition for the above matrix X c yields The system matrices in the companion form are given by