Refinement Methods for State Estimation via Sylvester-Observer Equation

We present new iterative methods based on refinement process for solving large sparse Sylvesterobserver equations applied in state estimation of a continuous-time system. These methods use projection methods to produce low-dimensional Sylvester-observer matrix equations that are solved by the direct methods. Moreover, the refinement process described in this paper has the capability of improving the results obtained by any other methods. Some numerical results will be reported to illustrate the efficiency of the proposed methods.


Introduction
where A ∈ R n×n , B ∈ R n×m , and C ∈ R r×n .We well know that all the state-feedback problems, such as the feedback stabilization, the LQR, and the state-feedback H ∞ -control problems, require the state vector x t explicitly 1 .In most practical situations, the states are not fully accessible and the designer only knows the output and the input vectors.The unavailable states somehow need to be estimated accurately from the knowledge of the matrices A, B, and C, the output vector y t , and the input vector u t .There are two common procedures for state estimation: one via eigenvalue assignment EVA and the other via solution of the Sylvester-observer equation.

Some Fundamental Results
The two basic concepts in control theory are controllability and observability of a control system.In this section, we will state some well-known facts about controllability and observability for convenient use later in the paper.For an excellent account of controllability and observability and related results, see 1 .
Definition 2.1.The system 1.1 is said to be controllable if, starting from any initial state x 0 , the system can be driven to any final state x f in some finite time t f , choosing the input vector u t , 0 ≤ t ≤ t f , appropriately.
Observability is the dual concept of controllability.
Definition 2.2.The system 1.1 is said to be observable if there exists t 1 > 0 such that the initial state x 0 can be uniquely determined from the knowledge of u t and y t for all 0 ≤ t ≤ t 1 .
Remark 2.3.Since matrix C does not have any role in the definition of controllability, the controllability of the system 1.1 see 1 , is often referred to as the controllability of the pair A, B .Similarly, since matrix B does not have any role in the definition of observability, the observability of the system 1.1 is often referred to as the observability of the pair A, C .Some well-known criteria of controllability and observability are now stated in the next two theorems.The proofs of Theorems 2.4 and 2.5 can be found in 1 .In the following,

Theorem 2.4. The pair A, B is controllable if and only if any one of the following equivalent conditions holds.
1 The controllability matrix 2 Rank A − λI, B n for every eigenvalue λ of A.
3 Let λ, x be an eigenpair of A, that is, Ay λy; then Cy / 0.

Theorem 2.5. The pair A, C is observable if and only if any one of the following equivalent conditions holds
1 The observability matrix has rank n.
2 The matrix λI−A C has rank n for every eigenvalue λ of A.
Definition 2.6.A matrix A is called a stable matrix if all of the eigenvalues of A have negative real parts.

State Estimation
In this section, we discuss how the states of a continuous-time linear system 1.1 can be estimated.The discussions here apply equally to the discrete-time systems possibly with some minor changes.So we concentrate on the continuous-time case only.We describe two common procedures for state estimation: one via eigenvalue assignment EVA and the other,via solution of the Sylvester-observer equation.

State Estimation via Eigenvalue Assignment
Now consider the linear time-invariant system 1.1 .Let x t be an estimate of the state vector x t .Obviously, we would like to construct the vector x t in such a way that the error e t x t − x t approaches zero as fast as possible, for all initial states x 0 and for every input u t .The following theorem shows that the problem of state estimation can be solved by finding a matrix K such that the matrix A − KC has a suitable desired spectrum.See 1 .Proof.See 1 .

State Estimation via Sylvester Equation
There is another approach for state estimation.Knowing A, B, C, u t , and y t , let us construct the system ż t Fz t Gy t Pu t , 3.
where F is n×n, G is n×r, and P is n×m, in such a way that for some constant n×n nonsingular matrix X the error vector e t z t − Xx t → 0 for all x 0 , z 0 , and for every input u t .The vector z t will then be an estimate of Xx t .The system 3.3 is then said to be the state observer for the system 1.1 .D. Luenberger originated the idea and is hence referred to in control theory as the Luenberger observer; see 1 .
We now show that the system 3.3 will be a state observer if the matrices X, F, G, and P satisfy certain requirements.

Inputs
The system matrices A, B, and C of order n × n, n × m, and r × n, respectively.

