列主元高斯消去 LU分解 迭代法_高斯列主元消去法步骤

列主元高斯消去 LU分解 迭代法由刀豆文库小编整理，希望给你工作、学习、生活带来方便，猜你可能喜欢“高斯列主元消去法步骤”。

数学与信息科学学院

实 验 报 告

课程名称： 《 数值计算方法 》 实验名称： 数值积分 实验类型： 验证性■综合性□设计性□ 实验室名称： 数学实验室 班级学号： 090721 学生姓名：

任课教师（教师签名）： 成绩：

实验日期： 2012-4-27

一、实验目的及题目

实验目的：

南昌航空大学数学与信息科学学院实验报告

熟悉线性方程组求解原理，运用列主消元高斯消去法，LU分解法及Jacobi迭代法与Gau-Seidel迭代法丢出线性方程的解 实验题目：

1、用列主元消去法解方程组：

x1x23x442xxxx1123

4

3xxx3x32341x12x23x3x442、用LU分解法界方程组Axb，其中

48240124122412124

A,b

062022626216

3、分别用jacobi迭代法和gau-seidel迭代法求解方程组：

10x1x22x3118xx3x1123

4

2xx10x6231x13x2x311x425

二、实验原理、程序框图、程序代码等

实验原理：

1、列主元高斯消去原理：

每次消去中，选取每列的绝对值最大的元素作为消去对象并作为主元素。然后换行使之变到主元位子上，在进行消元计算。设A(k)xb(k)，确定第k列主元所在位置ik，在交换ik行和k行后，在进行消元。根据矩阵理论，交换ik和k两个方程的位置，列主元素的消去过程相当于对交换后的新矩阵进行消元，即

LkIk,ikAkA(k1)

同时，右端向量b(k)变化为 LkIk,ikbkb(k1)

2、LU分解原理：

第 页

南昌航空大学数学与信息科学学院实验报告

直接三角分解：先将A分解为上三角矩阵U和下三角矩阵L（A=LU），则原为题等价于求解两个方程组：Lyb,Uxy3、迭代法

（1）jacobi迭代法（又称简单迭代法）考虑n阶线性代数方程组

a11x1a12x2.....a1nxnb1axax.....axb2112222nn

2

an1x1an2x2.....annxnbn其矩阵形式为

Axb

设该方程组的系数矩阵A非奇异且aii0(i1,2,.....,n)，可将A分解为：

ADLU

其中Ddiag[a11,a22,.....,ann]；

0a21La31an1 0a31an10an10a120，U=－0a13a230a1na2n

an1,n0然后化为如下等价形式

xD1(LU)xD1b

简记为

xBJxfJ

选取初始向量x(0)x1,x2(0)(0)xn(0)T,通过上式可得到线性代数方程组的迭代格式

x(k1)BJx(k)fJ,(k0,1,2)

其分量形式为

xi即(k1)n1k(biaijxj),(i1,2,n)(k0,1,2)aiii1j1第页

南昌航空大学数学与信息科学学院实验报告

(k1)1(k)(k)(k)x(baxaxax)11221331nn1a111(k1)(k)(k)(k)x(baxaxax)222112332nna 221(k1)(k)(k)(k)x(baxaxax)nnn1n2n,n1n112ann

以上计算过程称迭代法，矩阵BJ陈为jacobi迭代法的迭代矩阵。

（2）Gau-Seidel迭代法 将jacobi迭代格式改为如下形式

xi(k1)i1n1(k)(k)(biaijxjaijxj),(i1,2,n)(k0,1,2)aiii1ji1其矩阵形式为

x(k1)(DL)1Ux(k)(DL)1b

于是得Gau-Seidel迭代公式为

x(k)BGx(k)fG,(k0,1,2)

程序代码： 第一题程序代码： function GEpiv(A,b)[m,n]=size(A);nb=n+1;Ab=[A b];for i=1:m-1 [pivot,p]=max(abs(Ab(i:n,j)));ip=p+i-1;if ip~=i Ab([i ip],:)=Ab([ip i],:);end

第 页 称BG为Gau-Seidel迭代法的迭代矩阵。

