Fork me on GitHub

Euler-Maruyama discretization("欧拉-丸山"数值解法)

欧拉法的来源

  在数学和计算机科学中,欧拉方法(Euler method)命名自它的发明者莱昂哈德·欧拉,是一种一阶数值方法,用以对给定初值的常微分方程(即初值问题)求解。它是一种解决常微分方程数值积分的最基本的一类显型方法(Explicit method)。

什么是欧拉法

  欧拉法是以流体质点流经流场中各空间点的运动即以流场作为描述对象研究流动的方法。——流场法

  它不直接追究质点的运动过程,而是以充满运动液体质点的空间——流场为对象。研究各时刻质点在流场中的变化规律。将个别流体质点运动过程置之不理,而固守于流场各空间点。通过观察在流动空间中的每一个空间点上运动要素随时间的变化,把足够多的空间点综合起来而得出的整个流体的运动情况。

  常微分方程的数值解法的一种。基本思想是迭代。其中分为前进的EULER法、后退的EULER 法、改进的EULER法。所谓迭代,就是逐次替代,最后求出所要求的解,并达到一定的精度。误差可以很容易的计算出来。

欧拉算法

  微分方程的本质特征是方程中含有导数项,数值解法的第一步就是设法消除其导数值,这个过程称为离散化。实现离散化的基本途径是用向前差商来近似代替导数,这就是欧拉算法实现的依据。欧拉(Euler)算法是数值求解中最基本、最简单的方法,但其求解精度较低,一般不在工程中单独进行运算。所谓数值求解,就是求问题的解y(x)在一系列点上的值y(xi)的近似值yi。对于常微分方程:

  \frac{dy}{dx}=f(x,y),x∈[a,b]

  y(a) = y0

  可以将区间[a,b]分成n段,那么方程在第xi点有y'(xi) = f(xi,y(xi)),再用向前差商近似代替导数则为:\frac{(y(x_i+1)-y(x_i))}{h}= f(x_i,y(x_i)),在这里,h是步长,即相邻两个结点间的距离。因此可以根据xi点和yi点的数值计算出yi+1来:

  y_{i+1}= y_i+ h\times f(x_i ,y_i),i=0,1,2,L

  这就是欧拉格式,若初值yi + 1是已知的,则可依据上式逐步算出数值解y1,y2,L。

  为简化分析,人们常在yi准确yi = y(xi)的前提下估计误差 y(xi + 1) − yi + 1,这种误差称为局部截断误差。

  如果一种数值方法的局部截断误差为O(hp + 1),则称它的精度是p阶的,或称之为p阶方法。欧拉格式的局部截断误差为O(h2),由此可知欧拉格式仅为一阶方法。

欧拉公式

  y(x_i+1)=yi+h \times f(x_i,y_i)

  且x_i=x_0+i\times h (i=0,1,2,…,n-1)

  局部截断误差是O(h2)

非形式化的几何诠释

考虑计算下面一个未知曲线的形状:它具有给定的起点并且满足一个给定的微分方程。这里,所谓“微分方程”可以看作能够通过曲线上任意点的位置而计算出这一点的切线斜率的公式。

思路是,一开始只知道曲线的起点,曲线其他部份是未知的,不过通过微分方程,的斜率可以被计算出来,也就得到了切线。

顺着切线向前走一小步到点。如果我们假设是曲线上的一点(实际上通常不是),那么同样的道理就可以确定下一条切线,依此类推。在经过几步之后,一条折线就被计算出来了。大部分的情况下,这条折线与原先的未知曲线偏离不远,并且任意小的误差都可以通过减少步长来得到(虽然对于刚性方程而言会比较复杂)。

欧拉方法的推导

以以下微分方程为例 y'(t) = f(t,y(t)), \qquad \qquad y(t_0)=y_0,

希望用y 在点(t0,y(t0)) 附近的线性近似来得到其近似解。利用时间 tn 时的数值,若用单步的欧拉方法,可得到时间tn+1 =&nbsp ;tn + h 时的近似值如下:

y_{n+1} = y_n + hf(t_n,y_n). \qquad \qquad

欧拉方法是一种显型方法,也就是说yn + 1 的解是yi, i \leq n 的显函数。

欧拉方法可以求解一阶的微分方程,一个N阶的微分方程可以用导入N − 1 个变数来表示y' , y'', ..., y(N) 的方式,分解N个一阶的微分方程。因此可以用欧拉方法求解以下的向量\mathbf{y}(t)=(y(t),y'(t),y''(t),...,y^{(N )}(t)) 来得到高阶微分方程的解。

改进的欧拉算法

  先用欧拉法求得一个初步的近似值,称为预报值,然后用它替代梯形法右端的yi+1再直接计算 fi+1,得到校正值yi+1,这样建立的预报-校正系统称为改进的欧拉格式:

  预报值y~{i+1}=y_i+1 + h \times f(x_i,y_i)

  校正值 y_{i+1} =y_i+(\frac{h}{2})\times [f(x_i,y_i)+f(x_i+1,y~i+1)]

  它有下列平均化形式:

  y_p=y_i+h \times f(x_i,y_i)

  且 y_c=y_i+h \times f(x_i+1,y_p)

  且 y_{i+1}=\frac{(x_p+y_c)}{2}

  它的局部截断误差为O(h3),可见,改进欧拉格式较欧拉格式提高了精度,其截断误差比欧拉格式提高了一阶。

  注:欧拉法用差商 \frac{y(x_{i+1})-y(x_i)}{h} 近似代替y(xi)的导数,局部截断误差较大;改进欧拉法先用欧拉法求出预报值,再利用梯形公式求出校正值,局部截断误差比欧拉法低了一阶,较大程度地提高了计算精度。

 

posted @ 2016-12-08 13:44  stardsd  阅读(7335)  评论(0编辑  收藏  举报