# The LMI Approach for Stabilizing of Linear Stochastic Systems

Authors: Ivan Ivanov
Journal: International Journal of Stochastic Analysis
 (2013)

---

## Abstract

Stochastic linear systems subjected both to Markov jumps and to multiplicative white noise are considered. In order to stabilize such type of stochastic systems, the so-called set of generalized discrete-time algebraic Riccati equations has to be solved. The LMI approach for computing the stabilizing symmetric solution (which is in  fact the equilibrium point) of this system is studied. We construct a new modification of the standard LMI approach, and we show how to apply the new modification. Computer realizations of all modifications are compared. Numerical experiments are given where the LMI modifications are numerically compared. Based on the experiments the main conclusion is that the new LMI modification is faster than the standard LMI approach.


---

## Body

## 1. Introduction

In this paper we investigate general stochastic algebraic Riccati equations which are related to LQ control models for stochastic linear systems with multiplicative white noise and markovian jumping. We consider a stochastic system described by(1)dy(t)=(A0(η(t))y(t)+B0(η(t))u(t))dt+∑j=1k‍(Aj(η(t))y(t)+Bj(η(t))u(t))dWj(t)y(0)=y0,z(t)=C(η(t))x(t)+D(η(t))u(t),
whereW(t)=(W1(t),…,Wk(t))Tis a k-dimensional standard Brownian motion witht∈[0,+∞)andW(0)=0, defined by a filtered probability space(Ω,ℱ,ℱt,P). In addition, η(t), t≥0is a right continuous homogeneous Markov chain with the state space the set𝒟={1,…,N}and the probability transition matrixP(t)=(pij(t))=eQt, t≥0, andQ=(qij)with∑j=1k‍qij=0, i∈𝒟, and qij≥0if i=j. It is assumed that the W(t) andη(t)are independent stochastic processes andP{η(0)=i}>0for alli∈𝒟. The state vector x is an×1real vector, u denotes the vector of m control variables, and z is the regulated output vector with p components. The matrix coefficients Aj,Bj,0≤j≤k, C(i), D(i), i∈𝒟are constant matrices of appropriate dimensions with real elements.The stochastic systems with multiplicative white noise naturally arise in control problems of linear uncertain systems with stochastic uncertainty. It is important for applications to find a stabilizing controller for the above stochastic system (for more details see [1, 2]). For this purpose it is enough to compute the stabilizing solution to the following stochastic generalized Riccati algebraic equations:
(2)𝒯i(X)∶=A0(i)TX(i)+X(i)A0(i)+∑j=1k‍Aj(i)TX(i)Aj(i)+Q(i)+∑j=1N‍qijX(j)-S(i,X)T[R(i,X)]-1S(i,X)=0withR(i,X)>0,i∈𝒟,(3)R(i,X)=D(i)TD(i)+∑j=1k‍Bj(i)TX(i)Bj(i),Q(i)=CT(i)C(i),S(i,X)=B0(i)TX(i)+∑j=1k‍Bj(i)TX(i)Aj(i)+D(i)TC(i).The concepts of stabilizing solution of the Riccati-type equation (2) and of stabilizability for the triple(A,B,Q), where as usualA=(A0,A1,…,Ak),B=(B0,B1,…,bk), are defined in a standard way (see [1]).Applying Theorem 4.9 of Dragan and Morozan [2] we deduce that system (2) has a unique stabilizing solution(X~(1),X~(2),…,X~(N))with [D(i)TD(i)+∑j=1k‍Bj(i)TX~(i)Bj(i)]>0. The control
(4)u(t)=-[R(i,X~)]-1S(i,X~)x(t)
stabilizes system (1). An e1ective iterative convergent algorithm to compute these stabilizing solutions is presented in [3].Let us consider the special case of (1) where𝒟={1}. The stochastic linear quadratic model studied by Yao et al. in [4] is obtained. In this model the special functionalE∫0∞‍[x(t)TQx(t)+u(t)TRu(t)]dtis minimized (see equation (8) from [4]). The matricesQandRare so-called cost weighting matrices. It is allowed that the cost weighting matrices for the state and the control are singular, that is, Q=CTCandR=DTDare singular matrices. Such type models belong to a wide class of indefinite SLQ models. Applications of the indefinite LQ problems can be found in [5] for pollution control, [4, 6–9] in which problem appears in field of the mathematical finance. In this case the stochastic LQ problem can be solved via the following stochastic algebraic Riccati equation (the symmetric matrixPbeing the unknown one):
(5)𝒯(P):=ATP+PA+Q+∑j=1k‍CjTPCj-(BTP+∑j=1k‍DjTPCj)T×[R+∑j=1k‍DjTPDj]-1(BTP+∑j=1k‍DjTPCj)=0withR+∑j=1k‍DjTPDj>0.The coefficient matrices of (5)A, B, Q, R, Cj, andDjforj=1,…,kare given ones of sizesn×n,n×p,n×n,m×m,n×n, andn×m, respectively, andQdenotes the state matrix, andRis the control matrix. For a deterministic case (Cj=Dj=0,j=1…,k) it is assumed that the matrixQ is positive semidefinite andRis a positive definite one.Here the computation of stabilizing solution of system (1) is explicitly expressed in terms of the solutions of some linear matrix inequalities (LMI). The paper is devoted to the LMI approach and its modifications. The LMI approach is very important for the practice and real-world problems. Very often the LMI approach is an only method for solving a given class of problems. The application of the LMI approach to the solution of the optimal control problems is studied in [10–14]. We introduce a modification set of nonlinear equations equivalent to (1) which lead us to the new convex optimization problems. The LMI approach applied to the new optimization problem gives a fast way to find the stabilizing solution to (1). We will compare the numerical effectiveness of the introduced LMI solvers. Numerical simulations are used to demonstrate the performance of the considered solvers.The notations used in this note are standard. Here,X≥0 (orX>0); it is denoted thatX=XTis positive semidefinite (or positive definite) and∥·∥2denotes the spectral matrix norm.


