求解正定線性方程組的外推的PSS迭代方法

吳思婷， 鮑 亮， 黃景宣

（ 華東理工大學(xué)數(shù)學(xué)學(xué)院，上海 200237）

在科學(xué)與工程的很多領(lǐng)域，如固體力學(xué)、動力學(xué)、非線性規(guī)劃等方向[1-7]都需要求解如式（1）大型稀疏線性方程組

其中A∈Cn×n為正定矩陣，x∈Cn×n為未知向量，b∈Cn為已知向量。目前關(guān)于此類問題的求解主要分為直接法和迭代法。當(dāng)式(1)的系數(shù)矩陣規(guī)模較小時，直接法更加方便。但當(dāng)式(1)的系數(shù)矩陣規(guī)模較大時，直接法在計算過程中可能會產(chǎn)生大量的非零元，破壞矩陣的稀疏性，大大增加了計算機的儲存量和計算量。所以對于較大規(guī)模的稀疏線性方程組，通常使用迭代法進行求解。

2003年，Bai等[8]基于Hermitian和反Hermitian(HS)分裂，提出了一種求解非Hermitian正定線性方程組的Hermitian和反Hermitian(HSS)迭代方法。該方法由于形式簡單和無條件收斂的優(yōu)點，受到了許多學(xué)者的關(guān)注并進行了深入地研究，提出了一系列有效求解大型稀疏線性方程組的迭代方法[9-11]。后來Bai等[12]將系數(shù)矩陣A進行了正定和反Hermitian(PS)分裂，提出了一種求解正定線性方程組的PSS迭代方法。PSS迭代方法不僅可以求解非Hermitian正定線性方程組，還可以有效地求解Hermitian正定線性方程組。Cao等[13]提出了一種有效的廣義PSS(GPSS)迭代方法，它繼承了HSS迭代方法、PSS迭代方法和GHSS迭代方法的優(yōu)點。2019年，潘春萍等[14]運用外推的技術(shù)對HSS迭代方法進行加速，提出了一種求解非Hermitian正定線性方程組的外推的HSS(EHSS)迭代方法，大大加快了迭代方法的收斂速度。

HS分裂僅適用于非Hermitian正定矩陣，而PS分裂對于任意正定矩陣都有著很好的效果，同時外推的技術(shù)能夠?qū)⒌椒ǖ氖諗克俣却蟠蠹涌欤摷夹g(shù)目前僅用于非Hermitian正定線性方程組的求解，還未有效地用于Hermitian正定線性方程組的求解?；谏鲜霭l(fā)現(xiàn)，如果將對系數(shù)矩陣進行PS分裂，將外推的技術(shù)有效地推廣到大型稀疏正定線性方程組的求解中，這種方法不僅適用于非Hermitian正定線性方程組的求解，還可用于Hermitian正定線性方程組的求解，大大加快原有迭代方法的收斂速度。

本文將線性方程組的系數(shù)矩陣分解成正定矩陣和反Hermitian矩陣，再進行非對稱二步迭代，得到了一種求解正定線性方程組的外推的正定和反Hermitian(EPSS)迭代方法，并且給出了該方法收斂的充要條件；同時理論分析了滿足一定條件時，EPSS迭代方法的譜半徑比PSS迭代方法的譜半徑更小；最后數(shù)值實驗表明，通過參數(shù)值的選擇，EPSS迭代法比PSS迭代法的收斂速度更快和迭代次數(shù)更小，且選擇合適的參數(shù)值時可以大大提高新方法的收斂效率。

1 EPSS迭代方法

正定和PS分裂：任意的正定矩陣A∈Cn×n都有如下形式分解：

A=P+S

其中，P∈Cn×n為正定矩陣；S∈Cn×n為反Hermitian矩陣。關(guān)于矩陣P和S的選取，有如下兩種常見的PS分裂格式：

基于上述PS分裂，Bai等[12]給出如下求解正定線性方程組的PSS迭代方法。

1.1 PSS迭代方法

設(shè)A∈Cn×n為正定矩陣，給一個定初始向量x(0)∈Cn，對k=于0,1,直2···,到K迭代序{列x(k)}收斂，計算：

