求解非埃尔米特正定方程组的广义LHSS迭代法

初鲁， 鲍亮， 董贝贝

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

对于大规模正定线性方程组的求解问题，迭代法是目前的主要方法. 2001年，BAI等[1]提出了一类HSS迭代法，该方法因具有形式简单、易程序化等优点，成为当前的研究热点，得到了许多新的推广. 基于HSS迭代法，BAI等[2]给出了PSS算法, CAO等[3]对PSS算法做了改进，给出了广义PSS(GPSS)算法，LI等[4]给出了修正的GPSS算法；HUANG等[5]将系数矩阵分解成正定和半正定两部分，在这两部分之间进行迭代；2007年，LI等[6]提出了一类LHSS(Lopsided HSS)迭代算法，在系数矩阵的埃尔米特和非埃尔米特之间做非对称迭代，此算法格式简单，在参数较为松弛的约束条件下即可达到收敛；后来，POUR等[7]在LHSS法基础上提出了一类新的HSS方法，围绕系数矩阵的埃尔米特部分做二步非对称迭代.

本文对LHSS迭代法做进一步研究，将方程组的系数矩阵分解成2个正定部分，然后在这2个正定部分间做非对称迭代，给出了一类广义LHSS迭代法.数值结果表明，只要恰当地将系数矩阵分解为2个正定部分，广义LHSS迭代法在处理某些问题时能取得比HSS迭代法更好的效果.

1 背景知识

在多学科领域常出现Ax=b形式的大规模线性方程组，其中A∈Cn×n为非埃尔米特正定矩阵.为了求解该方程，BAI等[1]在矩阵的HS分解(埃尔米特和反埃尔米特分解)基础上提出了HSS迭代算法,该算法基于以下分解：

A=H+S，

(1)

然后给出一个初始向量x(0)∈Cn，通过以下2步计算x(k+1)，直到迭代序列{x(k)}收敛：

(2)

其中α为大于0的常数.对于HSS迭代法，BAI等[1]证明了其迭代矩阵的谱半径小于1，即HSS算法收敛.

2007年，LI等[6]提出了一类LHSS迭代法，该迭代法在系数矩阵A的埃尔米特部分H和反埃尔米特部分S之间做非对称二步迭代. LHSS迭代法的迭代过程如下：

给定初始向量x(0)∈Cn，通过以下2步计算x(k+1)，直到序列{x(k)}满足停止条件：

(3)

其中α为大于0的常数.

2 广义LHSS迭代法

本文对LHSS迭代法做进一步研究，给出了广义LHSS(GLHSS)迭代算法以及相应迭代法的格式.

首先，任意非埃尔米特正定矩阵A有以下分裂形式[3]：

A=G+K+S，

(4)

其中，S表示反埃尔米特矩阵，G和K表示埃尔米特正定矩阵.

其次，任意非埃尔米特正定矩阵A又可分裂为以下形式：

A=P1+P2，

(5)

其中，P1和P2为正定矩阵.

P1和P2的2种可供选择的格式如下：

最后，在系数矩阵A的分裂矩阵P1和P2之间做非对称的二步迭代，得到GLHSS迭代法.

算法1GLHSS迭代算法.

给定初始向量x(0)∈Cn，通过以下2步计算x(k+1)，直到迭代序列{x(k)}满足停止条件：

(6)

其中α为大于0的常数.

3 收敛性证明

为得到GLHSS迭代法的收敛性质，需要以下引理：

引理1[8]设A∈Cn×n，A=Mi-Ni(i=1,2)为矩阵A的2种分裂格式，x(0)∈Cn为一个给定的初始向量，由矩阵A的2种分裂格式得到的二步迭代格式为：

(7)

由上述二步迭代生成的迭代序列{x(k)}为

(8)

(9)

则对于任意初向量x(0)∈Cn，矩阵序列{x(k)}收敛于线性方程Ax=b的唯一解x*∈Cn.

