复线性方程组的预处理MCG算法张迎春;吕全义;肖曼玉【摘要】复线性方程组在科学与工程计算的诸多领域中有着重要的应用价值,如何高效的求解复线性方程组,一直是人们所关心的问题.目前对于复线性方程组,常用的处理方式有以下两种:一种是直接对方程组迭代求解,另外一种是将其转化为实线性方程组后进行求解.本文主要从两种处理方式讨论了共轭梯度法(CG法),并理论上证明了两种处理方式下的CG法具有相同的收敛性.之后基于变形共轭梯度法(MCG 法)收敛速度的本质与CG法类似,只需将MCG法推广到复线性方程组进行研究,并且为了提高MCG法的收敛速度,提出了一种预处理MCG法.最后,通过数值算例验证了算法与理论分析的一致性,以及预处理算法的有效性.%Complex linear equations have a wide application in science and engineering, and an important issue is how to solve it with high efficiency. Until now, complex linear equations are usually solved by either iteration methods or the solution of the real equations transformed from the original equations. Conjugate gradient method (CG method) is discussed from two different viewpoints, and it is proved theoretically that these two kinds of CG methods have the same convergence. Because the convergence speed of the modified conjugate gradient method (MCG method) and conjugate gradient method are essentially similar, MCG method is extended to solve complex linear equations. Besides, a preconditioned MCG method is proposed in order to improve the convergence speed. Finally, the consistency of algorithms and theoretical analysis and effectiveness of the proposed precondition algorithm are validated by numerical examples.【期刊名称】《工程数学学报》【年(卷),期】2018(035)003【总页数】11页(P308-318)【关键词】复线性方程组;变形共轭梯度法(MCG法);预处理方法;收敛性【作者】张迎春;吕全义;肖曼玉【作者单位】西北工业大学应用数学系,西安 710072;西北工业大学应用数学系,西安 710072;西北工业大学应用数学系,西安 710072【正文语种】中文【中图分类】O241.61 引言在科学计算和工程应用等诸多领域中，会产生大量的复线性方程组求解问题，包括时薛定谔方程[1]、Helmholtz方程[2]、声散射问题[3]等．因此，对于复线性方程组的求解就吸引了越来越多研究者的关注．在文献[4,5]中，讨论了系数矩阵为一般复矩阵的情况下，将CG法推广，采用复双共轭梯度法直接求解复线性方程组，并从理论上证明了该算法在复数域上是收敛的，该算法克服了CG法只能求解对称正定的情况，但是在求解系数矩阵为对称正定的情况下，复双共轭梯度法的收敛性却不如通常的CG法好．另外，对于复对称线性方程组，目前研究比较热门的一种算法是HSS迭代算法[6]，如在文献[7–9]中，借助复系数方程组的实系数形式，采用HSS迭代算法求解复对称线性方程组，但是在每步迭代的过程中涉及到求解对称复线性方程组．为了避免这种操作，之后人们又对HSS算法进行了改进，得到了求解复对称方程组的MHSS迭代算法，虽然收敛速度加快，但是这类方法如何选取最优收敛参数却比较困难．本文探讨复线性方程组的求解问题，我们首先讨论转化为两种形式求解的关系，形式一是直接从复线性方程组入手进行求解；形式二是将复线性方程组的系数矩阵实部与虚部分离进行重组，转换为实线性方程组进行求解．另外，在求解的过程中，分别讨论了当系数矩阵为Hermite正定矩阵时，采用CG法，当系数矩阵为非Hermite正定矩阵时，采用MCG法，并在此基础上进行预处理．文中欧氏空间Rn上的内积定义：[x,y]=xTy，其中x,y∈Rn，酉空间Cn上的内积定义：[x,y]=xHy，其中x,y∈Cn，σ(A)表示A的特征值集合，表示A的共轭矩阵，∥x∥2表示向量x的2-范数．2 Hermite正定线性方程组的CG算法考虑复线性方程组其中A∈Cn×n,x∈Cn,b∈Cn，且A为n阶Hermite正定矩阵．形式一为(1)式．形式二：对复线性方程组(1)进行实部与虚部分离，可得其中A1,x1,b1分别表示A,x,b的实部，A2,x2,b2分别表示A,x,b的虚部，且2.1 Hermite正定线性方程组的CG算法首先，我们给出形式一中复线性方程组(1)的复CG算法[4](记为CCG法)：第1步：任意给定初始向量x0∈Cn，置k:=0，计算第2步：若rk=0，计算停止；否则，计算第3步：计算第4步：置k:=k+1，转第2步．其次，我们再给出形式二中线性方程组(2)的CG法[10]：第1步：任意给定初始向量∈R2n，置k:=0，计算第2步：若rk=0，计算停止；否则，计算第3步：计算第4步：置k:=k+1，转第2步．2.2 两种算法的收敛性分析由于A为n阶Hermite正定矩阵，则有如下性质：性质1[11]A为n阶Hermite正定矩阵，则A的特征值均为实数．性质2 A为n阶Hermite正定矩阵，则：1)为n阶Hermite正定矩阵；2) 为n阶对称正定矩阵．证明 1) 由于为n阶Hermite矩阵．又由性质1知，A的特征值均为正实数，则A的特征值也为正实数，故为n阶Hermite正定矩阵．2) 因为故可得即为2n阶对称矩阵．令可得则可知矩阵C与矩阵相似．又由于A和均为正定矩阵，故为正定矩阵．综上所述，为2n阶对称正定矩阵．证毕与建立实系数线性方程组中的CG法的收敛性定理[10]类似，对于CCG法，我们也有相应的结论．定理1 设A为n阶Hermite正定矩阵，则由CCG法中得到的向量序列{rk}和{zk}满足定理2 设线性方程组(1)相容，则对于任意的初始向量x0∈Cn，CCG法最多经过n次迭代求得线性方程(1)的精确解．定理3[12] 共轭梯度法产生的近似解xk满足其中x∗为方程组(2)的精确解，Ck是k阶Chebyshev多项式，且满足因为定理3中不等式右端项满足由此可知系数矩阵的条件数κ越小，CG法收敛就越快，即CG法的收敛性与系数矩阵的谱分布有关．如果系数矩阵的谱均匀地分布在一个很长的区间上，则CG法收敛很慢，如果其谱除少数几个零散的特征值外，其余大部分凝聚在少数几个点附近，则CG法收敛非常快[12]．这表明CG法的收敛速度快慢与特征值的集中度相关，类似地，由于CG法与CCG法从算法本质上来说是相同的算法，故CCG法的收敛速度的快慢也同样由系数矩阵的特征值集中度决定．当然，对于形式二，即方程(2)，采用CG算法，最多经过2n次迭代求得线性方程(1)的精确解．而这并不能真正了解它们的收敛速度的关系．下面讨论方程(1)与(2)的系数矩阵的特征值的关系．定理4 设A为2n阶Hermite正定矩阵，A的特征值集合为σ(A)={λi|λi>0,1≤i≤n}⊂R，则矩阵的全部特征值的集合为σ()=σ(A)∪σ(A)．证明设由性质2的证明可知，故的全部特征值集合为σ(=σ(A)∪σ(A)．证毕由定理4知，CG法的特征值分布与CCG法的特征值分布相同，故两种算法具有相同的收敛速度．对于复线性方程组(1)的求解，由定理2可知，CCG法最多经过n次迭代即可求得线性方程的精确解，可见用CG法求解方程(2)也最多经过n次迭代就可求得方程的精确解．此外，分析CCG法与CG法在每步迭代过程中的计算量，可以发现在两个算法中，其计算量关键在于内积的运算和矩阵与向量相乘的运算，因此，我们可以通过分析两个算法的步骤2来说明这一问题．对于CCG法，计算αk所需的乘法运算量为4n2+4n，加法运算量为4n2+2n−3，而对于CG法，计算αk所需的乘法运算量为4n2+4n，加法运算量为4n2+2n−2，类似地，可以计算出其他步骤的运算量，最后通过理论上分析，发现CCG法与CG法的运算量几乎相等，即两个算法每步具有相同的运算复杂度．2.3 数值算例对于复线性方程组(1)，分别采用CG法、CCG法进行对比分析，检验是否与理论分析一致．算例中，初始向量x0=0,rk=b−Axk，其中xk表示第k步的迭代向量，rk表示第k步的迭代误差，终止条件为满足∥rk∥2/∥r0∥2≤ 10−6．误差用r表示，迭代次数用k表示，时间用t表示，其中r= ∥rk∥2/∥r0∥2．算例使用Matlab R2010a(7.10)软件，CPU2.40GHz微机计算．例1 对于复线性方程组(1)，A为n阶随机Hermite正定矩阵，b为n维随机复列向量，其中A和b的元素均按均匀分布选取．两种算法的误差、迭代次数和时间的结果对比见表1．表1: 例1的计算结果n=1000 n=2000 n=3000 CCG法CG法CCG法CG法CCG法CG法9.02e−7 k r 8.26e−7 8.26e−7 8.91e−7 8.91e−7 9.02e−7 224 t 153 153 211 211 224 0.96 0.97 4.85 5.09 11.96 13.10从表1的数据可以看出，CG算法与CCG算法的迭代次数相同，这与理论分析是一致的．而所用时间几乎一样，故可以发现CG算法与CCG算法每步有相同的运算复杂度．故两种方式的计算效果是相同的．3 复线性方程组的MCG算法对于复线性方程组(1)，考虑A为n阶非Hermite正定矩阵的情况．由于MCG算法的收敛速度本质上与CG法类似，由上面的讨论可知，我们只需采用形式一研究MCG算法的效果即可．3.1 复线性方程组的MCG算法将文献[10]中的MCG算法推广到求解复线性方程组(1)，我们得到如下的复MCG 算法(记为CMCG算法)：第1步：任意给定初始向量x0∈Cn，置k:=0，计算第2步：若rk=0，或者r k̸=0而zk=0，计算停止；否则，计算第3步：计算第4步：置k:=k+1，转第2步．3.2 MCG算法的预处理算法采用与实线性方程组加速求解的类似方法，对复线性方程组(1)进行预处理，通过构造预处理条件子M−1，使其近似等于系数矩阵A的逆矩阵，而后在复线性方程组(1)的两边同乘以预处理矩阵，即M−1Ax=M−1b．通过这种预处理技术，使得系数矩阵的奇异值分布更集中，从而加快CMCG算法的收敛速度．在这里，取M−1为其中如下给出复线性方程组(1)的预处理复MCG算法(记为PCMCG算法)：第1步：任意给定初始向量x0∈Cn，置k:=0，计算第2步：若rk=0，或者r k̸=0而zk=0，计算停止；否则，计算第3步：计算第4步：置k:=k+1，转第2步．从PCMCG算法的迭代过程可以看出，它要比CMCG算法多求解两个复线性方程组=rk与MHpk=，下面给出这两个复线性方程组求解的迭代格式．首先，给出求解Merk=rk的迭代格式：类似地，求解MHpk=的迭代格式为：3.3 收敛性分析与文献[10]中对求解实线性方程组的MCG算法的收敛性证明类似，我们同样有下面的结论：定理5 设A∈Cn×n，则由CMCG算法中得到的向量序列{rk}和{zk}满足定理6 设线性方程组(1)相容，则对于任意的初始向量x0∈Cn，CMCG算法最多经过n次迭代求得线性方程组(1)的精确解．对于PCMCG算法中所涉及的迭代次数q，在一定的范围内，系数矩阵特征值的密集程度会随着q的增大而更集中，从而加快CMCG算法的收敛速度．下面通过数值算例来说明PCMCG算法将如何改善CMCG算法的计算结果．3.4 数值算例对于复线性方程组(1)，分别采用CMCG算法、PCMCG算法、MHSS算法[7]进行对比分析．例2、例3中的初始向量、终止准则的选取以及符号的表示与例1保持一致，a表示MHSS算法[7]中的最优收敛参数．算例使用MatlabR2010a(7.10)软件，CPU 2.40GHz微机计算．例2[7] 考虑复线性方程组该问题来源于抛物方程[13]在单位正方形区域[0,1]×[0,1]上以空间步长上R22-Padé近似的中心差分离散，R22-Padé近似满足算子方程对该算子方程分解可得到如下形式的两个方程可以看出复线性方程组(3)即为离散形式的算子方程(6)．对于复线性方程组(3)，其中τ表示时间步长，且τ=h，L可以用张量积的形式表示为L=I⊗Vm+Vm⊗I,m=右端项b的定义如下另外，我们在系数矩阵的左右两侧同乘以h2进行规范化．三种算法的误差、迭代次数和时间的结果对比见表2和表3．表2:例2的计算结果(n=210,a=0.75)算法r k t 2.41 MHSS算法[7]CMCG算法9.67e−7 370 1.84 PCMCG算法q=1 9.61e−7 54 4.03 PCMCG算法q=29.49e−7 370 1.86 PCMCG算法q=4 7.77e−7 98 8.77e−7 51 1.77表3:例2的计算结果(n=212,a=0.54)算法r k t 93.45 MHSS算法[7]CMCG算法9.78e−7 906 100.88 PCMCG算法q=1 9.41e−7 73 160.10 PCMCG算法q=2 9.68e−7 906 66.63 PCMCG算法q=4 9.64e−7 228 8.88e−7 115 53.12例3[7] 考虑复线性方程组其中M和L分别是惯性矩阵和刚度矩阵，CV和CH分别是黏性矩阵和滞后阻尼矩阵，ω是角频率，µ是阻尼系数，取M=I,CV=10I,CH=µL,µ=0.02,ω=π，其中矩阵L与例2中相同，右端项b=(1+i)A1，其中1表示所有元素都为1的向量，并在系数矩阵的左右两侧同乘以h2进行规范化．三种算法的误差、迭代次数和时间的结果对比见表4和表5．表4:例3的计算结果(n=210,a=0.08)算法r k t 1.67 MHSS算法[7]CMCG算法7.49e−7 175 1.65 PCMCG算法q=1 7.42e−7 38 1.82 PCMCG算法q=26.62e−7 175 1.66 PCMCG算法q=4 3.44e−7 91 9.40e−7 51 1.63表5:例3的计算结果(n=212,a=0.04)算法r k t 79.46 MHSS算法[7]CMCG算法9.98e−7 632 94.74 PCMCG算法q=1 8.49e−7 50 113.11 PCMCG算法q=2 7.81e−7 632 78.79 PCMCG算法q=4 4.35e−7 310 6.26e−7 166 76.94从表4和表5中的数据可知，当内迭代次数q=1时，PCMCG算法与CMCG算法的迭代次数几乎接近，这是由于此时的预处理条件子M−1=D−1为对角阵，PCMCG算法相对于CMCG算法影响很小，因此，两种算法在这种情况下迭代次数几乎相同．从表4和表5中的迭代时间数据可知，对于同一阶数的复线性方程组求解，当内迭代次数q在一定的范围内时，PCMCG算法的收敛速度随着内迭代次数的增加而相应的加快，所用时间也相应减少．特别是，当PCMCG算法的内迭代次数q=4时，收敛速度最快，迭代次数最少，并且运行时间均少于CMCG算法和MHSS算法的运行时间．结果表明：CMCG算法经预处理后，无论从运行时间还是迭代次数上，PCMCG 算法的效果都要比CMCG算法的效果好．4 结论通过理论分析，我们表明两种形式下CG法是一致的，这为之后将MCG算法推广到求解复线性方程组的研究奠定了基础；之后，我们采用与实线性方程组加速求解的类似方法，对CMCG算法进行预处理，得到了求解非Hermite正定线性方程组的PCMCG算法，数值算例表明，当内迭代次数为q=4时，PCMCG算法的收敛速度优于CMCG算法和MHSS算法．参考文献：[1] Van D K,Toyama F M.Accurate numerical solutions of the time-dependent Schr¨odinger equation[J].Physical Review E:Statistical Nonlinear&Soft Matter Physics,2007,75(2):036707[2] Benzi M,Bertaccini D.Block preconditioning of real-valued iterative algorithms for complex linear systems[J].IMA Journal of Numerical Analysis,2008,28(3):598-618[3] Rao S M.An iterative method to solve acoustic scattering problems using a boundary integral equation[J].Journal of the Acoustical of America,2011,130(4):1792-1798[4] Jacobs D A H.A generalization of the conjugate-gradient method to solve complex systems[J].IMA Journal of Numerical Analysis,1986,6(4):447-452[5] Joly P,Meurant plex conjugate gradient methods[J].Numerical Algorithms,1993,4(3):379-406[6] 蒋尔雄.矩阵计算[M].北京：科学出版社，2008 Jiang E X.Matrix Computations[M].Beijing:Science Press,2008[7] Bai Z Z,Benzi M,Fang C.Modif i ed HSS iteration methods for a class of complex symmetric linear systems[J].Computing,2010,87(3-4)：93-111 [8] Bai Z Z,Benzi M,Fang C．On preconditioned MHSS iteration methods of complex symmetric linear systems[J].Numerical Algorithms,2011,56(2):297-317[9] Bai Z Z.On preconditioned iteration methods for complex linear systems[J]．Journal of Engineering Mathematics,2015,93(1):41-60[10] 张凯院，徐仲.数值代数(第2版)[M].北京：科学出版社，2010 Zhang K Y,Xu Z.Numerical Algebra(2nd Edition)[M].Beijing:Science Press,2010 [11] 吴建平，王正华，李晓梅.稀疏线性方程组的高效求解与并行计算[M].长沙：湖南科学技术出版社，2004 Wu J P,Wang Z H,Li X M.Efficient Solution and Parallel Computation of Sparse Linear Equations[M].Changsha:Hunan Science and Technology Press,2004[12]徐树方，钱江.矩阵计算六讲[M].北京：高等教育出版社，2011 Xu S F,Qian J.Six Lectures on Matrix Computations[M].Beijing:Higher Education Press,2011[13] Axelsson O,Kucherov A.Real valued iterative methods for solving complex symmetric linear systems[J].Numerical Linear Algebra with Applications,2000,7(4):197-218。