## 2. The LMI Approach to the Generalized Stochastic Riccati Equations

Rami and Zhou [15] have investigated the stochastic algebraic Riccati equation (5) in casek=1, where the numerical method to compute the maximal solution to (5) withk=1is derived. The method is based on the solution of a convex optimization problem over linear matrix inequalities. The LMI approach is considered as powerful tool in optimization. Further on, the authors in [4] have developed a computational approach to such SLQ models(k>1)using an LMI formulation. The LMI optimization is a successful method for solving (5) in this case. Although we cannot solve the stochastic algebraic Riccati equation we can still find the optimal control law via the LMI approach. It is well known [16] that if the stochastic system is stabilizable, in the mean square sense, the LMI optimization method always yields the maximal positive semidefinite solution to (5). The problem of stability and optimality of SLQ model (5) withk=1(the case of one-dimensional Brownian motion) is treated by [7]. Moreover, this study is continued in terms of the multidimensional model (k>1) by [4]. The optimization problem associated with (5) is
(6)max〈In,P〉×(Q+ATP+PA+∑j=1k‍CjTPCj(BTP+∑j=1k‍DjTPCj)TBTP+∑j=1k‍DjTPCjR+∑j=1k‍DjTPDj)≥0R+∑j=1k‍DjTPDj>0P=PT,
with respect to the variablePwhich is a symmetric matrix.Under notations(7)W(i,X)=(A0(i)+0.5qiiI)TX(i)+X(i)(A0(i)+0.5qiiI)=A~0(i)TX(i)+X(i)A~0(i),Π1(i,X)=∑j=1k‍Aj(i)TX(i)Aj(i)+∑j=i‍qijX(j),
the optimization problem associated with the set of equations (2) is
(8)max∑i=1N‍〈In,X(i)〉subjecttoi∈𝒟×(W(i,X)+Π1(i,X)+C(i)TC(i)S(i,X)TS(i,X)R(i,X))≥0R(i,X)>0X(i)=X(i)T.In this paper we investigate the numerical solvability of the semidefinite programming problem (8) for different types of matricesD(i)TD(i),  i∈𝒟. However, the numerical experiments for finding the maximal solution of (2) show that the LMI method (8) is slowly working for different types of matricesD(i)TD(i)in the case k=1. Here we introduce a new modification to accelerate the LMI method for solving the optimization problem (8). In many applications of control system theory the following fact is exploited. If any matrixR=RTis singular or zero and the matrixBhas the full rank, then there exists a symmetric matrix Z such thatR+BTZBis a positive definite one. In our investigation the matricesB(i)=(B1(i)B2(i)⋮Bk(i))have the full rank. Thus, we replaceX(i)=Z(i)+Y(i)inR(i,X). It is obtained
(9)R(i,Z+Y)=D(i)TD(i)+∑j=1k‍Bj(i)T(Z(i)+Y(i))Bj(i)=R~(i)+∑j=1k‍Bj(i)TY(i)Bj(i),
whereY(i)is the new unknown matrix. We apply this conclusion to set of equations (2). We construct the matricesZ(1),…,Z(N)such that all matrices
(10)R~(i)=D(i)TD(i)+∑j=1k‍Bj(i)TZ(i)Bj(i)
are positive definite ones fori∈𝒟. Then, the following set of Riccati equations is obtained regarding Y(1),…,Y(N):
(11)𝒲i(Y)∶=W(i,Y)+Q~(i)-S~(i,Y)T[R~(i,Y)]-1×S~(i,Y)=0,i∈𝒟,
where
(12)Q~(i)=C(i)TC(i)+A~0(i)TZ(i)+Z(i)A~0(i)+∑j=1k‍Aj(i)TZ(i)Aj(i)+∑j=i‍qijZ(j)R~(i,Y)=R~(i)+∑j=1k‍Bj(i)TY(i)Bj(i)S~(i,Y)=B0(i)TY(i)+∑j=1k‍Bj(i)TY(i)Aj(i)+L~(i),L~(i)=B0(i)TZ(i)+∑j=1k‍Bj(i)TZ(i)Aj(i)+D(i)TC(i).
Thus, the new optimization problem over LMI conditions related to (11) is derived:
(13)max∑i=1N‍〈In,Y(i)〉subject  toi∈𝒟×(W(i,Y)+Π1(i,Y)+Q~(i)S~(i,Y)TS~(i,Y)R~(i,Y))≥0R~(i,Y)>0Y(i)=Y(i)T.Our experience on computations with LMI approach shows that there are examples whereY(1),…,Y(N)are not positive semidefinite, but the matricesR~(i,Y)are still positive definite. In this reason we can consider a practical implementation of (13) where the constrainR~(i,Y)>0is omitted. This leads us to the following problem:
(14)max∑i=1N〈In,Y(i)〉‍subject  toi∈𝒟×(W(i,Y)+Π1(i,Y)+Q~(i)S~(i,Y)TS~(i,Y)R~(i,Y))≥0Y(i)=Y(i)T.


