# The Relaxed Gradient Based Iterative Algorithm for the Symmetric (Skew Symmetric) Solution of the Sylvester EquationAX+XB=C

Authors: Xiaodan Zhang; Xingping Sheng
Journal: Mathematical Problems in Engineering
 (2017)

---

## Abstract

In this paper, we present two different relaxed gradient based iterative (RGI) algorithms for solving the symmetric and skew symmetric solution of the Sylvester matrix equationAX+XB=C. By using these two iterative methods, it is proved that the iterative solution converges to its true symmetric (skew symmetric) solution under some appropriate assumptions when any initial symmetric (skew symmetric) matrix X0 is taken. Finally, two numerical examples are given to illustrate that the introduced iterative algorithms are more efficient.


---

## Body

## 1. Introduction

For the convenience of our statements, the following notations will be used throughout the paper:Rm×n represents the set of m×n real matrices. For A∈Rm×n, we write AT, R(A), tr(A), ρ(A), λmax(A), λmin(A), A2, and AF to denote the transpose, the range space, the trace, the spectral radius, the maximal eigenvalue, the minimal eigenvalue, the spectral norm, and the Frobenius norm of a matrix A, respectively; that is, A2=λmax(ATA) and AF=[tr(ATA)]1/2. σmax(A), σmin(A) are the maximal singular value and the minimal nonzero singular value of A. The symbol of In represents an identity matrix of order n, 1m×n is an m×n matrix whose elements are 1, and cond(A)=σmax(A)/σmin(A) is the condition number of matrix A. The inner product in space Rm×n is defined as A,B=tr(ATB); particularly, AF2=tr(ATA)=tr(AAT)=A,A. A⊗B denotes the Kronecker product, defined as A⊗B=(aijB), i=1,…,m, j=1,…,n. For any matrix X=(x1,x2,…,xn)∈Rm×n, the vector operator is defined as vec(X)=(x1T,x2T,…,xnT)T∈Rmn. Using vector operator and Kronecker product, we have vec(AXB)=(BT⊗A)vec(X).Considering the symmetric (skew symmetric) solution of the Sylvester matrix equation(1)AX+XB=C,where A∈Rn×n, B∈Rn×n, C∈Rn×n, and the unknown matrix X∈Rn×n.The Sylvester matrix equation (1) has many applications in linear system theory, for example, pole/eigenstructure assignment [1–4], robust pole assignment [5–8], robust partial pole assignment [9], observer design [10], model matching problem [11], regularization of descriptor systems [12, 13], disturbance decoupling problem [14], and noninteracting control [15].As is well known, (1) has a unique solution if and only if A and -B possess no common eigenvalues [16] and the solution can be computed by solving a linear system (I⊗A+BT⊗I)vec(X)=vec(C). Using this method, it will increase the computational cost and storage requirements, so that this approach is only applicable for small sized Sylvester equations.Due to these drawbacks, many other methods for the solution have appeared in the literature. The idea of transforming the coefficient matrix into a Schur or Hessenberg form to compute (1) have been presented in [16, 17]. When the linear matrix (1) is inconsistent, a finite iterative method to solving its Hermitian minimum norm solutions has been presented in [18]. An efficient iterative method based on Hermitian and skew Hermitian splitting has been proposed in [19]. Krylov subspace based methods have been presented in [20–26] for solving Sylvester equations and generalized Sylvester equations. Recently based on the idea of a hierarchical identification principle [27–29], some efficient gradient based iterative algorithms for solving generalized Sylvester equations and coupled (general coupled) Sylvester equations have been proposed in  [27, 30–32].  Particularly, for Sylvester equations of form (1), it is illustrated in [33] that the unknown matrix to be identified can be computed by a gradient based iterative algorithm. The convergence properties of the methods are investigated in  [27, 32].  Niu et al. [34] proposed a relaxed gradient based iterative algorithm for solving Sylvester equations. Wang et al. [35] proposed a modified gradient based iterative algorithm for solving Sylvester equations (1). More recently, Xie and Ma [36] gave an accelerated gradient based iterative algorithm for solving (1). In [37, 38] Xie et al. studied the special structure solution of matrix equation (1) by using iterative method.In this paper, inspired by [28, 34–36], we first derive a relaxed gradient based iterative (RGI) algorithm for solving the symmetric solution of matrix equation (1). Theoretical analysis shows that our method converges to the exact symmetric solution for any initial value with some appropriate assumptions. Then the proposed algorithm can be also applied to the skew symmetric solution of matrix equation (1). Numerical results illustrate that the proposed method is correct and feasible. We must point out that the ideas in this paper have some differences comparing with that in [28, 34–36].The rest of the paper is organized as follows. In Section2, some main preliminaries are provided. In Section 3, the relaxed gradient based iterative methods are studied. Finally, a numerical example is included to verify the superior convergence for the algorithms.


