非Hermitian正定线性方程组的外推的广义HSS方法*

吴思婷，鲍 亮

(华东理工大学数学学院,上海 200237)

1 引言

大型稀疏线性方程组在许多科学与工程领域都有着非常重要的作用，比如在最优化问题、流体力学和电磁学等领域都需要对此类问题进行求解[1-7]。因此，如何快速有效地求解大型稀疏线性方程组已经成为当下非常重要的课题之一。

本文考虑如式(1)所示的大型稀疏正定线性方程组的求解：

Ax=b

(1)

其中,A∈Cn×n为非Hermitian正定矩阵，x∈Cn为未知向量，b∈Cn为已知向量。当方程组(1)的系数矩阵规模较小时，使用直接法更加方便。而对于较大规模的稀疏非Hermitian正定线性方程组，研究人员一般采用迭代法进行求解。

在求解某些特定问题时，求解系数矩阵αI+S的线性方程组可能会出现一些问题，甚至可能会同求解原始的线性方程组Ax=b一样困难。而选择恰当的矩阵K会使得矩阵αI+K+S比矩阵αI+S更对角占优，因此求解系数矩阵αI+K+S的线性方程组会比求解系数矩阵αI+S的线性方程组更加方便。基于上述发现，本文将系数矩阵进行广义的HS分裂，再进行非对称的二步迭代。这可以将二步迭代中求解系数矩阵αI+S的线性方程组简化为求解系数矩阵αI+K+S的线性方程组，大大加快了迭代方法的收敛速度。

本文将系数矩阵进行广义的Hermitian和反Hermitian分裂，再通过引入新的变量同时结合外推的技术，给出一种外推的广义Hermitian和反Hermitian迭代方法EGHSS(Extrapolated and Generalized HSS)。还理论分析了本文方法的收敛性，给出了该方法收敛的充要条件，并将该方法迭代矩阵的谱半径与GHSS方法迭代矩阵的谱半径和EHSS方法迭代矩阵的谱半径进行比较。数值实验结果表明，在处理某些实际问题中，EGHSS迭代方法比GHSS迭代方法和EHSS迭代方法更有效，且在合适的参数值下新方法的收敛效率可以大大提高。

2 EGHSS迭代方法

首先简要回顾EHSS迭代方法[13]。

方法1EHSS迭代方法。

(2)

其中,α>0和ω>0均为给定的常数。

事实上，可以将式(2)改写成如式(3)～式(5)所示的等价形式：

x(k+1)=M1(α,ω)x(k)+N1(α,ω)b,

k=0,1,2,…

(3)

[HS-(1-ω)α(H+S)+α2I]

(4)

(5)

这里的M1(α,ω)是EHSS迭代方法的迭代矩阵。

考虑到上述迭代方法中求解系数矩阵为αI+S的线性方程组仍有一定的难度，为了将其简化为求解系数矩阵αI+K+S的线性方程组，本文将HS分裂中的Hermitian矩阵H作如式(6)所示的分裂：

H=G+K

(6)

其中，G和K为Hermitian半正定矩阵，从而得到系数矩阵A∈Cn×n有如式(7)所示广义的HS分裂：

A=G+K+S

(7)

其中，G和K为Hermitian半正定矩阵，S为反Hermitian矩阵。一种常见的矩阵分裂格式如式(8)所示：

(8)

基于上述系数矩阵的广义HS分裂，下面介绍EGHSS迭代方法。

方法2EGHSS迭代方法。

设A∈Cn×n为非Hermitian正定矩阵，给定一个初始向量x(0)∈Cn，对于k=0,1,2,…, 直到迭代序列{x(k)}收敛，计算如式(9)所示：

(9)

其中,α>0和ω>0均为给定的常数。

当ω=0时，EGHSS迭代方法就成为GHSS迭代方法；当K=0时，EGHSS迭代方法就成为EHSS迭代方法；当ω=0且K=0时，EGHSS迭代方法就成为HSS迭代方法。

观察EHSS方法的迭代格式，不难发现式(9)等价式(10)～式(12)所示的矩阵-向量形式：

x(k+1)=M(α,ω)x(k)+N(α,ω)b,k=0,1,2,…

(10)

