1、计算机、Microsoft Windows 7Matlab 7.05. 实习内容5.1 用欧拉方法,改进欧拉方法,4阶龙格库塔方法分别求下面微分方程的初值dy/dx=y*cos(x+2) y(-2)=1 x-2,05.1.1求精确解变量分离方程情形:形如的方程,这里分别是的连续函数.如果,我们可将方程改写成,这样,变量就”分离”开来了,两边同时积分即可:为任意常数.常数变易法:一阶线性微分方程,其中在考虑区间上是的连续函数.可先解出方程的解,这是属于变量分离方程情形,可解得:,这里是任意常数.然后将变易为的待定函数,令,将其代入原方程可得:所以可解得,这里是任意常数.将代入可得原方程的通解:恰
2、当微分方程情形:的一阶微分方程,这里假设在某矩形域内是的连续函数,且具有连续的一阶偏导数.若,则为恰当微分方程.判断为恰当微分方程后,则可用如下解法:设是原方程的解,则,所以设,则,所以,由此,由此可解得,所以原方程的通解为为任意常数。首先可以求得其精确解为:y=exp(sin(x+2) x=-2:0.1:2; y= exp(sin(x+2) plot(x,y,r.-); Data=x,yy = Columns 1 through 4 1.00000000000000 1.10498683033169 1.21977855600062 1.34382524373165 Columns 5 th
3、rough 8 1.47612194644573 1.61514629644208 1.75881884576699 1.90449653438673 Columns 9 through 12 2.04900865016427 2.188*461 2.31977682471585 2.43807150515633 Columns 13 through 16 2.53968253238078 2.62100592628670 2.67901644757271 2.71148101768216 Columns 17 through 20 2.71712300843128 2.69571859920
4、382 2.64811384739078 2.57616043684702 Columns 21 through 24 2.48257772801500 2.37075712617031 2.24453058057755 2.10792744704554 Columns 25 through 28 1.96494288855677 1.81933699108106 1.67447782737116 1.53323499677732 Columns 29 through 32 1.39792381994500 1.27029521881156 1.151*453 1.04245724559826
5、 Columns 33 through 36 0.94329695305041 0.85406694858102 0.77449730249739 0.70413637458179 Columns 37 through 40 0.64241520775037 0.58870142585634 0.54234231959371 0.50269776253737 Column 41 0.46916418587400Data = -2.00000000000000 1.00000000000000 -1.90000000000000 1.10498683033169 -1.8000000000000
6、0 1.21977855600062 -1.70000000000000 1.34382524373165 -1.60000000000000 1.47612194644573 -1.50000000000000 1.61514629644208 -1.40000000000000 1.75881884576699 -1.30000000000000 1.90449653438673 -1.20000000000000 2.04900865016427 -1.10000000000000 2.188*461 -1.00000000000000 2.31977682471585 -0.90000
7、000000000 2.43807150515633 -0.80000000000000 2.53968253238078 -0.70000000000000 2.62100592628670 -0.60000000000000 2.67901644757271 -0.50000000000000 2.71148101768216 -0.40000000000000 2.71712300843128 -0.30000000000000 2.69571859920382 -0.20000000000000 2.64811384739078 -0.10000000000000 2.57616043
8、684702 0 2.48257772801500 0.10000000000000 2.37075712617031 0.20000000000000 2.24453058057755 0.30000000000000 2.10792744704554 0.40000000000000 1.96494288855677 0.50000000000000 1.81933699108106 0.60000000000000 1.67447782737116 0.70000000000000 1.53323499677732 0.80000000000000 1.39792381994500 0.
9、90000000000000 1.27029521881156 1.00000000000000 1.151*453 1.10000000000000 1.04245724559826 1.20000000000000 0.94329695305041 1.30000000000000 0.85406694858102 1.40000000000000 0.77449730249739 1.50000000000000 0.70413637458179 1.60000000000000 0.64241520775037 1.70000000000000 0.58870142585634 1.8
10、0000000000000 0.54234231959371 1.90000000000000 0.50269776253737 2.00000000000000 0.469164185874005.1.2 用欧拉法求解设常微分方程的初始问题有唯一解。则由欧拉法求初值问题(1),(2)的数值解的差分方程为:程序如下:建立函数文f1.mfunction x,y=f1(fun,x_span,y0,h)x=x_span(1):h:x_span(2);y(1)=y0;for n=1:length(x)-1 y(n+1)=y(n)+h*feval(fun,x(n),y(n);endx=x;y=y在MAT
11、LAB输入以下程序: clear all fun=inline( y*cos(x+2) x,y1=f1(fun,-2,2,1,0.1); x,y1 plot(x,y1,g*-)结果及其图象:ans = -1.90000000000000 1.10000000000000 -1.80000000000000 1.20945045818058 -1.70000000000000 1.32798465534234 -1.60000000000000 1.45485187516708 -1.50000000000000 1.58885260659392 -1.40000000000000 1.7282
12、8754069001 -1.30000000000000 1.87092926670362 -1.20000000000000 2.01402582996363 -1.10000000000000 2.15434436081705 -1.00000000000000 2.28826055379421 -0.90000000000000 2.41189579915842 -0.80000000000000 2.52129845713651 -0.70000000000000 2.61265966186586 -0.60000000000000 2.68254800178024 -0.500000
13、00000000 2.72814250373577 -0.40000000000000 2.74744062038227 -0.30000000000000 2.73941822501564 -0.20000000000000 2.70412232942903 -0.10000000000000 2.64268410367377 0 2.55724888375039 0.10000000000000 2.45082978042675 0.20000000000000 2.32710059365817 0.30000000000000 2.19015046372483 0.40000000000
14、000 2.04422599002736 0.50000000000000 1.89348605020813 0.60000000000000 1.74179062418299 0.70000000000000 1.59253854452440 0.80000000000000 1.44856157120511 0.90000000000000 1.31207486378276 1.00000000000000 1.184*502 1.10000000000000 1.06739566199422 1.20000000000000 0.96074840947946 1.300000000000
15、00 0.86483739767581 1.40000000000000 0.77943645422926 1.50000000000000 0.70408067871132 1.60000000000000 0.63814657271418 1.70000000000000 0.58092024172055 1.80000000000000 0.53165239417811 1.90000000000000 0.48960040640242 2.00000000000000 0.454058731286725.1.3用改进欧拉法求解:计算公式为:即先用欧拉法得,进而由(3)的第一式得初始近似
16、值,然后再用(3)的第二式进行迭代,反复改进这个近似值,直到(为所允许的误差)为止,并把取作为的近似值,这个方法就称为改进欧拉法。通常称(3)为预报校正公式,其中第一式称为预报公式,第二式称为校正公式。这个公式还可以写为:(下文改进欧拉法计算就以下面为准)建立函数文件f2.mfunction x,y=f2(fun,x_span,y0,h)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; y*cos(x+2) x,y2=f2(fun,
17、-2,2,1,0.1); x,y2 plot(x,y2,b+- -1.90000000000000 1.10472522909029 -1.80000000000000 1.21920722936399 -1.70000000000000 1.34289777810139 -1.60000000000000 1.47479651303944 -1.50000000000000 1.61338861747877 -1.40000000000000 1.75660494566073 -1.30000000000000 1.90181495275894 -1.20000000000000 2.045
18、86183721764 -1.10000000000000 2.185*243 -1.00000000000000 2.31576355572445 -0.90000000000000 2.43368296895904 -0.80000000000000 2.53497167173468 -0.70000000000000 2.61603367901013 -0.60000000000000 2.67384966784500 -0.50000000000000 2.70619076790723 -0.40000000000000 2.71178326404468 -0.300000000000
19、00 2.69040521940878 -0.20000000000000 2.64290353044130 -0.10000000000000 2.57112934628408 0 2.47779956085378 0.10000000000000 2.36630057150132 0.20000000000000 2.24045633274254 0.30000000000000 2.10428512493058 0.40000000000000 1.96176831551518 0.50000000000000 1.81665028035196 0.60000000000000 1.67
20、228260188095 0.70000000000000 1.53151888526290 0.80000000000000 1.39666016389536 0.90000000000000 1.26944574589307 1.00000000000000 1.151*251 1.10000000000000 1.04229147409852 1.20000000000000 0.94339433622144 1.30000000000000 0.85437588457241 1.40000000000000 0.77496982261477 1.50000000000000 0.704
21、72971912662 1.60000000000000 0.64309273433755 1.70000000000000 0.58943293626443 1.80000000000000 0.54310392687123 1.90000000000000 0.50347143064050 2.00000000000000 0.469937065943505.1.4 用4阶龙格库塔求解标准的四阶龙格-库塔公式(亦称为经典的四阶龙格-库塔公式):公式的截断误差阶为。建立函数文件f3.mfunction x,y=f3(fun,x_span,y0,h)x=x_span (1): k2=feval
22、(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; clear all; x,y3=f3(fun,-2,2,1,0.1); x,y3 plot(x,y3, y*- -1.90000000000000 1.10498674569681 -1.80000000000000 1.21977837351828 -1.70000000000000 1.34382495362539 -1.60000
23、000000000 1.47612154334742 -1.50000000000000 1.61514577980337 -1.40000000000000 1.75881821961767 -1.30000000000000 1.90449580650198 -1.20000000000000 2.04900783085566 -1.10000000000000 2.188*744 -1.00000000000000 2.31977585752433 -0.90000000000000 2.43807048140885 -0.80000000000000 2.53968146311619
24、-0.70000000000000 2.62100482231058 -0.60000000000000 2.67901531971109 -0.50000000000000 2.71147987684658 -0.40000000000000 2.71712186540483 -0.30000000000000 2.69571746425094 -0.20000000000000 2.64811272993642 -0.10000000000000 2.57615934548225 0 2.48257667095154 0.10000000000000 2.37075611204203 0.
25、20000000000000 2.24452961927934 0.30000000000000 2.10792655020966 0.40000000000000 1.96494206934318 0.50000000000000 1.81933626317413 0.60000000000000 1.67447720334541 0.70000000000000 1.53323448621294 0.80000000000000 1.39792342776452 0.90000000000000 1.27029494424773 1.00000000000000 1.151*231 1.1
26、0000000000000 1.04245718123925 1.20000000000000 0.94329697235934 1.30000000000000 0.85406703400224 1.40000000000000 0.77449743625256 1.50000000000000 0.70413654020399 1.60000000000000 0.64241539118308 1.70000000000000 0.58870161605161 1.80000000000000 0.54234250864301 1.90000000000000 0.50269794543536 2.00000000000000 0.469164360045035.1.5 问题讨论与分析由以上数值分析结果绘制表格:精确解欧拉方法改进的欧拉方法四阶龙格-库塔方法xiyi误差-2.00 1.0000000.000000-1.90 1.1049871.1000000.0049871.1047250.000262-1.80 1.2197791.2094500.0103281.2192070.0005711.219778-1.70 1.3438251.3279850.0158411.3428980.0009271.34