常微分实习报告常微分方程数值求解问题.docx
《常微分实习报告常微分方程数值求解问题.docx》由会员分享，可在线阅读，更多相关《常微分实习报告常微分方程数值求解问题.docx（31页珍藏版）》请在冰豆网上搜索。
常微分实习报告常微分方程数值求解问题
课程实习报告
课程名称：
常微分方程课程实习
实习题目：
常微分方程数值求解问题的实习
姓名：
系：
专业：
年级：
学号：
指导教师：
职称：
年月日
课程实习报告结果评定
评语：
成绩：
指导教师签字：
评定日期：
5.1用不同格式对同一个初值问题的数值求解及其分析……………………..1
5.1.5问题讨论与分析……………………………………………………………………15
5.2一个算法不同不长求解同一个初值问题及其分析…………………………………..18
常微分方程课程实习
1.实习的目的和任务
目的：
通过课程实习能够应用MATLAB软来计算微分方程（组）的数值解；了解常微分方程数值解。
任务：
通过具体的问题，利用MATLAB软件来计算问题的结果，分析问题的结论。
2.实习要求
能够从案例的自然语言描述中，抽象出其中的数学模型；能够熟练应用所学的数值解计算方法；能够熟练使用MATLAB软件；对常微分方程数值解有所认识，包括对不同算法有所认识和对步长有所认识。
3.实习地点
数学实验室、学生宿舍
4.主要仪器设备
计算机、
MicrosoftWindows7
Matlab7.0
5.实习内容
5.1用欧拉方法，改进欧拉方法，4阶龙格—库塔方法分别求下面微分方程的初值dy/dx=y*cos（x+2）y（-2）=1x∈[-2，0]
5.1.1求精确解
①变量分离方程情形:
形如
的方程,这里
分别是
的连续函数.如果
我们可将方程改写成
这样,变量就”分离”开来了,两边同时积分即可:
为任意常数.
②常数变易法:
一阶线性微分方程
其中
在考虑区
间上是的连续函数.可先解出方程
的解,这是属于变量分
离方程情形,可解得:
这里
是任意常数.然后将变
易为
的待定函数
令
将其代入原方程可得:
所以可解得
这里
是任意常数.将
代入
可得原方程的通解:
为任意常数.
③恰当微分方程情形:
形如
的一阶微分方程,这里
假设
在某矩形域内是的连续函数,且具有连续的一阶偏导数.
若
则为恰当微分方程.判断为恰当微分方程后,则可用如下解法：
设
是原方程的解，则
，所以设
，
则
，所以
，由此
由此可解得
，所以原方程的通解为
为任意常数。
首先可以求得其精确解为：
y=exp（sin（x+2））
>>x=-2:
0.1:
2;
>>y=exp（sin（x+2））
>>plot（x,y,'r.-'）;
>>Data=[x',y']
y=
Columns1through4
1.000000000000001.104986830331691.219778556000621.34382524373165
Columns5through8
1.476121946445731.615146296442081.758818845766991.90449653438673
Columns9through12
2.049008650164272.188********4612.319776824715852.43807150515633
Columns13through16
2.539682532380782.621005926286702.679016447572712.71148101768216
Columns17through20
2.717123008431282.695718599203822.648113847390782.57616043684702
Columns21through24
2.482577728015002.370757126170312.244530580577552.10792744704554
Columns25through28
1.964942888556771.819336991081061.674477827371161.53323499677732
Columns29through32
1.397923819945001.270295218811561.151********4531.04245724559826
Columns33through36
0.943296953050410.854066948581020.774497302497390.70413637458179
Columns37through40
0.642415207750370.588701425856340.542342319593710.50269776253737
Column41
0.46916418587400
Data=
-2.000000000000001.00000000000000
-1.900000000000001.10498683033169
-1.800000000000001.21977855600062
-1.700000000000001.34382524373165
-1.600000000000001.47612194644573
-1.500000000000001.61514629644208
-1.400000000000001.75881884576699
-1.300000000000001.90449653438673
-1.200000000000002.04900865016427
-1.100000000000002.188********461
-1.000000000000002.31977682471585
-0.900000000000002.43807150515633
-0.800000000000002.53968253238078
-0.700000000000002.62100592628670
-0.600000000000002.67901644757271
-0.500000000000002.71148101768216
-0.400000000000002.71712300843128
-0.300000000000002.69571859920382
-0.200000000000002.64811384739078
-0.100000000000002.57616043684702
02.48257772801500
0.100000000000002.37075712617031
0.200000000000002.24453058057755
0.300000000000002.10792744704554
0.400000000000001.96494288855677
0.500000000000001.81933699108106
0.600000000000001.67447782737116
0.700000000000001.53323499677732
0.800000000000001.39792381994500
0.900000000000001.27029521881156
1.000000000000001.151********453
1.100000000000001.04245724559826
1.200000000000000.94329695305041
1.300000000000000.85406694858102
1.400000000000000.77449730249739
1.500000000000000.70413637458179
1.600000000000000.64241520775037
1.700000000000000.58870142585634
1.800000000000000.54234231959371
1.900000000000000.50269776253737
2.000000000000000.46916418587400
5.1.2用欧拉法求解
设常微分方程的初始问题
有唯一解。
则由欧拉法求初值问题
（1），
（2）的数值解的差分方程为：
程序如下：
建立函数文f1.m
function[x,y]=f1（fun,x_span,y0,h）
x=x_span
（1）:
h:
x_span
（2）;
y
（1）=y0;
forn=1:
length（x）-1
y（n+1）=y（n）+h*feval（fun,x（n）,y（n））;
end
x=x';y=y';
在MATLAB输入以下程序：
>>clearall
>>fun=inline（'y*cos（x+2）'）;
>>[x,y1]=f1（fun,[-2,2],1,0.1）;
>>[x,y1]
>>plot（x,y1,'g*-'）
结果及其图象：
ans=
-2.000000000000001.00000000000000
-1.900000000000001.10000000000000
-1.800000000000001.20945045818058
-1.700000000000001.32798465534234
-1.600000000000001.45485187516708
-1.500000000000001.58885260659392
-1.400000000000001.72828754069001
-1.300000000000001.87092926670362
-1.200000000000002.01402582996363
-1.100000000000002.15434436081705
-1.000000000000002.28826055379421
-0.900000000000002.41189579915842
-0.800000000000002.52129845713651
-0.700000000000002.61265966186586
-0.600000000000002.68254800178024
-0.500000000000002.72814250373577
-0.400000000000002.74744062038227
-0.300000000000002.73941822501564
-0.200000000000002.70412232942903
-0.100000000000002.64268410367377
02.55724888375039
0.100000000000002.45082978042675
0.200000000000002.32710059365817
0.300000000000002.19015046372483
0.400000000000002.04422599002736
0.500000000000001.89348605020813
0.600000000000001.74179062418299
0.700000000000001.59253854452440
0.800000000000001.44856157120511
0.900000000000001.31207486378276
1.000000000000001.184********502
1.100000000000001.06739566199422
1.200000000000000.96074840947946
1.300000000000000.86483739767581
1.400000000000000.77943645422926
1.500000000000000.70408067871132
1.600000000000000.63814657271418
1.700000000000000.58092024172055
1.800000000000000.53165239417811
1.900000000000000.48960040640242
2.000000000000000.45405873128672
5.1.3用改进欧拉法求解：
计算公式为：
即先用欧拉法得
，进而由（3）的第一式得初始近似值
，然后再用（3）的第二式进行迭代，反复改进这个近似值，直到
（
为所允许的误差）为止，并把
取作为
的近似值
，这个方法就称为改进欧拉法。
通常称（3）为预报校正公式，其中第一式称为预报公式，第二式称为校正公式。
这个公式还可以写为：
（下文改进欧拉法计算就以下面为准）
程序如下：
建立函数文件f2.m
function[x,y]=f2（fun,x_span,y0,h）
x=x_span
（1）:
h:
x_span
（2）;
y
（1）=y0;
forn=1:
length（x）-1
k1=feval（fun,x（n）,y（n））;
y（n+1）=y（n）+h*k1;
k2=feval（fun,x（n+1）,y（n+1））;
y（n+1）=y（n）+h*（k1+k2）/2;
end
x=x';y=y';
在MATLAB输入以下程序：
>>clearall
>>fun=inline（'y*cos（x+2）'）;
>>[x,y2]=f2（fun,[-2,2],1,0.1）;
>>[x,y2]
>>plot（x,y2,'b+-'）
结果及其图象：
ans=
-2.000000000000001.00000000000000
-1.900000000000001.10472522909029
-1.800000000000001.21920722936399
-1.700000000000001.34289777810139
-1.600000000000001.47479651303944
-1.500000000000001.61338861747877
-1.400000000000001.75660494566073
-1.300000000000001.90181495275894
-1.200000000000002.04586183721764
-1.100000000000002.185********243
-1.000000000000002.31576355572445
-0.900000000000002.43368296895904
-0.800000000000002.53497167173468
-0.700000000000002.61603367901013
-0.600000000000002.67384966784500
-0.500000000000002.70619076790723
-0.400000000000002.71178326404468
-0.300000000000002.69040521940878
-0.200000000000002.64290353044130
-0.100000000000002.57112934628408
02.47779956085378
0.100000000000002.36630057150132
0.200000000000002.24045633274254
0.300000000000002.10428512493058
0.400000000000001.96176831551518
0.500000000000001.81665028035196
0.600000000000001.67228260188095
0.700000000000001.53151888526290
0.800000000000001.39666016389536
0.900000000000001.26944574589307
1.000000000000001.151********251
1.100000000000001.04229147409852
1.200000000000000.94339433622144
1.300000000000000.85437588457241
1.400000000000000.77496982261477
1.500000000000000.70472971912662
1.600000000000000.64309273433755
1.700000000000000.58943293626443
1.800000000000000.54310392687123
1.900000000000000.50347143064050
2.000000000000000.46993706594350
5.1.4用4阶龙格—库塔求解
标准的四阶龙格-库塔公式（亦称为经典的四阶龙格-库塔公式）：
公式的截断误差阶为
。
程序如下：
建立函数文件f3.m
function[x,y]=f3（fun,x_span,y0,h）
x=x_span
（1）:
h:
x_span
（2）;
y
（1）=y0;
forn=1:
length（x）-1
k1=feval（fun,x（n）,y（n））;
k2=feval（fun,x（n）+h/2,y（n）+h/2*k1）;
k3=feval（fun,x（n）+h/2,y（n）+h/2*k2）;
k4=feval（fun,x（n+1）,y（n）+h*k3）;
y（n+1）=y（n）+h*（k1+2*k2+2*k3+k4）/6;
end
x=x';y=y';
在MATLAB输入以下程序：
>>clearall;
>>fun=inline（'y*cos（x+2）'）;
>>[x,y3]=f3（fun,[-2,2],1,0.1）;
>>[x,y3]
>>plot（x,y3,'y*-'）
结果及其图象：
ans=
-2.000000000000001.00000000000000
-1.900000000000001.10498674569681
-1.800000000000001.21977837351828
-1.700000000000001.34382495362539
-1.600000000000001.47612154334742
-1.500000000000001.61514577980337
-1.400000000000001.75881821961767
-1.300000000000001.90449580650198
-1.200000000000002.04900783085566
-1.100000000000002.188********744
-1.000000000000002.31977585752433
-0.900000000000002.43807048140885
-0.800000000000002.53968146311619
-0.700000000000002.62100482231058
-0.600000000000002.67901531971109
-0.500000000000002.71147987684658
-0.400000000000002.71712186540483
-0.300000000000002.69571746425094
-0.200000000000002.64811272993642
-0.100000000000002.57615934548225
02.48257667095154
0.100000000000002.37075611204203
0.200000000000002.24452961927934
0.300000000000002.10792655020966
0.400000000000001.96494206934318
0.500000000000001.81933626317413
0.600000000000001.67447720334541
0.700000000000001.53323448621294
0.800000000000001.39792342776452
0.900000000000001.27029494424773
1.000000000000001.151********231
1.100000000000001.04245718123925
1.200000000000000.94329697235934
1.300000000000000.85406703400224
1.400000000000000.77449743625256
1.500000000000000.70413654020399
1.600000000000000.64241539118308
1.700000000000000.58870161605161
1.800000000000000.54234250864301
1.900000000000000.50269794543536
2.000000000000000.46916436004503
5.1.5问题讨论与分析
由以上数值分析结果绘制表格：
精确解
欧拉方法
改进的欧拉方法
四阶龙格-库塔方法
xi
yi
yi
误差
yi
误差
yi
误差
-2.00
1.000000
1.000000
0.000000
1.000000
0.000000
1.000000
0.000000
-1.90
1.104987
1.100000
0.004987
1.104725
0.000262
1.104987
0.000000
-1.80
1.219779
1.209450
0.010328
1.219207
0.000571
1.219778
0.000000
-1.70
1.343825
1.327985
0.015841
1.342898
0.000927
1.34
升级会员 