陳興團，馬昌鳳

(福建師范大學(xué) 數(shù)學(xué)與信息學(xué)院，福建 福州，350117)

矩陣方程常常在多個領(lǐng)域出現(xiàn)，如控制理論[1-2]、系統(tǒng)理論[3-4]、穩(wěn)定性理論以及一些純應(yīng)用數(shù)學(xué)理論。因為矩陣方程[5-13]的應(yīng)用比較廣泛，所以，許多學(xué)者都致力于對它進(jìn)行研究。有很多方法可以用來求解矩陣方程，如共軛梯度法[14-15]、矩陣分解法等，本文提出了一種參數(shù)迭代法[16-19]來求解矩陣方程。還有一些其他方法，見參考文獻(xiàn)[20-24]。

設(shè)A，B，C，D∈Cn×n，考慮矩陣方程

AXB+CXD=F

(1)

方程(1)僅有一解的充要條件是BT?A+DT?C的特征值不為零。若B和C均為可逆矩陣，則矩陣方程(1)等價于

(C-1A)X+X(DB-1)=C-1FB-1

(2)

方程(2)僅有一解的充要條件是C-1A與DB-1沒有互為反號的特征值。此外，如果A，B，C和D都是埃爾米特正定矩陣時，那么，BT?A+DT?C也是埃爾米特正定矩陣，從而方程(1)有唯一解。

1 參數(shù)迭代法的建立與收斂性分析

選取實數(shù)α≠0，使αB+D和αC+A均為可逆矩陣，令

(3)

即

(4)

將方程(1)兩端乘以α，得

αAXB+αCXD=αF

X=MXN+Y0

(5)

易見，方程(5)與方程(1)等價。下面將給出相應(yīng)的參數(shù)迭代格式：

X(k+1)=MX(k)N+Y0，k=0，1，…

(6)

利用拉直算子將式(6)轉(zhuǎn)化為列向量形式

vec(X(k+1))=(NT?M)vec(X(k))+vec(Y0)，k=0，1，…

(7)

式(7)是一階線性定常迭代格式，它的迭代矩陣為NT?M，且有

ρ(NT?M)=ρ(NT)ρ(M)=ρ(M)ρ(N)。

有如下的收斂性定理。

定理1對任意給定的初始矩陣X(0)，格式(6)收斂的充要條件是ρ(M)ρ(N)<1。

證明當(dāng)矩陣C可逆時，設(shè)M的任一特征值為λM，對應(yīng)的特征向量為x，則有Mx=λMx，利用式(3)可求得

(8)

(9)

即M的特征值可由C-1A的特征值表示。

當(dāng)矩陣B可逆時，設(shè)N的任一特征值為λN，則λN也是NT的一個特征值，記對應(yīng)的特征向量為y，則有NTy=λNy，利用式(3)可求得

(10)

(11)

即N的特征值可由DB-1的特征值表示。

定理2設(shè)B和C都可逆，對任意給定的初始矩陣X(0)，格式(6)收斂的充要條件是

(12)

證明類似于定理1。

推論1若B和C是可逆矩陣，并且DB-1和C-1A的特征值實部都大(小)于零，則可得對任意正(負(fù))數(shù)α，格式(6)收斂。

證明由上述可知，若矩陣B可逆，且DB-1的特征值的實部都大于零，可得μj>0；若矩陣C可逆，且C-1A的特征值的實部都大于零，可得ηi>0；又由于α為正數(shù)，則有

即有格式(6)收斂。同理可證，當(dāng)B和C都可逆，并且DB-1和C-1A的特征值的實部全都小于零，這樣對任意負(fù)數(shù)α，格式(6)收斂。證畢。

從定理1可以得到，僅需選取適當(dāng)實數(shù)α，使得ρ(M)ρ(N)<1，那么格式(6)就能夠收斂，這樣子問題就轉(zhuǎn)化為怎么樣選取該實數(shù)α，才會使得格式(6)收斂最快，即ρ(M)ρ(N)最小，因此，將使得ρ(M)ρ(N)達(dá)到最小的實數(shù)α=αopt，稱為格式(6)的最優(yōu)參數(shù)。

2 最優(yōu)參數(shù)的選取

在這一小節(jié)中，將討論如何選取格式(6)的最優(yōu)參數(shù)，若矩陣A，B，C和D都是埃爾米特正定的。由于

故C-1A和DB-1的特征值全為正數(shù)。由上述推論1知道，對于任意的正數(shù)α，迭代格式(6)是收斂的。設(shè)C-1A和DB-1的特征值排序為

η1≥η2≥…≥ηn，μ1≥μ2≥…≥μn

(13)

(14)

(15)

由式(14)和(15)可得

(16)

當(dāng)?shù)玫降袷?6)的最優(yōu)參數(shù)即知道了αM和αN的取值時，便可計算C-1A和DB-1的最大和最小特征值。為了避免計算C-1和B-1，給出了最優(yōu)參數(shù)αM和αN的近似取法。