## 2. Preliminaries

In this section, we reviewed the ideas and principles of the gradient based iterative (GI) method, the relaxed gradient based iterative (RGI) method, and the modified gradient based iterative (MGI) method.Letμ>0 be the convergence factor or step factor. The gradient based iterative method for Ax=b is as follows: (2)xk=xk-1-μATAxk-1-b.The convergence of the gradient based iterative is stated as follows.Lemma 1 (see [32]).
Assume matrixA is full column-rank and 0<μ<2/A22; then the gradient based iterative sequences {x(k)} in (2) converge to x∗; that is, limk→∞x(k)=x∗ or the error x~(k)=x(k)-x∗ converges to zero for any initial value x(0). Moreover, the maximal 2-convergence rate is given by γmax=1/βμopt, where(3)βμopt=cond2A-1cond2A+1,μopt=2σmax2A+σmin2A.In this case, the error vector x~(k) satisfies (4)x~k≤βkμoptx~0,k≥0.In [28], Ding and Chen presented the following algorithm based on gradient for solving (1).Algorithm 2 (see [28] (the gradient based iterative (GI) algorithm)).
Step 1. Input matrices A,B,C∈Rn×n, given any small positive number ϵ. Choose the initial matrices X1(0) and X2(0). Compute X(0)=(X1(0)+X2(0))/2. Set k≔1.
Step 2. If δk=AX(k)+X(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (5)X1k=Xk-1-μATAXk-1+Xk-1B-C,X2k=Xk-1-μAXk-1+Xk-1B-CBT,Xk=X1k+X2k2.
Step 4. Set k≔k+1; return to Step 2.The authors [28] also pointed out that if the convergence factor μ is chosen in (0,2/(σ2A+σ2B)), Algorithm 2 will converge to the exact solution of (1).Niu et al. [34] gave a relaxed gradient based iterative algorithm for solving (1). When μ is in 0,2/ω1-ωσ2A+σ2B+σBAT the following algorithm has been proven to be convergent.Algorithm 3 (see [34] (the relaxed gradient based iterative (RGI) algorithm)).
Step 1. Input matrices A,B,C∈Rn×n, given any small positive number ϵ and appropriate positive number ω. Choose the initial matrices X1(0) and X2(0). Compute X(0)=ωX1(0)+(1-ω)X2(0). Set k ≔ 1.
Step 2. If δk=AX(k)+X(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (6)X1k=Xk-1-1-ωμATAXk-1+Xk-1B-C,X2k=Xk-1-ωμAXk-1+Xk-1B-CBT,Xk=ωX1k+1-ωX2k.
Step 4. Set k≔k+1; return to Step 2.Recently, in [35] Wang et al. proposed a modified gradient based iterative (MGI) algorithm to solve (1). The main difference is that, in the step of computing X2(k), the last approximate solution X1(k) has been considered fully to update X(k-1).Algorithm 4 (see [35] (the modified gradient based iterative (MGI) algorithm)).
Step 1. Input matrix A,B,C∈Rn×n, given any small positive number ϵ and appropriate positive number ω. Choose the initial matrices X1(0) and X2(0). Compute X(0)=(X1(0)+X2(0))/2. Set k≔1.
Step 2. If δk=AX(k)+X(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (7)X1k=Xk-1-μATAXk-1+Xk-1B-C.
Step 4. Compute (8)Xk-1=X1k+X2k-12.
Step 5. Update the sequences (9)X2k=Xk-1-μAXk-1+Xk-1B-CBT.
Step 6. Compute(10)Xk=X1k+X2k2.
Step 7. Set k≔k+1; return to Step 2.More recently, Xie and Ma [36] presented the following AGBI algorithm for solving (1) based on the idea of MGI.Algorithm 5 (see [36] (the accelerated gradient based iterative (AGBI) algorithm)).
Step 1. Input matrix A,B,C∈Rn×n, given any small positive number ϵ and appropriate positive number ω. Choose the initial matrices X1(0) and X2(0). Compute X(0)=ωX1(0)+(1-ω)X2(0). Set k≔1.
Step 2. If δk=AX(k)+X(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (11)X1k=Xk-1-1-ωμATAXk-1+Xk-1B-C.
Step 4. Compute (12)Xk-1=ωX1k+1-ωX2k-1.
Step 5. Update the sequences (13)X2k=Xk-1-ωμAXk-1+Xk-1B-CBT.
Step 6. Compute (14)Xk=ωX1k+1-ωX2k.
Step 7. Set k≔k+1; return to Step 2.


## 3. Main Results

In this section, we first study the necessary and sufficient conditions of the symmetric solution for (1). Then the relaxed gradient based iterative algorithm for the symmetric solution of equation (1) is proposed. Following the same line, the relaxed gradient based iterative algorithm for the skew symmetric solution of equation (1) is also presented.Theorem 6.
The matrix equation (1) has a unique symmetric solution Xs if and only if the following pair of the matrix equations (15)AX+XB=C,BTX+XAT=CThas a unique common solution X∗ and Xs=(X∗+(X∗)T)/2.Proof.
IfXs is a unique symmetric solution of (1), then (Xs)T=Xs and AXs+XsB=C; further we have(16)BTXs+XsAT=AXs+XsBT=CT.This shows that Xs is also the solution of the pair matrix equations (15).
Conversely, if the system of matrix equations (15) has a common solution X∗, let us denote Xs=(X∗+(X∗)T)/2; then we can check that (17)AXs+XsB=AX∗+X∗T2+X∗+X∗T2B=AX∗+X∗B2+AX∗T+X∗TB2=AX∗+X∗B2+BTX∗+X∗ATT2=C+CTT2=C.This implies that Xs is the unique symmetric solution of (1).According to Theorem6, if the unique common solution X∗ of equations (15) can be obtained, then the unique symmetric solution of (1) is Xs=(X∗+(X∗)T)/2.According to Theorem6, we construct a relaxed gradient based iterative algorithm to solve the symmetric solution of (1).Algorithm 7 (the relaxed gradient based iterative (RGI) algorithm for symmetric solution of (1)).
Step 1. Input matrices A,B,C∈Rn×n, given any small positive number ϵ and appropriative positive number ω such that 0<ω<1. Choose any initial matrix X(0). Set k≔1.
Step 2. If δk=AXs(k)+Xs(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (18)X1k=Xk-1-μ1ATAXk-1+Xk-1B-C+AXk-1+Xk-1B-CBT,X2k=Xk-1-μ2BBTXk-1+Xk-1AT-CT+BTXk-1+Xk-1AT-CTA,Xk=ωX1k+1-ωX2k,Xsk=Xk+XTk2.
Step 4. Set k≔k+1; return to Step 2.In the following paragraph, we will investigate the convergence of Algorithm7.Theorem 8.
Assume the matrix equations (15) have a unique solution X∗; then the iterative sequence {X(k)} generated by the Algorithm 7 converges to X∗, if (19)0<μ1<1ωA22+B22;0<μ2<11-ωA22+B22that is, limk→∞X(k)=X∗ or the error X(k)-X∗ converges to zero for any initial value X(0).Further the sequence{Xs(k)} converges to Xs, where Xs is the unique symmetric solution of (1).Proof.
Define the error matrix(20)X~k=Xk-X∗,θ1k-1≔AX~k-1+X~k-1B,θ2k-1≔BTX~k-1+X~k-1AT.We have (21)X~k=Xk-X∗=ωX1k+1-ωX2k-X∗=ωX1k-X∗+1-ωX2k-X∗=ωXk-1-μ1ATAXk-1+Xk-1B-C+AXk-1+Xk-1B-CBT-X∗+1-ωXk-1-μ2BBTXk-1+Xk-1AT-CT+BTXk-1+Xk-1AT-CTA-X∗=X~k-1-ωμ1ATAX~k-1+X~k-1B+AX~k-1+X~k-1BBT+1-ωμ2BBTX~k-1+X~k-1AT+BTX~k-1+X~k-1ATA=X~k-1-ωμ1ATθ1k-1+θ1k-1BT-1-ωμ2Bθ2k-1+θ2k-1A.The following two equalities are easily derived.(22)ATθ1k-1+θ1k-1BTF2≤A22+B22θ1k-1F2,Bθ2k-1+θ2k-1AF2≤A22+B22θ2k-1F2.Moreover (23)X~kF2=X~k-1-ωμ1ATθ1k-1+θ1k-1BT-1-ωμ2Bθ2k-1+θ2k-1AF2=X~k-1F2+ω2μ12ATθ1k-1+θ1k-1BTF2+1-ω2μ22Bθ2k-1+θ2k-1AF2-2ωμ1trX~k-1TATθ1k-1+θ1k-1BT-21-ωμ2trX~k-1TBθ2k-1+θ2k-1A+2ω1-ωμ1μ2trATθ1k-1+θ1k-1BTTBθ2k-1+θ2k-1A≤X~k-1F2+ω2μ12A22+B22θ1k-1F2+1-ω2μ22A22+B22θ2k-1F2-2ωμ1θ1k-1F2-21-ωμ2θ2k-1F2+2ω1-ωμ1μ2ATθ1k-1+θ1k-1BT,Bθ2k-1+θ2k-1A≤X~k-1F2+ω2μ12A22+B22θ1k-1F2+1-ω2μ22A22+B22θ2k-1F2-2ωμ1θ1k-1F2-21-ωμ2θ2k-1F2+2ωμ1ATθ1k-1+θ1k-1BTF1-ωμ2Bθ2k-1+θ2k-1AF≤X~k-1F2+ω2μ12A22+B22θ1k-1F2+1-ω2μ22A22+B22θ2k-1F2-2ωμ1θ1k-1F2-21-ωμ2θ2k-1F2+ω2μ12ATθ1k-1+θ1k-1BTF2+1-ω2μ22Bθ2k-1+θ2k-1AF2=X~k-1F2-2ωμ11-ωμ1A22+B22θ1k-1F2-21-ωμ21-1-ωμ2A22+B22θ2k-1F2.Therefore, the above equality implies that (24)X~kF2≤X~0F2-2ωμ11-ωμ1A22+B22∑i=1kθ1i-1F2-21-ωμ21-1-ωμ2A22+B22∑i=1kθ2i-1F2.It follows from 0<μ1<1/ωA22+B22 and 0<μ2<1/1-ωA22+B22 that(25)∑i=1∞θ1i-1F2<∞,∑i=1∞θ2i-1F2<∞.The above two inequalities (25) mean that(26)limk→∞θ1k-1=0,limk→∞θ2k-1=0.In other words(27)limk→∞AXk-1+Xk-1B=C,limk→∞BTXk-1+Xk-1AT=CT.Equation (27) implies that (28)limk→∞Xk=X∗.From Theorem6 and the above limited equation (28), we have (29)limk→∞Xsk=limk→∞Xk+XTk2=X∗+X∗T2=Xs.Following the same line, the idea of Algorithm7 can be extended to solve the skew symmetric solution of (1). First, we need the following theorem.Theorem 9.
The matrix equation (1) has a unique skew symmetric solution Xss if and only if the following pair of the matrix equations (30)AX+XB=C,BTX+XAT=-CThas a unique common solution X∗ and Xss=(X∗-(X∗)T)/2.The relaxed gradient based iterative algorithm for solving the skew symmetric solution of (1) can be stated as follows.Algorithm 10 (the relaxed gradient based iterative (RGI) algorithm for skew symmetric solution of (1)).
Step 1. Input matrices A,B,C∈Rn×n, given any small positive number ϵ and appropriative positive number ω such that 0<ω<1. Choose any initial matrix X(0). Set k≔1.
Step 2. If δk=AXss(k)+Xss(k)B-CF/CF<ϵ, stop; otherwise, go to Step 3.
Step 3. Update the sequences (31)X1k=Xk-1-μ1ATAXk-1+Xk-1B-C+AXk-1+Xk-1B-CBT,X2k=Xk-1-μ2BBTXk-1+Xk-1AT+CT+BTXk-1+Xk-1AT+CTA,Xk=ωX1k+1-ωX2k,Xssk=Xk-XTk2.
Step 4. Set k≔k+1; return to Step 2.Similarly, we have the following theorem to ensure the convergence of the Algorithm10.Theorem 11.
Assume the matrix equations (30) have a unique solution X∗; then the iterative sequence {X(k)} generated by the Algorithm 10 converges to X∗, if (32)0<μ1<1ωA22+B22,0<μ2<11-ωA22+B22;that is, limk→∞X(k)=X∗ or the error X(k)-X∗ converges to zero for any initial value X(0).Furthermore, the sequence{Xss(k)} converges to Xss, where Xss is the unique skew symmetric solution of (1).


## 4. Numerical Examples

In this section, two numerical examples are used to show the efficiency of the RGI Method. All the computations were performed on Intel® Core™ i7-4500U CPU@ 1.80 GHZ 2.40 GHZ system by using MATLAB 7.0. EER is the Frobenius norm of absolute error matrices which is defined to be EER=X(k)-X∗F, where X(k) is the kth iteration result for the RGI Method.Example 1.
In matrix equation (1), we choose (33)A=4.53.53.0-5.5-5.0-6.0-3.02.04.02.05.0-3.01.51.50.0-0.5,B=-1.5-6.5-3.0-2.51.06.03.01.0-4.0-5.0-4.00.0-1.5-1.50.0-0.5,C=-5.512.511.04.53.5-3.52.0-4.5-8.5-10.5-1.0-1.52.00.00.0-1.0.It is easy to show that the matrix equation (1) is consistent and has unique symmetric solution. By computing, the unique symmetric solution Xs is given as (34)Xs=-110110-100-11-110-10.TakingX(0)=10-612×3 and applying the RGI Method (Algorithm 7) to compute the symmetric solution of AX+XB=C, the sum of A22+B22 is 962.2175. When taking ω=0.4, μ1=0.0026, and μ2=0.0017, the iterative errors δk versus k are shown in Figure 1, where δk=Xs(k)-XsF2.Figure 1
Convergence cure of symmetric solution (ω=0.4, μ1=0.0026, and μ2=0.0017).Example 2.
In matrix equation (1), if we choose (35)A=4.58.53.0-5.5-5.0-9.0-3.07.08.014.05.0-12.01.51.50.0-0.5,B=-1.5-6.5-3.0-2.51.06.03.01.0-4.0-10.0-4.00.0-1.5-1.50.0-0.5,C=-11.08.028.50.57.5-2.5-27.06.5-19.5-19.534.0-3.54.025.013.53.5,it is easy to check that the matrix equation (1) is consistent and has unique skew symmetric solution. The unique skew symmetric solution Xss is computed as(36)Xss=0201-202-10-201-11-10.Also takingX(0)=10-612×3 and ω=0.5, μ1=0.0021, and μ2=0.0021, then using Algorithm 10, the convergence cure of skew symmetric solution is shown in Figure 2.Figure 2
Convergence cure of skew symmetric solution (ω=0.5, μ1=0.0021, and μ2=0.0021).



---

*Source: 1624969-2017-04-03.xml*