数值积分算法误差分析.doc

1
实验名称：数值积分算法误差分析
1.实验原理
1)欧拉法原理
在数学和计算机科学中，欧拉方法（Eulermethod）命名自它的发明者HYPERLINK"http://wiki.mbalib.com/w/index.php?title=0%23%20*)n%252%11%145'0%1D&action=edit"莱昂哈德·欧拉，是一种一阶数值方法，用以对给定初值的常微分方程(即初值问题)求解。它是一种解决常微分方程数值积分的最基本的一类显型方法(Explicitmethod)。
微分方程的本质特征是方程中含有导数项，数值解法的第一步就是设法消除其导数值，这个过程称为离散化。实现离散化的基本途径是用向前差商来近似代替导数，这就是欧拉算法实现的依据。欧拉(Euler)算法是数值求解中最基本、最简单的方法，但其求解精度较低，一般不在工程中单独进行运算。所谓数值求解，就是求问题的解y(x)在一系列点上的值y(xi)的近似值yi。对于常微分方程：

可以将区间[a,b]分成n段，那么方程在第xi点有y'(xi)=f(xi,y(xi))，再用向前差商近似代替导数则为：，在这里，h是步长，即相邻两个结点间的距离。因此可以根据xi点和yi点的数值计算出yi+1来：
，i=0,1,2,L
这就是欧拉格式，若初值yi+1是已知的，则可依据上式逐步算出数值解y1,y2,L。
龙哥库塔法原理
HYPERLINK"http://zh.wikipedia.org/wiki/%CA%FD%D6%B5%B7%D6%CE%F6"数值分析中，龙格－库塔法（Runge-Kutta）是用于模拟HYPERLINK"http://zh.wikipedia.org/wiki/%B3%A3%CE%A2%B7%D6%B7%BD%B3%CC"常微分方程的解的重要的一类隐式或显式迭代法。这些技术由数学家HYPERLINK"http://zh.wikipedia.org/wiki/%BF%A8%A0%96%A1%A4%FD%88%B8%F1"卡尔·龙格和HYPERLINK"http://zh.wikipedia.org/w/index.php?title=%C2%ED%B6%A1%A1%A4%CD%FE%B6%FB%BA%A3%C4%B7%A1%A4%BF%E2%CB%FE&action=edit&redlink=1"马丁·威尔海姆·库塔于1900年左右发明。
龙格库塔法的家族中的一个成员如此常用，以至于经常被称为“RK4”或者就是“龙格库塔法”。该方法主要是在已知方程导数和初值信息，利用计算机仿真时应用，省去求解微分方程的复杂过程。
令HYPERLINK"http://zh.wikipedia.org/wiki/%B3%F5%D6%B5%86%96%EE%7D"初值问题表述如下。

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

其中

这样，下一个值(yn+1)由现在的值(yn)加上时间间隔(h)和一个估算的斜率的乘积决定。该斜率是以下斜率的加权平均：
k1是时间段开始时的斜率；
k2是时间段中点的斜率，通过HYPERLINK"http://zh.wikipedia.org/wiki/%C5%B7%C0%AD%B7%A8"欧拉法采用斜率k1来决定y在点tn+h/2的值；
k3也是中点的斜率，但是这次采用斜率k2决定y值；
k4是时间段终点的斜率，其y值用k3决定。
当四个斜率取平均时，中点的斜率有更大的权值：

RK4法是四阶方法，也就是说每步的误差是h5HYPERLINK"http://zh.wikipedia.org/w/index.php?title=%B4%F3O%BC%C7%B7%A8&action=edit&redlink=1"阶，而总积累误差为h4阶。
2实验方案
目标方程为

对该方程进行分析可知，该方程一阶微分方程，为y对t的微分。
子函数程序：
subf.m文件
functiondy=subf(t,y)
dy=(1-t)*y;

目标方程解析解
>>H=dsolve('Dy=(1-t)*y','y(0)=1','t')

得到：H=exp(1)^(1/2)/exp((t-1)^2/2)

3.实验程序
欧拉法仿真脚本程序
clear%欧拉法
h=0.1;
tN=10;
t=0:h:tN;
N=length(t);
y(1)=1;
j=1;
forj=1:N
dataY(j)=y(1);
yn=y+h*subf(t(j),y);
y=yn;
end
%解析解法
ezplot('exp(1)^(1/2)/exp((t-1)^2/2)',[010])
holdon
plot(t,dataY,'r+')