Output
An estimate x t of the state vector x t .
Assumptions 1. A, C is observable.
Step 1. Find a nonsingular solution of the Sylvester-observer equation 1.2 by choosing F as a stable matrix and choosing G in such a way that the resulting solution X is nonsingular.
Step 3. Construct the observer z t by solving the system of differential equations z t Fz t Gy t Pu t , z 0 z 0 .

3.5
Step 4. Find an estimate x t of x t : x t X −1 z t .

Characterization of Nonsingular Solutions of the Sylvester Equation
In this section, we describe some necessary conditions for a unique solution of the Sylvester equation to have such properties.The following theorem was proved by Bhattacharyya and de Souza.The proof here has been taken from 1 .
Theorem 3.5.Let A, F, G, and C, respectively, be of order n × n, n × n, n × r, and r × n.Let X be a unique solution of the Sylvester-observer equation 1.2 .Then, necessary conditions for X to be nonsingular are that A, C is observable and F, G is controllable.
Proof.See 1 .Example 3.8.In this example we show how the Sylvester-observer equation 1.2 can be applied in state estimation of a continuous time system 1.1 .In this sense, at first we use the MATLAB function ode23 for directly solving 1.1 with u t as the unit step function and x 0 6, 0 T .Then, we apply Algorithm 3.4 for computing the estimate x t .Also, the differential equation 3.12 was solved with z 0 0 and MATLAB function ode23.The comparison of the state x t and estimate x t for this example is shown in Figure 1.The solid line corresponds to the exact state and the dotted line corresponds to the estimate state.Let According to criteria 1 of Theorem 2.5, A, C is observable.Thus, we can use Algorithm 3.4 for state estimation of 1.1 .
Step 5. Choose According to criteria 1 of Theorem Remark 3.9.According to Algorithm 3.4 the most important step is Step 1.In the case that n is small, there are many reliable algorithm s for solving the Sylvester-observer equation and the states of a continuous-time system can be estimated.However, for large and sparse systems solving the Sylvester-observer equation by the available methods can be costly.In Sections 4 and 5, we introduce two iterative refinement methods for solving large sparse Sylvesterobserver equations.

Block Refinement Method
As we already mentioned, so far many numerical methods have been developed by different authors; see 1, 13, 17 .For example, the Hessenberg-Schur method is now widely used as an effective computational method for the Sylvester-observer equation.But numerical stability of this method has not been investigated.As the iterative methods are very efficient for the solution of computational problems, we thought it will be good idea to create an iterative method for solving the Sylvester-observer equation XA − FX GC where A, F, X ∈ R n×n , G ∈ R n×r , and C ∈ R r×n .In this section we propose to show that the obtained approximate solution of the Sylvester-observer equation by any method can be improved, in other words the accuracy can be increased.If this idea is applicable, then we have found an iterative method for solving of the Sylvester-observer equation.
Proof.Let V m v 1 , . . ., v m and W m w 1 , . . ., w m be two orthogonal bases constructed by the block Arnoldi process for the matrices A T and F, respectively.Thus, we have Also, the square block Hessenberg matrices A m and F m m r , where r and are the dimensions of blocks whose nonzero entries are the scalars a ij and f ij , constructed by the Block Arnoldi process, can be expressed as If we set

4.6
Since Y m is the solution of 4.1 and by using 4.3 , we have According to Theorem 4.1, we can develop an iterative method for solving the Sylvester-observer equation when the matrices A, F, G, and C are large and sparse.For achieving this idea, if we choose m < n, then instead of solving XA − FX GC we can solve 4.1 .In other words, in this method, first we transform the initial Sylvester-obsever equation to another Sylvester equation with less dimensions and then in each iteration step solve this matrix equation and extend the obtained solution to the solution of the initial equation by 4.5 .The algorithm is as follows.
Algorithm 4.2 block refinement method . 1 Start: choose an initial approximate solution X 0 , and a tolerance ε.
2 Select two numbers r and for dimensions of block and set m r m < n .
3 Compute R X 0 GC − X 0 A − FX 0 .4 Construct the orthonormal bases V m and W m ∈ R n×m by the block Arnoldi process, such that R X 1 and go to step 3 .Remark 4.3.By choosing m < n, Algorithm 4.2 reduces the original large sparse Sylvesterobserver equation 1.2 to a low-dimensional Sylvester-observer equation 4.1 .In step 5 , we solve this low-dimensional matrix equation by any direct method such as the Hessenberg-Schur method.Then, in step 6 by using relation 4.5 , we extend the obtained solution to the solution of the original matrix equation.Also, according to Theorem 4.1, Algorithm 4.2 is the convergence for any initial matrix X 0 .

