计算方法综合实践.docx
- 文档编号:27507205
- 上传时间:2023-07-02
- 格式:DOCX
- 页数:18
- 大小:99.78KB
计算方法综合实践.docx
《计算方法综合实践.docx》由会员分享,可在线阅读,更多相关《计算方法综合实践.docx(18页珍藏版)》请在冰豆网上搜索。
计算方法综合实践
1.应用自己熟悉的算法语言编写程序,使之尽可能具有通用性。
2.上机前充分准备,复习有关算法,写出计算步骤,反复检查,调试程序。
(注:
在练习本上写,不上交)
3.完成计算后写出实验报告,容包括:
算法步骤叙述,变量说明,程序清单,输出计算结果,结构分析和小结等。
(注:
具体题目具体分析,并不是所有的题目的实验报告都包含上述容!
)
4.独立完成,如有雷同,一律判为零分!
5.上机期间不允许做其他任何与课程设计无关的事情,否则被发现一次扣10分,被发现三次判为不及格!
非特殊情况,不能请假。
旷课3个半天及以上者,直接判为不及格。
一、基本技能训练
1、误差分析
实验1.2误差传播与算法稳定性
实验目的:
体会稳定性在选择算法中的地位。
误差扩的算法是不稳定的,是我们所不期望的;误差衰减的算法是稳定的,是我们努力寻求的,这是贯穿本课程的目标。
问题提出:
考虑一个简单的由积分定义的序列
显然
。
当
时,
而对于
时,利用分部积分易得
另一方面,我们有
实验容:
由以上递推关系,我们可得到计算序列
的两种方法。
(I)
(II)
symsnIn5In6In7;
In5=vpa((exp(-1)),5);
In6=vpa((exp(-1)),6);
In7=vpa((exp(-1)),7);
fprintf('%.5f%.6f%.7f\n',eval(In5),eval(In6),eval(In7));
forn=2:
10
In5=vpa((1-n*In5),5);
In6=vpa((1-n*In6),6);
In7=vpa((1-n*In7),7);
fprintf('%.5f%.6f%.7f\n',eval(In5),eval(In6),eval(In7));
end
五位六位七位
0.367880.3678790.3678794
0.264240.2642410.2642411
0.207280.2072770.2072766
0.170890.1708930.1708934
0.145530.1455330.1455329
0.126800.1268020.1268024
0.112380.1123840.1123835
0.100930.1009320.1009320
0.091610.0916120.0916123
0.083880.0838770.0838771
symsnEn;
En5=vpa(0,5);
En6=vpa(0,6);
En7=vpa(0,7);
fprintf('%.5f%.6f%.7f\n',eval(En5),eval(En6),eval(En7));
forn=10:
-1:
2
En5=vpa(((1-En5)/n),5);
En6=vpa(((1-En6)/n),6);
En7=vpa(((1-En7)/n),7);
fprintf('%.5f%.6f%.7f\n',eval(En5),eval(En6),eval(En7));
end
五位六位七位
0.000000.0000000.0000000
0.100000.1000000.1000000
0.100000.1000000.1000000
0.112500.1125000.1125000
0.126790.1267860.1267857
0.145540.1455360.1455357
0.170890.1708930.1708929
0.207280.2072770.2072768
0.264240.2642410.2642411
0.367880.3678790.3678795
(1)由所得数据可以了解,两种算法随着n的增大,误差越来越大,而第二种算法随着n的减小,数据越来越精确,且三种有效数字结果大致相同,所以第二种更精确。
(2)算法一E1的误差e1,由于En=1-n*En-1,n=2,3,4…..,通过推算可得误差
=
,n越大,其误差
就越大,那么最后
算出的结果是e1的n!
倍;
算法二
的误差
,由于
,推倒可得误差
=
,
比
缩小了(N-n)!
倍,则n减少,所以
就越小。
所以算法一越往后算,一步步的误差会使误差越大,而算法二由于是从后面递推回来,其误差会被缩小。
所以算法二更优。
(3)算法一中,n增大,误差增大,算法二中,n减小,误差减小。
所以当某一步发生误差
后随着n变大,算法一的误差越来越大,而算法二由于是往回推,所以误差变小。
(4)通过以上可一推出,算法一随着n变大,容易使误差越来越大,而算法二会使误差变小,所以算法二更加稳定。
2、求解非线性方程
3、插值
实验目的:
掌握Lagrange插值法和Newton插值法
问题提出:
,已知
的函数值表如下
x
0
0.1
0.2
0.3
0.4
0.5000
0.5398
0.5793
0.6179
0.7554
用插值法求
和
的近似值。
实验容:
(1)分别用Lagrange插值法和Newton插值法编程求解;
(2)求出插值多项式系数,对比计算结果。
拉格朗日差值
functionyh=lagrange(x,y,xh)
tic
n=length(x);
m=length(xh);
x=x(:
);
y=y(:
);
xh=xh(:
);
yh=zeros(m,1);
c1=ones(1,n-1);
c2=ones(m,1);
fori=1:
n,
xp=x([1:
i-1i+1:
n]);
yh=yh+y(i)*prod((xh*c1-c2*xp')./(c2*(x(i)*c1-xp')),2);
end
toc
x=[00.10.20.30.4];
y=[0.50000.53980.57930.61790.7554];
xh=[0.130.36];
lagrange(x,y,xh)
时间已过0.026678秒。
ans=
0.5537
0.6780
多项式为
((5*xh-1)*((5*xh)/2-1)*(10*xh-1)*((10*xh)/3-1))/2+(5793*xh*(5*xh-2)*(10*xh-1)*(10*xh-3))/2000-(6179*xh*(5*xh-1/2)*(10*xh-2)*(10*xh-4))/3000+(3777*xh*(5*xh-1)*(10*xh-3)*((10*xh)/3-1/3))/2000-(2699*xh*(5*xh-3/2)*(10*xh-2)*((10*xh)/3-4/3))/500
作图
x0=[00.10.20.30.4];
y0=[0.50000.53980.57930.61790.7554];
xh=[0:
0.01:
0.4];
y=lagrange(x0,y0,xh);
时间已过0.000158秒。
>>plot(xh,y);holdon;
牛顿差值
functionf=Newton(x,y,x0,x1)
tic
symst;
if(length(x)==length(y))
n=length(x);
c(1:
n)=0.0;
else
disp('x和y的维数不相等!
');
return;
end
f=y
(1);
y1=0;
l=1;
for(i=1:
n-1)
for(j=i+1:
n)
y1(j)=(y(j)-y(i))/(x(j)-x(i));
end
c(i)=y1(i+1);
l=l*(t-x(i));
f=f+c(i)*l;
y=y1;
end
f=simplify(f);
g=subs(f,'t',x0)
g1=subs(f,'t',x1)
A=zeros(n,n-1);
A=[y',A];
forj=2:
n
fori=j:
n
A(i,j)=(A(i,j-1)-A(i-1,j-1))/(x(i)-x(i+1-j));
end
end
disp('差商表为');
disp(A);
toc
x0=[00.10.20.30.4];
y0=[0.50000.53980.57930.61790.7554];
Newton(x0,y0,0.13,0.36)
g=
9228481/0000000
g1=
91741/00000
差商表为
1.0e+04*
00000
0.00000.0004000
-0.0000-0.0004-0.004100
-0.0000-0.00010.00160.01900
0.00420.04190.21010.69481.6896
时间已过0.368255秒。
ans=
(251*t^4)/6-(97*t^3)/+(501*t^2)/00+(7*t)/00+1/2
作图
>>t=[0:
0.01:
0.4];
>>f=[];
>>f=(251.*t.^4)./6.-(97.*t.^3)./.+(501.*t.^2)./00.+(7.*t)./00.+1/2;
>>plot(t,f);holdon;
两种方法作图曲线为
拉格朗日插值法时间复杂度为O(n),牛顿插值法时间复杂度为O(n^2),实际上拉格朗日算法计算时间为0.000158,牛顿算法计算时间0.368255。
说明拉格朗日算法时间复杂度低。
4、数值积分
实验目的:
掌握Romberg积分法
实验容:
用Romberg积分法计算下列定积分:
function[R,y]=Romberg(a,b,n)
k=1;
while1
T=zeros(k+1,k+1);
T(1,1)=1/2*(b-a)*(f(a)+f(b));
fori=1:
k
h=(b-a)/2^i;
s=0;
forj=1:
2^(i-1)
s=s+f(a+(2*j-1)*h);
end
T(i+1,1)=T(i,1)/2+h*s;
end
forj=1:
k
c=1/(4^j-1);
form=j:
k
T(m+1,j+1)=T(m+1,j)+c*(T(m+1,j)-T(m,j));
end
end
ifabs(T(k,k)-T(k+1,k+1)) break; end k=k+1; end y=T(k,k); R=T; 输入积分公式 (1) functiony=f(x) y=1/(1+x^4); [R,y]=Romberg(0,1,10) R= 0.75000 0.84560.8775 y= 0.7500 (2) functiony=f(x) y=cos(sin(x)); [R,y]=Romberg(0,pi/4,10) R= 0.69120 0.70990.7161 y= 0.6912 (3) functiony=f(x) y=exp(x)/x; [R,y]=Romberg(1,2,10) R= 3.20640 3.09713.0607 y= 3.2064 (4) functiony=f(x) y=x.^x*log(3+x); [R,y]=Romberg(2,3,10) R= 27.40760 22.127120.3669 y= 27.4076 龙贝格计算量较小 二、提高技能训练 1、 (2)实验目的与要求 (1)了解Kepler方程; (2)能够用牛顿迭代方法计算Kepler方程; (3)通过调节轨道偏心率e查看运算收敛情况。 实验容及数据来源 编写解kepler方程的方法,给定偏心率为0.01、平近点角为32度时计算出偏近点角 。 clear; x0=32; symsE; f=inline('E-0.01*sin(E)-32'); g=inline('1-0.01*cos(E)'); x1=x0-f(x0)/g(x0); y=abs(x1-x0); x0=x1;i=i+1; whiley>10^(-5) x1=x0-f(x0)/g(x0); y=abs(x1-x0); x0=x1; i=i+1; end E=vpa(x1) 结果: E=32.83 (3)通过调节轨道偏心率e查看运算收敛情况; e0.00111.82.210 E32.0005532.9999133.5355930.8838131.35097 导0.999161.013181.93911-0.89582-8.97892 可以看出,e在0~1之间,收敛;在1~1.8之间,发散;在1.8~2.2左右,收敛;2.2之后,发散. 2、 (日照时间分布)已知某地区在不同月份的平均日照时间的观测数据如下表所示 月份 1 2 3 4 5 6 7 8 9 10 11 12 日照/(h/月) 80.9 67.2 67.1 50.5 32.0 33.6 35.6 46.8 52.3 62.0 64.1 71.2 试分析日照时间的变化规律。 x=1: 1: 12; y=[80.967.267.150.532.033.636.646.852.362.064.171.2]; x1=1: 0.01: 12; y1=interp1(x,y,x1); y2=interp1(x,y,x1,'spline'); plot(x,y,'go',x1,y1,'r-',x1,y2,'b') 其中红线是分段插值得到的,因此为折线,绿色是三次样条曲线,‘o‘处是插值节点。 从图中可以看出,日照时间从一月开始减少至5月,之后又开始增加至12月。 (水道测量问题)水深数据如下表所示,设船的吃水深度为1.5m,问在 的区域,哪些地方船只避免进入,画船只避免进入的区域。 x 129.0 140.0 108.5 88.0 185.5 195.5 105.5 y 7.5 141.5 28.0 147.0 22.5 137.5 85.5 z 1.22 2.44 1.83 2.44 1.83 2.44 2.44 x 157.5 107.5 77.0 81.0 162.0 162.0 117.5 y -6.5 -81.5 3.0 56.5 -66.5 84 -38.5 z 2.74 2.74 2.44 2.44 2.74 1.22 2.74 x=[129.0140.0108.588.0185.5195.5105.5157.5107.577.081.0162.0162.0117.5]; y=[7.5141.528.0147.022.5137.585.5-6.5-81.53.056.5-66.584-38.5]; z=-[1.222.441.832.441.832.442.442.742.742.442.442.741.222.74]; x1=linspace(min(x),max(x),40); y1=linspace(min(y),max(y),40); [Xi,Yi]=meshgrid(x1,y1); [Xi,Yi,Zi]=griddata(x,y,z,Xi,Yi,'cubic'); mesh(Xi,Yi,Zi); xlabel('X'); ylabel('Y'); zlabel('Z'); figure; [c,h]=contour(Xi,Yi,Zi,[-1.2: -0.1: -2.8]); clabel(c,h); xlabel('X'); ylabel('Y'); 三、本课程设计的心得体会(500字左右) 通过这次计算方法的课程设计,我比过去得到了更多的收获。 不仅仅是知识的学习,更明白了巩固旧知识的重要性,明白了什么是温故而知新。 也懂得了人际关系的重要性,保持好的人际关系,遇到问题,别人才乐于伸出援手。 在试验的进行中,我们使用matlab对问题进行求解。 刚开始的时候,使用matlab解决此类问题有一些难度,包括如何操作,具体算法编排都不是很了解,但是通过互相探讨,查找资料,渐渐的掌握了一些技巧,能够解决一些基础问题。 同时我们也明白了maatlab的强大,它对数学数据和图形的处理能力是c比不上的。 在本次实验中,我的计算方法知识得到了很好的巩固与更新。 在不断的出错,不断的查找,不断的求索与修改之中,过去知识中的漏洞与迷茫渐渐不复存在,也明白了自己掌握的东西还是太少,太浅薄,需要努力寻找新的东西。 而寻找东西的能力恰恰也在本次实验中得到了培养。 由于一次次的求解失误,只能努力的寻找解决方法,书本,网络,同学,这些都在本次试验中运用到。 课程设计已经结束,但是我收获了软件操作知识,数学知识,以及查找资料等能力,这些都是我对专业,对课程有了一个更加深入的了解。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 计算方法 综合 实践
![提示](https://static.bdocx.com/images/bang_tan.gif)