# Perturbation Analysis of the Nonlinear Matrix EquationX-∑i=1m‍Ai*XpiAi=Q

Authors: Jing Li
Journal: Abstract and Applied Analysis
 (2013)

---

## Abstract

Consider the nonlinear matrix equationX-∑i=1m‍Ai*XpiAi=Q with 0<pi<1. Two perturbation bounds and the backward error of an approximate solution to the equation are derived. Explicit expressions of the condition number for the equation are obtained. The theoretical results are illustrated by numerical examples.


---

## Body

## 1. Introduction

In this paper we consider the Hermitian positive definite solution of the nonlinear matrix equation(1)X-∑i=1mAi*XpiAi=Q,
where 0<pi<1(i=1,2,…,m), A1,A2,…,Am are n×n complex matrices, m is a positive integer, and Q is a positive definite matrix. Here, Ai* denotes the conjugate transpose of the matrix Ai.Whenm>1, (1) is recognized as playing an important role in solving a system of linear equations. For example, in many physical calculations, one must solve the system of linear equation
(2)Mx=f,
where
(3)M=(I0⋯0A10I⋯0A2⋮⋮⋱⋮⋮00⋯IAmA1*A2*⋯Am*-Q)
arises in a finite difference approximation to an elliptic partial differential equation (for more information, refer to [1]). We can rewrite M as M=M~+D, where
(4)M~=(X-p10⋯0A10X-p2⋯0A2⋮⋮⋱⋮⋮00⋯X-pmAmA1*A2*⋯Am*-Q),D=(I-X-p10⋯000I-X-p2⋯00⋮⋮⋱⋮⋮00⋯I-X-pm000⋯00).M~ can be factored as
(5)M~=(-I0⋯000-I⋯00⋮⋮⋱⋮⋮00⋯-I0-A1*Xp1-A2*Xp2⋯-Am*Xpm-I)×(-X-p10⋯0-A10-X-p2⋯0-A2⋮⋮⋱⋮⋮00⋯-X-pm-Am00⋯0X)
if and only if X is a solution of equation X-∑i=1mAi*XpiAi=Q. When   m=1, this type of nonlinear matrix equations arises in ladder networks, dynamic programming, control theory, stochastic filtering, statistics, and so forth [2–7].For the similar equationsX±A*X-pA=Q, Xs±A*X-tA=Q, X+∑i=1mAi*X-1Ai=I, and X=Q-A*X-1A+B*X-1B, there were many contributions in the literature to the theory, numerical solutions, and perturbation analysis [8–32]. Jia and Gao [33] derived two perturbation estimates for the solution of the equation X-A*XqA=Q with 0<q<1. In addition, Duan et al. [34] proved that the equation X-∑i=1mAi*XδiAi=Q(0<|δi|<1) has a unique positive definite solution. They also proposed an iterative method for obtaining the unique positive definite solution. However, to our best knowledge, there has been no perturbation analysis for (1) with m>1 in the known literature.As a continuation of the previous results, the rest of the paper is organized as follows. In Section2, some preliminary lemmas are given. In Section 3, two perturbation bounds for the unique solution to (1) are derived. Furthermore, in Section 4, we obtain the backward error of an approximate solution to (1). In Section 5, we also discuss the condition number of the unique solution to (1). Finally, several numerical examples are presented in Section 6.We denote by𝒞n×n the set of n×n complex matrices, by ℋn×n the set of n×n Hermitian matrices, by I the identity matrix, by i the imaginary unit, by ∥·∥ the spectral norm, by ∥·∥F the Frobenius norm, and by λmax(M) and λmin(M) the maximal and minimal eigenvalues of M, respectively. For A=(a1,…,an)=(aij)∈𝒞n×n and a matrix B, A⊗B=(aijB) is a Kronecker product, and vecA is a vector defined by vecA=(a1T,…,anT)T. For X,Y∈ℋn×n, we write X≥Y (X>Y,resp.) if X-Y is Hermitian positive semidefinite (definite, resp.).


## 2. Preliminaries

Lemma 1 (see [35]).
IfA≥B>0 and 0≤γ≤1, then Aγ≥Bγ.Lemma 2 (see [33]).
For any Hermitian positive definite matrixX and Hermitian matrix ΔX, one has  (i)
Xq=(sinqπ/π)∫0∞‍X1/2(λI+X)-1X1/2λq-1dλ,  0<q<1;(ii)
Xq=(sinqπ/(1-q)π)∫0∞‍X1/2(λI+X)-1X(λI+X)-1X1/2λq-1dλ,0<q<1.In addition, if X+ΔX≥(1/ν)X>0 and 0<q<1, then  (iii)
∥X-1/2A*((X+ΔX)q-Xq)AX-1/2∥≤q(∥X-1/2ΔXX-1/2∥)+ν∥X-1/2ΔXX-1/2∥2)∥Xq/2AX-1/2∥2.Lemma 3 (see [34]).
The matrix equationX-∑i=1mAi*XδiAi=Q(0<|δi|<1) always has a unique positive definite solution X. The matrix sequence Xk:
(6)Xs+m+1=Q+∑i=1mAi*Xs+iδiAi,s=0,1,2,…,
converges to the unique positive definite solution X for arbitrary initial positive definite matrices X1,X2,…,Xm.