Weighted Block Refinement Method
In this section we discuss a new iterative method based upon a modified block refinement method.The new process uses instead of the Euclidean scalar product another one, denoted by •, • D where D is a chosen diagonal matrix.The idea of changing the inner product is to accelerate the convergence of the components of the residual which are far away from zero.To achieve this, an appropriate weight is associated to each term of the inner product.A natural choice of these weights is the entries of the first residual.The following method is based on reduction of A and F to the Hessenberg matrix with the use of weighted block Arnoldi process.Before giving a complete description of the new algorithm, let us define the D-scalar product.
If u and v are two vectors of R n , their D-scalar product is where D diag d 1 , d 2 , . . ., d n is a diagonal matrix.This inner product is well defined if and only if the matrix D is positive definite, that is, d i > 0, for all i ∈ {1, . . ., n}.
In this case, we can define the D-norm • D associated with this inner product by

5.2
Theorem 5.1.Let X 0 be the approximate solution obtained by an arbitrary method for the matrix equation 1.2 , and let R X 0 GC − X 0 A − FX 0 be the corresponding residual.

If m m n steps of the weighted block Arnoldi process by the diagonal matrices D and D, respectively, for matrices A and F have been run and Y m ∈ R m×m is the solution of the low-dimensional Sylvester equation
Proof.By using the weighted Arnoldi process, we generate the bases U m u 1 , . . ., u m and V m v 1 , . . ., v m that are, respectively, D-orthonormal and D-orthonormal; thus it holds that where U m , V m ∈ R n×m and D, D ∈ R n×n are two diagonal matrices.Moreover, the square Hessenberg matrices A m and F m whose nonzero entries are the scalars a ij and f ij , constructed by the weighted Arnoldi process, can be expressed in the form Now, we set where is the solution of the Sylvester-obsever equation 5.3 .Thus, the new residual matrix becomes

5.8
Multiplying the above relation on the left by V T m and on the right by U m , we have Now, by using 5.3 , 5.5 , and 5.6 we get In order to get Y m ∈ R m×m , we need to solve the low-dimensional Sylvester equation 5.3 .
According to the results, we can develop an iterative method for solving of the Sylvesterobserver equation.The algorithm is as follows.
Algorithm 5.2 weighted block refinement WBR method . 1 Start: choose an initial solution X 0 , new dimension m lesser than n and a tolerance .
3 Construct the D-orthonormal basis U m ∈ R n×m and D-orthonormal basis V m ∈ R n×m by the weighted Arnoldi process, such that R X 1 and go to step 2 .
Remark 5.3.By choosing m < n, Algorithm 5.2 reduces the original large sparse Sylvesterobserver equation 1.2 to a low-dimensional Sylvester-observer equation 5.3 .In step 4 , we solve this low-dimensional matrix equation by any direct method such as the Hessenberg-Schur method.Also, according to Theorem 5.1, Algorithm 5.2 is the convergence for any initial matrix X 0 .

Numerical Experiments
In this section, we present some numerical examples to illustrate the effectiveness of our new iterative methods for solving large sparse Sylvester-obsever equation.In Examples 6.1 and 6.2, we apply Algorithms 2 and 3 for solving matrix equation 1.2 .In Example 6.3, we compare the Hessenberg-Schur method described in 1 with our new algorithms for solving large sparse Sylvester-obsever equation.In order to show the efficiency of our algorithms, we choose the matrices A and C arbitrary in these three examples.But in Example 6.4, we use four matrices from MATLAB matrix collection with the large estimation of condition numbers.
The initial approximate solution is X 0 0 n×n .The error is monitored by means of the test  The desired accuracy has been chosen as 10 −6 , but the model works well with any choice of 10 −t : 1, 0, 0, . . ., 0 .

