黄金比例分割法确定对称逐次超松弛迭代法的最佳松弛因子//.paper.edu -1- 黄金比例分割法确定对称逐次超松弛迭代法的最佳松弛因子 张德宣，杜成斌，孙立国 河海大学工程力学系，南京 (210098) 摘 要：本文提出了将黄金分割法确定松弛因子与对称逐次超松弛法的改进迭代格式相结合的迭代算法算法应用黄金比例分割法确定最佳松弛因子，成功的将其与运行速度和效率很高的对称逐次超松弛法的改进迭代格式相结合，并给出了迭代收敛性证明，编写了相应的程序，对一实际结构进行的算例计算表明，与大型商业软件的计算结果相比较，本文所提算法具有精度高，收敛快的优点 关键词：黄金分割；对称逐次超松弛法；松弛因子 1．引言 线性代数方程组 bAx = (其中A为非奇异系数矩阵,b为列向量)的迭代解法通常有Jacobi迭代法,Gauss Seidel迭代法和对称逐次超松弛迭代（SSOR）法三种，此外还有USAOR[1]迭代法，USSOR[2]迭代法等 但是当结构的自由度达到较大的数目（如10000个自由度以上），Jacobi法和Gauss Seidel法存在着运算速率很慢这一明显缺点，而且这些迭代方法的结果的精度不能保证。
因此本文针对运算速率很高的对称逐次超松弛法的改进迭代格式提出了用黄金分割法确定松弛因子的思路，编写了黄金分割计算松弛因子与相应的迭代法相结合的程序，实际算例表明，本文方法可提高求解大型线性方程组的求解效率和精度 2．SSOR-PCG法改进的迭代格式 SSOR法在每步迭代中，需计算系数矩阵与方向向量的乘积，下面引用了该法的改进迭代格式，避免了该乘积的计算，比原迭代格式节省计算量8％—50％ SSOR法的改进迭代格式[3]为： 置初值 0,,,,, 0000010000 ==??==??= ???? kzWdVyzgWyRKg Tδδ ， (1) R： ( )kk Vyy ,=α 式中 ( ， )为内积表示，下同 如果 εα ≤ ，则停止，否则 //.paper.edu -2- ??????????????????????????????+==+??==??++=+=??=+??+++++??++RkkzWdzVyzVyyVyyVdzWdyydVdzdVyykTkkkkkkkkkkkkkkkkkkkkkkkkkk转到,1,,),,/(),()),((,),2,/(),(111111111ββττδδτ (2) 从迭代格式可见，由于V是对角阵，计算Vy和Vd为2n次乘法运算。
改进的迭代格式由于避免了计算Ad，比原迭代格式可省 ( )nra 2?? 次乘法运算， ar 为A的各行非零元素的平均个数因此， 3=ar 是改进迭代格式节省计算量约8％； 62 >>ar 时，可节省计算量约50％，各迭代步计算量接近CG（共轭梯度）法的计算量 3．松弛因子的选取 本文采用上述对称逐次超松弛法的改进迭代格式求解大型线性方程组，将黄金分割法确定松弛因子与之相结合，松弛因子的取值对迭代公式的收敛速度影响很大，它的好坏直接影响到收敛的快慢为了保证迭代过程的收敛，取 21 << ω ，但是在1和2之间仍然有很多的取值，究竟如何取值，目前并没有统一的规定，一般是将ω的取值区间(1,2)进行M等分，松弛因子ω分别取 M/11+ , M/21+ ,… …, MM /)1(1 ??+ ,然后分别代入迭代公式求解的同时比较出最少的迭代次数，这样就得到了最佳的ω值但是这种方法如果要得到比较精确的ω值，M的取值就要变大，如果要在程序中实现，就导致计算量增加很多，运算时间就会延长许多，是很不经济的 黄金分割比例法是一种经典的算法，本文即采用此种方法来确定最佳松弛因子，其具体步骤和主要思想如下： (1)利用优选法思想,在(1,2)之间选取四个点，即： dot1=1，dot2=dot4-0.618(dot4-dot1)，dot3=dot1+0.618(dot4-dot1)，dot4=2 。
(2)分别取dot2与 dot3作为松弛因子代入迭代程序，比较出最少的迭代次数，如果dot2对应的迭代次数少，则选取(dot1, dot2)作为收敛区间，如果是dot3对应的迭代次数少，则选取(dot3, dot4)作为收敛区间 (3)在所选取的收敛区间里循环进行上述的两个步骤，直到选取出迭代次数最少时所对应的最佳松弛因子 通过本文中所采用的算例验证，本文方法确定松弛因子的加速效果是明显的，计算的时间相对较短在SSOR法改进迭代格式中采用黄金比例分割的综合算法就可以根据不同的具体工程问题确定最佳松弛因子，该算法简单明了，而且效果显著，大大节省了选取最佳松弛因子的时间，提高了大型有限元方程组计算的收敛速度 4．收敛性证明 定义1[4] 设A nnC ??∈ ， nλλλ ,,, 21 L 为A的特征值，我们称 //.paper.edu -3- iiA λρ max)( = (3) 为A的谱半径 定义2[5] 由 1)(0 << Bρ ，且误差向量 ????= xx kk )()(ε 满足 )0(εε kk B= ，推出欲使初始误差缩小 s??10 所需最小迭代次数为)(ln10lnBsk ρ??≥ ，称 )(ln)( BBR ρ??= 为迭代法收敛速度。
定义3[6] nnijaA ×= )( ，如果A的元素满足 ∑≠=≥nijjijij aa1),2,1,( nji L= ，称A为对角占优矩阵，记为 0DA∈ ；而当 ∑≠=>nijjijij aa1),2,1,( nji L= 时，称A为严格对角占优矩阵，记为DA∈ 引理1[7] 设有等价的方程组 fBxxbAx +=??= (4) 构造的一阶定常迭代法 fBxx kk +=+ )()1( (5) 对于任意选取的初始向量 )0(x ，迭代法收敛的充要条件是迭代矩阵B的谱半径 1)( <Bρ 定理1 设 bAx = ，本文采用改进的SSOR迭代法收敛的充分条件为： (1) nnRA ×∈ 为严格对角占优矩阵，且 ULDA ????= （D为单位对角阵，L为单位下三角阵，U为单位上三角阵）； (2) 21 << ω ； (3) 取 [ ] [ ]IkD = ，当 )21,0(∈nk 时； 则本文采用改进的SSOR迭代法收敛。
证明：本文所采用的SSOR-PCG法改进的迭代格式所采用的迭代矩阵为： ))1(()()( 11 UILIACDIL L ωωωωωω +????=????= ???? (6) 采用反证法，假设λ为 ωL 的任一特征值，且 1≥λ ，则 =?? ωλ LI ωωωλ LLILI ?????? ?? )()( 1 [ ]UILILI ωωωλω ??????????= ?? )1()()( 1 [ ]ULILI ωλωωλω ????+????= ?? )1()( 1 0)1(1 =????+????= ?? ULILI ωλωωλω (7) 由于 ULDA ????= 为严格对角占优矩阵，所以 0≠????= ULDA ，由于 LI ω?? 也为对角阵，从而 01 ≠?? ??LI ω 分析式（7）后半部分，并由 1≥λ 可得 //.paper.edu -4- ∑∑∑ ∑+=??=??= +=????+??≥????+??nijijiijijiiiiijnijijiijiiii aaIaaI11111 1)1()1( ωωλωλωωλωλ )()1()()1(1,11 111 1∑ ∑ ∑∑ ∑=??= +=??= +=++??+??>+??+??>njiijnijijiijiijiiiiijnijijiijiiii aaaIaaI ωωωλωλωωλωλ∑∑∑=== ???????????????? ??+??=???????????????? ??+??=??+??>niiininiiniiininiiiiiii akkakaI111)1()1(1)1(ωωλωλωλωλωλ由 21 << ω ，进而 01 >??iω ，又由 )21,0(∈nk ，推得 01 >?? ink ω ，于是可以得出 ∑=≥???????????????? ??+?? niiiniinakk10)1()1( ωωλ，此式表明 ( ) ULI ωλωωλ ????+??1 也是严格对角占优矩阵，其行列式 0)1( ≠????+?? ULI ωλωωλ ，这与式（7）矛盾，故 ωL 的特征值均有 1<λ ，结合定义1可知谱半径 1)( <ωρ L ，则该迭代方法收敛。
本文论证出了改进的SSOR迭代法收敛的一个充分条件，而由引理1并结合定义2可知，迭代次数与收敛速度成反比，即 )(Bρ 越小，迭代收敛速度R(B)越快而一般来说迭代矩阵B的谱半径 )(Bρ 对ω极为敏感，但是究竟是ω取值大能减少迭代次数，还是ω取值小可以减少迭代次数，这并无规律，完全取决于计算当时的情况，如结构划分单元的多少及计算精度等因素都有一定的关系 5．算例 。