雅可比迭代法和高斯超松弛迭代

雅可比迭代法和高斯-赛德尔迭代法以及超松弛迭代对于给定的方程用下式逐步代入求近似解的方法称为迭代法。如xk(当)的极限存在，此极限即方程组的真正解，此迭代法收敛，否则称迭代法收敛。1、雅可比（Jacobi）迭代法设有方程组Ax=b

(56)其展开形式为(57)系数矩阵A为非奇异阵，且(i=1-n)A可分解为A＝D＋A0（58）其中：改写线性方程组（57)式，将第i个方程（i=1~n）表示为xi的表达式：（59）改写后的（57）式的矩阵表达式为：（60）其中很容易看出（56）式和（60）式间系数矩阵A与及自由项列阵b与之间存在如下关系（61）用迭代法解（60）式，其迭代公式为（62）初始向量为x0。此法称雅可比迭代法，称为雅可比迭代法的迭代矩阵。k为迭代次数。雅可比迭代法的分量形式为（63）雅可比迭代法分量形式（63）式也可改写为（64）（64）式更方便于编程求解。雅可比迭代法公式简单，迭代思路明确。每迭代一次只需计算n个方程的向量乘法，程序编制时需设二个数组分别存放xk和xk+1便可实现此迭代求解。2、高斯－赛德尔（Gauss-seidel）迭代法由雅可比迭代法可知，在计算xk+1的过程中，采用的都是上一迭代步的结果xk。考察其计算过程，显然在计算新分量xik+1时，已经计算得到了新的分量，。有理由认为新计算出来的分量可能比上次迭代得到的分量有所改善。希望充分利用新计算出来的分量以提高迭代解法的效率，这就是高斯－赛德尔迭代法（简称G－S迭代法）对（64）式进行改变可以得到G－S迭代法的分量形式（65）G－S迭代法的分量形式亦可表示为（66）

(65)式也可写成矩阵形式。方程组的系数A在（58）式基础上还可进一步分解，若将A0继续分解为一个下三角阵A0L和一个上三角阵A0U，则系数矩阵的分解则可表达为（67）其中：（68）（65）式可写成矩阵形式按向量的迭代次数归并可得若（D＋）的逆存在，则有高斯－赛德尔迭代的矩阵形式可表达为（69）高斯－赛德尔迭代法每步迭代的计算量与雅可比迭代相当，但在计算机进行计算时，只需存放x一个数组。用迭代法求解时要注意迭代法的收敛性①。无论是用雅可比迭代法或高斯-赛德尔迭代法都存在解的收敛问题，甚至有的线性方程组，用雅可比迭代法解是收敛的，而用高斯-赛德尔迭代法解却是发散的。反之亦然。当线性方程组(56)中的系数矩阵A(n阶矩阵)为严格对角优势矩阵，即A的每一行对角元素的绝对值都严格大于同行其它元素之和:则可证明上述二种迭代法收敛。定理如果线性方程组的系数矩阵A为严格对角优势矩阵，则对求解线性方程组Ax=b的雅可比迭代法和高斯-赛德尔迭代法均收敛。有限单元法的求解方程Ka=P中，系数矩阵K具有主元占优的特点，但却不能保证是严格对角优势矩阵，因此采用雅可比迭代法时要注意解的收敛问题。高斯-赛德尔迭代法的收敛性要求在下一节超松弛迭代法时一并讨论。2、逐次超松弛迭代法逐次超松弛迭代法(SuccessiveOverRelaxationMerhod,简称SOR法)是高斯-赛德尔迭代法的一种加速收敛的方法。是大型稀疏矩阵线性方程组的有效解法之一。首先用另一种变换形式讨论一下迭代过程。对于(56)式的线性方程组Ax=b

