李貝貝, 崔靜靜, 黃政閣, 謝曉鳳

(廣西民族大學(xué)數(shù)學(xué)與物理學(xué)院, 廣西 南寧 530006)

1 引言

考慮求解如下大型稀疏線性方程組

其中A∈Cn×n是一個復(fù)對稱的非Hermitian 正定矩陣, 其形式為

其中W,T∈Rn×n分別是對稱正定和對稱半正定矩陣. 這里是虛數(shù)單位, 假設(shè)T≠0, 則A是非Hermitian 正定矩陣.

在科學(xué)與工程計算的領(lǐng)域中, 如固體力學(xué)、科學(xué)計算、動力學(xué)、工程計算、非線性規(guī)劃等[1-4]都需要求解復(fù)對稱線性系統(tǒng)(1). 目前, 求解復(fù)對稱線性系統(tǒng)(1) 常用的方法有直接法和迭代法兩大類. 當(dāng)系數(shù)矩陣階數(shù)較小時通?？刹捎弥苯臃ㄇ蠼? 但當(dāng)用直接法求解大型稀疏矩陣方程組時會破壞系數(shù)矩陣的稀疏性, 從而增加計算機的存儲量,降低計算效率. 而迭代法具有可充分利用和保持系數(shù)矩陣的稀疏性, 節(jié)省計算機存儲空間等優(yōu)點, 因此在求解大型稀疏線性系統(tǒng)時, 更傾向于使用迭代法. 近年來, 數(shù)值迭代方法在許多數(shù)學(xué)領(lǐng)域中都有極為廣泛的應(yīng)用, 例如: 四元數(shù)方向[5], 微分方程方向[6-7]等. 鑒于復(fù)對稱線性系統(tǒng)來源的廣泛性, 本文針對此類問題構(gòu)造一種有效的迭代算法,促進(jìn)相關(guān)領(lǐng)域的發(fā)展. 在求解復(fù)對稱線性系統(tǒng)的已有迭代法中, 其中最著名的迭代法是基于Hermitian 和反Hermitian 方法的迭代法. 首先對系數(shù)矩陣A進(jìn)行Hermitian 和反Hermitian 分裂:

這里H和S分別是矩陣A的Hermitian 和反Hermitian 部分. 根據(jù)A的Hermitian 和反Hermitian 分裂,白中治、Colub 和Ng 提出了求解復(fù)對稱線性系統(tǒng)(1)的Hermitian和反Hermitian 分裂(HSS) 迭代法[8]; 基于HSS 迭代方法, 白中治、Benzi 和陳芳提出了改良的HSS(MHSS) 迭代方法[9]; 潘春平、王紅玉和曹文芳提出了外推的HSS(EHSS)迭代方法[10]. 此外, 基于HSS 迭代法, 許多改進(jìn)及推廣的相應(yīng)迭代法被提出. 例如, LHSS 迭代方法[11], PHSS 迭代方法[12], PMHSS 迭代方法[13], LPMHSS迭代方法[14], AHSS 迭代方法[15], APMHSS 迭代方法[16]和QHSS 迭代方法[17]等.

潘春平等對HSS 方法使用外推技術(shù)得到EHSS 方法. EHSS 方法在迭代次數(shù)和計算時間方面都優(yōu)于HSS 迭代方法, 具體可參考文獻(xiàn)[10]. 下面簡單介紹一下EHSS 迭代方法.

EHSS 迭代方法:設(shè)A∈Cn×n是正定矩陣, 給定一個初始向量x(0)∈Cn, 計算

直到迭代序列{x(k)} 收斂, 這里α是給定的正常數(shù),ω是非負(fù)常數(shù),I是單位矩陣.

EHSS 迭代法可等價地表示為

其中這里M(α,ω) 是EHSS 迭代方法的迭代矩陣. 實際上, 迭代格式(2) 可由系數(shù)矩陣A進(jìn)行如下分裂得到