其中， α>0 為給定的常數(shù)，I為單位矩陣。將式(2)寫成如下矩陣-向量形式：

其中：M(α)=(αI+S)?1(αI?P)(αI+P)?1(αI?S) 和N(α)=2α(αI+S)?1(αI+P)?1，M(α) 為PSS方法的迭代矩陣。

為了加快迭代速度，本文將對稱的二步迭代格式進行變形，推導(dǎo)出一種新的非對稱的二步迭代格式。首先，將式(2)中的(a)式進行變形，得到

再將式(4)代入式(2b)中消去矩陣P，得

接著利用

(αI+S)x=(S?αI)x+2αx=(S+αI)x

可將式(5)變成如下形式

最后在式(6)中引入一個新的松弛變量 ω ，即：

由此可得求解大型稀疏正定線性方程組的外推的PSS迭代方法，記為EPSS迭代方法。下面給出EPSS迭代方法具體的迭代格式。

1.2 EPSS方法的迭代過程

設(shè)A∈Cn×n為正定矩陣，給定一個初始向量x(0)∈Cn，對k于=0,1直,2,··到·迭代序{列x(k)}收斂，計算：

其中， α>0 ， ω>0 均為給定的常數(shù)。

當(dāng) ω=0 時，EPSS迭代方法就成為了PSS迭代方法；當(dāng) ω=1 時，EPSS迭代方法就成為了式(6)的方法。

可以將式(7)改寫成式(8)的等價的矩陣-向量形式：

其中：M(α,ω)=(αI+S)?1[S?(1?ω)αI+(2?ω)α(αI+P)?1(αI?S)]=(αI+S)?1(αI+P)?1{(αI+P)[S?(1?ω)αI] +(2?ω)α(αI?S)}= (αI+S)?1(αI+P)?1{PS?(1?ω)α(P+S)+α2I}，N(α,ω)=(2?ω)α(αI+S)?1(αI+P)?1，這里的M(α,ω)是EPSS方法的迭代矩陣。

通過計算發(fā)現(xiàn)EPSS迭代方法也是系數(shù)矩陣A經(jīng)過以下分裂得到：

A=B(α,ω)?C(α,ω),

2 收斂性分析

本文根據(jù)得到的迭代格式對EPSS迭代方法的收斂性進行了詳細的分析，給出了迭代格式收斂的充要條件。下面先給出PSS迭代方法收斂的相關(guān)引理。

引理1[12]設(shè)A∈Cn×n為正定矩陣，令A(yù)=P+S，P和S分別為正定矩陣和反Hermitian矩陣，M(α)為PSS方法的迭代矩陣，V(α)=(αI?P)(αI+P)?1，則迭代矩陣M(α) 的譜半徑 ρ(M(α)) 以 //V(α)//2為上界，因此有

ρ(M(α))≤//V(α)//2＜1,?α>0,

即PSS迭代方法收斂到線性方程組式(1)的精確解x?∈Cn。

根據(jù)上述EPSS方法的迭代矩陣M(α,ω) 的形式，不難發(fā)現(xiàn)EPSS迭代方法和PSS迭代方法迭代矩陣之間存在一定的聯(lián)系。接下來本文進一步詳細地分析EPSS方法的迭代矩陣M(α,ω) 與PSS方法的迭代矩陣M(α) 之間的關(guān)系。

定理1設(shè)A∈Cn×n為正定矩陣，P和S分別為矩陣A的正定部分和反Hermitian部分， α>0 和ω≥0為兩個參數(shù)M，(α)為PSS方法的迭代矩陣，則EPSS方法的迭代矩陣M(α,ω)滿足：

證明：

ωI+(2?ω)M(α)=ωI+(2?ω)(αI+S)?1(αI?P)(αI+P)?1·(αI?S)=(αI+S)?1(αI+P)?1[ω(αI+P)(αI+S)+(2?ω)(αI+P)(αI?P)(αI+P)?1(αI?S)]=(αI+S)?1·(αI+P)?1[ω(αI+P)(αI+S)+(2?ω)(αI?P)(αI+P)·(αI+P)?1(αI?S)]=(αI+S)?1(αI+P)?1[ω(α2I+αP+αS+PS)+(2?ω)(αI?P)(αI?S)]=(αI+S)?1(αI+P)?1·[2PS+2α2I?2(1?ω)α(P+S)]

