数值分析实验报告二.doc

《数值分析实验报告二.doc》由会员上传分享，免费在线阅读，更多相关内容在工程资料-天天文库。

1、数值实验报告二一、实验名称解线性方程组的列主元素高斯消去法和LU分解法二、实验目的通过数值实验，从屮体会解线性方程组选主元的必耍性和LU分解法的优点,以及方程组系数矩阵和右端向量的微小变化对解向量的影响。三、实验内容解下列两个线性方程纽'3.016.031.99、(2)四、算法描述1、列主元素高斯消去法记:aij}=atj(i,j=1,2,3…斤)=

2、bj(i=1,2,3)消元过程：对J*k=1,2,3⑴选行号L'使

3、哦卜max

4、。絆

5、。(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

6、=yo先由Ly=b解出向量y,再由Ux二y解出向量x,即为原方程组Ax二b的解。五、程序流程图1、列主元素高斯消去法的M文件主程序:function[RA,RB,n,X]=1iezhu(A,b)B=[Ab];n=length(b);RA=rank(A);RB=rank(B);zhica=RB-RA;ifzhica>0,dispC请注意：因为RA〜=RB,所以此方程组无解.)returnendifRA==RBifRA==ndisp(请注意：因为RA=RB=n,所以此方程组有唯一解・')X=zeros(n,1);C=zeros(1zn+1);fo

7、rp=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;fork=p+l:nm=B(kzp)/B(p,p)；B(kzp:n+l)=B(kzp:n+1)-m*B(p,p:n+l);endendb=B(1:n,n+l);A=B(1:n,1:n);X(n)=b(n)/A(n,n);forq=n-l:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(qrq)；endelsedisp「请注意：因为RA=RB

8、穷多解.')endend在MATLABT作窗口输入的程序为：»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=

9、1;1;1],det(A)JRA,RB,n,X]=liezhu(A,b)2、LU分解法的M文件程序：functionhl=zhijieLU(A)[n

10、n]=size(A);RA=rank(A);ifRA〜=ndisp(谛注意：因为A的n阶行列式hl等丁•零，所以A不能进行LU分解.A的秩RA如下门RA,hl=det(A);returnendifRA==nforp=l:nh(p)=det(A(l:p,1:p));endhl=h(1:n);fori=l:nifh(lzi)==0disp「请注意:因为A的i•阶主子式等丁•零，所以A不能进行LU分解.A的秩RA和各阶顺序主了式值hl依次如下门returnendendifh(l,i)~=0dispCffl注意:因为A的各阶主子式都不等丁•家所以A能

11、进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下门forj=l:nU(l,j)=A(lzj)；endfork=2:nfori=2:nforj=2:nL(l,l)=l;L(i,i)=l;ifi>jL(1Z1)=1；L(2Z1)=A(2,1)/U(1,1)；L(izl)=A(izl)/U(l,l)；L(i,k)=(A(izk)-L(i,l:k-l)*U(1:k-l,k))/U(kzk)；elseU(k,j)=A(k,j)-L(k,l:k-l)*U(l:k-l,j);endendendendhl;RA,U,Lendend在MATLABT作窗口

12、输入的程序为：»A=ll0,-7,0,l；-3,2.099999,6,2;5,-l,5,-l;2,l,0,2J,b=[8;5.900001;5;lJ,det(A),