其中

HSS 迭代方法在求解線性系統(tǒng)時每一步迭代都需要求解一個位移反Hermitian 線性系統(tǒng), 這無疑是困難和耗費時間的. 為解決這個問題, 白中治、Benzi 和陳芳提出了改良的HSS(MHSS) 迭代方法, 文獻(xiàn)[14] 中的數(shù)值實驗結(jié)果表明MHSS 迭代法無論是在迭代次數(shù)還是計算時間上都優(yōu)于HSS 迭代方法. 下面簡單介紹一下MHSS 迭代方法.

MHSS 迭代方法:給定一個初始向量x(0)∈Cn, 計算

直到迭代序列{x(k)} 收斂, 這里α是給定的正常數(shù),I是單位矩陣.

迭代格式(3) 可以改寫為

其中

M(α) 是MHSS 迭代方法的迭代矩陣. MHSS 迭代方法也可以看做是由矩陣A進(jìn)行如下分裂所得到的

其中

然而, EHSS 迭代方法在求解線性系統(tǒng)時每一步迭代也都需要求解一個位移反Hermitian 線性系統(tǒng), 且MHSS 迭代法在求解某些復(fù)對稱線性系統(tǒng)時收斂速度比較慢, 為了克服這些缺點, 受EHSS 迭代法思想的啟發(fā), 考慮對MHSS 迭代方法使用外推技術(shù), 期望能得到比EHSS 和MHSS 方法更好的結(jié)果.

本文的具體布局為: 在第二節(jié)中構(gòu)造了外推的MHSS(EMHSS) 迭代方法; 第三節(jié)研究了EMHSS 迭代方法的迭代矩陣與MHSS 迭代方法的迭代矩陣之間的關(guān)系, 并討論了EMHSS 迭代方法的收斂條件; 第四節(jié)給出的數(shù)值實驗說明了本文方法的有效性。

2 EMHSS 迭代方法

受EHSS 迭代法思想的啟發(fā), 考慮對MHSS 迭代方法使用外推技術(shù), 構(gòu)造了外推的MHSS 迭代法, 簡稱為EMHSS 迭代法.

根據(jù)EHSS 方法, 考慮迭代格式(3) 中的第一步迭代

在上述等式兩邊同時乘以i, 并移項可得

將上式帶入格式(3) 中的迭代格式

可得

再結(jié)合

可得如下迭代格式

聯(lián)立MHSS 迭代法中的第一步迭代, 有

在(5) 式中引入松弛參數(shù)ω, 可得

結(jié)合迭代格式(3) 中的第一步迭代和(7) 式可得到如下的外推MHSS(EMHSS) 迭代法.

外推的MHSS(EMHSS) 迭代法:任給一個初始向量x(0)∈Cn, 計算

直到迭代序列{x(k)} 收斂, 這里α是給定的正常數(shù),ω是非負(fù)常數(shù),I是單位矩陣.當(dāng)ω=1 時, 該方法即為迭代格式(6).

經(jīng)過簡單的計算, 外推的MHSS(EMHSS) 迭代方法可表示為

其中

和

這里L(fēng)(α,ω) 是EMHSS 迭代法的迭代矩陣.

如果引入矩陣

和

則有下列等式成立

則EMHSS 迭代法(8) 可看作是對A進(jìn)行上述分裂(10) 所得到.

3 EMHSS 迭代方法的收斂性分析

定理3.1設(shè)A=W+iT∈Cn×n, 其中W∈Rn×n,T∈Rn×n分別是對稱正定矩陣和對稱半正定矩陣. 令α,ω分別是正實數(shù)和非負(fù)實數(shù),M(α) 和L(α,ω) 分別是MHSS迭代方法和EMHSS迭代方法的迭代矩陣, 則有

其中η=(1+i)ω-2i.

證明由(4) 式可得

則