## 3. Perturbation Bounds

Here the perturbed equation(7)X~-∑i=1m‍A~i*X~piA~i=Q~.
is considered, where 0<pi<1  and  A~i and Q~ are small perturbations of Ai and Q in (1), respectively. We assume that X and X~ are the solutions of (1) and (7), respectively. Let ΔX=X~-X, ΔQ=Q~-Q and ΔAi=A~i-Ai.By Lemma3, we know that (1) always has a unique positive definite solution X; then in this section two perturbation bounds for the unique positive definite solution of (1) are developed. The relative perturbation bound in Theorem 5 does not depend on any knowledge of the actual solution X of (1). Furthermore, a sharper perturbation bound in Theorem 8 is derived.To prove the next theorem, we first verify the following lemma.Lemma 4.
IfX is a solution of (1), then
(8)X≥(λmin(Q)+∑i=1m‍λmin(Ai*Ai)λminpi(Q))I=βI.Proof.
By Lemma3, (1) with 0<pi<1 always has a unique positive definite solution X. Then X>0, and it follows that Xpi>0. Therefore X≥Q. By Lemma 1 and (1), we have X≥Q+∑i=1mAi*QpiAi≥  (λmin(Q)+∑i=1mλmin(Ai*Ai)λminpi(Q))I=βI.The next theorem generalizes [33, Theorem 4] with m=1,  ∥ΔQ∥=0 to arbitrary integer m≥1,  ∥ΔQ∥>0.Theorem 5.
Letb=β+∥ΔQ∥-∑i=1mpiβpi∥Ai∥2,s=∑i=1mβpi∥ΔAi∥(2∥Ai∥+∥ΔAi∥). If
(9)0<b<2(β-s),b2-4(β-s)(s+∥ΔQ∥)≥0,
then
(10)∥X~-X∥∥X∥≤ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥≜ξ1,
where
(11)ϱ=2s∑i=1m∥ΔAi∥(b+b2-4(β-s)(s+∥ΔQ∥)),ω=2b+b2-4(β-s)(s+∥ΔQ∥).Proof.
Let(12)Ω={∑i=1m∥ΔAi∥ΔX∈ℋn×n:∥X-1/2ΔXX-1/2∥≤ϱ∑i=1m∥ΔAi∥+ω∥ΔQ∥}.
Obviously, Ω is a nonempty bounded convex closed set. Let
(13)f(ΔX)=∑i=1m‍(A~i*(X+ΔX)piA~i-Ai*XpiAi)+ΔQ,ΔX∈Ω.
Evidently, f:Ω↦ℋn×n is continuous. We will prove that f(Ω)⊆Ω.
For everyΔX∈Ω, it follows that ∥X-1/2ΔXX-1/2∥≤ϱ∑i=1m∥ΔAi∥+ω∥ΔQ∥. Thus
(14)(ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥)I≥X-1/2ΔXX-1/2≥(-ϱ∑i=1m‍∥ΔAi∥-ω∥ΔQ∥)I,(1+ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥)X≥X+ΔX≥(1-ϱ∑i=1m‍∥ΔAi∥-ω∥ΔQ∥)X.
According to (9), we have
(15)ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥=2(∥ΔQ∥+s)b+b2-4(β-s)(s+∥ΔQ∥)≤2(∥ΔQ∥+s)b≤b2(β-s)<1.
Therefore
(16)(1-ϱ∑i=1m‍∥ΔAi∥-ω∥ΔQ∥)X>0.
From Lemmas 2 and 4, it follows that
(17)∥X-1/2[∑i=1m‍Ai*((X+ΔX)pi-Xpi)Ai]X-1/2∥≤(∥X-1/2ΔXX-1/2∥21-ϱ∑i=1m∥ΔAi∥-ω∥ΔQ∥∥X-1/2ΔXX-1/2∥+∥X-1/2ΔXX-1/2∥21-ϱ∑i=1m∥ΔAi∥-ω∥ΔQ∥)×(∑i=1m‍pi∥Xpi/2AiX-1/2∥2)≤(∥X-1/2ΔXX-1/2∥+∥X-1/2ΔXX-1/2∥21-ϱ∑i=1m∥ΔAi∥-ω∥ΔQ∥)×(∑i=1m‍piβ1-pi∥Ai∥2).
Therefore
(18)∥X-1/2f(ΔX)X-1/2∥-=∥X-1/2[∑i=1m‍A~i*(X+ΔX)piA~i-Ai*XpiAi]-----×X-1/2+X-1/2ΔQX-1/2∑i=1m‍∥-≤∥∑i=1m‍X-1/2Ai*((X+ΔX)pi-Xpi)AiX-1/2∥--+∥X-1/2ΔQX-1/2∥--+∥∑i=1m‍X-1/2[ΔAi*(X+ΔX)pi(Ai+ΔAi)-----+Ai*(X+ΔX)piΔAi]X-1/2∑i=1m‍X-1/2∑i=1m‍∥-≤(∥X-1/2ΔXX-1/2∥+∥X-1/2ΔXX-1/2∥21-ϱ∑i=1m∥ΔAi∥-ω∥ΔQ∥)--×(∑i=1m‍piβ1-pi∥Ai∥2)--+∑i=1m∥ΔAi∥(2∥Ai∥+∥ΔAi∥)β1-pi‍------×(1+ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥)+∥ΔQ∥β-≤(ξ1+ξ121-ξ1)(∑i=1m‍piβ1-pi∥Ai∥2)--+sβ(1+ξ1)+∥ΔQ∥β-=ξ1.
That is, f(Ω)⊆Ω. By Brouwer's fixed point theorem, there exists a ΔX∈Ω such that f(ΔX)=ΔX. Moreover, by Lemma 3, we know that X and X~ are the unique solutions to (1) and (7), respectively. Then
(19)∥X~-X∥∥X∥=∥ΔX∥∥X∥=∥X1/2(X-1/2ΔXX-1/2)X1/2∥∥X∥≤∥X-1/2ΔXX-1/2∥≤ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥.Remark 6.
According to(20)ϱ∑i=1m‍∥ΔAi∥+ω∥ΔQ∥=2(∑i=1m∥ΔAi∥(2∥Ai∥+∥ΔAi∥)+∥ΔQ∥)b+b2-4(β-s)(s+∥ΔQ∥),
we get  ϱ∑i=1m∥ΔAi∥+ω∥ΔQ∥→0  for  ∥ΔQ∥→0 and ∥ΔAi∥→0(i=1,2,…,m). Therefore (1) is well posed.Next, a sharper perturbation estimate is derived.Subtracting (1) from (7), we have
(21)ΔX+∑i=1m‍sinpiππ∫0∞[(λI+X)-1X1/2Ai]*ΔX+∑i=1m‍sinpiππ×ΔX[(λI+X)-1X1/2Ai]λpi-1dλ=E+h(ΔX),
where
(22)Bi=XpiAi,i=1,2,…,m,E=∑i=1m‍(Bi*ΔAi+ΔAi*Bi)+∑i=1m‍ΔAi*XpiΔAi+ΔQ,Zi(ΔX)=sinpiππ∫0∞X1/2(λI+X)-1ΔX(λI+X+ΔX)-1----------×ΔX(λI+X)-1X1/2λpi-1dλ,Vi(ΔX)=sinpiππ∫0∞X1/2(λI+X)-1ΔX---------×(λI+X+ΔX)-1X1/2λpi-1dλ,h(ΔX)=∑i=1m‍[Ai*Zi(ΔX)Ai-A~i*Vi(ΔX)ΔAih(ΔX)-ΔAi*Vi(ΔX)AiA~i*A~i*Vi].Lemma 7.
If∑i=1m(∥Ai∥2/β1-pi)<1, then the linear operator L:ℋn×n→ℋn×n defined by
(23)LW=W+∑i=1m‍sinpiππ×∫0∞[(λI+X)-1X1/2Ai]*‍×W[(λI+X)-1X1/2Ai]×λpi-1dλ,W∈ℋn×n
is invertible.Proof.
It suffices to show that the following equation:(24)LW=V
has a unique solution for every  V∈ℋn×n. Define the operator M:ℋn×n→ℋn×n by
(25)MZ=∑i=1m‍sinpiππ×∫0∞‍X-1/2Ai*X1/2×(λI+X)-1X1/2ZX1/2×(λI+X)-1X1/2AiX-1/2×λpi-1dλ,Z∈ℋn×n.
Let Y=X-1/2WX-1/2. Thus (21) is equivalent to
(26)Y+MY=X-1/2VX-1/2.
According to Lemma 2, we have
(27)∥MY∥≤∑i=1m∥sinpiππ‍×∫0∞X-1/2Ai*X1/2×(λI+X)-1X(λI+X)-1×X1/2AiX-1/2λpi-1dλsinpiππ∥∥Y∥≤∑i=1m‍(1-pi)∥Xpi/2AiX-1/2∥2∥Y∥≤∑i=1m‍∥Ai∥2β1-pi∥Y∥<∥Y∥,
which implies that ||M||<1 and I+M is invertible. Therefore, the operator L is invertible.Furthermore, we define operatorsPi:𝒞n×n→ℋn×n by
(28)PiZi=L-1(Bi*Zi+Zi*Bi),Zi∈𝒞n×n,i=1,2,…,m.
Thus, we can rewrite (21) as
(29)ΔX=L-1ΔQ+∑i=1m‍PiΔAi+L-1(∑i=1m‍ΔAi*XpiΔAi)+L-1(h(ΔX)).
Define
(30)∥L-1∥=maxW∈ℋn×n∥W∥=1∥L-1W∥,∥Pi∥=maxZ∈𝒞n×n∥Z∥=1∥PiZ∥,i=1,2,…,m.
Now we denote that
(31)l=∥L-1∥-1,ζ=∥X-1∥,ξi=∥Xpi∥,ni=∥Pi∥,θ=ζ2l∑i=1m‍ξi∥Ai∥2,i=1,2,…,m,ϵ=1l∥ΔQ∥+∑i=1m‍(ni∥ΔAi∥+ξil∥ΔAi∥2),σ=ζl∑i=1m‍ξi(2∥Ai∥+∥ΔAi∥)∥ΔAi∥.Theorem 8.
Suppose thatX and X~ are the solutions of (1) and (7), respectively. If ∑i=1m∥Ai∥2/β1-pi<1,
(32)σ<1,ϵ<(1-σ)2ζ+σζ+2θ+2((ζ+θ)(σζ+θ)),
then
(33)∥X~-X∥≤2ϵ1+ϵζ-σ+(1+ζϵ-σ)2-4ϵ(ζ+θ)≡ν.Proof.
Let(34)f(ΔX)=L-1ΔQ+∑i=1m‍PiΔAi+L-1(∑i=1m‍ΔAi*XpiΔAi)+L-1(h(ΔX)).
Obviously, f:ℋn×n→ℋn×n is continuous. The condition (32) ensures that the quadratic equation (ζ+θ)x2-(1+ζϵ-σ)x+ϵ=0 with respect to the variable x has two positive real roots. The smaller one is
(35)ν=2ϵ1+ϵζ-σ+(1+ζϵ-σ)2-4ϵ(ζ+θ).
Define Ω={ΔX∈ℋn×n:∥ΔX∥≤ν}. Then for any ΔX∈Ω, by (32), we have
(36)∥X-1ΔX∥≤∥X-1∥∥ΔX∥≤ζν≤ζ·2ϵ1+ζϵ-σ=1+ζϵ+σ-11+ζϵ-σ≤1+-2(1-σ)(ζσ+θ)(ζ+σζ+2θ)(1+ζϵ-σ)<1.
It follows that I-X-1ΔX is nonsingular and
(37)∥I-X-1ΔX∥≤11-∥X-1ΔX∥≤11-ζ∥ΔX∥.
Using (22) and Lemma 2, we have
(38)∥Zi(ΔX)∥≤(1-pi)∥ΔX∥2∥X-1∥2×∥(I+X-1ΔX)-1∥∥Xpi∥≤ξiζ2∥ΔX∥21-ζ∥ΔX∥,∥Vi(ΔX)∥≤∥Xpi∥∥ΔX∥∥X-1∥∥(I+X-1ΔX)-1∥≤ξiζ∥ΔX∥1-ζ∥ΔX∥,(39)∥h(ΔX)∥≤∑i=1m‍(∥Ai∥2∥Zi(ΔX)∥+(2∥Ai∥+∥ΔAi∥)×∥ΔAi∥∥Vi(ΔX)∥∥Ai∥2)≤∑i=1m‍(∥ΔX∥1-ζ∥ΔX∥ξiζ2∥Ai∥2∥ΔX∥21-ζ∥ΔX∥+(2∥Ai∥+∥ΔAi∥)×∥ΔAi∥ξiζ∥ΔX∥1-ζ∥ΔX∥).
Noting (31) and (34), it follows that
(40)∥f(ΔX)∥≤1l∥ΔQ∥+∑i=1m‍(ni∥ΔAi∥+ζil∥ΔAi∥2)+1l∥h(ΔX)∥≤ϵ+σ∥ΔX∥1-ζ∥ΔX∥+θ∥ΔX∥21-ζ∥ΔX∥≤ϵ+σν1-ζν+θν21-ζν=ν,
for ΔX∈Ω. That is, f(Ω)⊆Ω. According to Schauder fixed point theorem, there exists ΔX*∈Ω such that f(ΔX*)=X*. It follows that X+ΔX* is a Hermitian solution of (7). By Lemma 3, we know that the solution of (7) is unique. Then ∥ΔX*∥=∥X~-X∥≤ν.Remark 9.
From Theorem8, we get the first order perturbation bound for the solution as follows:
(41)∥X~-X∥≤1l∥ΔQ∥+∑i=1m‍ni∥ΔAi∥+O(∥(ΔA1,ΔA2,…,ΔAm,ΔQ)∥F2),
as(ΔA1,ΔA2,…,ΔAm,ΔQ)→0.
Combining this with (29) gives
(42)ΔX=L-1ΔQ+L-1∑i=1m‍(Bi*ΔAi+ΔAi*Bi)+O(∥(ΔA1,ΔA2,…,ΔAm,ΔQ)∥F2).
as (ΔA1,ΔA2,…,ΔAm,ΔQ)→0.