可得：

由引理1和定理1可得EPSS迭代方法收斂的充要條件。

定理2設(shè)A∈Cn×n為正定矩陣，P和S分別為矩陣A的正定部分和反Hermitian部分，對任意的 α>0 ，若 ω 滿足 0≤ω＜2 ，則EPSS迭代方法收斂[15-16]。

證明：根據(jù)定理1可知PSS方法的迭代矩陣M(α) 和EPSS方法的迭代矩陣M(α,ω) 存在如下關(guān)系：

設(shè) λ 和 μ 分別為PSS方法的迭代矩陣M(α) 和EPSS方法的迭代矩陣M(α,ω) 的特征值，則

若 λ 為實數(shù)，則

若 λ 為復(fù)數(shù)，設(shè) λ=a+bi，i為虛數(shù)單位，則

綜上所述可得：

根據(jù)引理1可知 ρ(M(α))＜1,?α>0 ，因此當(dāng) 0≤ω＜2時，即0＜2?ω≤2時，得ρ(M(α,ω))≤[ω+(2?ω)ρ(M(α))]＜[ω+(2?ω)]=1，故EPSS迭代方法收斂。

由定理2可知，EPSS迭代方法是收斂的，并且EPSS迭代方法的收斂與參數(shù) α 和 ω 的選取有關(guān)，不是無條件收斂的。除此之外，在定理2的證明中發(fā)現(xiàn)，EPSS方法的迭代矩陣M(α,ω) 的特征值與PSS方法的迭代矩陣M(α) 特征值存在一定的聯(lián)系。下面本文繼續(xù)通過理論分析EPSS方法的迭代矩陣M(α,ω) 的譜半徑與PSS方法的迭代矩陣M(α)譜半徑之間的大小關(guān)系，進一步判斷兩種迭代方法收斂速度的快慢。

證明：根據(jù)定理2可得

其中， μ 為EPSS方法的迭代矩陣M(α,ω) 的特征值。為了便于求解，令

f(ω)=ω2+2ω(2?ω)a+(2?ω)2(a2+b2)?4(a2+b2)

當(dāng)f(ω)＜0時，可得

ω[(1?a)2+b2]＜(2a?1)2+4b2?1

根據(jù)定理3可知，當(dāng)PSS方法的迭代矩陣M(α)的特征值以及松弛變量 ω 滿足一定條件時，EPSS迭代方法不僅收斂，而且迭代矩陣的譜半徑小于PSS迭代方法的譜半徑。由于迭代矩陣的譜半徑越小，則迭代速度越快，迭代次數(shù)越少。因此本文通過理論證明，在處理某些問題時，選取合適的參數(shù)值，EPSS迭代方法會比PSS迭代方法有更好的效果。

3 數(shù)值實驗

本文對EPSS迭代方法進行了數(shù)值實驗，并將其結(jié)果與PSS迭代方法進行比較。本次數(shù)值實驗是在英特爾四核處理器（2.40 GHz，8 GB RAM）的環(huán)境下應(yīng)用Matlab編程語言進行的。

例1 考慮區(qū)域為 ?=[0,1] 的一類微分方程[9]：

?u′′+qu′=f

邊界條件為齊次邊界條件。利用中心差分離散方程得到系數(shù)矩陣A有如下格式：

A=tridiag(?1+qh/2,2,?1?qh/2)∈Rn×n

其中，h=為離散網(wǎng)格值，N為矩陣規(guī)模；q為常數(shù)； tridiag(·,·,·) 表示三對角矩陣。

例2 考慮區(qū)域為 ?=[0,1]×[0,1]×[0,1] 的一類三維對流擴散方程[8]：

?(uxx+uyy+uzz)+q(ux+uy+uz)=f(x,y,z)

其中，q為常數(shù)，并且區(qū)域邊界滿足Dirichlet邊界條件。下面采用向前差分離散方程得到的系數(shù)矩陣為如下形式的線性方程組：

A=Tx?I?I+I?Ty?I+I?I?Tz∈Rn3×n3,

