大连理工大学矩阵与数值分析上机作业.docx

1、矩阵与数值分析上机作业 学校： 大连理工大学 学院： 班级： 姓名： 学号： 授课老师： 注：编程语言Matlab程序：Norm.m函数function s=Norm(x,m)%求向量x的范数%m取1,2,inf分别 表示1,2，无穷范数n=length(x);s=0;switch m case 1 %1-范数 for i=1:n s=s+abs(x(i); end case 2 %2-范数 for i=1:n s=s+x(i)2; end s=sqrt(s); case inf %无穷-范数 s=max(abs(x);end 计算向量x，y的范数Test1.mclear all;clc;n1

2、=10;n2=100;n3=1000;x1=1./1:n1;x2=1./1:n2;x3=1./1:n3;y1=1:n1;y2=1:n2;y3=1:n3;disp(n=10时);disp(x的1-范数:);disp(Norm(x1,1);disp(x的2-范数:);disp(Norm(x1,2);disp(x的无穷-范数:);disp(Norm(x1,inf);disp(y的1-范数:);disp(Norm(y1,1);disp(y的2-范数:);disp(Norm(y1,2);disp(y的无穷-范数:);disp(Norm(y1,inf);disp(n=100时);disp(x的1-范数:)

3、;disp(Norm(x2,1);disp(x的2-范数:);disp(Norm(x2,2);disp(x的无穷-范数:);disp(Norm(x2,inf);disp(y的1-范数:);disp(Norm(y2,1);disp(y的2-范数:);disp(Norm(y2,2);disp(y的无穷-范数:);disp(Norm(y2,inf);disp(n=1000时);disp(x的1-范数:);disp(Norm(x3,1);disp(x的2-范数:);disp(Norm(x3,2);disp(x的无穷-范数:);disp(Norm(x3,inf);disp(y的1-范数:);disp(N

4、orm(y3,1);disp(y的2-范数:);disp(Norm(y3,2);disp(y的无穷-范数:);disp(Norm(y3,inf);运行结果：n=10时x的1-范数:2.9290；x的2-范数:1.2449； x的无穷-范数:1y的1-范数:55； y的2-范数:19.6214； y的无穷-范数:10n=100时x的1-范数:5.1874；x的2-范数: 1.2787； x的无穷-范数:1y的1-范数:5050； y的2-范数:581.6786； y的无穷-范数:100n=1000时x的1-范数:7.4855； x的2-范数:1.2822； x的无穷-范数:1y的1-范数: 50

5、0500； y的2-范数:1.8271e+004；y的无穷-范数:1000程序Test2.mclear all;clc;n=100;%区间h=2*10(-15)/n;%步长x=-10(-15):h:10(-15);%第一种原函数f1=zeros(1,n+1);for k=1:n+1 if x(k)=0 f1(k)=log(1+x(k)/x(k); else f1(k)=1; endendsubplot(2,1,1);plot(x,f1,-r);axis(-10(-15),10(-15),-1,2);legend(原图);%第二种算法f2=zeros(1,n+1);for k=1:n+1 d=1

6、+x(k); if(d=1) f2(k)=log(d)/(d-1); else f2(k)=1; endendsubplot(2,1,2);plot(x,f2,-r);axis(-10(-15),10(-15),-1,2);legend(第二种算法);运行结果：显然第二种算法结果不准确，是因为计算机中的舍入误差造成的，当时，计算机进行舍入造成恒等于1，结果函数值恒为1。程序：秦九韶算法:QinJS.mfunction y=QinJS(a,x)%y输出函数值%a多项式系数，由高次到零次%x给定点n=length(a);s=a(1);for i=2:n s=s*x+a(i);endy=s;计算p(

7、x):test3.mclear all;clc;x=1.6:0.2:2.4;%x=2的邻域disp(x=2的邻域：);xa=1 -18 144 -672 2016 -4032 5376 -4608 2304 -512;p=zeros(1,5);for i=1:5 p(i)=QinJS(a,x(i);enddisp(相应多项式p值：);pxk=1.95:0.01:20.5;nk=length(xk);pk=zeros(1,nk);k=1;for k=1:nk pk(k)=QinJS(a,xk(k);endplot(xk,pk,-r);xlabel(x);ylabel(p(x);运行结果：x=2的

8、邻域：x =1.6000 1.8000 2.0000 2.2000 2.4000相应多项式p值：p = 1.0e-003 * -0.2621 -0.0005 0 0.0005 0.2621p(x)在1.95,20.5上的图像程序：LU分解，LUDecom.mfunction L,U=LUDecom(A)%不带列主元的LU分解N = size(A);n = N(1);L=eye(n);U=zeros(n);for i=1:n U(1,i)=A(1,i);L(i,1)=A(i,1)/U(1,1);endfor i=2:n for j=i:n z=0; for k=1:i-1 z=z+L(i,k)*

9、U(k,j); end U(i,j)=A(i,j)-z; end for j=i+1:n z=0; for k=1:i-1 z=z+L(j,k)*U(k,i); end L(j,i)=(A(j,i)-z)/U(i,i); endendPLU分解，PLUDecom.mfunction P,L,U =PLUDecom(A)%带列主元的LU分解m,m=size(A);U=A;P=eye(m);L=eye(m);for i=1:m for j=i:m t(j)=U(j,i); for k=1:i-1 t(j)=t(j)-U(j,k)*U(k,i); end end a=i;b=abs(t(i); fo

10、r j=i+1:m if babs(t(j) b=abs(t(j); a=j; end end if a=i for j=1:m c=U(i,j); U(i,j)=U(a,j); U(a,j)=c; end for j=1:m c=P(i,j); P(i,j)=P(a,j); P(a,j)=c; end c=t(a); t(a)=t(i); t(i)=c; end U(i,i)=t(i); for j=i+1:m U(j,i)=t(j)/t(i); end for j=i+1:m for k=1:i-1 U(i,j)=U(i,j)-U(i,k)*U(k,j); end endendL=tril

11、(U,-1)+eye(m);U=triu(U,0);（1） （2）程序：Test4.mclear all;clc;for n=5:30 x=zeros(n,1); A=-ones(n); A(:,n)=ones(n,1); for i=1:n A(i,i)=1; for j=(i+1):(n-1) A(i,j)=0; end x(i)=1/i; end disp(当n=);disp(n); disp(方程精确解:); x b=A*x; %系数b disp(利用LU分解方程组的解:); L,U=LUDecom(A); %LU分解 xLU=U(Lb) disp(利用PLU分解方程组的解:); P,L,U =PLUDecom(A); %PLU分解 xPLU=U(L(Pb) %求解A的逆矩阵 disp(A的准确逆矩阵:); InvA=inv(A) InvAL=zeros(n); %利用LU分解求A的逆矩阵 I=eye(n); for i=1:n InvAL(:,i)=U(LI(:,i); end disp(利用LU分解的A的逆矩阵:); InvALEnd运行结果：（1） 只列出n=5,6,7的结果