数值分析第二章上机

《数值分析第二章上机》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。

1、用列主元消元法解线性方程组AX=b的MATLAB程序function[RA,RB,n,X]=liezhu(A,b)B=[Ab];n=length(b);RA=rank(A);RB=rank(B);zhica=RB-RA;ifzhica>0,disp('请注意：因为RA~=RB，所以此方程组无解.')returnendifRA==RBifRA==ndisp('请注意：因为RA=RB=n，所以此方程组有唯一解.')X=zeros(n,1);C=zeros(1,n+1);forp=1:n-1[Y,j]=max(abs(B(p:n,p)))

2、;C=B(p,:);B(p,:)=B(j+p-1,:);B(j+p-1,:)=C;fork=p+1:nm=B(k,p)/B(p,p);B(k,p:n+1)=B(k,p:n+1)-m*B(p,p:n+1);endendb=B(1:n,n+1);A=B(1:n,1:n);X(n)=b(n)/A(n,n);forq=n-1:-1:1X(q)=(b(q)-sum(A(q,q+1:n)*X(q+1:n)))/A(q,q);endelsedisp('请注意：因为RA=RB

3、工作窗口输入程序>>A=[111;139;1749];b=[6;5;2];[RA,RB,n,X]=liezhu(A,b)运行后输出结果请注意：因为RA=RB=n，所以此方程组有唯一解.RA=3RB=3n=3X=6.3750-0.3333-0.0417将矩阵A进行直接LU分解的MATLAB程序functionhl=zhjLU(A)[nn]=size(A);RA=rank(A);ifRA~=ndisp('请注意：因为A的n阶行列式hl等于零，所以A不能进行LU分解.A的秩RA如下:'),RA,hl=det(A);returnendifR

4、A==nforp=1:nh(p)=det(A(1:p,1:p));endhl=h(1:n);fori=1:nifh(1,i)==0disp('请注意：因为A的r阶主子式等于零，所以A不能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下：'),hl;RAreturnendendifh(1,i)~=0disp('请注意：因为A的各阶主子式都不等于零，所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下：')forj=1:nU(1,j)=A(1,j);endfork=2:nfori=2:nforj=2:nL(1,1)=1;

5、L(i,i)=1;ifi>jL(1,1)=1;L(2,1)=A(2,1)/U(1,1);L(i,1)=A(i,1)/U(1,1);L(i,k)=(A(i,k)-L(i,1:k-1)*U(1:k-1,k))/U(k,k);elseU(k,j)=A(k,j)-L(k,1:k-1)*U(1:k-1,j);endendendendhl;RA,U,Lendend习题3.41.(1)在MATLAB工作窗口输入程序>>A=[24-6;153;132];hl=zhjLU(A)运行后输出结果请注意：因为A的各阶主子式都不等于零，所以A能进行LU分解.

6、A的秩RA和各阶顺序主子式值hl依次如下：RA=3U=2.00004.0000-6.000005.00006.0000003.8000L=1.0000000.50001.000000.50000.20001.0000hl=2618(2)在MATLAB工作窗口输入程序>>A=[116;-129;1-23];hl=zhjLU(A)运行后输出结果请注意：因为A的各阶主子式都不等于零，所以A能进行LU分解.A的秩RA和各阶顺序主子式值hl依次如下：RA=3U=1.00001.00006.000002.000015.00000019.5000

7、L=1.000000-1.00001.000001.0000-1.50001.0000hl=1336用P范数讨论AX=b解和A的性态的MATLAB程序functionAcp=zpjwc(A,jA,b,jb,p)Acp=cond(A,p);dA=det(A);X=Ab;dertaA=A-jA;PndA=norm(dertaA,p);dertab=b-jb;Pndb=norm(dertab,p);ifPndb>0jX=Ajb;Pnb=norm(b,p);PnjX=norm(jX,p);dertaX=X-jX;PnjdX=norm(d

8、ertaX,p);jxX=PnjdX/PnjX;PnjX=norm(jX,p);PnX=norm(X,p);jxX=PnjdX/PnjX;xX=PnjdX/PnX;Pndb=norm(dertab,p);xAb=Pndb/Pnb;Pn