下面通過改變系數(shù)矩陣的矩陣規(guī)模N、參數(shù) α和 ω ，得到EPSS迭代方法和PSS迭代方的譜半徑、迭代次數(shù)和迭代時間，并將得到的數(shù)據(jù)進行比較分析，最終得到3種迭代方法的收斂速度之間的比較。下文數(shù)值實驗中的參數(shù) α 和 ω 均為數(shù)值實驗得出的近似最優(yōu)值。

3.1 收斂速度分析

表1給出了例1中qh分別為 10、100、1000時，PSS迭代方法和EPSS迭代方法的譜半徑、收斂迭代次數(shù)和時間的比較。由表1可以發(fā)現(xiàn)，當(dāng)矩陣規(guī)模相同時，EPSS迭代方法的譜半徑、迭代次數(shù)和迭代時間基本上都小于PSS迭代方法的譜半徑、迭代次數(shù)和迭代時間，尤其是當(dāng)矩陣規(guī)模N不斷變大時，EPSS迭代方法的譜半徑、迭代次數(shù)和迭代時間都明顯小于PSS迭代方法的譜半徑、迭代次數(shù)和迭代時間。因此在同一精度范圍內(nèi)，選取合適的參數(shù)后，EPSS迭代方法比PSS迭代方法會有更好的效果。

表1 例1中PSS和EPSS迭代方法譜半徑，收斂迭代次數(shù)和時間Table 1 Spectral radius, IT and CPU of PSS and EPSS iterative method for example 1

表2給出了例2取q=1000 且 ω=0.1 時，PSS迭代方法、EPSS迭代方法和EHSS迭代方法的譜半徑、收斂迭代次數(shù)和時間的比較。本次實驗選取ω=0.1，但由于PSS迭代方法的迭代格式中沒有變量，ω故該方法沒有標(biāo)注ω的取值。由表2可以發(fā)現(xiàn)，當(dāng)矩陣規(guī)模相同時，EPSS迭代方法的譜半徑、迭代次數(shù)和迭代時間基本上都小于PSS迭代方法和EHSS迭代方法的譜半徑、迭代次數(shù)和迭代時間。因此在同一精度范圍內(nèi)，選取合適的參數(shù)后，EPSS迭代方法比PSS迭代方法和EHSS迭代方法會有更好的效果。

表2 例2中PSS、EPSS和EHSS迭代方法譜半徑，收斂迭代次數(shù)以及時間（ q =1 000）Table 2 Spectral radius, IT and CPU of PSS, EPSS and EHSS iterative method for example 2 when q=1000

為了直觀地展示新方法與其他方法的迭代速度的比較，下面給出EPSS迭代方法與其他相關(guān)迭代方法的殘量下降速度的比較。圖1展示了例1中系數(shù)矩陣規(guī)模為 512×512 階且qh為100和1000時兩種迭代方法的殘量隨著迭代步數(shù)的變化趨勢。由圖1可以看出，兩種EPSS迭代方法的殘量下降曲線位于PSS迭代方法的殘量下降曲線下方。這也可以說明，當(dāng)系數(shù)矩陣規(guī)模為 512×512 階且qh為100、1000時，EPSS迭代方法的收斂速度更快，并且當(dāng)qh為1000時，EPSS迭代方法的迭代次數(shù)遠小于PSS迭代方法的收迭代次數(shù)。因此在求解某些問題時，EPSS迭代方法相較于PSS迭代方法具有迭代次數(shù)更少，迭代時間更短等優(yōu)點。

圖1 例1中EPSS和PSS迭代方法的殘量速度比較Fig. 1 Comparison of the residuals speed in EPSS and PSS for example 1

圖2展示了例2中當(dāng)參數(shù)q=1000 且 ω=0.1 時，系數(shù)矩陣規(guī)模分別為 1728×1728 、 2744×2744 階的EPSS、PSS和EHSS 3種迭代方法的殘量隨著迭代步數(shù)的變化趨勢，顯然EPSS迭代方法的殘量下降曲線位于PSS迭代方法和EHSS迭代方法的殘量下降曲線之下。這也可以說明，EPSS迭代方法的收斂速度更快。因此在求解某些問題時，EPSS迭代方法相較于PSS迭代方法以及EHSS迭代方法具有迭代次數(shù)更少，迭代時間更短等優(yōu)點。