引理2对∀α>0,若矩阵P为复数域上的正定矩阵，则‖P(αI+P)-1‖2<1，其中‖‖2表示矩阵的二范数.

证明因P为正定矩阵，P*也为正定矩阵，α为一个大于0的常数，则对任意非零向量y有

〈(α2I+P+P*)y,y〉>0，

(10)

其中〈,〉表示向量内积，式(10)两端同时加上〈P*Py,y〉有

〈(α2I+P+P*+P*P)y,y〉>〈P*Py,y〉，

(11)

整理后即得

〈(αI+P)y,(αI+P)y〉>〈Py,Py〉.

(12)

令x=(αI+P)y，由于矩阵αI+P非奇异，所以对于任意非零向量x，都可找到一个非零向量y与之对应，经变换，式(12)可化为

〈x,x〉>〈P(αI+P)-1x,P(αI+P)-1x〉，

(13)

于是可得

(14)

由向量x的任意性可得

‖P(αI+P)-1‖2<1.

(15)

证毕！

‖(αI-P)P-1‖2<1.

〈αy,y〉<〈Hy,y〉.

代入H=(P+P*)，得到

〈αy,y〉<〈(P+P*)y,y〉.

(16)

式(16)两边做进一步整理，得

〈(α2I-αP-αP*+P*P)y,y〉<〈P*Py,y〉.

(17)

于是有

〈(αI-P)y,(αI-P)y〉<〈Py,Py〉.

(18)

令x=Py，由于y为任意非零向量且矩阵P为非奇异矩阵，因此，x可选取为任意非零向量，经变换，式(18)可化为

〈(αI-P)P-1x,(αI-P)P-1x〉<〈x,x〉，

(19)

变形后得到

(20)

由于x的选取是任意的，因此，对于任意大于0的常数α有

‖(αI-P)P-1‖2<1 .

(21)

证毕！

证明由引理1，GLHSS迭代法对应的迭代矩阵为

(22)

式(22)相似于：

(23)

对M′(α)谱半径ρ(M′(α))，有：

(24)

综上,即可得到M′(α)的谱半径ρ(M′(α))<1；再由矩阵的相似性质，即可得到ρ(M(α))<1.

证毕！

4 数值实验

数值实验中用Matlab编程求解以下大型稀疏矩阵的线性方程组：

Ax=b，A=I⊗B+BT⊗I，

其中，

M=tridiag(-1,2,-1)N=tridiag(0.5,0,-0.5)，

tridiag表示三对角矩阵，⊗表示克罗内克积.当n=8和16时，系数矩阵的维度实际上为64×64阶和256×256阶，记达到截断误差所需的迭代时间为CPU，迭代次数为IT，实验中截断误差选为10-6，使用英特尔四核处理器(2.70 GHZ，8 GB RAM).表1～表2，图1～图4为2种算法的比较.

表1 n=8时HSS和GLHSS法的迭代次数和迭代时间比较

表2 n=16时HSS和GLHSS法的迭代次数和迭代时间比较

图1 n=8时HSS与GLHSS法谱半径比较Fig.1 Comparison of spectral radius between HSS and GLHSS when n=8

图2 n=16时HSS与GLHSS法谱半径比较Fig.2 Comparison of spectral radius between HSS and GLHSS when n=16

图3 n=8时HSS与GLHSS法残量下降速度比较Fig.3 Comparison of residual decline between HSS and GLHSS when n=8

图4 n=16时HSS与GLHSS法残量下降速度比较Fig.4 Comparison of residual decline between HSS and GLHSS when n=16

从图1～图4、表1和表2中可以看出，对于本文中的算例，GLHSS法到达截断误差所需的迭代步数及迭代时间均优于HSS，GLHSS迭代矩阵的最小谱半径优于HSS，且残量下降速度亦优于HSS.

5 总 结

提出了一种广义的LHSS(GLHSS)迭代法用于求解大规模正定线性方程组，数值结果表明，该方法在求解某些问题时较经典的HSS迭代法效果更好.将来的研究可以围绕如何选取更佳的正定矩阵来开展.