数值分析课程设计报告文档格式.docx
- 文档编号:17276046
- 上传时间:2022-11-30
- 格式:DOCX
- 页数:33
- 大小:393.39KB
数值分析课程设计报告文档格式.docx
《数值分析课程设计报告文档格式.docx》由会员分享,可在线阅读,更多相关《数值分析课程设计报告文档格式.docx(33页珍藏版)》请在冰豆网上搜索。
I0=int(1/(5+x),x,0,1);
I0=vpa(I0);
fori=1:
20
I0=-5*I0+1/i
End
经过20次迭代,最终结果为I0=0.007997523028232163831452112940177
②从
较粗糙的估计值出发,用递推公式:
计算从
先求出
的值近似精确值,再用递推公式
I0=int(0.2*x^30,x,0,1);
fori=30:
-1:
2
I0=-0.2*I0+0.2/i
经过30次迭代,最终结果为I0=.0883********
③分析所得结果的可靠性以及出现这种现象的原因。
方法2比较可靠,因为方法1每次都将误差放大5倍,而方法2则将误差缩小到原来的0.2倍。
实验二
2.1用高斯消元法的消元过程作矩阵分解。
设
消元过程可将矩阵A化为上三角矩阵U,试求出消元过程所用的乘数
、
并以如下格式构造下三角矩阵L和上三角矩阵U
验证:
矩阵A可以分解为L和U的乘积,即A=LU。
算法分析:
由于matlab中含有将A分解成一个下三角矩阵和一个上三角矩阵,所以直接调用。
程序及其运行结果:
A=[2023;
181;
2-315]
[L,U]=lu(A)
A=
2023
181
2-315
L=
1.000000
0.05001.00000
0.1000-0.40511.0000
U=
20.00002.00003.0000
07.90000.8500
0015.0443
矩阵A可以分解为L和U的乘积
2.3验证希尔伯特矩阵的病态性:
对于三阶矩阵
取右端向量
,验证:
(1)向量
是方程组
的准确解;
(2)取右端向量b的三位有效数字得
,求方程组的准确解
,并与X的数据
作比较。
说明矩阵的病态性。
(1)直接验证即可
A=[11/21/3;
1/21/31/4;
1/31/41/5];
B=[1;
1;
1];
b=sym(A*B)
结果与上述结果相同
(2)先求出解
,与数据
b=[1.83;
1.08;
0.783];
B=A\b
结果为
B=
.0799********
0.540000000000004
1.439999999999997
结果误差变化较大,说明数值的微小变化导致结果的巨大误差,矩阵A为病态的
实验三
3.1用泰勒级数的有限项逼近正弦函数
用计算机绘出上面四个函数的图形。
程序:
x=0:
pi/100:
pi;
m=0:
pi/2;
y0=sin(x);
y1=m;
y2=m-m.^3/6;
y3=m-m.^3/6+t.^5/120;
plot(x,y0,m,y1,m,y2,m,y3)
3.2绘制飞机的降落曲线
一架飞机飞临北京国际机场上空时,其水平速度为540km/h,飞行高度为1000m。
飞机从距机场指挥塔的横向距离12000m处开始降落。
根据经验,一架水平飞行的飞机其降落曲线是一条三次曲线。
建立直角坐标系,设飞机着陆点为原点O,降落的飞机为动点
,则
表示飞机距指挥塔的距离,
表示飞机的飞行高度,降落曲线为
该函数满足条件:
(1)试利用
满足的条件确定三次多项式中的四个系数;
(2)用所求出的三次多项式函数绘制出飞机降落曲线。
(1)程序:
A=[1000;
11200012000^212000^3;
0100;
01240003*(12000^2)];
b=[0,1000,0,0]'
;
inv(A)*b
ans=
1.0e-004*
0
0.2083
-0.0000
(2)程序:
50:
12000;
y=0.207e-004*(x.^2)+-0.00001158e-004*(x.^3);
plot(x,y)
实验四
4.1曾任英特尔公司董事长的摩尔先生早在1965年时,就观察到一件很有趣的现象:
集成电路上可容纳的零件数量,每隔一年半左右就会增长一倍,性能也提升一倍。
因而发表论文,提出了大大有名的摩尔定律(Moore’sLaw),并预测未来这种增长仍会延续下去。
下面数据中,第二行数据为晶片上晶体数目在不同年代与1959年时数目比较的倍数。
这些数据是推出摩尔定律的依据:
年代
1959
1962
1963
1964
1965
增加倍数
1
3
4
5
6
试从表中数据出发,推导线性拟合的函数表达式。
y=[13456];
x2=[19591962196319641965];
P=polyfit(x2,y,1);
disp(P);
运行结果为:
0.0008-1.6255
结果分析:
由运行结果可知线性拟合的函数表达式为y=0.0008x-1.6255,由于是近似值,所以不免有误差。
实验五
5.1用几种不同的方法求积分
的值。
(1)牛顿-莱布尼茨公式;
(2)梯形公式;
(3)辛卜生公式;
(4)复合梯形公式。
(1)牛顿-莱布尼茨公式
输入:
y=4/(1+x^2);
int_y=int(y,x,0,1)
输出:
int_y=
pi
(2)梯形公式
functiony=dada(x)
y=4./(1+x.^2);
以'
dada.m'
保存
T=(1-0)/2*(dada
(1)+dada(0))
结果:
T=3
(3)辛卜生公式
S=quad('
dada'
0,1)
S=
3.141592682924567
(4)复合梯形公式
输入d=1/100;
d:
f=dada(x);
FT=trapz(f)*d
FT=
3.141575986923129
5.2设X为标准正态随机变量,即X~N(0,1)。
现分别取
,试设计算法计算30个不同的概率值;
,并将计算结果与概率论教科书中的标准正态分布函数表作比较。
(提示:
m=0.1:
0.1:
3;
P=zeros(30,1);
fori=1:
30
P(i)=1/sqrt(2*pi)*int(sym('
exp(-x^2/2)'
),m(i),4);
end
P
),-4,m(i));
P=
0.4601
0.4207
0.3821
0.3445
0.3085
0.2742
0.2419
0.2118
0.1840
0.1586
0.1356
0.1150
0.0968
0.0807
0.0668
0.0548
0.0445
0.0359
0.0287
0.0227
0.0178
0.0139
0.0107
0.0082
0.0062
0.0046
0.0034
0.0025
0.0018
0.0013
0.5398
0.5792
0.6179
0.6554
0.6914
0.7257
0.7580
0.7881
0.8159
0.8413
0.8643
0.8849
0.9032
0.9192
0.9332
0.9452
0.9554
0.9640
0.9713
0.9772
0.9821
0.9861
0.9892
0.9918
0.9938
0.9953
0.9965
0.9974
0.9981
0.9986
5.5求空间曲线L:
弧长公式为
运行程序为:
symst;
a=diff(cos(t));
b=diff(sin(t));
c=diff(2-cos(t)-sin(t));
y=int((a^2+b^2+c^2)^0.5,t,0,2*pi);
vpa(y)
ans=8.737752570984804741619351957708
实验六
6.1用欧拉公式和四阶龙格-库塔法分别求解下列初值问题;
欧拉公式:
function[t,x]=Euler(fun,t0,tt,x0,N)
h=(tt-t0)/N;
t=t0+[0:
N]'
*h;
x(1,:
)=x0'
fork=1:
N
f=feval(fun,t(k),x(k,:
));
f=f'
x(k+1,:
)=x(k,:
)+h*f;
Euler.m'
functionf=Euler_fun(t,x)
f=0.9*x./(1+2*t);
Euler_fun.m'
functionmain_Euler
[t,x]=Euler('
Euler_fun'
0,1,1,20);
fh=dsolve('
Dx=0.9*x/(1+2*t)'
'
x(0)=1'
21
ft(k)=t(k);
fx(k)=subs(fh,ft(k));
[t,x]
main_Euler.m'
保存
main_Euler
01.0000
0.05001.0450
0.10001.0878
0.15001.1285
0.20001.1676
0.25001.2051
0.30001.2413
0.35001.2762
0.40001.3100
0.45001.3427
0.50001.3745
0.55001.4055
0.60001.4356
0.65001.4649
0.70001.4936
0.75001.5216
0.80001.5490
0.85001.5758
0.90001.6021
0.95001.6278
1.00001.6531
四阶龙格-库塔法
functionR=rk4(f,a,b,ya,N)
h=(b-a)/N;
T=zeros(1,N+1);
Y=zeros(1,N+1);
T=a:
h:
b;
Y
(1)=ya;
forj=1:
k1=h*feval(f,T(j),Y(j));
k2=h*feval(f,T(j)+h/2,Y(j)+k1/2);
k3=h*feval(f,T(j)+h/2,Y(j)+k2/2);
k4=h*feval(f,Y(j)+h,Y(j)+k3);
Y(j+1)=Y(j)+(k1+2*k2+2*k3+k4)/6;
ans=[T'
Y'
]
rk4.m'
functionz=f(x,y)
z=0.9*y./(1+2*x);
f.m'
rk4('
f'
0,1,1,20)
0.05001.0392
0.10001.0765
0.15001.1121
0.20001.1463
0.25001.1791
0.30001.2108
0.35001.2415
0.40001.2712
0.45001.3000
0.50001.3281
0.55001.3554
0.60001.3821
0.65001.4082
0.70001.4337
0.75001.4586
0.80001.4831
0.85001.5071
0.90001.5306
0.95001.5538
1.00001.5765
欧拉公式
function[t,x]=Eulers(fun,t0,tt,x0,N)
Eulers.m'
functionf=Euler_funs(t,x)
f=-t*x/(1+t^2);
Euler_funs.m'
functionmains_Euler
Euler_funs'
0,1,2,20);
mains_Euler.m'
mains_Euler
02.0000
0.05002.0000
0.10001.9950
0.15001.9851
0.20001.9706
0.25001.9516
0.30001.9287
0.35001.9021
0.40001.8725
0.45001.8402
0.50001.8058
0.55001.7696
0.60001.7323
0.65001.6941
0.70001.6554
0.75001.6165
0.80001.5777
0.85001.5392
0.90001.5012
0.95001.4639
1.00001.4274
四阶龙格-库塔法
functionR=rk4s(f,a,b,ya,N)
rk4s.m'
functionz=fs(x,y)
z=-x*y/(1+x^2);
fs.m'
rk4s('
fs'
0,1,2,20)
0.05001.9918
0.10001.9795
0.15001.9632
0.20001.9433
0.25001.9200
0.30001.8936
0.35001.8645
0.40001.8331
0.45001.7998
0.50001.7650
0.55001.7291
0.60001.6923
0.65001.6551
0.70001.6176
0.75001.5801
0.80001.5429
0.85001.5060
0.90001.4697
0.95001.4340
1.00001.3990
6.2用求解常微分方程初值问题的方法计算积分上限函数
的值,实际上是将上面表达式两端求导化为常微分方程形式,并用初值条件
试用MATLAB中的指令ode23解决定积分
的计算问题。
odefun=@(x,y)2*pi*exp(x)/x;
%定义函数
tspan=[45];
%求解区间
y0=0;
%初值
[x,y]=ode23(odefun,tspan,y0)
x=
4.0000
4.0001
4.0007
4.0036
4.0182
4.0911
4.1911
4.2911
4.3911
4.4911
4.5911
4.6911
4.7911
4.8911
5.0000
y=
0.0001
0.0005
0.0125
0.0625
0.3129
1.5732
8.0862
17.6281
27.9249
39.0424
51.0526
64.0337
78.0708
93.2571
109.6940
129.1469
当x=5时,y=129.1469
实验七
7.1用二分法求下列方程在指定区间[a,b]上的实根近似值:
(要求误差不超过0.01)
(1)x-sinx-1=0,[a,b]=[1,3];
(2)xsinx=1,[a,b]=[0,1.5].
解;
(1)运行程序:
a=1;
b=3;
k=0;
y1=a-sin(a)-1;
while(b-a)>
0.01
k=k+1;
x0=(a+b)/2;
x(k)=x0;
y0=x0-sin(x0)-1;
y(k)=y0;
ify0*y1<
b=x0;
else
a=x0;
y1=y0;
end
disp([x'
y'
])
[2,408488074749405/4503599627370496]
[3/2,-4481036472577419/9007199254740992]
[7/4,-1053779023151395/4503599627370496]
[15/8,-712341393175443/9007199254740992]
[31/16,8975041611503/2251799813685248]
[61/32,-85593294866561/2251799813685248]
[123/64,-154268929443443/9007199254740992]
[247/128,-7430218095237/1125899906842624]
[495/256,-11835035444387/9007199254740992]
[991/512,93879030099/70368744177664]
[1981/1024,10840787111/1125899906842624]
[3961/2048,-5875158237047/9007199254740992]
[7923/4096,-723616715653/2251799813685248]
[15847/8192,-701966501541/4503599627370496]
[31695/16384,-164654758197/2251799813685248]
(2)运行程序:
a=0.1;
b=1.5;
y1=a*sin(a)-1;
while(b-a)>
k=k+1;
x0=(a+b)/2;
y0=x0*sin(x0)-1;
ify0*y1<
b=x0;
else
a=x0;
y1=y0;
disp([x'
[4/5,-3838103856873717/9007199254740992]
[23/20,55933053762297/1125899906842624]
[39/40,-1738305320249353/9007199254740992]
[17/16,-323478390340705/4503599627370496]
[177/160,-98943562990273/9007199254740992]
[361/320,21826383383843/1125899906842624]
[143/128,18951300402293/4503599627370496]
[1423/1280,-7626403754673/2251799813685248]
[1981/1024,10840787111/1125899906842
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数值 分析 课程设计 报告