## 4. Backward Error

In this section, a backward error of an approximate solution for the unique solution to (1) is obtained.Theorem 10.
LetX~>0 be an approximation to the solution X of (1). If Σ=∑i=1mpi∥X~pi/2AiX~-1/2∥2<1 and the residual R(X~)≡Q+∑i=1mAi*X~piAi-X~ satisfies
(43)∥R(X~)∥<θ12∥X~-1∥min{1,θ12},whereθ1≡1+∥X~-1∥∥R(X~)∥-Σ>0,
then
(44)∥X~-X∥≤μ∥R(X~)∥,whereμ=2∥X~∥∥X~-1∥θ1+θ12-4∥X~-1∥∥R(X~)∥.Proof.
Let(45)Ψ={ΔX∈ℋn×n:∥X~-1/2ΔXX~-1/2∥≤θ2∥R(X~)∥},
where θ2=μ/∥X~∥. Obviously, Ψ is a nonempty bounded convex closed set. Let
(46)g(ΔX)=∑i=1m‍Ai*[(X~+ΔX)pi-X~pi]Ai+R(X~).
Evidently g:Ψ↦ℋn×n is continuous. We will prove that g(Ψ)⊆Ψ. For every ΔX∈Ψ, we have
(47)∥X~-1/2ΔXX~-1/2∥≤θ2∥R(X~)∥.
Hence
(48)X~-1/2ΔXX~-1/2≥-θ2∥R(X~)∥I.
That is,
(49)X~+ΔX≥(1-θ2∥R(X~)∥)X~.
Using (43), one sees that
(50)θ2∥R(X~)∥=2∥X~-1∥∥R(X~)∥θ1+θ12-4∥X~-1∥∥R(X~)∥<2∥X~-1∥∥R(X~)∥θ1<1.
Therefore (1-θ2∥R(X~)∥)X~>0.
According to (17), we obtain
(51)∥X~-1/2g(ΔX)X~-1/2∥≤(∥X~-1/2ΔXX~-1/2∥+∥X~-1/2ΔXX~-1/2∥21-θ2∥R(X~)∥)Σ+∥X~-1/2R(X~)X~-1/2∥≤(θ2∥R(X~)∥+(θ2∥R(X~)∥)21-θ2∥R(X~)∥)Σ+∥X~-1∥∥R(X~)∥=θ2∥R(X~)∥,
for ΔX∈Ψ. That is, g(Ψ)⊆Ψ. By Brouwer’s fixed point theorem, there exists a ΔX∈Ψ such that g(ΔX)=ΔX. Hence X~+ΔX is a solution of (1). Moreover, by Lemma 3, we know that the solution X of (1) is unique. Then
(52)∥X~-X∥=∥ΔX∥≤∥X~∥∥X~-1/2ΔXX~-1/2∥=θ2∥X~∥∥R(X~)∥=μ∥R(X~)∥.