因此, 根據(jù)上式和(9) 式可得

定理3.2設(shè)A=W+iT∈Cn×n, 這里W∈Rn×n,T∈Rn×n分別是對稱正定矩陣和對稱半正定矩陣. 假設(shè)A是正規(guī)矩陣, 則W和T滿足WT=TW. 設(shè)λj和μj(j=1,2,··· ,n) 分別是矩陣W和T的特征值. 設(shè)λmax和λmin分別為矩陣W的最大特征值和最小特征值,μmax為矩陣T的最大特征值. 如果0 ≤ω<2, 則有

(1) 若0<α≤1,μmax-λmin≤, 則ρ(L(α,ω))<1;

證明如果W∈Rn×n,T∈Rn×n滿足WT=TW, 則存在正交矩陣Q, 使得

其中

由(9) 式可知

ρ(L(α,ω))<1 成立需滿足不等式而

若2αλj(α+λj)(α+μj)>2α2(μ2j+λ2j) 成立, 即

則必有不等式(11) 成立.

(1) 當(dāng)0 <α≤1 時, 如果μmax-λmin≤, 必然滿足μj-λj≤λ2j, 則有不等式2αλj(α+λj-αλj)+2μj(αλj+λ2j-αμj)>0, 即ρ(L(α,ω))<1.

(2) 當(dāng)α> 1, 如果, 必然滿足, 則2αλj(α+λj-αλj)+2μj(αλj+λ2j-αμj)>0, 即ρ(L(α,ω))<1.

4 數(shù)值例子

本節(jié)通過數(shù)值實驗驗證EMHSS 迭代方法求解復(fù)對稱線性系統(tǒng)的有效性, 并在迭代次數(shù)(IT), 計算時間(CPU) 和相對殘差(RES) 方面將其結(jié)果與EHSS 和MHSS 方法進(jìn)行了比較. 所有數(shù)值實驗均在Intel(R)Core(TM)i7-8700 CPU @3.20GHz 和RMA 16.0GB,Win7 系統(tǒng)的電腦上用MATLAB R2018b 進(jìn)行的.

在本次實驗中，所有測試迭代法選取均以零向量x(0)=0作為初始迭代向量. 所有測試迭代法的終止準(zhǔn)則是當(dāng)k步的相對殘差滿足

時, 計算停止.

例4.1[9]Ax=(W+iT)x=b滿足線性系統(tǒng)(1), 其中

e1和em分別是第一個元素和第m個元素為1 的單位向量. 選取右端向b=(1+i)A1,這里1是所有元素都為1 的向量.

表1 列出了在不同問題規(guī)模下MHSS 迭代方法和EMHSS 迭代方法求解例4.1 時的實驗最優(yōu)參數(shù). 表2 給出了在表1 的參數(shù)下MHSS 迭代方法和EMHSS 迭代方法求解例4.1 時的迭代次數(shù)(IT), 計算時間(CPU) 和相對誤差(RES). 另外, 由于EHSS 迭代方法求解例4.1 時是不收斂的, 因此在表1 和表2 中未列出EHSS 迭代方法的實驗結(jié)果. 觀察表2 可知, 當(dāng)矩陣規(guī)模相同時, EMHSS 迭代方法的迭代次數(shù)(IT) 和計算時間(CPU) 都少于MHSS 迭代方法的計算時間和迭代次數(shù). 除了m= 64 外, EMHSS迭代方法的相對誤差(RES) 都要優(yōu)于MHSS 迭代方法. 因此, 從數(shù)值實驗的結(jié)果可知EMHSS 迭代法的計算效率要優(yōu)于MHSS 迭代法.

表1 EMHSS 和MHSS 迭代方法求解例4.1 時的最優(yōu)參數(shù)取值

表2 EMHSS 和MHSS 迭代方法求解例4.1 時的最優(yōu)參數(shù)取值

例4.2[11,18]考慮復(fù)Helmholtz 方程

