实验名称微分方程数值解..docVIP

探索实验8 常微分方程初值问题数值解 实验目的 了解常微分方程初值问题的数值解概念，掌握解常微分方程初值问题的Euler方法及改进的Euler方法和Runge-Kutta方法解常微分方程初值问题的算法构造和计算。能用程序实现解常微分方程初值问题的Euler方法、改进的Euler方法和经典的Runge-Kutta方法，学习用计算机求常微分方程初值问题数值解的一些科学计算方法和简单的编程技术。 二、概念与结论 1. 常微分方程初值问题： 常微分方程特解问题称为初值问题，通常其形式为 2．常微分方程初值问题数值解： 常微分方程初值问题的解y(x)在[a,b]上的有限个值y(xk)的近似值yk称为常微分方程初值问题数值解，其中 xk= xk-1 + hk-1 ，k=1,2,3,…,N。xk称为节点， hk 称为步长。 通常，步长h不变，取为等距步长 hk=h=（b-a）/N，N为等分区间[a,b]分割数。 3.常微分方程初值问题数值解法： 求常微分方程初值问题数值解yk的方法称为微分方程初值问题数值解法。 4.单步法： 在计算yk+1之时只用到yk 的方法，其计算公式有： 显式单步法计算公式: yk+1=yk+hΦ(xk,yk;h） 隐式单步法计算公式: yk+1 = yk +hΦ(xk ,yk,yk+1,h） 式中函数Φ是连续函数，称为增量函数。 5.多步法： 在计算yk+1之时不仅用到yk ，还要用yk-1，yk-2，…,一般m 步法要用到yk,yk-1,yk-2，… yk-m+1， 多步法也有显式方法和隐式方法之分。 6.数值解法的局部截断误差： 假设某常微分方程初值问题数值解法在x= xk没有误差，即y(xk)= yk，称 Tk+1 =y(xk+1)- yk+1 为该数值方法的局部截断误差。 显式单步法有其局部截断误差为： Tk+1 =y(xk+1)- y(xk)- hΦ(xk,y(xk),h） 7.数值解法的阶： 如果某常微分方程初值问题数值解法的局部截断误差为： Tk+1 =?(h p+1) 则称该数值方法的阶为P,或称该数值方法为P阶方法。 数值方法的阶越高，方法越好。 三、程序中Mathematica语句解释： 1.k++,++k,k--,--k k++ 表示赋值关系 k = k+1 , 如: k=1;Table[++k,{5}]获得表{2,3,4,5,6} ++k 表示先处理k的值,再做赋值 k=k+1, 如: k=1;Table[k++,{5}]获得表{1,2,3,4,5} k-- 表示赋值关系 k = k-1, 如: k=1;Table[k--,{5}]获得表{1, 0, -1, -2, -3} --k 表示先处理k的值,再做赋值 k=k-1,如:k=1;Table[--k,{4}]获得表{0,-1,-2,-3} 2. x+=k, x*=k x+=k 表示 x = x + k x*=k 表示 x = x * k 3.For[stat,test,incr,body] 以stat为初值,重复计算incr和body直到test为False终止 。这里start为初始值,test为条件,incr为循环变量修正式,body为循环体,通常由incr项控制test的变化。 注意: 上述命令形式中的start可以是由复合表达式提供的多个初值,如果循环体生成 Break[ ] 语句,则退出For循环; 如果循环体生成Continue[ ] 语句,则由incr的增量进入For语句的下一次循环。 四、方法与程序 在实际问题中，常常会遇到一些常微分方程初值问题。对这些问题如果只求助于高等数学中介绍的用求通解加确定边界条件的方法去求解往往是无能为力的。因为一方面通解可能求不出来，另一方面所求的解可能是不能用初等函数表达的形式。人们处理这类问题主要采用常微分方程初值问题数值解的方法来求解。这类方法有很多，这里主要介绍解常微分方程初值问题的Euler方法、改进的Euler方法和Runge-Kutta方法的构造过程和算法程序。 Euler方法 Euler方法是最简单的求微分方程数值解的方法。这个方法由于精度不高，实用中较少使用。人们常用Euler方法来说明求微分方程数值解涉及到的一些问题。 Euler方法的构造过程： Euler方法涉及的计算公式有很多方法，这里介绍在求微分方程数值解用的较多的Taylor展开法。 设 f(x,y) 充分光滑, xn+1=xn+h,将y(