数值分析实验综述.docx
- 文档编号:153180
- 上传时间:2022-10-04
- 格式:DOCX
- 页数:26
- 大小:813.98KB
数值分析实验综述.docx
《数值分析实验综述.docx》由会员分享,可在线阅读,更多相关《数值分析实验综述.docx(26页珍藏版)》请在冰豆网上搜索。
数值实验
目录
数值实验 1
第一章实验 1
1.根号5的极限 1
2.求pi的近似值 2
3.画图y=sinx与其泰勒展开 3
4.三维图mesh与surf实现。
4
第二章实验 6
1.列主元三角分解A. 6
2.追赶法求方程,电路的电流 8
3.方程组的性态与矩阵条件数的实验 9
4.Wilson矩阵。
求特征值,条件数。
误差 12
第三章实验 15
1.分别用Jacobi,Gauss-Seild,共轭梯度法解方程 15
2.利用共轭梯度法计算矩阵--10^5阶 19
3.利用cgs,bicg,bicgstab,等计算矩阵解 20
第五章实验 22
1.Newton插值f(x)=1/(1+4x^2).图形,误差 22
2.f(x)=1/(1+4x^2)插值 24
3.飞机的外形轮廓 24
第九章四阶龙格库塔方法(自选) 25
1. 解算微分方程组 25
2. Matlab的四阶Rk方法 25
3. 作图比较,一定误差 27
第一章实验
1.根号5的极限
代码如下:
x0=5;
fori=0:
1:
100
xi=sqrt(x0);
ei=x0-xi;
fprintf('NO:
%d,xi=%0.8f,Theerroris:
%0.8f\n',i,xi,ei);
x0=xi;
ifei<10^(-8)
break;
end
end
可以看出极限是1.
2.求pi的近似值
代码如下:
sum0=0;
forn=1:
1:
50000
y=(-1)^(n+1)/(2*n-1);
sum1=y+sum1;
pi_1=4*sum0;
pi_2=4*sum1;
error=pi_2-pi_1;
sum0=sum1;
fprintf('NO:
%d,pi_2=%0.8f,Theerror:
%0.8f\n',n,pi_2,error);
ifabs(error)<10^(-4)
break;
end
end
求得结果是:
3.画图y=sinx与其泰勒展开
x=0:
pi/100:
2*pi;
y=sin(x);
y1=0;
y2=0;
y3=0;
fori=0:
2
y1=y1+(-1)^(i)*x.^(2*i+1)/factorial(2*i+1);
end
fori=0:
5
y2=y2+(-1)^(i)*x.^(2*i+1)/factorial(2*i+1);
end
fori=0:
10
y3=y3+(-1)^(i)*x.^(2*i+1)/factorial(2*i+1);
end
plot(x,y,'*r',x,y1,'b',x,y2,'-g',x,y3,'k')
axis([02*pi-1.51.5])
看出n=2发散,n=10几乎与y=sinx重合。
4.三维图mesh与surf实现。
首先用mesh函数.
fori=0:
1:
2;
k=[0.2,0.1,0.05];
[xi,yi]=meshgrid(-10:
k(1,i+1):
10);
z=exp(-abs(xi))+cos(xi+yi)+1./(xi.^2+yi.^2+1);
subplot(2,2,i+1)
title('mesh')
mesh(xi,yi,z)
end
再用surf函数.
fori=0:
1:
2;
k=[0.2,0.1,0.05];
[xi,yi]=meshgrid(-10:
k(1,i+1):
10);
z=exp(-abs(xi))+cos(xi+yi)+1./(xi.^2+yi.^2+1);
subplot(2,2,i+1)
title('mesh')
mesh(xi,yi,z)
end
第二章实验
1.列主元三角分解A.
clc
A=[11111;12345;1361015;14102035;15153570]
E=eye(5)
[m,n]=size(A);
ifm~=n
error('不是方阵')
return
end
ifdet(A)==0
error('不能分解')
end
u=A;p=eye(m);l=eye(m);
fori=1:
m
forj=i:
m
t(j)=u(j,i);
fork=1:
i-1
t(j)=t(j)-u(j,k)*u(k,i);
end
end
a=i;b=abs(t(i));
forj=i+1:
m
ifb b=abs(t(j)); a=j; end end ifa~=i forj=1: m c=u(i,j); u(i,j)=u(a,j); u(a,j)=c; end forj=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); forj=i+1: m u(j,i)=t(j)/t(i); end forj=i+1: m fork=1: i-1 u(i,j)=u(i,j)-u(i,k)*u(k,j); end end end l=tril(u,-1)+eye(m) u=triu(u,0) B=A/E 2.追赶法求方程,电路的电流 A=[2-2000000;-25-200000;0-25-20000;00-25-2000; 000-25200;0000-25-20;00000-25-2;000000-25] bb=[220/270000000] n=size(A,1); s=zeros(n,1); %-----È¡³öÈý¶Ô½Ç-------- b=diag(A); a=diag(A,-1); c=diag(A,1); d=zeros(n,1); u=zeros(n-1,1); fori=1: n-1 d (1)=b (1); u(i)=c(i)/d(i); d(i+1)=b(i+1)-a(i)*u(i); end %-----×·µÄ¹ý³Ì------------ y=zeros(n,1); y (1)=bb (1)/d (1); fori=2: n y(i)=(bb(i)-a(i-1)*y(i-1))/d(i); end %-----¸ÏµÄ¹ý³Ì--------------- s(n)=y(n); fori=n-1: -1: 1 s(i)=y(i)-u(i)*s(i+1); end s 3.方程组的性态与矩阵条件数的实验 clc n=5; A=zeros(n); C=zeros(n); b=zeros(1,n); fori=1: n x(i)=1+0.1*i; forj=1: n A(i,j)=x(i)^(j-1); C(i,j)=1/(i+j-1); b(i)=b(i)+A(i,j); end end A C b d=cond(A,2) DD=A\b' (1)当n=5,10,20时,cond(A)=5.3615e+05,8.6823e+11,1.5428e+22。 看出越来越病态了。 (2)当n=5,10,20时,解分别如下: 说明条件数越大,越病态,发散了。 当n=10时就一个值有一定误差。 (3)n=10,解方程(A+delta)x=b. clc n=10; A=zeros(n); C=zeros(n); b=zeros(1,n); fori=1: n x(i)=1+0.1*i; forj=1: n A(i,j)=x(i)^(j-1); %C(i,j)=1/(i+j-1); b(i)=b(i)+A(i,j); end end A(2,2)=A(2,2)+10^(-8) A(10,10)=A(10,10)+10^(-8) A d=cond(A,2) x=A\b' 那么结果如下: 可以看出很小的扰动导致解误差变得较大了。 4.Wilson矩阵。 求特征值,条件数。 误差 (1)行列式,条件数,特征值 (2)(A+A0)(x+x0)=b,求此时的x0与||x0||2. 代码如下: clc A=[10787;7565;86109;75910] A1=[107.28.16.9;7.085
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 实验 综述