M(α,ω)=(αI+K+S)-1[K+S-(1-ω)αI+

(2-ω)α(αI+G)-1(αI-K-S)]=

(αI+K+S)-1(αI+G)-1{(αI+G)[K+

S-(1-ω)αI]+(2-ω)α(αI-K-S)}=

(αI+K+S)-1(αI+G)-1{G(K+S)-

(1-ω)α(G+K+S)+α2I}

(11)

(12)

这里的M(α,ω)为EGHSS迭代方法的迭代矩阵。

当然EGHSS迭代方法也可由系数矩阵A经过如式(13)～式(15)所示的分裂得到：

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

(13)

(14)

(1-ω)α(G+K+S)+α2I]=

(15)

3 收敛性分析

本节讨论了EGHSS迭代方法的收敛性，并给出了迭代方法收敛的充要条件。

引理1[12]设A∈Cn×n为正定矩阵，令A=G+K+S,G和K为Hermitian半正定矩阵，S为反Hermitian矩阵，M(α)为GHSS迭代方法的迭代矩阵。如果G或K为正定矩阵，且α>0，则迭代矩阵M(α)的谱半径ρ(M(α))满足式(16)：

(16)

即GHSS迭代方法收敛到线性方程组式(1)的精确解x*∈Cn。

从上述EGHSS方法的迭代格式推导过程不难发现，EGHSS方法的迭代矩阵与GHSS方法的迭代矩阵存在一定的关系。为了证明EGHSS迭代方法的收敛性，本节根据EGHSS迭代方法的迭代矩阵M(α,ω)的形式，进一步分析了EGHSS迭代方法的迭代矩阵M(α,ω)与GHSS迭代方法的迭代矩阵M(α)之间的关系。

定理1设A∈Cn×n为非Hermitian正定矩阵，令A=G+K+S,G和K为Hermitian半正定矩阵，S为反Hermitian矩阵，对任意的α>0和ω≥0，GHSS迭代方法的迭代矩阵如式(17)所示：

(17)

那么EGHSS迭代方法的迭代矩阵M(α,ω)满足式(18)：

(18)

证明

ωI+(2-ω)M(α)=

ωI+(2-ω)(αI+K+S)-1(αI-G)

(αI+G)-1(αI-K-S)=

(αI+K+S)-1(αI+G)-1[ω(αI+G)

(αI+K+S)+(2-ω)(αI+G)(αI-G)(αI+G)-1

(αI-K-S)]=(αI+K+S)-1(αI+G)-1

[ω(αI+G)(αI+K+S)+

(2-ω)(αI-G)(αI+G)(αI+G)-1

(αI-K-S)]=(αI+K+S)-1(αI+G)-1

[ω(α2I+αG+αK+αS+GK+GS)+(2-ω)

(αI-G)(αI-K-S)]=(αI+K+S)-1(αI+G)-1

[2GK+2GS+2α2I-2(1-ω)α(G+K+S)]

从而可以得到式(19)：

[G(K+S)-(1-ω)α(G+K+S)+α2I]=

(19)

□

由引理1和定理1可以得出EGHSS迭代方法收敛的充要条件。

定理2[14,15]设A∈Cn×n为非Hermitian正定矩阵，令A=G+K+S,G和K为Hermitian半正定矩阵，S为反Hermitian矩阵，如果0≤ω<2，那么对任意的α>0，EGHSS方法的迭代格式式(9)收敛。

证明由定理1知EGHSS迭代方法的迭代矩阵M(α,ω)与GHSS迭代方法的迭代矩阵M(α)之间存在如式(20)所示的关系：

(20)

设λ和μ分别为GHSS迭代方法的迭代矩阵M(α)和EGHSS迭代方法的迭代矩阵M(α,ω)的特征值，则有式(21):

(21)

若λ为实数，则式(22)成立:

(22)

若λ为复数，设λ=a+bi,i为虚数单位，则式(23)成立:

(23)

综上所述,式(24)成立：

(24)

又由引理1知∀α>0，ρ(M(α))<1,所以当0≤ω<2时，即0<2-ω≤2时，有式(25):

(25)

因此,EGHSS迭代方法收敛。

□

