《数值分析》第五讲：常微分方程数值解.ppt
常微分方程数值解,理学院 陈丽娟,数值分析第五讲,第五章：常微分方程数值解,5.1 引言,1、常微分方程与解,为n阶常微分方程。,如果函数 在区间a,b内n阶可导，称方程,为方程满足定解条件的解。,第五章：常微分方程数值解,解的图示,第五章：常微分方程数值解,本教材重点讨论定解问题(初值问题）,定解条件（初始条件）,是否能够找到定解问题的解取决于,仅有极少数的方程可以通过“常数变易法”、“可分离变量法”等特殊方法求得初等函数形式的解，绝大部分方程至今无法理论求解。,2、数值解的思想,第五章：常微分方程数值解,（1）将连续变量 离散为,（2）用代数的方法求出解函数 在 点的近似值,如果找不到解函数数学界还关注：解的存在性解的唯一性解的光滑性解的振动性解的周期性解的稳定性解的混沌性,5.2 Euler方法,第五章：常微分方程数值解,第一步：连续变量离散化,第二步：用直线步进,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卷。,第五章：常微分方程数值解,例 P106,第五章：常微分方程数值解,初值问题,Bernoulli方程,由Bernoulli方程的求解方法可得解析解,Euler格式为,步进计算结果见P106表5.1,第五章：常微分方程数值解,Euler值,Euler格式的误差分析,事实上Euler格式的每一步都存在误差，为了方便讨论算法的好坏，假定第 n 步准确的前提下分析第n+1步的误差，称为局部截断误差。,第五章：常微分方程数值解,由Taylor公式,第五章：常微分方程数值解,Euler格式的误差为,2、后退Euler格式,令,得,令,Euler格式,第五章：常微分方程数值解,后退Euler格式,令,得,令,后退Euler格式的值,Euler格式的值,3、梯形格式,4、改进的Euler格式,第五章：常微分方程数值解,为方便计算，一般用以下改进格式计算,用改进格式计算例5.1的结果见P110表5.2,第五章：常微分方程数值解,5、两步Euler格式,第五章：常微分方程数值解,已知Euler格式,即Euler格式具有一阶精度,如果令,并假定,第五章：常微分方程数值解,则有,记,则,其中,第五章：常微分方程数值解,特别要注意的是：一般,两点预测格式具有二阶精度。,当 时,所以,第五章：常微分方程数值解,考察两点校正格式的精度,则,第五章：常微分方程数值解,比较,得,即梯形格式具有二阶精度，因此两步格式从预测到校正均达到二阶精度。,因此得具有二阶精度的两步Euler格式,预测,校正,第五章：常微分方程数值解,5.3 Lunge-Kutta方法,1、二阶Lunge-Kutta方法（P113-P115),第五章：常微分方程数值解,依据精度要求的待定系数法,令,确定 使,具有二阶精度,第五章：常微分方程数值解,对照,第五章：常微分方程数值解,可解得,得,改进的Euler格式,3、三阶Lunge-Kutta方法,第五章：常微分方程数值解,补充,确定参数 使,具有二阶精度,第五章：常微分方程数值解,使,具有三阶精度,即,分别将,在点 Taylor展开，代入（P116),与 的Taylor展开比较,第五章：常微分方程数值解,得,可解得,第五章：常微分方程数值解,得一个三阶精度的Runge-Kutta格式,4、四阶Lunge-Kutta方法见P117,5.4 几种方法的数值计算,例5.1 P106,第五章：常微分方程数值解,四阶经典Lunge-Kutta方法,例5.1 P106,第五章：常微分方程数值解,几种方法的结果与误差,第五章：常微分方程数值解,参考程序-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;for k=1:10 y0(k+1)=y0(k)+h*(y0(k)-2*x0(k)/y0(k);endfor i=1:10 y1(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);hold on;plot(x0,y0,or);hold on;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;for k=1:10 k1=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;endhold on;plot(x,y,b);hold on;plot(x0,y0,or);,第五章：常微分方程数值解,5.5 线性多步方法,1、Adams显式格式,第五章：常微分方程数值解,P27(2.5.12),第五章：常微分方程数值解,第五章：常微分方程数值解,第五章：常微分方程数值解,r=3 时见P125（5.5.6）,第五章：常微分方程数值解,（5.5.6）式的误差分析,第五章：常微分方程数值解,第五章：常微分方程数值解,例5.1 P106,第五章：常微分方程数值解,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;for k=4:10 y0(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;endhold on;plot(x,y,b);hold on;plot(x0,y0,or);,第五章：常微分方程数值解,几种方法的结果与误差,第五章：常微分方程数值解,2、Adams隐式格式,隐式格式见P125,第五章：常微分方程数值解,3、Adams改进格式,校正,预测,第五章：常微分方程数值解,4、线性多步的一般格式,1、由（1）获得一个2阶精度显示格式,（1）,第五章：常微分方程数值解,只需令,解得,Euler两步格式,第五章：常微分方程数值解,四阶格式见P128,Milne格式,Hamming格式,见P130,第五章：常微分方程数值解,
