Iterations for a General Class of Discrete-Time Riccati-Type Equations: A Survey and Comparison

The discrete-time linear control systems have been applied in the wide area of applications such as engineering, economics, biology. Such type systems have been intensively considered in the control literature in both the deterministic and the stochastic framework. The stability and optimal control of stochastic differential equations withMarkovian switching has recently received a lot of attention, see Freiling and Hochhaus [8], Costa, Fragoso, and Marques [2], Dragan and Morozan [4, 5]. The equilibrium in these discrete-time stochastic systems can be found via the maximal solution of the corresponding set of discrete-time Riccati equations.

The matrix F involved in the above definition is called stabilizing feedback gain.
We will investigate the computation of the maximal solution of a set of equations (1).A solution X of (1) is called maximal if X ≥ X for any solution X.
We will consider three cases.In the first case the weighting matrices R(i)=D T (i) D(i), i = 1, . . ., N are assumed to be positive definite and Q(i)=C T (i) C(i), i = 1, . . ., N are positive semidefinite.Thus, the matrices R(i, X)=R(i)+∑ r l=0 B l (i) T E i (X)B l (i), i = 1,...,N are positive definite.We present an overview of several computational methods [6,10] to compute the maximal and stabilizing solutions of a considered class of discrete-time Riccati equations.In addition, we apply a new approach, where the variables are changed and an equivalent set of the same type Riccati equations is obtained.A new iteration for the maximal solution to this equivalent set of nonlinear equations is proposed.This is the subject of section 2.
We investigate the applicability of the existing methods for the maximal solution to (1) where the weighting matrices R(i), Q(i) are indefinite in the second case.These weighting matrices are indefinite, but the matrices R(i, X), i = 1, . . ., N are still positive definite.Similar investigations have been executed by Rami and Zhou [12,14] in case in the infinite time horizon.The important tool for finding the maximal solution is a semidefinite programming associated with the linear matrix inequalities (LMI).The method for the maximal solution used an LMI optimization problem is called the LMI approach or the LMI method.Rami and Zhou [14] have described a technics for applying the LMI method to indefinite linear quadratic problem in infinite time horizon.Here we will extend their findings and will modify their technics to indefinite linear quadratic problem of Markovian jump linear systems.We propose a new optimization problem suitable to the occasion.The investigation is accompanied by comparisons of the LMI approach on different numerical examples.This is the subject of section 3.
The third case is considered in section 4. Here, the solution of (1) under the assumption that at least one of matrices R(i, X), i = 1,...,N is positive semidefinite is analysed.In this case set of equations ( 1) can be written as for i = 1,...,N with the additional conditions R(i, X) ≥ 0, and Such type generalized Riccati equation is introduced in [15].The notation Z † stands for the Moore-Penrose inverse of a matrix Z.We derive a suitable iteration formula for computing the maximal solution of (3) -( 4) and the convergence properties of the induced matrix sequence are proved.In addition, the LMI approach is modified and applied to the case of semidefinite matrices R(i, X)(i = 1, . . ., N). Numerical simulations for comparison the derived methods are presented in the section.
We are executing some numerical experiments in this investigation.Based on the results from experiments the considered methods are compared in all cases.In the examples we consider a MJLS with three operation modes describing an economic system, adapted from [17] which studies a time-variant macroeconomic model where some of the parameters are allowed to fluctuate in an exogenous form, according to a Markov chain.The operation modes are interpreted as the general situation: "neutral", "bad" or "good" (N = 3).See [17] and references therein for more details.Our experiments are carried out in the MATLAB on a 1,7GHz PENTIUM computer.In order to execute our experiments the suitable MATLAB procedures are used.
The notation H n stands for the linear space of symmetric matrices of size n over the field of real numbers.For any X, Y ∈H n , we write ,...,X(N)) ∈H n and the inequality Y ≥ Z mean that for i = 1,...,N, X(i) ∈H n and Y(i) ≥ Z(i), respectively.The 149 Iterations for a General Class of Discrete-Time Riccati-Type Equations: A Survey and Comparison linear space H n is a Hilbert space with the Frobenius inner product < X, Y >= trace(XY).Let .denote the spectral matrix norm.