6.2
Example 6.2.Consider Example 6.1 again.We apply the weighted Block refinement method for solving XA − FX GC and take 10 −6 .In Table 2, we report the results for different values of m.Example 6.3.According to the results in Tables 1 and 2, we see that the weighted block refinement method in comparison with block refinement method works better.Now, consider that A and C are the same matrices used in Example 6.1.We apply our two iterative methods with 2 iterations and the Hessenberg-Schur method to solve the Sylvester-observer equation when the dimensions of the matrices are large.Results are shown in Table 3. Example 6.4.In this example we show that the convergence of our proposed algorithms independent of the matrices structure.In this sense, we use four matrices from MATLAB collection for the matrix A. The first matrix is a sparse, random finite element matrix with the condition number 1.84E 03.The second matrix is a symmetric, positive semidefinite SPD Toeplitz matrix that is composed of the sum of 800 rank 2 SPD Toeplitz matrices with the condition number 6.18E 04.The third matrix is a row diagonally dominant matrix with the condition number 1.35E 010.The last matrix is a sparse singular, row diagonally dominant matrix resulting from discrediting the Neumann problem with the usual five-point operator on a regular mesh.The estimated condition number is 5.56E 017.For all of these examples, the matrix C is C sprand 1, n, d , where is a random, sparse matrix with approximately d • n • n uniformly distributed nonzero entries with d 0.5.We choose the matrices F and G completely satisfying the controllability and observability requirement of the pairs F, G and G, C .We apply the Hessenberg-Schur method and weighted block refinement method Algorithm 5.2 with 3 iterations for solving the Sylvester-observer equation XA − FX GC.The results are shown in Table 4.
It is also obvious from Table 4 that the performance of the weighted block refinement method is much better than that of the Hessenberg-Schur method, specifically for the illconditioned matrices.

Comments and Conclusion
In this paper, we propose two new iterative algorithm s for solving the large sparse Sylvesterobsever matrix equations.The existing projection methods use the Arnoldi process, but the methods described in this paper are based on the weighted block Arnoldi process.Moreover,

Theorem 4 . 1 .Figure 1 :
Figure 1: The a first and b second variables of the state x t and estimate x t for Example 3.8.
Theorem 3.1.If A, C is observable, then the states x t of the system 1.1 can be estimated by

Theorem 3.2. The
system 3.3 is a state observer of the system 1.1 , that is, z t is an estimate of Xx t in the sense that the error e t z t − Xx t → 0 as t → ∞ for any initial conditions x 0 , ∈ R n×n , G ∈ R n×r , and C ∈ R r×n , is called the Sylvester-observer equation.
Remark 3.7.According to Theorem 3.5 and Corollary 3.6, the controllability of F, G and observability of A, C guarantee the existence of a nonsingular solution of the Sylvesterobserver equation 1.2 in Step 1 of Algorithm 3.4.Moreover, there are other choices for F and G in Step 1 of Algorithm 3.4 provided F, G is controllable.Also, we can use Theorems 2.4 and 2.5 for analyzing the controllability of F, G and the observability of A, C .
then necessary and sufficient conditions for the unique solution X of 1.2 to be nonsingular are that F, G is controllable and A, C is observable.Proof.See 1 .
2.4, F, G is controllable.Thus, by Corollary 3.6, the nonsingular solution X of XA − FX GC is 1with the value of depending on the examples.The time is given in seconds for all examples.All numerical tests are performed in MATLAB software on a PC with 2.20 GHz with main memory 2 GB.For the first test we use two arbitrary matrices A and C. We choose the matrices F and G completely satisfying the controllability requirement of the pair G, C .Now We apply block refinement method for solving Sylvester-observer equation XA − FX GC with n 200.Also, we take 10 −6 .In Table1, we report the results for different values of m.In Table1, the results show that by increasing the values of m and l, the number of iterations decreases.The last column of Table 1 also shows the decreasing of time consumption.Note that the fourth and fifth columns of this table are the errors of the orthogonalization method.

Table 1 :
Implementation of block refinement method to solve the Sylvester equation with different values of m.

Table 2 :
Implementation of weighted block refinement method to solve the Sylvester-observer equation with different values of m.

Table 3 :
Implementation of new Iterative methods with 2 iterations and the Hessenberg-Schur method for solving the Sylvester equation.

Table 4 :
Effectiveness of the Hessenberg-Schur method and WBR algorithm with 3 iterations for randomly generated matrices.refinement process presented in Sections 4 and 5 has the capability of improving the results obtained by an arbitrary method.The numerical examples show the efficiency of the proposed schemes. the