2023年LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.pdf

《2023年LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.pdf》由会员分享，可在线阅读，更多相关《2023年LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.pdf（11页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、 LU 分 解 法、列 主 元 高 斯 法、J a c o b i 迭 代 法、Ga u s s-S e i d e l法 的 原 理 及 M a t l a b程 序 精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 一、实验目的及题目 1.1 实验目的：（1）学会用高斯列主元消去法，LU 分解法，Jacobi 迭代法和 Gauss-Seidel 迭代法解线性方程组。（2）学会用 Matlab 编写各种方法求解线性方程组的程序。1.2 实验题目：1.用列主元消去法解方程组：1 2 41 2 3 41 2 3 41 2 3 43 42 12 3 32 3 4x x xx

2、 x x xx x x xx x x x 2.用 LU 分解法解方程组,Ax b 其中 48 24 0 1224 24 12 120 6 20 26 6 2 16A，4422b 3.分别用 Jacobi 迭代法和 Gauss-Seidel 迭代法求解方程组：1 2 32 3 41 2 31 2 3 410 2 118 3 112 10 63 11 25x x xx x xx x xx x x x 二、实验原理、程序框图、程序代码等 2.1 实验原理 2.1.1 高斯列主元消去法的原理 Gauss 消去法的基本思想是一次用前面的方程消去后面的未知数，从而将方程组化为等价形式：解方程组用分解法解方

3、程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 11 1 12 2 1 122 2 2 2n nn nnn n nb x b x b x gb x b x gb x g 这个过程就是消元，然后再回代就好了。具体过程如下：对于 1,2,1 k n，若()0,kkka 依次计算()()(

4、1)()()(1)()()/,1,k kik ik kkk k kij ij ik kjk k ki i ik km a aa a m ab b m b i j k n 然后将其回代得到：()()()()()1/()/,1,2,1n nn n nnnk k kk k kj j kkj kx b ax b a x a k n n 以上是 高斯消去。但是高斯消去法在消元的过程中有可能会出现()0kkka 的情况，这时消元就无法进行了，即使主元数()0,kkka 但是很小时，其做除数，也会导致其他元素数量级的严重增长和舍入误差的扩散。因此，为了减少误差，每次消元选取系数矩阵的某列中绝对值最大的元素作

5、为主元素。然后换行使之变到主元位置上，再进行销元计算。即 高斯列主元消去法。2.1.2 直接三角分解法（LU分解）的原理 先将矩阵 A直接分解为 A LU 则求解方程组的问题就等价于求解两个三角形方程组。直接利用矩阵乘法，得到矩阵的三角分解计算公式为：1 11 1 111111,1,2,/,2,1,2,3,()/,1,2,i ii ikkj kj km mjmkik ik im mk kkmu a i nl a u i nu a l u j k k nk nl a l u u i k k n k n 且 由上面的式子得到矩阵 A的 LU分解后，求解 Ux=y 的计算公式为 解方程组用分解法解方

6、程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 1 1111,2,3,/()/,1,2,1ii i ij jjn n nnni i ij j iij iy by b l y i nx y ux y u x u i n n 以上为 LU分解法。2.1.3Jacobi 迭代法和 Gaus

7、s-Seidel 迭代法的原理（1）Jcaobi 迭代 设线性方程组 b Ax(1)的系数矩阵 A可逆且主对角元素nna,.,a,a22 11均不为零,令 nna,.,a,a diag D22 11 并将 A分解成 D D A A(2)从而(1)可写成 b x A D Dx 令 1 1f x B x 其中b D f,A D I B 1111.(3)以1B为迭代矩阵的迭代法(公式)1 11 f x B x k k(4)称为 雅可比(Jacobi)迭代法,其分量形式为 解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除