The positive definite case
Let us assume that the weighting matrices R(i), i = 1,...,N are positive definite and Q(i), i = 1, . . ., N are positive semidefinite.Thus the matrices R(i, X), i = 1,...,N are positive definite.In this section, we consider set of equations (1) where a matrix X belongs to the domain: Note that X ∈ Dom P implies Y ∈ Dom P for all Y ≥ X and that Dom P is open and convex.We consider the map Dom P→H n .We investigate some iterations for finding the maximal solution to (1).For the matrix function P (i, X) we introduce notations and we present set of equations (1) as follows: Then, for the matrix function P (i, X) we rewrite We will study the system X(i)=P (i, X) for i = 1, . . ., N. We start by some useful properties to P (i, X).For briefly we use Ãl (i, Z)=A l (i)+B l (i)F(i, Z), l = 0, 1, . . ., r for some Z ∈H n .Lemma 2.1.[10] Assuming Y ∈H n and Z ∈H n are symmetric matrices, then the following properties of P (i, X) , i = 1, . . ., N Dragan, Morozan and Stoica [6] have been proposed an iterative procedure for computing the maximal solution of a set of nonlinear equations (1).The proposed iteration [6, iteration (4.7)] is : k = 1, 2, 3 . . ., and ε is a small positive number.Note that iteration ( 5) is a special case of the general iterative method given in [6,Theorem 3.3].Based on the Gauss-Seidel technique the following modification is observed by Ivanov [10]: where The convergence properties of matrix sequences defined by ( 5) and ( 6) are derived in the corresponding papers.
In this section, we are proving that iteration (5) has a linear rate of convergence.Theorem 2.2.Assume that conditions a),b),c) of theorem 2.1 are fulfilled for a symmetric solution X ∈ Dom P of set of equations (1).Then, the sequence {X (k) } ∞ k=1 defined by ( 5) converges to the maximal solution Proof.Following the course of the proof of theorem 2.1 it has been proved that X (k) ≥ X for all k.Therefore, for i = 1, . . ., N we conclude R(i, and then, the limit F(i, X)=lim k→∞ F(i, X (k) ) exists and Based on the proof of theorem 2.1 and the properties of lemma 2.1 the following equality is established Consider the difference Then Thus where In addition, using X (0) ≥ X (k) we in fact have 153 Iterations for a General Class of Discrete-Time Riccati-Type Equations: A Survey and Comparison Furthermore, we estimate for i = 1, . . ., N Note that Now, assuming that the inequality max 1≤i≤N Thus, the proof of the theorem is complete.
Let us consider the following example in order to compare iterations ( 5) and ( 6).

The LMI approach
There exists an increasing interest to consider a computational approach to stochastic algebraic Riccati equations via a semidefinite programming problem over linear matrix inequalities.Similar studies can be found in [12][13][14].The main result from such type studies is the equivalence between the feasibility of the LMI and the solvability of the corresponding stochastic Riccati equation.Moreover, the maximal solution of a given stochastic algebraic Riccati equation can be obtained by solving a corresponding convex optimization problem (an LMI approach).
Further on, following the classical linear quadratic theory [13,14] we know that the optimization problem is associated with (1) has the form (for example see [1,7]): However, we can apply the same approach to equivalent system (7).As a result we formulate a new optimization problem assigned to (7) and we will use it to find the maximal solution to (7).
The corresponding optimization problem, associated to the maximal solution to (7), is given by: max ∑ N i=1 I, Y(i) It is very important to analyze a case where the weighting matrices R(i), Q(i), i = 1,...,N are indefinite in the field of linear quadratic stochastic models.This case has a practical importance.There are studies where the cost matrices are allowed to be indefinite (see [3,16] and reference there in).In this paper we will investigate this special case to considered general discrete-time Riccati equations (1).We will interpret iterations ( 6) and ( 9) in a case where matrices R(i), i = 1,...,N are indefinite and however, we will look for a maximal solution from Dom P. Based on the next example, we compare the LMI approach (through optimization problems (10) and ( 11)) for solving the maximal solution to set of nonlinear equations (1).
We take the n × n matrices Q(1), Q(2) and Q(3) as follows: Example 3.1.The coefficient matrices A 0 (i), A 1 (i), B 0 (i), B 1 (i), L(i), i = 1, 2, 3 for system (1) are given through formulas (using the MATLAB notations):  The MATLAB function mincx is applied with the relative accuracy equals to 1.e − 10 for solving the corresponding optimization problems.Our numerical experiments confirm the effectiveness of the LMI approach applied to the optimization problems (10) and (11).We have compared the results from these experiments in regard of number of iterations and time

The positive semidefinite case
We will investigate new iterations for computing the maximal solution to a set of Riccati equations ( 1) where the matrices R(i, X), i = 1, . . ., N are positive semidefinite.It is well known the application of a special linear quadratic stochastic model in the finance [18] where the cost matrix R is zero and the corresponding matrix R + B T XB is singular.So, this special case where it is necessary to invert a singular matrix is important to the financial modelling process.Without loss of generality we assume that all matrices R(i, X), i = 1, . . ., N in (1) are positive semidefinite.Thus, we will investigate set of equations ( 3)-( 4) for existence a maximal solution.Investigations on similar type of generalized equations have been done by many authors (see [8,9] and literature therein).
We introduce the following new iteration: We will prove that the matrix sequence defined by (12) converges to the maximal solution X of (3)-( 4) and R(i, X), i = 1, . . ., N are positive semidefinite.Thus, the iteration (12) constructs a convergent matrix sequence.
159 Iterations for a General Class of Discrete-Time Riccati-Type Equations: A Survey and Comparison and Combining ( 13) and ( 14) we write down for X ∈H n and H ∈H n .Obviously W X (i, H) ≥ 0 and W X (i, X)=0.
Proof.Let us consider the difference 161 Iterations for a General Class of Discrete-Time Riccati-Type Equations: A Survey and Comparison Thus, it is received a nonincreasing matrix sequence {X i } ∞ i=1 of symmetric matrices bounded below by X which converges to X and X ≥ X.Hence, X ∈ dom P † by Lemma 4.3.
The theorem is proved.Thus, following the above theorem we could compute the maximal solution X of the set of equations ( 3)-( 4).We should apply iteration (12).The next question is: How to apply the LMI method in this case?
Let us consider the modified optimization problem: Theorem 4.2.Assume that (A, B) is stabilizable and there exists a solution to the inequalities P (i, X) − X(i) ≥ 0 for i = 1,...,N.Then there exists a maximal solution X + of ( 3)-( 4) if and only if there exists a solution X for the above convex programming problem ( 16) with X + ≡ X.
Proof.Note that the matrix X =( X(1),...,X(N)) satisfies the restrictions of optimization problem (16) if and only if The last statement follows immediately by lemma 4.2.
Let us consider set of equations ( 7)-( 8) under the assumption that R(i)+∑ r l=0 B l (i) T Y(i)B l (i) ≥ 0, i = 1, . . ., N. Thus, optimization problem (11) is transformed to the new optimization problem: max ∑ N i=1 I, Y(i) It is easy to verify that the solution of ( 18) is the maximal solution to ( 7)-( 8) with the positive semidefinite assumption to matrices R(i)+∑ r l=0 B l (i) T Y(i)B l (i) ≥ 0, i = 1,...,N.We investigate the numerical behavior of the LMI approach applied to the described optimization problems LMI: (16) and LMI(Y): (18) for finding the maximal solution to set of discrete-time generalized Riccati equations ( 3) -( 4).In addition, we compare these LMI solvers with derived recurrence equations (12) for the maximal solution to the same set of equations.We will carry out some experiments for this purpose.Example 4.1.We consider the case of r = 1, n = 6, 7, where the coefficient real matrices A 0 (i), A 1 (i), A 2 (i), B 0 (i), B 1 (i), B 2 (i), L(i), i = 1, 2, 3 are given as follows (using the MATLAB notations): In our definitions the functions randn(p,k) and sprand(q,m,0.3)return a p-by-k matrix of pseudorandom scalar values and a q-by-m sparse matrix respectively (for more information see the MATLAB description).
Results from experiments are given in table 4. The parameters mIt and avIt are the same as the previous tables.In addition, the CPU time in seconds is included.The optimization problems ( 16) and ( 18) need the equals iteration steps (the column avIt) for finding the maximal solution to set of equations ( 3) -( 4).However, the executed examples have demonstrated that LMI problem (18) faster than LMI problem (16).Moreover, iterative method (12)   This choice of the matrices B 0 (i), B 1 (i), i = 1, 2, 3 guaranteed that the matrices R(i, X) are positive semidefinite, i.e. there are symmetric matrices X which belongs to dom P † .The remain coefficient matrices are already in place.
We find the maximal solution to (3) -(4) for the constructed example with iterative method (12) and the LMI approach applied to optimization problems (16) and (18).The results are the following.Iteration (12)   The LMI approach for optimization problem (16) does not give a satisfactory result.After 32 iteration step the calculations stop with the computed maximal solution V.However, the norm of the difference between two solutions W and V is W(1) − V(1) = 1.0122e − 9, W(2) − V(2) = 2.8657e − 6, W(3) − V(3) = 3.632e − 6.
The LMI approach for optimization problem (18) needs 28 iteration steps to compute the maximal solution to (1).This solution Z has the same eigenvalues as in (19).The norm of the difference between two solutions W and Z is W The results from this example show that the LMI approach applied to optimization problem (18) gives the more accurate results than the LMI method for (16).Moreover, the results obtained for problem (16) are not applicability.A researcher has to be careful when applied the LMI approach for solving a set of general discrete time equations in positive semidefinite case.

Conclusion
This chapter presents a survey on the methods for numerical computation the maximal solution of a wide class of coupled Riccati-type equations arising in the optimal control of discrete-time stochastic systems corrupted with state-dependent noise and with Markovian jumping.In addition, computational procedures to compute this solution for a set of discrete-time generalized Riccati equations ( 7)-( 8) are derived.Moreover, the LMI solvers for this case are implemented and numerical simulations are executed.The results are compared and the usefulness of the proposed solvers are commented.

Table 1 .
(6)have executed hundred examples of each value of n.We report the maximal number of iteration steps mIt and the average number of iteration steps avIt of each size for all examples needed for achieving the accuracy.The used accuracy equals to 1.e − 10.The results are listed in Table1.The average number of iteration steps for method (6) smaller than the corresponding average number for method(5).The last column of Table1shows how much is the acceleration of method(6).Results for Example 2.1.Comparison between iterations for 100 runs.

Table 2 .
Comparison between iterations for Example 3.1.

Table 3 .
Comparison between iterations for Example 3.1.
(10)xecution.The executed four tests of examples have demonstrated that LMI problem(10)performance needs more computational work than LMI problem(11)and thus, LMI method(11)is faster than LMI method(10)(see the CPU times displayed in tables 2 and 3).

Table 4 .
is much faster than the LMI approaches and it achieves the same accuracy.Comparison between methods for the maximal solution in Example 4.1.
We introduce an additional example where the above optimization problems are compared.Example 4.2.The parameters of this system are presented as follows.The coefficient matrices are (r = 1, k = 3, n = 6):