南昌航空大学数学与信息科学学院实验报告

pivot=Ab(i,i);for k=i+1:m Ab(k,i:nb)=Ab(k,i:nb)-(Ab(k,i)/pivot)*Ab(i,i:nb);end end x=zeros(n,1);x(n)=Ab(n,nb)/Ab(n,n);for i=n-1:-1:1 x(i)(Ab(i,nb)-Ab(i,i+1:n)*x(i+1:n,1))/Ab(i,i);end for k=1:n fprintf('x[%d]=%fn',x(k));end 运行结果为 输入

A=[1 1 0 3;2 1-1 1;3-1-1 3;-1 2 3-1];b=[4;1;-3;4];得结果为 x =

-1.3333 2.3333-0.3333 1.0000

第二题代码：

function x=lufj(A,b)%lu分解求值函数 [l u p]=lu(A);y=lx(l,b);x=us(u,y);end

第页

南昌航空大学数学与信息科学学院实验报告

function y=lx(A,b)%求系数矩阵为下三角方程函数 [n,m]=size(A);y(1)=b(1)/A(1,1);for i=2:n s=0;for k=1:i-1 s=s+A(i,k)*y(k);end y(i)=(b(i)-s)/A(i,i);end end

function x=us(A,b)%求系数矩阵为上三角方程函数 [n,m]=size(A);x(n)=b(n)/A(n,n);for i=0:n-2 s=0;for k=0:i s=s+A(n-i-1,n-k)*x(n-k);end x(n-i-1)=(b(n-i-1)-s)/A(n-i-1,n-i-1);end end 输入矩阵

A=[48-24 0 12;-24 4 12 12;0 6 20 2;-6 6 2 16];b=[4;4;-2;-2];结果为 x =

第页

南昌航空大学数学与信息科学学院实验报告

0.5212 1.0055-0.3757-0.2597

第三题代码：（1）Jacobi迭代：

function x = agui_jacobi(a,b)n=length(b);N=100;e=1e-4;x0=zeros(n,1);x=x0;x0=x+2*e;k=0;d=diag(diag(a));l=-triu(a,-1);u=-triu(a,1);while norm(x0-x,inf)>e&k

（2）Gau-seidel迭代 function x=gau(A,b,x)[d u l]=fenjie(A);B=inv(d-l)*u;f=inv(d-l)*b;

第页

南昌航空大学数学与信息科学学院实验报告

for i=1:9 y=x;x=B*y+f;end

function [d u l]=fenjie(A)[n m]=size(A);for i=1:n for j=1:n d(i,j)=0;end end for i=1:n d(i,i)=A(i,i);end for i=1:n for j=1:i l(i,j)=-A(i,j);u(j,i)=-A(j,i);end for k=i:n l(i,k)=0;u(k,i)=0;end end（1）

输入矩阵：a=[10-1 2 0;0 8-1 3;2-1 10 0;-1 3-1 11];>> b=[-11;-11;6;25];>> x = agui_jacobi(a,b)结果为

第页

南昌航空大学数学与信息科学学院实验报告

x=-1.4674-2.3587 0.6576 2.8424

（2）A=[10-1 2 0;0 8-1 3;2-1 10 0;-1 3-1 11];b=[-11;-11;6;25];x=[-1.3;-2.1;0.55;2.66];gau(A,b,x)结果为 x=-1.4674-2.3587 0.6576 2.8424

四、实验中存在的问题及解决方案

命令需多次修改，如Gau迭代中，m文件中因end的不匹配使得程序无法运行处结果，则需细心调试，去除多余的end，才能使程序正常运行。

五、心得体会

在线性代数中，从理论上解决了线性方程组的求解问题，但常用的高斯消元法有时会导致结果产生严重的偏差，因此需要研究数值解法。Jacobi迭代与Gau-Seidel有一个共同特点：新的近似解是已知近似解的线性函数，并且新的近似解只与已知近似解有关。

编写程序是一个繁杂的过程，需要不断调试，并且善于发现程序中的细小问题，才能够正常运行处程序。

第页