8、谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8,.,k,n,.,ix a baxni jj)k(j j i iii)k(i2 1 0 2 1111(5)其中 Tnx,.x,x x 0 02010 为初始向量.（2）Gauss-Seidel 迭代 由雅可比迭代公式可知,在迭代的每一步计算过程中是用 k x的全部分量来计算 1 k x的所有分量,显然在计算第 i 个分量 1 kix时,已经计算出

9、的最新分量 1111 kikx,.,x没有被利用。把矩阵 A分解成 U L D A(6)其中 nna,.,a,a diag D22 11,U,L 分别为A的主对角元除外的下三角和上三角部分,于是,方程组(1)便可以写成 b Ux x L D 即 2 2f x B x 其中 b L D f,U L D B1212(7)以2B为迭代矩阵构成的迭代法(公式)2 21 f x B x k k(8)称为 高斯 塞德尔迭代法,用分量表示的形式为,.,k,n,ix a x a baxijni j)k(j ij)k(j ij iii)k(i2 1 0 2 1111 11 1 2.2 程序代码 2.2.1 高斯

10、列主元的代码 function Gauss(A,b)%A 为系数矩阵，b 为右端项矩阵 解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 m,n=size(A);n=length(b);for k=1:n-1 pt,p=max(abs(A(k:n,k);%找出列中绝

11、对值最大的数 p=p+k-1;if pk t=A(k,:);A(k,:)=A(p,:);A(p,:)=t;%交换行使之变到主元位置上 t=b(k);b(k)=b(p);b(p)=t;end m=A(k+1:n,k)/A(k,k);%开始消元 A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n);b(k+1:n)=b(k+1:n)-m*b(k);A(k+1:n,k)=zeros(n-k,1);if flag=0 Ab=A,b;end end x=zeros(n,1);%开始回代 x(n)=b(n)/A(n,n);for k=n-1:-1:1 x(k)=(b(k)-

12、A(k,k+1:n)*x(k+1:n)/A(k,k);end for k=1:n fprintf(x%d=%fn,k,x(k);end 2.2.2 LU 分解法的程序 function LU(A,b)%A 为系数矩阵，b 为右端项矩阵 m,n=size(A);%初始化矩阵 A，b，L 和 U 解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的

13、严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 n=length(b);L=eye(n,n);U=zeros(n,n);U(1,1:n)=A(1,1:n);%开始进行 LU 分解 L(2:n,1)=A(2:n,1)/U(1,1);for k=2:n U(k,k:n)=A(k,k:n)-L(k,1:k-1)*U(1:k-1,k:n);L(k+1:n,k)=(A(k+1:n,k)-L(k+1:n,1:k-1)*U(1:k-1,k)/U(k,k);end L%输出 L 矩阵 U%输出 U 矩阵 y=zeros(n,1);%开始解方程组 Ux=y y(1)=b(1);fo

14、r k=2:n y(k)=b(k)-L(k,1:k-1)*y(1:k-1);end x=zeros(n,1);x(n)=y(n)/U(n,n);for k=n-1:-1:1 x(k)=(y(k)-U(k,k+1:n)*x(k+1:n)/U(k,k);end for k=1:n fprintf(x%d=%fn,k,x(k);end 2.2.3 Jacobi 迭代法的程序 function Jacobi(A,b,eps)%A 为系数矩阵，b 为后端项矩阵，epe 为精度 m,n=size(A);解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主

15、交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 D=diag(diag(A);%求矩阵 D L=tril(A)-D;%求矩阵 L U=triu(A)-D;%求矩阵 U temp=1;x=zeros(m,1);k=0;while abs(max(x)-temp)eps temp=max(abs(x);k=k+1;%记录循环次数 x=-inv(D)*(L+U)*x+in

16、v(D)*b;%雅克比迭代公式 end for k=1:n fprintf(x%d=%fn,k,x(k);end 2.2.4 Gauss-Seidel 迭代程序 function Gauss_Seidel(A,b,eps)%A 为系数矩阵，b 为后端项矩阵，epe 为精度 m,n=size(A);D=diag(diag(A);%求矩阵 D L=D-tril(A);%求矩阵 L U=D-triu(A);%求矩阵 U temp=1;x=zeros(m,1);k=0;while abs(max(x)-temp)eps 解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代

17、码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 temp=max(abs(x);k=k+1;%记录循环次数 x=inv(D-L)*U*x+inv(D-L)*b;%Gauss_Seidel 的迭代公式 end for k=1:n fprintf(x%d=%fn,k,x(k);end 三、实验过程中需要记录的实验数据表格 3.1 第一题（高斯列主元

18、消去）的数据 A=1 1 0 3;2 1-1 1;3-1-1 3;-1 2 3-1;b=4;1;-3;4;Gauss(A,b)x1=-1.333333 x2=2.333333 x3=-0.333333 x4=1.000000 3.2 第二题（LU 分解法）数据 A=48-24 0-12;-24 24 12 12;0 6 20 2;-6 6 2 16;b=4;4;-2;-2;LU(A,b)L=1.0000 0 0 0-0.5000 1.0000 0 0 0 0.5000 1.0000 0-0.1250 0.2500-0.0714 1.0000 U=48.0000-24.0000 0-12.000

19、0 0 12.0000 12.0000 6.0000 0 0 14.0000-1.0000 0 0 0 12.9286 解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8 x1=0.521179 x2=1.005525 x3=-0.375691 x4=-0.25966

20、9 3.3 第三题 Jacobi 迭代法的数据 A=10-1 2 0;0 8-1 3;2-1 10 0;-1 3-1 11;b=-11;-11;6;25;Jacobi(A,b,0.00005)x1=-1.467396 x2=-2.358678 x3=0.657604 x4=2.842397 3.4 第三题用 Gauss_Seidel 迭代的数据 A=10-1 2 0;0 8-1 3;2-1 10 0;-1 3-1 11;b=-11;-11;6;25;Gauss_Seidel(A,b,0.00005)x1=-1.467357 x2=-2.358740 x3=0.657597 x4=2.84240

21、5 四、实验中存在的问题及解决方案 4.1 存在的问题（1）第一题中在 matlab 中输入 Gauss(A,b)（数据省略）得到 m=4 n=4?Undefined function or variable Ab.Error in=Gauss at 8ap,p=max(abs(Ab(k:n,k);没有得到想要的结果。（2）第二题中在 matlab 中输入 y=LU(A,b)得到 y=4.0000 6.0000-5.0000-3.3571不是方程组的解。解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文

22、档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重精品好文档，推荐学习交流 仅供学习与交流，如有侵权请联系网站删除 谢谢8（3）第三题中在用高斯赛德尔方法时在 matlab 中输入 Gauss-Seidel(A,b,eps)结果程序报错?Error using=Gauss Too many output arguments.得不到想要的结果。4.2 解决方案（1）针对第一题中由于程序的第二行漏了一个分号导致输出了 m 和 n 的值，第 8 行中将Ab 改为 A 问题就解决了。（

23、2）由于程序后面出现了矩阵 y 故输出的事矩阵 y 的值，但是我们要的事 x 的值，故只需要将 y 改成 x，或者直接把 y 去掉就解决了问题。（3）在 function 文件中命名不能出现“-”应该将其改为下划线“_”，所以将 M 文件名“Gauss-Seidel(A,b,eps)”改成“Gauss_Seidel(A,b,eps)”就解决问题了。五、心得体会 本次试验涉及到了用高斯列主元消去法，LU 分解法，Jacobi 迭代法以及 Gauss-seidel迭代法等四种方法。需要对这些方法的原理都要掌握才能写出程序，由于理论知识的欠缺，我花了很大一部分时间在看懂实验的原理上，看懂了实验原理之后就开始根据原理编写程序，程序中还是出现了很多的低级错误导致调试很久才能运行。通过这次试验使我深刻的体会到理论知识的重要性，没有理论知识的支撑是写不出程序来的。写程序时还会犯很多低级的错误，以后一定要加强理论知识的学习，减少编程时低级错误的产生。解方程组用分解法解方程组其中分别用迭代法和迭代法求解方程组二实验原理程序框图程序代码等实验原理高斯列主 交流如有侵权请联系网站删除谢谢精品好文档推荐学习交流仅供学习与交流如有侵权请联系网站删除谢谢这个过程就 的过程中有可能会出现的情况这时消元就无法进行了即使主元数但是很小时其做除数也会导致其他元素数量级的严重