多空指标公式解微分方程欧拉法,R-K法及其MATLAB实例

.
解微分方程的欧拉法，龙格-库塔法及其MATLAB简单
实例
欧拉方法（Euler method)用以对给定初值的常微分方程(即初值问题)求解
分为前进EULER法、后退EULER法、改进的EULER法。
缺点：
欧拉法 简单地取切线的端点作为下一步的起点进行计算，当步数增多时，误差会
因积累而越来越大。因此欧拉格 式一般不用于实际计算。
改进欧拉格式：
为提高精度，需要在欧拉格式的基础上进行改进。 采用区间两端的斜率的平均值
作为直线方程的斜率。改进欧拉法的精度为二阶。
算法为：
微分方程的本质特征是方程中含有导数项，数值解法的第一步就是设法消除其导
数值。
对于常微分方程：

x∈[a,b]
y(a) = y0
可以将区间[a,b]分成n段，那么方程在第xi点有y'(xi) = f(xi,y(xi))，
再用向前差商近似代替导数则为：

在这里，h是步长， 即相邻两个结点间的距离。因此可以根据xi点和yi点的数
值计算出yi+1来：

i=0,1,2,L
这就是向前欧拉格式。
改进的欧拉公式：
将向前欧拉公式中的导数f(xi,yi)改为微元两端导数的平均，即
..
.
上式便是梯形的欧拉公式。
可见，上式是隐式格式，需要迭代求解。为了便于求解，使用改进的欧拉公式：

数值分析中，龙格－库塔法（Runge- Kutta）是用于模拟常微分方程的解的重要
的一类隐式或显式迭代法。
实际上，龙格-库 塔法是欧拉方法的一种推广，向前欧拉公式将导数项简单取为
f(xn,yn),而改进的欧拉公式将导 数项取为两端导数的平均。
龙格-库塔方法的基本思想：
在区间[xn,xn+1]内多取几个点，将他们的斜率加权平均，作为导数的近似。
< br>龙格库塔法的家族中的一个成员如此常用，以至于经常被称为“RK4”或者就是
“龙格库塔法” 。

令初值问题表述如下。

则，对于该问题的RK4由如下方程给出：

其中




这样，下一个值(
y< br>n
+1
)由现在的值(
y
n
)加上时间间隔(
h)和一个估算的斜率的乘积
决定。该斜率是以下斜率的加权平均：
..
.
k1是时间段开始时的斜率；
k2是时间段中点的斜率，通过欧拉法采用斜率k1来决定y在点tn + h2
的值；
k3也是中点的斜率，但是这次采用斜率k2决定y值；
k4是时间段终点的斜率，其y值用k3决定。
当四个斜率取平均时，中点的斜率有更大的权值：

RK4法是四阶方法，也就是说每步的误差是
h
5
阶，而总积累误 差为
h
4
阶。
注意上述公式对于标量或者向量函数(
y
可以是向量)都适用。


例子：

下面给出了数值求解该微分方程的简单程序。
其中y1,y2,y3,y4分别为向前欧拉公式，改进的欧拉公式，4级4阶龙格- 库塔
公式及精确解。
h=0.1;
x=0:h:1;
y1=zeros(size(x));
y1(1)=1;
y2=zeros(size(x));
y2(1)=1;
y3=zeros(size(x));
y3(1)=1;

for i1=2:length(x)
y1(i1)=y1(i1-1)+h*(y1(i1-1)-2*x(i1-1)y1(i1-1));
k1=y2(i1-1)-2*x(i1-1)y2(i1-1);
k2=y2(i1-1)+h*k1-2*x(i1)(y2(i1-1)+h*k1);
y2(i1)=y2(i1-1)+h*(k1+k2)2;


k1=y2(i1-1)-2*x(i1-1)y2(i1-1);
k2=y2(i1-1)+h*k12-2*(x(i1-1)+h2)(y2(i1-1)+h*k12);
k3=y2(i1-1)+h*k22-2*(x(i1-1)+h2)(y2(i1-1)+h*k22);
..
.
k4=y2(i1-1)+h*k3-2*(x(i1-1)+h)(y2(i1-1)+h*k3);
y3(i1)=y3(i1-1)+(k1+2*k2+2*k3+k4)*h6;

end
y4=sqrt(1+2*x);
%plot(x,y1,x,y2,x,y3,x,y4)
%legend('y1','y2','y3','y4')

plot(x,y4-y1,x,y4-y2,x,y4-y3)
legend('y1','y2','y3')



..