## 3. Numerical Experiments

We carry out numerical simulations to present the numerical behaviour of introduced methods. In our experiments we apply three semidefinite programming problems (8), (13), and (14) for solving the stochastic algebraic Riccati equations (2). Our experiments are executed in MATLAB on an 2,16 GHz Intel(R) Dual CPU computer. The solutions of above optimization problems are obtained under the MATLAB lmi solvers which are executed with relative accuracy tol=1.0e-9.For all examples we takeQ(1)=Q(2)=Q(3)to be diagonal matrices with entries(1,1,…,1,0,0). We have executed a set of examples with different values ofnand constant weighting matricesR(1)=R(2)=R(3)=zeros(3,3). We compare all iterations introducing the following parameters: “m It”—the biggest number of iterations, “av It”— the average number of iterations. To determine the numbers “m It” and “av It” we count those examples of each size for which the corresponding iteration converges.We consider a family of examples in caseN=3, k=2, and n=8,…,12,15,20, where the coefficient real matrices are given as follows:A0(i), A1(i), A2(i), B0(i), B1(i), B2(i), L(i), i=1,2,3 were constructed using the MATLAB notations
(15)A0(1)=randn(n,n)10-eye(n,n);A0(2)=randn(n,n)10-eye(n,n);A0(3)=randn(n,n)10-eye(n,n);A1(1)=randn(n,n)10;A1(2)=randn(n,n)10;A1(3)=randn(n,n)10;A2(1)=randn(n,n)10;A2(2)=randn(n,n)10;A2(3)=randn(n,n)10;B0(1)=2*rand(n,3);B0(2)=2*rand(n,3);B0(3)=2*rand(n,3);B1(1)=randn(n,3)100;B1(2)=randn(n,3)100;B1(3)=randn(n,3)100;B2(1)=randn(n,3)100;B2(2)=randn(n,3)100;B2(3)=randn(n,3)100;L(1)=zeros(n,3);L(2)=zeros(n,3);L(3)=zeros(n,3).In our definitions the functions randn(p,k) and rand(p,k) return a p-by-k matrix of pseudorandom scalar values (for more information see the MATLAB description). The following transition probability matrix (see [3])
(16)(qij)=(0.3330.60.60.3330.30.10.3330.10.3)
is applied for all examples.For our purpose we have executed hundred examples of each value ofn for the test. The maximal number of iterations “m It” and average number of iterations “av It” of each size for all examples needed for achieving the relative accuracy are reported in the table. There are three columns where the maximal errors ET=max1,…,100maxi=1,2,3∥𝒯i(X~)∥2 and EW=max1,…,100maxi=1,2,3∥𝒲i(Y~)∥2 for each n for the test are presented. HereX~=(X(1),X(2),X(3))is a computed solution to (2) via (8), andY~=(Y(1),Y(2),Y(3))is a computed solution to (11) via (13) and (14). In addition, the time of execution for each method in casesn=15andn=20is reported. Results from experiments are given in Table 1.Table 1
Results from 100 runs for each value ofn. Z(1)=Z(2)=Z(3)=0.6In.
n
LMI: (8)
LMI: (13)
LMI: (14)
m It
av It
E
T
m It
av It
E
W
m It
av It
E
W
8
62
31.0
3.1274
e
-
08
62
31.5
6.9395
e
-
08
42
32.5
6.9395
e
-
08
9
54
34.5
2.7652
e
-
08
56
35.5
1.6774
e
-
08
41
31.3
6.7014
e
-
09
10
49
32.7
1.6131
e
-
08
51
34.2
2.0807
e
-
09
40
32.7
1.1040
e
-
09
11
57
34.5
1.6131
e
-
08
59
36.7
1.4413
e
-
09
46
35.7
8.7308
e
-
10
12
51
33.3
2.2023
e
-
08
53
34.9
4.2220
e
-
09
46
35.7
1.0796
e
-
09
15
47
37.7
2.4223
e
-
08
48
39.1
6.0639
e
-
09
47
34.5
6.0639
e
-
09
20
46
35.3
6.1281
e
-
08
47
35.5
4.9302
e
-
08
36
33.1
4.6979
e
-
08
CPU time for executing for 10 runs
15
59.8 s
61.2 s
39.5 s
20
170.19 s
172.6 s
117.5 s


## 4. Conclusion

We have made numerical experiments for computing this solution, and we have compared the numerical results. Our numerical experiments confirm the effectiveness of the proposed new transformations which lead us to the equivalent semidefinite programming problem. We have compared the results from the experiments in regard to number of iterations and CPU time for executing the above optimization problems forn=15, n=20. The solution of the optimization problems achieves the same accuracy for different number of iterations. The executed examples have demonstrated that the LMI problem performances (8) and (13) require the same average numbers of iterations (see the corresponding columns “av It” for all tests). In addition, the LMI performance (14) requires less number of iterations than the remaining approaches, and it is faster than the others—see the CPU time for execution from Table 1. This property of (14) follows from the structure of problem (14)—there is one inequality less than (13). However, this property is not valid in all cases.



---

*Source: 281473-2013-08-29.xml*