其中σ1,σ2是實系數(shù)函數(shù),u在[0,1]×[0,1] 上滿足Dirichlet 邊界條件. 使用步長為的五點差分格式去處理上述方程, 可得如下類似線性系統(tǒng)(1.1) 的方程

這里

H是一個階數(shù)為n的塊三對角矩陣, 且n=m2, 選取右端向量b= (i-1)A1, 這里1是所有元素都為1 的向量. 通常在上述等式的兩邊同時乘h2來正規(guī)化系數(shù)矩陣. 在實際計算中取σ1=σ2=10, 則矩陣H+σ1I和矩陣σ2I是對稱正定的.

表3 列出了在不同問題規(guī)模下MHSS, EHSS 和EMHSS 迭代方法求解例4.2 時當(dāng)?shù)螖?shù)達(dá)到最少時參數(shù)α和ω的取值范圍. 在表4 中, 給出了MHSS 和EHSS迭代方法的迭代次數(shù)和計算時間均最少的實驗最優(yōu)參數(shù)取值, 及在最優(yōu)參數(shù)下的迭代次數(shù)(IT), 計算時間(CPU) 和相對誤差(RES). EMHSS 迭代方法中參數(shù)α,ω選取為α=10, 及在α=10 的情況下ω的最優(yōu)取值, 及在此情況下的迭代次數(shù)(IT), 計算時間(CPU) 和相對誤差(RES). 從表4 中可以看到當(dāng)問題規(guī)模相同時, 雖然EMHSS迭代方法的相對誤差與MHSS 和EHSS 的相對誤差幾乎相同, 但EMHSS 迭代方法的迭代次數(shù)(IT) 和計算時間(CPU) 要少于MHSS 和EHSS 迭代方法迭代次數(shù)(IT) 和計算時間(CPU). 因此可知當(dāng)α= 10 時, EMHSS 迭代方法已經(jīng)優(yōu)于MHSS 和EHSS迭代方法. 若α和ω都選取為最優(yōu)參數(shù)時, EMHSS 迭代方法必然有更好的數(shù)值實驗結(jié)果. 總之, EMHSS 迭代方法求解復(fù)對稱線性系統(tǒng)時的計算效率比MHSS 和EHSS 迭代方法的計算效率高.

表3 EMHSS, EHSS 和MHSS 迭代方法求解例4.2 時的最優(yōu)參數(shù)取值范圍

表4 求解例4.2 時EMHSS, EHSS 和MHSS 方法的迭代次數(shù), 計算時間和相對殘差

例4.3[9]考慮形為如下的方程

這里M和K是慣性矩陣和剛度矩陣,CV和CH是粘性矩陣和滯后衰減矩陣,?是驅(qū)動循環(huán)頻率. 令CH=εK, 這里ε是一個衰減系數(shù),M=I,CV= 10I,K是在均勻網(wǎng)格[0,1]×[0,1] 且網(wǎng)格大小為上使用中心差分格式沿著負(fù)Laplacian 算符方向均勻的逼近Dirichlet 邊界條件形成的矩陣. 這里

另外, 選取右端向量b=(1+i)A1, 這里1是所有元素都為1 的向量. 通常在上述等式的兩邊同時乘h2來正規(guī)化系數(shù)矩陣. 在實際計算中分別取?=1 和ε=0.002.

表5 給出了MHSS 和EMHSS 迭代方法求解例4.3 時迭代次數(shù)達(dá)到最少的參數(shù)的取值范圍. 表6 在表5 給出的參數(shù)范圍的基礎(chǔ)上, 展示了迭代時間最少時的參數(shù)值, 迭代次數(shù), 迭代時間和相對殘差. 從表6 可知隨著問題規(guī)模的增大MHSS 和EMHSS 迭代方法的迭代次數(shù)逐漸增大, 且EMHSS 迭代方法的迭代次數(shù)遠(yuǎn)小于MHSS 迭代方法的迭代次數(shù). 同時, EMHSS 迭代方法的迭代時間也少于MHSS 迭代方法的迭代時間.另外, EHSS 迭代方法求解例4.3 時的時間花費非常大, 最優(yōu)參數(shù)的選擇是困難的, 因此在表5 和表6 中未列出EHSS 迭代方法的數(shù)值實驗結(jié)果.

