第五讲：常微分方程数值解

常微分方程数值解,数值分析第五讲,第五章：常微分方程数值解,5.1引言,1、常微分方程与解,为n阶常微分方程。,如果函数在区间a,b内n阶可导，称方程,为方程满足定解条件的解。,第五章：常微分方程数值解,解的图示,第五章：常微分方程数值解,本教材重点讨论定解问题(初值问题）,定解条件（初始条件）,是否能够找到定解问题的解取决于,仅有极少数的方程可以通过“常数变易法”、“可分离变量法”等特殊方法求得初等函数形式的解，绝大部分方程至今无法理论求解。,2、数值解的思想,第五章：常微分方程数值解,（1）将连续变量离散为,（2）用代数的方法求出解函数在点的近似值,如果找不到解函数数学界还关注：解的存在性解的唯一性解的光滑性解的振动性解的周期性解的稳定性解的混沌性,5.2Euler方法,第五章：常微分方程数值解,第一步：连续变量离散化,第二步：用直线步进,1、Euler格式,思想：向前差商近似代替导数,18世纪最杰出的数学家之一，13岁时入读巴塞尔大学，15岁大学毕业，16岁获得硕士学位。1727年-1741年（20岁-34岁）在彼得堡科学院从事研究工作，在分析学、数论、力学方面均有出色成就，并应俄国政府要求，解决了不少地图学、造船业等实际问题。24岁晋升物理学教授。1735年（28岁）右眼失明。,第五章：常微分方程数值解,1741年-1766（34岁-59岁）任德国科学院物理数学所所长，任职25年。在行星运动、刚体运动、热力学、弹道学、人口学、微分方程、曲面微分几何等研究领域均有开创性的工作。1766年应沙皇礼聘重回彼得堡，在1771年（64岁）左眼失明。Euler是数学史上最多产的数学家，平均以每年800页的速度写出创造性论文。他去世后，人们用35年整理出他的研究成果74卷。,第五章：常微分方程数值解,第五章：常微分方程数值解,后退Euler格式,令,得,令,后退Euler格式的值,Euler格式的值,2.后退的Euler格式,思想：向后差商代替导数,3、梯形格式,4、改进的Euler格式,第五章：常微分方程数值解,为方便计算，一般用以下改进格式计算,用改进格式计算例5.1的结果见P110表5.2,第五章：常微分方程数值解,例P106,第五章：常微分方程数值解,初值问题,Bernoulli方程,由Bernoulli方程的求解方法可得解析解,Euler格式为,步进计算结果见P106表5.1,第五章：常微分方程数值解,Euler值,Euler格式的误差分析,事实上Euler格式的每一步都存在误差，为了方便讨论算法的好坏，假定第n步准确的前提下分析第n+1步的误差，称为局部截断误差。,第五章：常微分方程数值解,由Taylor公式,第五章：常微分方程数值解,Euler格式的误差为,5、两步Euler格式,第五章：常微分方程数值解,如果令,则有,记,两步Euler格式公式,6.、两点预测校正公式,6、代数精度,第五章：常微分方程数值解,已知Euler格式,即Euler格式具有一阶精度,第五章：常微分方程数值解,记,则,其中,假定,第五章：常微分方程数值解,特别要注意的是：一般,两点预测格式具有二阶精度。,当时,所以,第五章：常微分方程数值解,考察两点校正格式的精度,则,第五章：常微分方程数值解,比较,得,即梯形格式具有二阶精度，因此两步格式从预测到校正均达到二阶精度。,因此得具有二阶精度的两步Euler格式,预测,校正,第五章：常微分方程数值解,5.3Lunge-Kutta方法,1、二阶Lunge-Kutta方法（P113-P115),第五章：常微分方程数值解,依据精度要求的待定系数法,令,确定使,具有二阶精度,第五章：常微分方程数值解,对照,第五章：常微分方程数值解,可解得,得,改进的Euler格式,3、三阶Lunge-Kutta方法,第五章：常微分方程数值解,补充,确定参数使,第五章：常微分方程数值解,使,具有三阶精度,即,分别将,在点Taylor展开，代入（P116),与的Taylor展开比较,第五章：常微分方程数值解,得,可解得,第五章：常微分方程数值解,得一个三阶精度的Runge-Kutta格式,4、四阶Lunge-Kutta方法见P117,5.4几种方法的数值计算,例5.1P106,第五章：常微分方程数值解,四阶经典Lunge-Kutta方法,例5.1P106,第五章：常微分方程数值解,几种方法的结果与误差,第五章：常微分方程数值解,参考程序-Euler,x=0:0.01:1;y=sqrt(1+2.*x);a=0.0;b=1.0;n=10;h=(b-a)/n;x0=a:h:b;y0(1)=1.0;fork=1:10y0(k+1)=y0(k)+h*(y0(k)-2*x0(k)/y0(k);endfori=1:10y1(1)=1.0;y1(i+1)=y1(i)+h*(y1(i)-2*x0(i)/y1(i);y1(i+1)=y1(i)+h*(y1(i)-2*x0(i)/y1(i)+y1(i+1)-2*x0(i+1)/y1(i+1)/2;endplot(x,y,b);holdon;plot(x0,y0,or);holdon;plot(x0,y1,*);,第五章：常微分方程数值解,参考程序-Lunge_Kutta,x=0:0.01:1;y=sqrt(1+2.*x);a=0.0;b=1.0;n=10;h=(b-a)/n;x0=a:h:b;y0(1)=1.0;fork=1:10k1=y0(k)-2*x0(k)/y0(k);k2=y0(k)+h*k1/2-(2*x0(k)+h)/(y0(k)+h*k1/2);k3=y0(k)+h*k2/2-(2*x0(k)+h)/(y0(k)+h*k2/2);k4=y0(k)+h*k3-2*(x0(k)+h)/(y0(k)+h*k3);y0(k+1)=y0(k)+h*(k1+2*k2+2*k3+k4)/6;endholdon;plot(x,y,b);holdon;plot(x0,y0,or);,第五章：常微分方程数值解,5.5线性多步方法,1、Adams显式格式,第五章：常微分方程数值解,P27(2.5.12),第五章：常微分方程数值解,第五章：常微分方程数值解,第五章：常微分方程数值解,r=3时四步显式Adams公式,第五章：常微分方程数值解,误差分析,第五章：常微分方程数值解,第五章：常微分方程数值解,例5.1P106,第五章：常微分方程数值解,Adams程序,x=0:0.01:1;y=sqrt(1+2.*x);a=0.0;b=1.0;n=10;h=(b-a)/n;x0=a:h:b;y0(1)=1.0;y0(2)=1.0954;y0(3)=1.1832;y0(4)=1.2649;fork=4:10y0(k+1)=y0(k)+h*(55*(y0(k)-2*x0(k)/y0(k)-59*(y0(k-1)-2*x0(k-1)/y0(k-1)+37*(y0(k-2)-2*x0(k-2)/y0(k-2)-9*(y0(k-3)-2*x0(k-3)/y0(k-3)/24;endholdon;plot(x,y,b);holdon;plot(x0,y0,or);,第五章：常微分方程数值解,几种方法的结果与误差,第五章：常微分方程数值解,2、Adams隐式格式,隐式格式见P125,第五章：常微分方程数值解,3、A