求解非埃尔米特正定方程组的广义LHSS迭代法
2018-11-26初鲁鲍亮董贝贝
初鲁, 鲍亮, 董贝贝
(华东理工大学 理学院, 上海 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迭代法效果更好.将来的研究可以围绕如何选取更佳的正定矩阵来开展.