表5 EMHSS 和MHSS 迭代方法求解例4.3 時的最優(yōu)參數(shù)取值范圍

表6 求解例4.3 時EMHSS 和MHSS 方法的迭代次數(shù), 計算時間和相對殘差

此外, 為了更直觀地觀察EMHSS, MHSS 迭代方法對參數(shù)α的敏感性, 在圖1 中給出了當(dāng)m=64 時三個例子的迭代次數(shù)(IT) 隨參數(shù)α的變化曲線. 另外, 因為EHSS迭代方法求解例4.2 時是收斂的, 因此中間圖還展示了EHSS 的迭代次數(shù)隨參數(shù)α的變化曲線. 從圖1 可以看到當(dāng)參數(shù)ω不變時, MHSS 和EMHSS 迭代方法求解三個例子時迭代次數(shù)隨著參數(shù)α的增大均是先減小后增大. 且也可看出不管參數(shù)α如何選取,EMHSS 迭代法的迭代次數(shù)均小于MHSS 迭代法的迭代次數(shù), 且迭代次數(shù)最小值是遠(yuǎn)小于MHSS 迭代方法的最小值.

圖1 m=64 時EMHSS, MHSS 和EHSS 方法的迭代次數(shù)隨參數(shù)α 的變化圖

圖2 直觀地展示了EMHSS 迭代方法當(dāng)m=64 且參數(shù)α固定時(例4.1:α=0.49;例4.2:α= 10; 例4.3:α= 25) 單一變量ω對迭代次數(shù)的影響. 從圖2 可以觀察到三個例子的迭代次數(shù)均是在ω接近于0 時達(dá)到最少.

圖2 m=64 時EMHSS 方法的迭代次數(shù)隨參數(shù)ω 的變化圖

同時, 圖3 畫出了當(dāng)m=64 時, EMHSS 迭代方法求解這三個例子時, 迭代次數(shù)隨參數(shù)α和ω變化的三維圖像. 從圖3 的左側(cè)圖可看到EMHSS 迭代方法求解例4.1 時,當(dāng)α在0.5 附近時迭代次數(shù)達(dá)到最少. 從圖3 的中間圖可以看出EMHSS 迭代方法求解例4.2 時, 當(dāng)α在0.2 附近,ω在10 附近時迭代次數(shù)達(dá)到最少. 且從圖1 - 圖3 可知,EMHSS 迭代法求解復(fù)對稱線性系統(tǒng)時具有較高的計算效率.

圖3 EMHSS 迭代方法的迭代次數(shù)隨參數(shù)α 和ω 變化的三維圖像

5 總結(jié)

本文基于EHSS 迭代法的思想對MHSS 迭代法采用外推技術(shù)提出了外推的MHSS(EMHSS) 迭代方法. 并給出了EMHSS 迭代法的迭代矩陣和MHSS 迭代法的迭代矩陣之間的關(guān)系, 且討論了EMHSS 迭代法收斂的條件. 最后通過數(shù)值例子, 從迭代次數(shù),計算時間和相對誤差方面說明了EMHSS 迭代法的有效性. 另外, 還給出了EMHSS 迭代方法的迭代次數(shù)隨著參數(shù)α變化曲線圖, 以及參數(shù)α和ω對EMHSS 方法的迭代次數(shù)影響的三維變化圖, 說明了求解復(fù)對稱線性系統(tǒng)時EMHSS 迭代方法優(yōu)于EHSS和MHSS 迭代法.