## 5. Condition Number

In this section, we apply the theory of condition number developed by Rice [36] to study condition numbers of the unique solution to (1).


### 5.1. The Complex Case

Suppose thatX and X~ are the solutions of the matrix equations (1) and (7), respectively. Let ΔAi=A~i-Ai, ΔQ=Q~-Q and ΔX=X~-X. Using Theorem 8 and Remark 9, we have
(53)ΔX=X~-X=L-1ΔQ+L-1∑i=1m‍(Bi*ΔAi+ΔAi*Bi)+O(∥(ΔA1,ΔA2,…,ΔAm,ΔQ)∥F2),
as (ΔA1,ΔA2,…,ΔAm,ΔQ)→0.By the theory of condition number developed by Rice [36], we define the condition number of the Hermitian positive definite solution X to (1) by
(54)c(X)=limδ→0sup∥(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)∥F≤δ∥ΔX∥Fξδ,
where ξ, ρ, and ηi,  i=1,2,…,m, are positive parameters. Taking ξ=ηi=ρ=1 in (54) gives the absolute condition number cabs(X), and taking ξ=∥X∥F, ηi=∥Ai∥F, and ρ=∥Q∥F in (54) gives the relative condition number crel(X).Substituting (53) into (54), we get (55)c(X)=1ξmax(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)≠0ΔAi∈𝒞n×n,ΔQ∈ℋn×n∥L-1(ΔQ+∑i=1m(Bi*ΔAi+ΔAi*Bi))∥F∥(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)∥F=1ξmax(E1,E2,…,Em,H)≠0Ei∈𝒞n×n,H∈ℋn×n∥L-1(ρH+∑i=1mηi(Bi*Ei+Ei*Bi))∥F∥(E1,E2,…,Em,H)∥F.LetL be the matrix representation of the linear operator L. Then it is easy to see that
(56)L=I+∑i=1m‍sinpiππ∫0∞‍[(λI+X)-1X1/2Ai]TI+∑i=1m‍sinpiππ⊗[(λI+X)-1X1/2Ai]*×λpi-1dλ.
Let
(57)L-1=S+iΣ,L-1(I⊗Bi*)=L-1(I⊗(XpiAi)*)=Ui1+iΩi1,L-1(BiT⊗I)Π=L-1((XpiAi)T⊗I)Π=Ui2+iΩi2,Sc=[S-ΣΣS],Ui=[Ui1+Ui2Ωi2-Ωi1Ωi1+Ωi2Ui1-Ui2],i=1,2,…,m,vecH=x+iy,  vecEi=ai+ibi,  g=(xT,yT,a1T,b1T,…,amT,bmT)T,  M=(E1,E2,…,Em,H), where x,y,ai,bi∈ℛn2,S,Σ,Ui1,Ui2,Ωi1,Ωi2∈ℛn2×n2,  i=1,2,…,m,  and  Π is the vec-permutation matrix, such that
(58)vecAT=ΠvecA.
Then we obtain that
(59)c(X)=1ξmaxM≠0∥L-1(ρH+∑i=1mηi(Bi*Ei+Ei*Bi))∥F∥(E1,E2,…,Em,H)∥F=1ξmaxM≠0(∥∑i=1m‍((I⊗Bi*)vecEiρL-1vecH-------+∑i=1m‍ηiL-1((I⊗Bi*)vecEi-------+(BiT⊗I)vecEi*)(∥∑i=1m‍((I⊗Bi*)vecEiρL-1vecHηiL-1∥)-------×(∥vec(E1,E2,…,Em,H)∥)-1=1ξmaxM≠0(∥∑i=1m‍((I⊗Bi*)vecEiρ(S+iΣ)(x+iy)-------+∑i=1m‍ηi[(Ui1+iΩi1)(ai+ibi)----------+(Ui2+iΩi2)(ai-ibi)]∑i=1m‍((I⊗Bi*)vecEi∥)-------×(∥vec(E1,E2,…,Em,H)∥)-1=1ξmaxg≠0∥(ρSc,η1U1,η2U2,…,ηmUm)g∥∥g∥=1ξ∥(ρSc,η1U1,η2U2,…,ηmUm)∥,Ei∈𝒞n×n,H∈ℋn×n.Then we have the following theorem.Theorem 11.
If∑i=1m∥Ai∥2/β1-pi<1, then the condition number c(X) defined by (54) has the explicit expression
(60)c(X)=1ξ∥(ρSc,η1U1,η2U2,…,ηmUm)∥,
where the matrices Sc and Ui are defined as in (57).Remark 12.
From (60) we have the relative condition number
(61)crel(X)=∥(∥Q∥FSc,∥A1∥FU1,∥A2∥FU2,…,∥Am∥FUm)∥∥X∥F.


### 5.2. The Real Case

In this subsection we consider the real case. That is, all the coefficient matricesAi, Q of (1) are real. In such a case the corresponding solution X is also real. Completely similar arguments as Theorem 11 give the following theorem.Theorem 13.
LetAi, Q be real and let c(X) be the condition number defined by (54). If ∑i=1m∥Ai∥2/β1-pi<1, then c(X) has the explicit expression
(62)c(X)=1ξ∥(ρSr,η1U1,η2U2,…,ηmUm)∥,
where
(63)Sr=([×X12(λI+X)-1X1/2Ai]TI+∑i=1msinpiππ‍∫0∞‍[(λI+X)-1X1/2Ai]T⊗[×X12(λI+X)-1X1/2Ai]Tλpi-1dλ)-1Ui=Sr[I⊗(AiTXpi)+((AiTXpi)⊗I)Π],i=1,2,…,m.Remark 14.
In the real case the relative condition number is given by(64)crel(X)=∥(∥Q∥FSr,∥A1∥FU1,∥A2∥FU2,…,∥Am∥FUm)∥∥X∥F.


## 5.1. The Complex Case

Suppose thatX and X~ are the solutions of the matrix equations (1) and (7), respectively. Let ΔAi=A~i-Ai, ΔQ=Q~-Q and ΔX=X~-X. Using Theorem 8 and Remark 9, we have
(53)ΔX=X~-X=L-1ΔQ+L-1∑i=1m‍(Bi*ΔAi+ΔAi*Bi)+O(∥(ΔA1,ΔA2,…,ΔAm,ΔQ)∥F2),
as (ΔA1,ΔA2,…,ΔAm,ΔQ)→0.By the theory of condition number developed by Rice [36], we define the condition number of the Hermitian positive definite solution X to (1) by
(54)c(X)=limδ→0sup∥(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)∥F≤δ∥ΔX∥Fξδ,
where ξ, ρ, and ηi,  i=1,2,…,m, are positive parameters. Taking ξ=ηi=ρ=1 in (54) gives the absolute condition number cabs(X), and taking ξ=∥X∥F, ηi=∥Ai∥F, and ρ=∥Q∥F in (54) gives the relative condition number crel(X).Substituting (53) into (54), we get (55)c(X)=1ξmax(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)≠0ΔAi∈𝒞n×n,ΔQ∈ℋn×n∥L-1(ΔQ+∑i=1m(Bi*ΔAi+ΔAi*Bi))∥F∥(ΔA1/η1,ΔA2/η2,…,ΔAm/ηm,ΔQ/ρ)∥F=1ξmax(E1,E2,…,Em,H)≠0Ei∈𝒞n×n,H∈ℋn×n∥L-1(ρH+∑i=1mηi(Bi*Ei+Ei*Bi))∥F∥(E1,E2,…,Em,H)∥F.LetL be the matrix representation of the linear operator L. Then it is easy to see that
(56)L=I+∑i=1m‍sinpiππ∫0∞‍[(λI+X)-1X1/2Ai]TI+∑i=1m‍sinpiππ⊗[(λI+X)-1X1/2Ai]*×λpi-1dλ.
Let
(57)L-1=S+iΣ,L-1(I⊗Bi*)=L-1(I⊗(XpiAi)*)=Ui1+iΩi1,L-1(BiT⊗I)Π=L-1((XpiAi)T⊗I)Π=Ui2+iΩi2,Sc=[S-ΣΣS],Ui=[Ui1+Ui2Ωi2-Ωi1Ωi1+Ωi2Ui1-Ui2],i=1,2,…,m,vecH=x+iy,  vecEi=ai+ibi,  g=(xT,yT,a1T,b1T,…,amT,bmT)T,  M=(E1,E2,…,Em,H), where x,y,ai,bi∈ℛn2,S,Σ,Ui1,Ui2,Ωi1,Ωi2∈ℛn2×n2,  i=1,2,…,m,  and  Π is the vec-permutation matrix, such that
(58)vecAT=ΠvecA.
Then we obtain that
(59)c(X)=1ξmaxM≠0∥L-1(ρH+∑i=1mηi(Bi*Ei+Ei*Bi))∥F∥(E1,E2,…,Em,H)∥F=1ξmaxM≠0(∥∑i=1m‍((I⊗Bi*)vecEiρL-1vecH-------+∑i=1m‍ηiL-1((I⊗Bi*)vecEi-------+(BiT⊗I)vecEi*)(∥∑i=1m‍((I⊗Bi*)vecEiρL-1vecHηiL-1∥)-------×(∥vec(E1,E2,…,Em,H)∥)-1=1ξmaxM≠0(∥∑i=1m‍((I⊗Bi*)vecEiρ(S+iΣ)(x+iy)-------+∑i=1m‍ηi[(Ui1+iΩi1)(ai+ibi)----------+(Ui2+iΩi2)(ai-ibi)]∑i=1m‍((I⊗Bi*)vecEi∥)-------×(∥vec(E1,E2,…,Em,H)∥)-1=1ξmaxg≠0∥(ρSc,η1U1,η2U2,…,ηmUm)g∥∥g∥=1ξ∥(ρSc,η1U1,η2U2,…,ηmUm)∥,Ei∈𝒞n×n,H∈ℋn×n.Then we have the following theorem.Theorem 11.
If∑i=1m∥Ai∥2/β1-pi<1, then the condition number c(X) defined by (54) has the explicit expression
(60)c(X)=1ξ∥(ρSc,η1U1,η2U2,…,ηmUm)∥,
where the matrices Sc and Ui are defined as in (57).Remark 12.
From (60) we have the relative condition number
(61)crel(X)=∥(∥Q∥FSc,∥A1∥FU1,∥A2∥FU2,…,∥Am∥FUm)∥∥X∥F.


## 5.2. The Real Case

In this subsection we consider the real case. That is, all the coefficient matricesAi, Q of (1) are real. In such a case the corresponding solution X is also real. Completely similar arguments as Theorem 11 give the following theorem.Theorem 13.
LetAi, Q be real and let c(X) be the condition number defined by (54). If ∑i=1m∥Ai∥2/β1-pi<1, then c(X) has the explicit expression
(62)c(X)=1ξ∥(ρSr,η1U1,η2U2,…,ηmUm)∥,
where
(63)Sr=([×X12(λI+X)-1X1/2Ai]TI+∑i=1msinpiππ‍∫0∞‍[(λI+X)-1X1/2Ai]T⊗[×X12(λI+X)-1X1/2Ai]Tλpi-1dλ)-1Ui=Sr[I⊗(AiTXpi)+((AiTXpi)⊗I)Π],i=1,2,…,m.Remark 14.
In the real case the relative condition number is given by(64)crel(X)=∥(∥Q∥FSr,∥A1∥FU1,∥A2∥FU2,…,∥Am∥FUm)∥∥X∥F.


## 6. Numerical Examples

To illustrate the results of the previous sections, in this section three simple examples are given, which were carried out using MATLAB 7.1. For the stopping criterion we takeεk+1(X)=∥Xk-∑i=1mAi*XkpiAi-Q∥<1.0e-10.Example 15.
We consider the matrix equation(65)X-A1*X1/2A1-A2*X1/3A2=I,
with
(66)A1=(1/3)+2×10-2∥A∥A,A2=(1/6)+3×10-2∥A∥A,A=(00.9501).Suppose that the coefficient matricesA1 and A2 are perturbed to A~i=Ai+ΔAi,i=1,2, where
(67)ΔA1=10-j∥CT+C∥(CT+C),ΔA2=3×10-j-1∥CT+C∥(CT+C),
and C is a random matrix generated by MATLAB function randn.We now consider the corresponding perturbation bounds for the solutionX in Theorems 5 and 8.The assumptions in Theorem5 are
(68)ass1=2(β-s)-b>0,ass2=β+∥ΔQ∥-∑i=1m‍piβpi∥Ai∥2>0,ass3=b2-4(β-s)(s+∥ΔQ∥)≥0.The assumptions in Theorem8 are
(69)ass4=1-σ>0,ass5=(1-σ)2ζ+σζ+2θ+2((ζ+θ)(σζ+θ))-ϵ>0.
By computation, we list them in Table 1.Table 1
Assumptions check for Example15 with different values of j.
j
4
5
6
7
ass
1
1.0749
1.0753
1.0753
1.0753
ass
2
0.9248
0.9247
0.9247
0.9247
ass
3
0.8543
0.8550
0.8550
0.8550
ass
4
0.9999
1.0000
1.0000
1.0000
ass
5
0.7645
0.7648
0.7648
0.7648The results listed in Table1 show that the assumptions of Theorems 5 and 8 are satisfied.By Theorems5 and 8, we can compute the relative perturbation bounds ξ1 and ξ2=ν/∥X∥, respectively. These results averaged as the geometric mean of 10 randomly perturbed runs. Some results are listed in Table 2.Table 2
Perturbation bounds for Example15 with different values of j.
j
4
5
6
7
∥
X
~
-
X
∥
/
∥
X
∥
2.5627
×
1
0
-
5
3.8447
×
1
0
-
6
5.1681
×
1
0
-
7
2.1776
×
1
0
-
8
ξ
1
2.1885
×
1
0
-
4
1.9891
×
1
0
-
5
2.4026
×
1
0
-
6
1.8251
×
1
0
-
7
ξ
2
8.0828
×
1
0
-
5
7.4741
×
1
0
-
6
8.8011
×
1
0
-
7
6.9496
×
1
0
-
8The results listed in Table2 show that the perturbation bound ξ2 given by Theorem 8 is fairly sharp, while the bound ξ1 given by Theorem 5 which does not depend on the exact solution is conservative.Example 16.
We consider the matrix equation(70)X-A1*X0.5A1-A2*X0.25A2=I,
with
(71)A1=(1/3)+2×10-2∥A∥A,A2=(1/6)+3×10-2∥A∥A,A=(2100012100012100012100012).
Choose X~1=A, X~2=2A. Let the approximate solution X~k of X be given with the iterative method (6), where k is the iteration number.The residualR(X~k)≡I+A1*X~k0.5A1+A2*X~k0.25A2-X~k satisfies the conditions in Theorem 10.By Theorem10, we can compute the backward error bound for X~k as follows:
(72)∥X~k-X∥≤μ∥R(X~k)∥,
where
(73)μ=2∥X~k∥∥X~k-1∥θ1+θ12-4∥X~k-1∥∥R(X~k)∥,θ1≡1+∥X~k-1∥∥R(X~k)∥-(0.5∥X~k1/4A1X~k-1/2∥2+0.25∥X~k1/8A2X~k-1/2∥2).
Some results are listed in Table 3.Table 3
Backward error bound for Example16 with different values of k.
k
8
10
12
14
∥
X
~
k
-
X
∥
6.1091
×
1
0
-
4
4.0865
×
1
0
-
5
2.6837
×
1
0
-
6
1.7372
×
1
0
-
7
μ
∥
R
(
X
~
k
)
∥
7.1435
×
1
0
-
4
4.7784
×
1
0
-
5
3.1381
×
1
0
-
6
2.0318
×
1
0
-
7The results listed in Table3 show that the error bound given by Theorem 10 is fairly sharp.Example 17.
We study the matrix equation(74)X-A1*X1/2A1-A2*X1/3A2=Q,
with
(75)A1=(00.55+10-k00),A2=12A1,Q=(1101).
By Remark 14, we can compute the relative condition number crel(X). Some results are listed in Table 4.
Table4 shows that the unique positive definite solution X is well conditioned.Table 4
Relative condition number for Example17 with different values of k.
k
1
3
5
7
9
c
rel
(
X
)
1.0717
1.0228
1.0225
1.0225
1.0225



---

*Source: 979832-2013-06-11.xml*