外推的MHSS迭代法求解一類大型稀疏的復對稱線性系統

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

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

1 引言

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

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

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

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

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

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

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

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

EHSS 迭代法可等價地表示為

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

其中

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

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

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

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

其中

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

其中

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

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

2 EMHSS 迭代方法

受EHSS 迭代法思想的啟發, 考慮對MHSS 迭代方法使用外推技術, 構造了外推的MHSS 迭代法, 簡稱為EMHSS 迭代法.

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

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

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

可得

再結合

可得如下迭代格式

聯立MHSS 迭代法中的第一步迭代, 有

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

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

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

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

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

其中

和

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

如果引入矩陣

和

則有下列等式成立

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

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

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

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

證明由(4) 式可得

則

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

定理3.2設A=W+iT∈Cn×n, 這里W∈Rn×n,T∈Rn×n分別是對稱正定矩陣和對稱半正定矩陣. 假設A是正規矩陣, 則W和T滿足WT=TW. 設λj和μj(j=1,2,··· ,n) 分別是矩陣W和T的特征值. 設λ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) 當0 <α≤1 時, 如果μmax-λmin≤, 必然滿足μj-λj≤λ2j, 則有不等式2αλj(α+λj-αλj)+2μj(αλj+λ2j-αμj)>0, 即ρ(L(α,ω))<1.

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

4 數值例子

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

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

時, 計算停止.

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

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

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

表1 EMHSS 和MHSS 迭代方法求解例4.1 時的最優參數取值

表2 EMHSS 和MHSS 迭代方法求解例4.1 時的最優參數取值

例4.2[11,18]考慮復Helmholtz 方程

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

這里

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

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

表3 EMHSS, EHSS 和MHSS 迭代方法求解例4.2 時的最優參數取值范圍

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

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

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

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

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

表5 EMHSS 和MHSS 迭代方法求解例4.3 時的最優參數取值范圍

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

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

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

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

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

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

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

5 總結

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