1、设计题目： 超松弛迭代法解线性方程组 2摘要本文是在matlab环境下熟悉的运用计算机编程语言并结合超松弛变量超松弛迭代法的理论基础对方程组求解。首先，本文以微分方程边值问题为例，导出了离散化后线性方程组即稀疏线性方程组，转化对稀疏线性方程组求解问题。其次，用超松弛( SOR) 迭代法编写matlab程序，对产生的稀疏线性方程组进行迭代法求解。然后，分别改变松弛因子和分段数n的值，分析其收敛性和收敛速度，做出各个方面的分析和比较得到相关结论。最后，将超松弛迭代算法在计算机上运用matlab语言实现, 得出了一组与精确解较接近的数值解,并画图比较，验证逐次超松弛( SOR) 迭代法的精确性。关键
2、词 ： 稀疏线性方程组 逐次超松弛迭代法 松弛因子 matlab 编程31、问题提出考虑两点边值问题 .1,0,0,2yadx容易知道它的精确解为 .1axeayx为了把微分方程离散，把 区间 等分，令 ， ， 得,0nnh1ih,1,2n到差分方程 ,2111 ahyhyyiiiii 简化为 ,211yyiii 从而离散后得到的线性方程组的系数矩阵为 hhhhA 2222 对 ， ， ，分别用 、 和 的超松弛迭代法14.0a2n15.0.1求解线性方程组，要求有 4 位有效数字，然后比较与精确解的误差，探讨使超松弛迭代法收敛较快的 取值，对结果进行分析。改变 ，讨论同样问题。n二、超松弛迭
3、代法产生的背景对从实际问题中得到维数相当大的线性代数方程组的求解仍然十分困难, 以至使人们不能在允许的时间内用直接方法得到解, 因此, 客观上要求用新的方法来解决大维数方程组的求解问题。现有大多数迭代法不是对各类线性方程组都有收敛性, 在解题时, 要对原方程组矩阵作一根本的变换, 从而可能使条件数变坏, 也可能破坏了变换前后方程组的等价性, 以及丧失使原方程组的对称性等。探求新的有效的解题方法依然是迫切的任务。逐次超松弛(Successive Over Relaxation)迭代法是在高斯-赛德4尔(GS)迭代法基础上为提高收敛速度，采用加权平均而得到的新算法。在求解过程中由于线性方程组的系数
4、矩阵维数较大, 采用计算机编写算法来求解, 从而实现了对解析模型的计算机数值逼近的计算方法#本论文以逐次超松弛迭代法为主要的求解方法。三、超松弛迭代法的理论基础（1）逐次超松弛迭代法逐次超松弛(Successive Over Relaxation)迭代法，简称 SOR 迭代法，它是在 GS 法基础上为提高收敛速度，采用加权平均而得到的新算法，设解方程(7.1.3)的 GS 法记为(1)再由 与 加权平均得这里 0 称为松弛参数，将(1)代入则得(2)该法称为 SOR 迭代法，WTBX0 称为松弛因子，当 =1 时(2)式即为高斯-赛德尔迭代法，简记 GS 法，将(2）写成矩阵形式，则得即 于是
5、得 SOR 迭代的矩阵表示(3)其中(4)分解后，有 .（2）逐次超松弛迭代法的收敛性根据迭代法收敛性定理，SOR 法收敛的充分必要条件为 ，收敛的充分条件为 ，但要计算 比较复杂，通常都不用此结论，而直接根5据方程组的系数矩阵 A 判断 SOR 迭代收敛性，下面先给出收敛必要条件.定理 1 设 ，则解方程 的 SOR 迭代法收敛的必要条件是 02.该定理为 SOR 迭代法收敛的必要条件。定理 2 若 对称正定，且 02，则解 Ax=b 的 SOR 迭代法对迭代收敛.对于 SOR 迭代法，松弛因子的选择对收敛速度影响较大，关于最优松弛因子研究较为复杂，且已有不少理论结果.下面只给出一种简单且便
6、于使用的结论。定理 3 设 为对称正定的三对角矩阵， 是解方程的 J 法迭代矩阵，若 ，记 ，则 SOR 法的最优松弛因子 为(5)且(6)根据定理， ,如图 1 所示.由(6)可知，当 =1，时，收敛速度为.说明 GS 法比 J 法快一倍.图 16定理 4 设 ,如果：（1）A 为严格对角占优矩阵；（2）0=1.则解 的 SOR 迭代法收敛。4、实验内容1.自定义函数 sor(A, b, nm, e, w)，以实现 SOR 方法求解线性方程组 AX=B，其中A系数矩阵；b常数列向量；w松弛因子；nm迭代的最大次数e 达到的精度上限(1)()kkX由离散后的差分方程： ,2211ahyhyhi
7、ii ,1,ni得到的线性方程组的系数矩阵为 hhhhA 2222 7常数列向量 b= yahyah20220)(其中 ， ， ， ，则有14.0nnh1。A 为（ aij） 200*200 型矩阵，b 为（bij）0.,51,5.2)( 2ahh200*1 型矩阵。在本次试验中，由于所提供数据较小，当最大迭代次数 nm 较小时，在 nm 迭代次数范围内，不能判断该超松弛迭代法是否收敛，此次取 nm=30000。迭代精度 e 也应取较小值才能使误差更小，此次取 e=0.00001。由定理 1 可知，本次试验中， 的取值范围为：0 2才能保证迭代法收敛。取 ，为 的矩阵。用 SOR 迭代公式得T
8、x)1,()020 .052/05.201. ;0./1;./. 5205201 ;./.1 )()1(9)(20)1(20 )(9899 )(5)(4)(3)(4)1(4 4333 )(3)()()(2)(2 2111 kkkk kkkkkk kkkkk xxx xxxx取不同值时，对应的迭代次数、与精确解的误差如下表 1。表 1 取不同值时对应的迭代次数与误差松弛因子 区间等分数 =200n迭代次数 误差0.5 18111 1.158181 10368 0.36931.5 4891 0.10531.6 3945 0.07311.7 3017 0.04601.8 2087 0.02611.9
9、 1113 0.0236满足误差 的迭代次数42*10xk图 1 计算值与精确值图形比较 从本组的实验中，可以看出 w 值的取定十分重要，它对求解的迭代次数影响十分明显。一个不好的 w 值甚至会导致迭代超过 10000 次仍未能求得需要精度的值。由表 1 可得，当 =1.9 时 SOR 迭代法收敛速度最快，误差最小。取 =1.9，nm=30000, 等各个因子相同时，当分Tx)6.0,.60,.()0段点 n 取不同值时，对应的迭代次数、与精确解的误差如下表 2。9区间等分数 n松弛因子 =1.9迭代次数 误差120 416 0.2477150 416 0.2477200 416 0.2477
10、250 416 0.2477满足误差 的迭代次数42*10xk从本组的实验中，当其他各个因子取适当值时，改变分段数时，对结果没有影响。图 3 精确图形由图 3 可得，当各个参数取值适当时，用 SOR 迭代法所得线性方程组的解与精确解误差极小，从而验证了 SOR 迭代法的准确性。10五、结论1. 通过本次的课程设计，可知逐次超松弛迭代法与 Jacobi 迭代法, Seidel 迭代法相比, 收敛速度较快。由逐次超松弛迭代法求出的方程组的数值解与该方程组的精确解十分接近, 离散化后线性方程组的逐次超松弛迭代法的精确性较高。逐次超松弛迭代法可以广泛地应用于实际。该算法不仅可以用来求解高阶稀疏线性方程
11、组, 还可以用来求解热传导问题这样可以大大减少计算量和计算机的内存储量, 从而提高计算效率。本次的课程设计，我们运用了matlab 语言来实现相关的计算，这样不仅对逐次超松弛迭代法有了更深层的了解掌握，还提高了对 matlab 的操作技术，深刻体会到了 MATLAB 功能的强大之处。通过本次试验，我掌握了用 Jacobi、Gauss-Seidel、SOR 迭代法求解线性方程组的方法；六、参考文献1 李庆扬，王能超，易大义.数值分析M, 清华大学出版社，2008.2 刘卫国. MATLAB 程序设计与应用M，高等教育出版社，2008.3 王诗然. 稀疏线性方程组求解的逐次超松弛迭代法J，沈阳师范
12、大学学报,4,407-409,2006.4 李建宇，黎燕. 牛顿一 SOR 迭代方法中最佳松弛因子的算法J，四川大学学报，4,381-382,1995.5 蔡大用.数值分析与实验学习指导M，清华大学出版社，2001.附录1.超松弛迭代法function n,x=cscdd(A,b,X,nm,e,w)n=1;D=diag(diag(A); %令 A=D-L-U,计算矩阵 DL=tril(-A)+D; %令 A=D-L-U,计算矩阵 LU=triu(-A)+D; %令 A=D-L-U,计算矩阵 UM=inv(D-w*L)*(1-w)*D+w*U); %计算迭代矩阵g=w*inv(D-w*L)*b;
13、 %计算迭代格式中的常数项%下面是迭代过程while n=nm11x=M*X+g; %用迭代格式进行迭代r=norm(x-X,inf);if rereturn;endX=x; n=n+1;enddisp(在最大迭代次数内不收敛!)2.输入初始值并调用 SOR 迭代法n0=200;m=1;a=0.4;h=1/n0;A=zeros(n0,n0);for i=1:n0A(i,i)=-(2*m+h);endfor i=2:n0-1A(i,i-1)=m; A(i,i+1)=m+h;endA(1,2)=m+h;A(n0,n0-1)=m;for i=1:n0-1b(i,1)=a*h2;endb(n0,1)=
14、a*h2-(m+h);for i=1:n0xi=i/n0;12y0(i,1)=(1-0.4)/(1-exp(-1)*(1-exp(-xi)+0.4*xi;x0(i,1)=1; endn1,x1=cscdd(A,b,x0,30000,0.00001,0.5);n1u1=norm(x1-y0,2)n2,x2=cscdd(A,b,x0,30000,0.00001,1.0);n2u2=norm(x2-y0,2)n3,x3=cscdd(A,b,x0,30000,0.00001,1.5);n3u3=norm(x3-y0,2)t=1/200:1/200:1;plot(t,y0)hold onplot(t,x1,g)hold onplot(t,x2,r)hold onplot(t,x3,k)legend(精确解，w=0.5,w=1,w=1.5);title(计算值与精确值图形比较)hold off