由定理2可知，EGHSS迭代方法是条件收敛的。由定理2的证明过程可知，EGHSS方法的迭代矩阵M(α,ω)的特征值与GHSS方法的迭代矩阵M(α)的特征值存在着一定的联系。接下来进一步分析EGHSS方法的迭代矩阵M(α,ω)的谱半径与GHSS方法的迭代矩阵M(α)谱半径之间的关系。

定理3设A∈Cn×n为非Hermitian正定矩阵，令A=G+K+S,G和K为Hermitian半正定矩阵，S为反Hermitian矩阵，λ=a+bi为GHSS方法的迭代矩阵M(α)的特征值，则有：

证明根据定理2可得式(26):

(26)

(27)

(28)

此时|μ|<|λ|，于是当a2+b2<1时，有0≤ω<2，则式(29)成立:

ρ(M(α,ω))≤ρ(M(α))<1

(29)

ρ(M(α))≤ρ(M(α,ω))

(30)

ρ(M(α))≤ρ(M(α,ω))<1

(31)

□

根据定理3可知，通过选择合适的参数，EGHSS方法的迭代矩阵M(α,ω)的谱半径会小于GHSS方法的迭代矩阵M(α)的谱半径。下一节将通过具体的数值实验对EGHSS迭代方法、GHSS迭代方法和EHSS迭代方法的收敛速度和迭代次数进行比较，并详细分析矩阵规模、α和ω对EGHSS迭代方法收敛速度的影响。

4 数值实验

本节对EGHSS迭代方法进行了数值实验，并将其结果与GHSS迭代方法和EHSS迭代方法进行比较。本次数值实验是在英特尔四核处理器(2.40 GHz，8 GB RAM)环境下应用Matlab编程语言实现的。

例1考虑文献[16]中的数值例子,如式(32)所示：

-u″+qu′=f

(32)

边界条件为齐次边界条件。根据中心差分方法对式(32)进行离散，可以得到如式(33)所示的系数矩阵：

A=tridiag(-1+qh/2,2,-1-qh/2)∈RN×N

(33)

其中，tridiag(·,·,·)表示三对角矩阵。

例2考虑区域为Ω=[0,1]×[0,1]×[0,1]的一类三维对流扩散方程[8],如式(34)所示：

-(uxx+uyy+uzz)+q(ux+uy+uz)=

f(x,y,z)

(34)

其中,q=1000为给定常数，且区域边界满足Dirichlet边界条件。下面采用向前差分离散方程得到如式(35)所示的系数矩阵：

A=Tx⊗I⊗I+I⊗Ty⊗I+

I⊗I⊗Tz∈RN3×N3

(35)

4.1 收敛速度分析

表1～表3分别给出了例1中GHSS迭代方法、EHSS迭代方法和EGHSS迭代方法的收敛迭代次数和收敛时间，并且给出了3种迭代方法近似最优变量的取值。其中N的取值表示的矩阵规模分别为256×256阶，512×512阶,1024×1024阶，2048×2048阶。从表1～表3可以发现，当数值例子相同且矩阵规模相同时，EGHSS迭代方法的谱半径、迭代次数和迭代时间普遍都小于GHSS迭代方法的和EHSS迭代方法的。且随着qh和N的不断增大，EGHSS迭代方法的迭代矩阵谱半径、迭代次数和迭代时间都远远小于GHSS迭代方法的和EHSS迭代方法的。因此，不难从表1～表3中发现，在合适参数下EGHSS迭代方法会相比GHSS迭代方法和EHSS迭代方法有更好的效果。

表4给出了例2中GHSS迭代方法、EHSS迭代方法和EGHSS迭代方法的收敛迭代次数和收敛时间，并且给出了3种迭代方法近似最优变量的取值。其中N的取值表示的矩阵规模为1728×1728阶和4096×4096阶。从表4可以发现，当矩阵规模相同时，EGHSS迭代方法的迭代速度比GHSS迭代方法的和EHSS迭代方法的更快。甚至当迭代次数相近时，EGHSS迭代方法的迭代时间明显小于GHSS迭代方法的。因此不难发现，在选择合适参数下EGHSS迭代方法会相比GHSS迭代方法和EHSS迭代方法有更好的效果。