当系数矩阵主元aii=1(I=1,2,..,n)时将式中系数矩阵分解为（70）则可（56）式的等价方程组（71）它的迭代公式为当进行了k次迭代得到xk后，xk一般地与真解x间仍存在差异。如何改进xk得到下一次迭代的结果xk＋1，可以引入xk的剩余向量rk（73）此时，迭代公式（72）可表示为（74）由此可见应用迭代法得到逐次改进的解xk＋1，实质上是用k次迭代后的剩余向量rk来改进解的第k次近似rk。因此可以引进一个加速迭代的模式来改进迭代法。令（75）其中ω称为松弛因子。式（75）是迭代公式（74）的一个改进，可以选择松弛因子ω加速迭代过程的收敛。式（75）的分量形式为（76）若对上述改进的迭代公式，按高斯－赛德尔迭代法尽量利用最新迭代得到的分量的原则，又可得到新的迭代公式（77）当线性方程组的系数矩阵A具有非零主元（aii≠0,i=1,2,3,…n）的特点时，可以得到主元为1的方程组形式（78）此时迭代公式（77）可改写为（79）迭代公式（79）式称为松弛因子迭代法。当=1，（79）式就是高斯-赛德尔迭代法；当ω<1时（79）式称为低松弛法；当ω＞1时（79）式称为超松弛法，即SOR方法。由于加速迭代收敛一般选取ω＞1，因此（79）式一般称为超松弛迭代法（SOR方法）。可以证明如果线性方程组的系数矩阵A为对称正定矩阵。ω的选取0＜ω＜2时，则解方程组的SOR方法一定收敛。由于有限单元法的求解方程其系数矩阵具有对称、正定的特点，因此SOR方法是迭代求解的一种常见方法。用超松弛迭代法求解时，应选择超松弛因子，适当的超松弛因子将会加快收敛速度。超松弛因子一般可取1.2左右。超松弛代法和高斯-赛德尔迭代要优于雅可比迭代。而由于选择了超松弛因子ω,一般超松弛迭代效率要高于高斯－赛德尔迭代。超松弛因子无法事先确定最优值，可在迭代过程中根据收敛速度进行调整。用超松弛迭代法求解思路明确，编程简单，存储迭代解xk只需一个数组。但需注意超松弛因子的选择。在用迭代法求解有限元方程时，仍需按系数矩阵(刚度矩阵)在计算机内的存储方式，与直接解法中的高斯消去法和三角分解法所讨论的相似，考虑系数矩阵特点，修改(79)式中元素的下标及迭代涉及的元素。如采用雅可比迭代或高斯-赛德尔迭代时，也应作此修正。读者可以进行修改并验证。附程序#include<iostream>#include<math.h>#include<iomanip>usingnamespacestd;#definekk50//定义最大方程元数#defineexp1e-5//定义结束门限值intn,ll,hh,i,j;doubleA[kk][kk],x[kk][kk],b[kk],y[kk],a[kk],z[kk],m,nn,d,e=1,w;voidmain(){cout<<"**********************************************************************"<<endl;cout<<"***本程序可以用雅可比迭代法,塞德尔迭代法,逐次超松弛法求解线性方程组***"<<endl;cout<<"***************************制作人：李小龙*****************************"<<endl;cout<<"***********************说明：方程最多个数为**************************"<<endl;cout<<"**********************************************************************"<<endl;//*********************************数据的输入********************************bb:cout<<"输入的方程元数"<<endl;cin>>n;cout<<"请输入方程系数矩阵:"<<"\n";for(i=0;i<n;i++)for(j=0;j<n;j++)cin>>A[i][j];cout<<"请输入右边向量:"<<"\n";for(i=0;i<n;i++)cin>>b[i];//*******************************判断是否对角占优*************************for(i=0;i<n;i++){for(j=0;j<n;j++){nn=0;if(i==j){d=fabs(A[i][i]);}elsenn=nn+fabs(A[i][j]);}if(nn>d){cout<<"该方程不对角占优，迭代不收敛"<<endl;cout<<"是否继续？是（），否（）"<<endl;cin>>hh;if(hh!=1)gotobb;elseexit(1);}}//********************************计算出迭代矩阵*************************for(i=0;i<n;i++){b[i]=b[i]/A[i][i];for(j=0;j<n;j++){if(i==j){x[i][i]=0;}else{x[i][j]=-A[i][j]/A[i][i];}}}//*******************************输出迭代矩阵*****************************cout<<"计算出迭代矩阵为："<<endl;for(i=0;i<n;i++){for(j=0;j<n;j++)cout<<x[i][j]<<"";cout<<b[i]<<"";cout<<endl;}//****************************迭代方法的选择*****************************cout<<"请你选择迭代方法！"<<endl;cout<<endl;cout<<endl;cout<<"选用雅可比迭代法，请输入（）！"<<endl;cout<<endl;cout<<endl;cout<<"选用塞德尔迭代法，请输入（）！"<<endl;cout<<endl;cout<<endl;cout<<"选用逐次超松弛法，请输入（）！"<<endl;cout<<endl;cout<<endl;cin>>ll;//*****************************赋迭代初值***********************************for(i=0;i<n;i++){y[i]=1;z[i]=1;}intf=1;switch(ll){case1:gotocc;break;case2:gotoaa;break;case3:gotodd;}//***************************雅可比迭代法************************************cc:while(e>exp){for(i=0;i<n;i++){z[i]=y[i];nn=0;for(j=0;j<n;j++){nn=nn+x[i][j]*z[j];y[i]=nn+b[i];}e=fabs(z[i]-y[i]);if(i==0){cout<<"第"<<f++<<"次迭代";}cout<<setw(8)<<setprecision(8)<<y[i]<<"";}cout<<endl;}cout<<endl;cout<<endl;cout<<"最后结果为："<<endl;cout<<endl;for(i=0;i<n;i++){cout<<"X"<<"["<<i<<"]"<<"="<<y[i];cout<<endl;}exit(1);//***********************************塞德尔迭代法********************************aa:while(e>exp){for(i=0;i<n;i++){z[i]=y[i];nn=0;for(j=0;j<n;j++){nn=nn+x[i][j]*y[j];y[i]=nn+b[i];}e=fabs(z[i]-y[i]);if(i==0){cout<<"第"<<f++<<"次迭代";}cout<<setw(8)<<setprecision(8)<<y[i]<<"";}cout<<endl;}cout<<endl;cout<<endl;cout<<"最后结果为："<<endl;cout<<endl;for(i=0;i<n;i++){cout<<"X"<<"["<<i<<"]"<<"="<<y[i];cout<<endl;}exit(1);//***********