定理4設(shè)A，B，C和D都是埃爾米特正定矩陣，且它們的特征值依次為

a1≥a2≥…≥an，b1≥b2≥…≥bn，

c1≥c2≥…≥cn，d1≥d2≥…≥dn，

證明根據(jù)廣義特征值問題的極值原理，由式(8)可得(下面的最大值和最小值是針對0≠x∈Cn來求取的)：

(17)

(18)

(19)

易見ρ(M)≤q(α)。要讓ρ(M)最小，只需使q(α)最小即可。因此，可得使q(α)達(dá)到最小的實數(shù)為

據(jù)式(13)，正數(shù)η1，ηn，μ1，μn之間的關(guān)系可歸結(jié)為以下6種：①ηn≤μn≤η1≤μ1；②ηn≤μn≤μ1≤η1；③ηn≤η1≤μn≤μ1；④μn≤ηn≤μ1≤η1；⑤μn≤ηn≤η1≤μ1；⑥μn≤μ1≤ηn≤η1。

定理5在引理1的條件下，若

(20)

則αopt=αM；否則，αopt=αN。

證明對關(guān)系式①有αM≤αN，令

可得

改寫式(20)為

(21)

或者

由引理1中的關(guān)系式②的證明過程知αopt=αN。

同理可證，對關(guān)系式③～⑥，定理也成立。證畢。

3 加速的參數(shù)迭代法

下面討論關(guān)于迭代格式(6)的加速方法。

若取初始向量X(0)=Y0時，則由數(shù)學(xué)歸納法可得，迭代格式(6)可寫為

(22)

由定理1可知，當(dāng)ρ(NT?M)<1時，上述迭代格式收斂，且有

(23)

其中：X*是原方程的精確解。

(24)

其中：S(0)=X(0)=Y0。

于是由式(24)可得部分和的迭代格式為

S(k+1)=S(k)+M2kS(k)N2k，S(0)=Y0，k=0，1，…

(25)

4 數(shù)值實驗

下面將給出幾個矩陣方程數(shù)值例子來說明參數(shù)迭代法的有效性和可行性。以下2個實驗都是在Windows 10系統(tǒng)下Matlab R2014b的環(huán)境下運行的。

例1考慮矩陣方程AXB+CXD=I，其中，

則C-1A和DB-1的特征值分別為

η1=6.411 5，η2=1.815 2，η3=0.773 3，

μ1=1.153 5，μ2=0.660 5，μ3=0.064 0。

根據(jù)定理3求得αM=2.226 7，αN=0.271 8。由式(16)求得

ρ(M)|αM=0.484 5，ρ(N)|αN=0.618 7。

由定理5知αopt=αM。

B和D的特征值分別為

b1=5.247 0，b2=3.555 0，b3=2.198 1，

d1=4.791 3，d2=2.000 0，d3=0.208 7。

根據(jù)給出的A，B，C和D，分別用迭代格式(6)和(24)求解，假設(shè)終止條件為10-10，若對于不同的參數(shù)α，初始矩陣為X(0)=Y0，則數(shù)值結(jié)果見表1。

表1 不同α對于例1的數(shù)值結(jié)果Table 1 Numerical results about different α for Case 1

從表1可以看出，對于迭代格式(6)而言，取不同的參數(shù)α，迭代次數(shù)也發(fā)生變化，選取最優(yōu)參數(shù)，則可在一定程度上降低迭代所需的次數(shù)，α越接近最優(yōu)參數(shù)，迭代次數(shù)越少。而對于加速的迭代格式(24)而言，在同樣的參數(shù)取值下，它所需要的迭代次數(shù)明顯少于迭代格式(6)的迭代次數(shù)，這說明了加速的參數(shù)迭代法是有效的，并且可以提高迭代格式(6)的收斂速度。

例2設(shè)A，B，C和D都是分塊的三對角矩陣，即

其中：

表2 不同α對于例2的數(shù)值結(jié)果Table 2 Numerical results about different α for Case 2

從表2中可以發(fā)現(xiàn)：若取α=0.02，則迭代格式(6)需要迭代479次，但迭代格式(24)只需要迭代10次即可滿足精度要求；若取α=1.076 0，可知迭代格式(6)要迭代11次，而迭代格式(24)僅僅需要迭代4次就能滿足精度要求。因而，計算最優(yōu)參數(shù)有助于選取合適的α，從而降低迭代次數(shù)和時間。加速的參數(shù)迭代法有較快的收斂速度，從而提升了實驗效率。

5 結(jié)論

文中提出了一種新的迭代方法即參數(shù)迭代法求解矩陣方程AXB+CXD=F，并且給出了在適當(dāng)條件下最優(yōu)參數(shù)和近似最優(yōu)參數(shù)的選取方法。此外，還對參數(shù)迭代法進(jìn)行了加速處理。最后，給出幾個數(shù)值結(jié)果來說明所提迭代法是有效和可行的。