整理LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.docx
《整理LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.docx》由会员分享，可在线阅读，更多相关《整理LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序.docx（13页珍藏版）》请在冰豆网上搜索。
整理LU分解法列主元高斯法Jacobi迭代法GaussSeidel法的原理及Matlab程序
精品文档
、实验目的及题目
1.1实验目的:
（1）学会用高斯列主元消去法，LU分解法，Jacobi迭代法和Gauss-Seide迭代法解线性
方程组。
（2）学会用Matlab编写各种方法求解线性方程组的程序。
1.2实验题目：
用列主元消去法解方程组：
iX”2+3x4=412为+X2-X3+&=1
12X|—x?
-X3+3x^=-3[-为+2X2+3X3-X4=4
1.
2.
用LU分解法解方程组
广48
-24
0-
12
2、
—24
24
12
，b=
4
0
6
20
2|
-2
<~6
6
2
16
、乂）
Ax=b,其中
A=
3.
分别用Jacobi迭代法和
Gauss-Seide迭代法求解方程组:
10%—X2+2x3=-11
8X2-X3+3x4=-11
2x^x2+10x3=6
一為+3x2—X3+11%=25
二、实验原理、程序框图、程序代码等
2.1实验原理
2.1.1高斯列主元消去法的原理
Gauss消去法的基本思想是一次用前面的方程消去后面的未知数，从而将方程组化为等
价形式:
》1必+Gx2+III+gXn=g
Ib22X^|（+b2nX^g2
bnnXn=gn
这个过程就是消元，然后再回代就好了。
具体过程如下:
对于k=1,2,川,n-1，若akJ工0,依次计算
精品文档
mik=aik）/akk）
ai（k—aj（k）-mikak：
）
b（5=b（k）-mkbkk）,i,j=k+1,川,n
然后将其回代得到:
Xn尽）/^
n
}xk=（bkk）-2akjk）Xj）/akk）,k=n—1,n—2,川,1
jM
以上是高斯消去。
但是高斯消去法在消元的过程中有可能会出现akk^0的情况，这时消元就无法进行了,
即使主元数akk^O,但是很小时，其做除数，也会导致其他元素数量级的严重增长和舍入误差的扩散。
因此，为了减少误差，每次消元选取系数矩阵的某列中绝对值最大的元素作为主元素。
然后换行使之变到主元位置上，再进行销元计算。
即高斯列主元消去法。
2.1.2直接三角分解法（LU分解）的原理
先将矩阵A直接分解为A=LU则求解方程组的问题就等价于求解两个三角形方程组。
直接利用矩阵乘法，得到矩阵的三角分解计算公式为:
|u1i=a1i,i二1,2,川,n
山=ai1/U11,i=2,ill,n
「kA
|Ukj=akj-送lkmUmj,,j=k,k+1,川，n
ml,k=2,3川In
lik=（ak—送limUmk）/Ukk,i=k+1,k+2」||,n且kHn
mzt
由上面的式子得到矩阵A的LU分解后，求解Ux=y的计算公式为
i-4
yi=D
IVi=bi—送lijyj,i=2,3，川n
[y
Xn=Vn/Unn
*n
Xi=（yi-2UijXj）/Uii,i=n—1,n—2,川，1j止十
以上为LU分解法。
精品文档
2.1.3Jacobi迭代法和Gauss-Seidel迭代法的原理
（1）Jcaobi迭代
设线性方程组
Ax=b
的系数矩阵A可逆且主对角元素a11,a22,...,ann均不为零，令
D=diag（an,a22,...,ann
并将A分解成
A=（A-D）中D
从而
（1）可写成
X=BiX+fi
-1
11
其中B=1-DA,f1=Db
以Bi为迭代矩阵的迭代法（公式）
称为雅可比Jacobi）迭代法,其分量形式为
（kdt）1r_J（k）
X=—Lb-壬aX八i^iiJj
a••J
u"J#
「=1,2,...n,k=01,2,...
其中x（J=（xFW）,..人（0斤为初始向量.
（2）Gauss-Seidel迭代
由雅可比迭代公式可知,在迭代的每一步计算过程中是用X（k）的全部分量来计算x（k+的
（k+）（k+）（k+）
所有分量,显然在计算第i个分量Xi时，已经计算出的最新分量X1,...,Xi~没有被利用。
把矩阵A分解成
A=D-L-U
其中D二diag（a11,a22,...,ann\~L"U分别为A的主对角元除外的下三角和上三角部分,
精品文档
于是,方程组
（1）便可以写成
（D-L风=Ux+b
X=B2X+f2
其中
以B2为迭代矩阵构成的迭代法（公式）
称为高斯一塞德尔迭代法，用分量表示的形式为
functionGauss（A,b）%A为系数矩阵，b为右端项矩阵
[m,n]=size（A）;
n=length（b）;
fork=1:
n-1
[pt,p]=max（abs（A（k:
n,k）））;%找出列中绝对值最大的数
p=p+k-1;
ifp>k
t=A（k,:
）;A（k,:
）=A（p,:
）;A（p,:
）=t;%交换行使之变到主元位置上
t=b（k）;b（k）=b（p）;b（p）=t;
A（k+1:
n,k+1:
n）=A（k+1:
n,k+1:
n）-m*A（k,k+1:
n）;
ifflag~=0
精品文档
精品文档
Ab=[A,b];
end
end
fork=n-1:
-1:
1
x（k）=（b（k）-A（k,k+1:
n）*x（k+1:
n））/A（k,k）;
end
fork=1:
n
fprintf（'x[%d]=%f\n',k,x（k））;
end
L（2:
n,1）=A（2:
n,1）/U（1,1）;
fork=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矩阵
%输出U矩阵
y=zeros（n,1）;%开始解方程组Ux=y
y
（1）=b
（1）;
fork=2:
n
y（k）=b（k）-L（k,1:
k-1）*y（1:
k-1）;
endx=zeros（n,1）;精品文档
精品文档x（n）=y（n）/U（n,n）;fork=n-1:
-1:
1
x（k）=（y（k）-U（k,k+1:
n）*x（k+1:
n））/U（k,k）;endfork=1:
n
fprintf（'x[%d]=%f\n',k,x（k））;
end
2.2.3Jacobi迭代法的程序
[m,n]=size（A）;
temp=1;
x=zeros（m,1）;
k=0;
whileabs（max（x）-temp）>eps
temp=max（abs（x））;
k=k+1;
%记录循环次数
%雅克比迭代公式
x=-inv（D）*（L+U）*x+inv（D）*b;
endfork=1:
n
fprintf（'x[%d]=%f\n',k,x（k））;
end
精品文档
2.2.4Gauss-Seidel迭代程序
functionGauss_Seidel（A,b,eps）%A为系数矩阵，b为后端项矩阵，epe为精度
[m,n]=size（A）;
temp=1;
x=zeros（m,1）;
k=0;
whileabs（max（x）-temp）>eps
temp=max（abs（x））;
endfork=1:
n
fprintf（'x[%d]=%f\n',k,x（k））;
end
三、实验过程中需要记录的实验数据表格
3.1第一题（高斯列主元消去）的数据>>A=[1103;21-11;3-1-13;-123-1];>>b=[4;1;-3;4];
>>Gauss（A,b）
x[1]=-1.333333
x[2]=2.333333
x[3]=-0.333333
x[4]=1.0000003.2第二题（LU分解法）数据
>>A=[48-240-12;-24241212;06202;-66216];
>>b=[4;4;-2;-2];
>>LU（A,b）
精品文档
1.0000
-0.5000
0
-0.1250
0
1.0000
0.5000
0.2500
0
0
1.0000
-0.0714
0
0
0
1.0000
-24.0000
12.0000
0
0
-12.0000
6.0000
-1.0000
12.9286
x[1]=0.521179
x[2]=1.005525
x[3]=-0.375691
x[4]=-0.259669
>>A=[10-120;08-13;2-1100;-13-111];b=[-11;-11;6;25];
Jacobi（A,b,0.00005）
x[1]=-1.467396
x[2]=-2.358678
x[3]=0.657604
x[4]=2.8423973.4第三题用Gauss_Seide迭代的数据>>A=[10-120;08-13;2-1100;-13-111];
>>b=[-11;-11;6;25];
>>Gauss_Seidel（A,b,0.00005）
x[1]=-1.467357
x[2]=-2.358740
x[3]=0.657597
x[4]=2.842405
四、实验中存在的问题及解决方案
4.1存在的问题
（1）第一题中在matlab中输入》Gauss（A,b）（数据省略）得到m=4n=4?
?
?
Undefinedfunctionorvariable"Ab".Errorin==>Gaussat8[ap,p]=max（abs（Ab（k:
n,k）））;没有得到想要的结
精品文档
果。
2）第二题中在matlab中输入>>y=LU（A,b）得到y=4.00006.0000-5.0000-3.3571不是方程组的解。
（3）第三题中在用高斯赛德尔方法时在matlab中输入>>Gauss-Seidel（A,b,eps结果程序报
错?
?
?
Errorusing==>GaussToomanyoutputargumentS得不至U想要的结果。
4.2解决方案
（1）针对第一题中由于程序的第二行漏了一个分号导致输出了m和n的值，第8行中将
Ab改为A问题就解决了。
（2）由于程序后面出现了矩阵y故输出的事矩阵y的值，但是我们要的事X的值，故只需
要将y改成X,或者直接把y去掉就解决了问题。
3）在function文件中命名不能出现“-”应该将其改为下划线“_”，所以将M文件名
Gauss-Seidel（A,b,eps”）改成“Gauss_Seidel（A,b,eps”）就解决问题了。
五、心得体会
本次试验涉及到了用高斯列主元消去法，LU分解法，Jacobi迭代法以及Gauss-seide迭
代法等四种方法。
需要对这些方法的原理都要掌握才能写出程序,由于理论知识的欠缺,我
花了很大一部分时间在看懂实验的原理上，看懂了实验原理之后就开始根据原理编写程序，程序中还是出现了很多的低级错误导致调试很久才能运行。
通过这次试验使我深刻的体会到理论知识的重要性,没有理论知识的支撑是写不出程序精品文档
精品文档
来的。
写程序时还会犯很多低级的错误，以后一定要加强理论知识的学习，减少编程时低级错误的产生。
升级会员 