圖2 例2中取不同矩陣規(guī)模時，EPSS、PSS和EHSS迭代方法的殘量速度比較Fig. 2 Comparison of the residuals speed in EPSS, PSS and EHSS for example 2 when matrix size is different

3.2 參數(shù) α 和 ω 的靈敏度分析

迭代矩陣的譜半徑對迭代方法收斂速度起著十分重要的影響，因此本節(jié)主要分析不同數(shù)值例子中，矩陣規(guī)模N、參數(shù) α 和 ω 對迭代矩陣譜半徑的影響。首先分析單一松弛變量 α 對PSS迭代方法和EPSS迭代方法中的譜半徑的影響，再通過分析變量α 和 ω 共同作用對EPSS迭代方法的譜半徑的影響，進而得到參數(shù) α 和 ω 的靈敏度的分析。參數(shù)的靈敏性分析例1中的參數(shù) α 和 ω 對EPSS迭代方法的譜半徑的影響。

圖3(a)示出了系數(shù)矩陣規(guī)模為 512×512 階，qh=100且ω=0.6時，PSS和EPSS迭代方法中的參數(shù)α與譜半徑的關(guān)系。圖3(b)示出了系數(shù)矩陣規(guī)模為階512×512，qh=1000且ω=0.7時，PSS和EPSS迭代方法中的參數(shù)α與譜半徑的關(guān)系。由圖3可以看出，EPSS迭代方法和PSS迭代方法的譜半徑均隨著 α 先減小后增大，并且隨著 α 不斷增大，EPSS迭代方法和PSS迭代方法的譜半徑曲線逐漸趨于重合。并且當(dāng) α 和 ω 處于一個適當(dāng)?shù)娜≈捣秶鷥?nèi)，EPSS迭代方法的譜半徑小于PSS迭代方法的譜半徑，且EPSS迭代方法的最優(yōu)譜半徑要小于PSS迭代方法的最優(yōu)譜半徑（小于1），這說明了EPSS迭代方法具有比PSS迭代方法更好的迭代效果。

圖3 例1中EPSS和PSS迭代方法的譜半徑與α的關(guān)系Fig. 3 Relationship between the parameter α and the spectral radius in EPSS and PSS for example 1

圖4示出了系數(shù)矩陣規(guī)模為 512×512 階且qh分別為100、1000時，EPSS迭代方法中的參數(shù) α 和ω 與譜半徑的關(guān)系。通過觀察圖4發(fā)現(xiàn)，當(dāng) α 和 ω 處于一個適當(dāng)?shù)娜≈捣秶鷥?nèi)，EPSS迭代方法的譜半徑小于1，該實驗結(jié)果符合定理2中收斂的證明。當(dāng)數(shù)值例子滿足qh=100 、且 α=3.9 和 ω= 0.6 時，EPSS迭代方法的譜半徑近似取最小值（小于1）。當(dāng)數(shù)值例子滿足qh=1000 、且 α=4.7 和 ω= 0.7 時，EPSS迭代方法的譜半徑近似取得最小值（小于1）。因此圖2分別選取 ω=0.6 和 ω=0.7 是合理的。

圖4 例1中EPSS迭代方法中參數(shù) α 和 ω 與譜半徑的關(guān)系Fig. 4 Relationship between the parameters α , ω and the spectral radius in the EGHSS iterative method for example 1

4 總 結(jié)

本文探討了求解正定線性方程組的EPSS迭代方法，給出了該迭代方法收斂的充要條件，并進行了數(shù)值實驗。通過與PSS迭代方法的譜半徑、迭代步數(shù)和迭代時間的比較，EPSS迭代方法具有迭代次數(shù)少，迭代時間短等優(yōu)點。本文還分析了參數(shù) α 和 ω對EPSS迭代方法譜半徑的影響，當(dāng)然理論上最優(yōu)參數(shù) α 和 ω 的確定還有待進一步的研究。總的來說，EPSS迭代方法是一種相較于PSS迭代方法更有競爭力的方法，特別是選取了合適的參數(shù) α 和 ω 后，EPSS迭代方法的收斂速度可以大大提高。