Table 1 IT and CPU of GHSS iterative method for example 1

Table 2 IT and CPU of EHSS iterative method for example 1

Table 3 IT and CPU of EGHSS iterative method for example 1

Table 4 IT and CPU of GHSS,EGHSS and EHSS iterative method for example 2

图1展示了例1取qh=10时，系数矩阵规模分别为256×256阶和512×512阶的情况下，EGHSS迭代方法、GHSS迭代方法和EHSS迭代方法的残量随着迭代步数的变化趋势。显然EGHSS迭代方法的残量下降曲线位于GHSS迭代方法和EHSS迭代方法的残量下降曲线之下。这说明对上述qh=10的数值例子，EGHSS迭代方法相较于GHSS迭代方法和EHSS迭代方法迭代次数更少，迭代时间更短。并且系数矩阵规模越大，EGHSS迭代方法相较于其他2种方法的优势越明显。因此，在处理某些实际问题中，EGHSS迭代方法比GHSS迭代方法和EHSS迭代方法具有迭代次数更少、迭代时间更短等优点。

4.2 参数α和ω的灵敏度分析

通常来说，迭代矩阵的谱半径越小，迭代方法的收敛速度越快。因此，本节主要讨论各个变量对EGHSS方法的迭代矩阵谱半径的影响。首先分析单一变量α对EGHSS迭代方法、GHSS迭代方法和EHSS迭代方法的迭代矩阵谱半径的影响；再分析2个变量α和ω共同对EGHSS迭代方法的谱半径的影响。本节参数灵敏度分析仅考虑例1中的参数α和ω对EHSS迭代方法谱半径的影响。

图2a表示例1选取矩阵规模为256×256阶，且qh=10和ω=0.6时，EGHSS迭代方法、GHSS迭代方法和EHSS迭代方法中参数α与迭代矩阵谱半径的关系。图2b表示例1选取矩阵规模为512×512阶，且qh=10,ω=0.5时，EGHSS迭代方法、GHSS迭代方法和EHSS迭代方法中参数α与迭代矩阵谱半径的关系。观察图2可以发现，3种方法迭代矩阵的谱半径均随着α增大先减小后增大。除此之外,EGHSS迭代方法的最小谱半径明显小于GHSS迭代方法和EHSS迭代方法的最小谱半径(都小于1)，这说明了EGHSS迭代方法相比GHSS迭代方法和EHSS迭代方法有更好的迭代效果。

图3a表示例1选取系数矩阵规模为256×256阶且qh=10时，EGHSS迭代方法中的参数α和ω与迭代矩阵谱半径的关系。图3b表示例1选取系数矩阵规模为512×512阶且qh=10时，EGHSS迭代方法中的参数α和ω与迭代矩阵谱半径的关系。观察图3可以发现，当α和ω处于一个适当的取值范围时，EGHSS迭代方法的谱半径小于1，符合定理2的收敛性证明。当系数矩阵规模为256×256阶且qh=10时，α=1.6和ω=0.6时,EGHSS迭代方法的谱半径近似取得最小值(小于1)。当系数矩阵规模为512×512阶且qh=10时，α=1.1和ω=0.5时,EGHSS迭代方法的谱半径近似取得最小值(小于1)。因此,图3a和图3b分别选择ω=0.6和ω=0.5是合理的。

5 结束语

本文提出了一种求解大型稀疏非Hermitian正定线性方程组的EGHSS迭代方法，给出了该方法收敛的充要条件，并进行了数值实验。在与GHSS迭代方法和EHSS迭代方法的谱半径、迭代步数及迭代时间的比较过程中不难发现，EGHSS迭代方法在处理某些问题时相较于其他2种迭代方法具有一定的优越性。除此之外，本文还分析了单一参数α对EGHSS方法迭代矩阵谱半径的影响，以及2个参数α和ω共同对EGHSS方法迭代矩阵谱半径的影响。当然，理论上的最优参数α和ω的确定还有待进一步的研究。总的来说，EGHSS迭代方法是一种较GHSS迭代方法和EHSS迭代方法更有竞争力的方法，特别是选取了合适的参数α和ω后，EGHSS方法的收敛速度大大提高。