用改进的欧拉方法和四阶龙格-库塔方法解初值问题

用改进的欧拉方法和四阶龙格-库塔方法解初值问题
一、题目：
取步长2.0=h ，分别用改进的欧拉方法和四阶龙格—库塔方法解初值问题
⎩
⎨⎧=≤≤+=.1)0(;10,'y x y x y 并比较结果。

二、基本思想：
1.改进的欧拉方法
改进的欧拉方法用梯形公式计算()()()()+1+1-=,.n n x n n x y x y x f x y x dx ⎰中的积分，并以
n y 和+1n y 分别表示()n y x 和()+1n y x 的近似值，就得到
()()()+1+1+1=+,+,.2
n n n n n n h y y f x y f x y （梯形公式） 梯形公式也是一个一步法公式。

由于公式右端也含有未知的+1n y ，故被称作是隐式的。

隐式格式实际上是关于+1n y 的一个函数方程。

为了避免解方程，可以采用欧拉方法计算初始值，再由梯形公式计算。

这样建立起来的计算格式称为改进的欧拉格式：
()()()()+1+1+1=+,,=+,+,.2
p n n n n n n n n n y y hf x y h y y f x y f x y ⎧⎪⎨⎪⎩ 梯形公式也可以采用迭代法求解。

如果仍然采用欧拉方法计算迭代初值，那么计算格式就是
[]()[]()[]()()
0+1+1+1+1+1=+,,=+,+,.2n n
n n k k n n n n n n y y hf x y h y y f x y f x y ⎧⎪⎨⎪⎩ 由于已假定(),f x y 满足里普希兹条件，所以有
[][][]()[]()
[][]+1-1-1+1+1+1,+1+1+1+1+1-=-,-.22
k k k k k k n n n n n n n n h hL y y f x y f x y y y ≤ 从而，迭代的收敛条件是0<<1.2hL 2.四阶龙格-库塔方法
龙格—库塔方法不是用求导数的办法，而是用计算不同点上()y x f ,的函数值，然后对这些函数值作线性拟合，构造近似公式。

组合的原则是使得近似公式与泰勒展开式有尽可能多的项吻合，以达到较高的精度。

改进的欧拉格式
()()()()+1+1+1=+,,=+,+,.2
p n n n n n n n n n y y hf x y h y y f x y f x y ⎧⎪⎨⎪⎩ 可以改写为
()2112
1K K y y n n ++
=+, 其中 ()().,,
,1121K y x hf K y x hf K n n n n +==+
可见改进的欧拉格式可这样理解：首先计算()n n y x ,和()11,K y x n n ++两点的函数值()n n y x f ,和()11,K y x f n n ++，然后以它们的算术平均值作为()y x f ,的近似值，而点()+11,+n n x y K 中的1K 则是通过已知值n y 来预报的。

这个处理过程告诉我们：如果设法多计算几个点处的函数值，然后以这些值得某种加权平均值作为(),f x y 的近似值，则有可能构造出更精确的计算格式。

一般地龙格—库塔格式可以写成
+1=1=+,N
n n i i i y y C K ∑
其中
()1-1=1,,
=+,+, =2,3,
,.n n i i n i n ij j j K hf x y K hf x a h y b K i N =⎛⎫ ⎪⎝⎭∑
选择参数,,i i ij C a b 的原则是：要求+1=1
=+
,N n n i i i y y C K ∑的右端在(),n n x y 处的泰勒展开式按h 的幂重新排列后得到的 23+1123=++++
n n y y d h d h d
h
与微分方程的解()+1n y x 在点n x 的展开式 ()()()()()'2''3'''+111=++++2!3!n n n n n y x y x hy x h y x h y
x
又尽可能多的项重合。

经典的四阶龙格—库塔方法的计算公式是：
()()()+11234121324
31=+22,6,,11+,+,2211+,+,22+,+.n n n n n n n n n n y y K K K K K hf x y K hf x h y K K hf x h y K K hf x h y K ⎧+++⎪⎪⎪=⎪⎪⎛⎫=⎨ ⎪⎝⎭⎪⎪⎛⎫=⎪ ⎪⎝⎭⎪⎪=⎩ 三、程序设计：
1.改进的欧拉方法
eulerpro.m 文件
function [x,y]=eulerpro(h,x0,y0);
x(1)=x0;y(1)=y0;
for i=1:10;
x(i+1)=x(i)+h;
y1=y(i)+h*(x(i)+y(i));
y2=y(i)+h*(x(i+1)+y1);
y(i+1)=(y1+y2)/2;
end
计算：[x,y]=eulerpro(0.1,0,1)
2.四阶龙格-库塔方法
lgkt4j.m 文件
function [x,y]=lgkt4j(x0,y0,h)
x(1)=x0;y(1)=y0;
for i=1:10
x(i+1)=x(i)+h;
K1=h*(x(i)+y(i));
K2=h*(x(i)+1/2*h+y(i)+1/2*K1);
K3=h*(x(i)+1/2*h+y(i)+1/2*K2);
K4=h*(x(i)+h+y(i)+K3); y(i+1)=y(i)+1/6*(K1+2*K2+2*K3+K4);
end
计算：[x,y]=lgkt4j(0,1,0.1)
四、结果分析：
改进的欧拉方法运行结果：
x = 0 0.1000 0.2000 0.3000
0.4000 0.5000 0.6000 0.7000
0.8000 0.9000 1.0000
y = 1.0000 1.1100 1.2421 1.3985
1.5818 1.7949
2.0409 2.3231
2.6456
3.0124 3.4282
四阶龙格-库塔方法运行结果：
x = 0 0.1000 0.2000 0.3000
0.4000 0.5000 0.6000 0.7000
0.8000 0.9000 1.0000
y =1.0000 1.1103 1.2428 1.3997
1.5836 1.7974
2.0442 2.3275
2.6511
3.0192 3.4366
准确解：
x = 0 0.1000 0.2000 0.3000
0.4000 0.5000 0.6000 0.7000
0.8000 0.9000 1.0000
y= 1.0000 1.1103 1.2428 1.3997
1.5836 1.7974
2.0442 2.3275
2.6511
3.0192 3.4366
比较计算结果可得四阶龙格-库塔方法更接近准确解。

改进的欧拉方法的局部误差为()3hο，整体截断误差为()2hο，而四阶龙格-库塔方法的局部截断误差为()5hο，整体截断误差为()4hο，误差明显小于改进的欧拉方法的误差，所以龙格库塔方法更接近于准确值。