数值分析实验报告二.doc

数值实验报告二 一、 实验名称 解线性方程组的列主元素高斯消去法和LU分解法 二、 实验目的 通过数值实验，从屮体会解线性方程组选主元的必耍性和LU分解法的优点, 以及方程组系数矩阵和右端向量的微小变化对解向量的影响。 三、实验内容 解下列两个线性方程纽 3.01 6.03 1.99、 r (1) 1.27 4.16 -1.23 兀2 — 1 、0.987 -4.81 9.34 丿 1 ‘10 -7 0 1、 / \ 8 、 -3 2.099999 6 2 兀2 5.900001 5 -1 5 -1 兀3 5 、2 \ I 0 2丿 卫4丿 1 (2) 四、算法描述 1、列主元素高斯消去法 记: a\ij} = atj (i, j = 1,2,3 …斤) =bj (i = 1,2,3 ) 消元过程：对J * k = 1,2,3 ⑴选行号L使|哦卜max|。絆|。 (2)交换硝)与4? (j = k,k+l,k+2??F)以及附与曙所含的数值。 (3)对于 i = k,k+l,k+2…〃,计算 回代过程： 心=（4）一工砧％）/血） （k = n-l,n-2, n-3 …1） j=k+\ 在此算法中的4?称为第k个列主元素，它的数值总耍被交换到第k个主对 角线元素的位置上。 2、LU分解法 通过MATLAB 口有的函数，把系数矩阵A分解成A=LU,其中：L是下三 角矩阵，U是上三和矩阵，这时方程组Ax=b就可以分解成两个容易求解的三角 形方程组Ly=b, Ux=yo先由Ly=b解出向量y,再由Ux二y解出向量x,即为原 方程组Ax二b的解。 五、程序流程图 1、列主元素高斯消去法的M文件主程序: function [RA,RB,n,X]=1iezhu(A,b) B=[A b]; n=length(b); RA=rank(A); RB=rank(B);zhica=RB-RA; if zhica0, dispC请注意：因为RA?=RB,所以此方程组无解.) return end if RA==RB i f RA==n disp (请注意：因为RA=RB=n,所以此方程组有唯一解?) X=zeros(n,1); C=zeros(1z n+1); for p= 1:n-1 [Y,j]=max(abs(B(p:n,p))); C=B (p,:); B (p, ：)= B(j+p-1, :); B ( j+p-1, :)=C; for k=p+l:n m= B(kzp)/ B(p,p)； B(kz p:n + l)= B(kz p:n + 1)-m* B (p,p:n + l); end end b=B(1:n,n+l);A=B(1:n,1:n); X (n) =b (n) /A (n, n); for q=n-l:-1:1 X(q) = (b(q)-sum(A(q,q+1:n) *X(q+1:n) ) )/A(qr q)； end else disp「请注意：因为RA=RBn,所以此方程组有无穷多解.) end end 在MATLAB T作窗口输入的程序为： ?A=[3.01,6.03,1.99;1.27,4.16,-1.23;0.987,-4.81,9.34],b=[l;l;l],det(A),[RA,RB, n,X]=liezhu(A,b) ★变化系数矩阵在MATLAB工作窗口中输入： ?A=[3.00,6.03,1.99; 1.27,4」6,-1,23;0.987,-4.81,9.34],b=| 1; 1; 1 ],det( A)JR A,RB, n,X]=liezhu(A,b) 2、LU分解法的M文件程序： function hl = zhi jieLU(A) [n n]=size(A);RA=rank(A); if RA?=n disp (谛注意：因为A的n阶行列式hl等丁?零，所以A不能进行LU分 解.A的秩RA如下门 RA,hl=det (A); return end if RA==n for p=l:n h(p)=det(A(l:p,1:p)); end hl=h (1:n); for i=l:n if h(lzi)==0 disp「请注意:因为A的i?阶主子式等丁?零，所以A不能进行LU分 解.A的秩RA和各阶顺序主了式值hl依次如下门 return end end if h(l,i)~=0 dispCffl注意:因为A的各阶主子式都不等丁?家所以A能进行LU分 解.A的秩RA和各阶顺序主子式值hl依次如下门 for j=l:n U(l, j)=A(lz j)； end for k=2:n for i=2:n for j=2:n L(l,l)=l;L(i,i)=l; if i j L(1Z1)=1；L(2Z1)=A(