第9章常微分方程初值问题数值解法教学课件

数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法1第九章第九章 常微分方程初值问题常微分方程初值问题数值解法数值解法/*Numerical Method for Ordinary Differential Equations*/许多实际问题的数学模型是微分方程或微分许多实际问题的数学模型是微分方程或微分方程的初值问题方程的初值问题,如物体运动如物体运动,电路震荡电路震荡,化学反化学反映及生物群体的变化等。映及生物群体的变化等。能用解析方法求出精确解的微分方程为数不能用解析方法求出精确解的微分方程为数不多多,而且有的方程即使有解析解而且有的方程即使有解析解,也可能由于解的也可能由于解的表达式非常复杂而不易计算表达式非常复杂而不易计算,因此有必要研究微因此有必要研究微分方程的数值解法。分方程的数值解法。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法2常微分方程中介绍的微分方程主要有常微分方程中介绍的微分方程主要有:(1)(1)变量可分离的方程变量可分离的方程 (2)(2)一阶线性微分方程一阶线性微分方程(贝努利方程贝努利方程)(3)(3)可降阶的一类高阶方程可降阶的一类高阶方程 (4)(4)二阶常系数齐次微分方程二阶常系数齐次微分方程 (5)(5)二阶常系数非齐次微分方程二阶常系数非齐次微分方程 (6)(6)全微分方程全微分方程本章主要介绍本章主要介绍一阶一阶常微分方程初值问题的常微分方程初值问题的数值解法数值解法。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法3图形解 xyo简单的微分方程复杂、大型的微分方程解析解 y=f(x)数值解(xi,yi)欧拉方法改进欧拉方法 梯形法龙格-库塔法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法4 初值问题及其初值问题及其数值解数值解的概念的概念1 引言引言常用的一些常用的一些解析解法解析解法：常数常数变易变易法、法、Lapalace变换等变换等分离变量分离变量法、变量法、变量代换代换、一阶一阶常微分方程初值问题：常微分方程初值问题：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法5对于初值问题对于初值问题 ，如果，如果 在下列区域内连续：在下列区域内连续：（解的（解的存在唯一存在唯一性）性）且关于且关于 满足满足Lipschitz条件，即存在常数条件，即存在常数 ，使，使则初值问题则初值问题 存在唯一解，且解是存在唯一解，且解是连续可微连续可微的。的。所谓所谓数值解数值解是指：在解的是指：在解的存在区间存在区间上取一系列点上取一系列点逐个求出逐个求出 的近似值的近似值等距等距节点：节点：步长步长数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法6 初值问题初值问题 的解析解及其数值解的的解析解及其数值解的几何几何意义：意义：初值问题初值问题 的解表示过点的解表示过点 的一条的一条曲线曲线初值问题初值问题 的数值解表示一组的数值解表示一组离散点列离散点列可用可用拟合拟合方法求该组数据方法求该组数据 的的近似曲线近似曲线积分积分曲线曲线数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法7建立微分方程数值解法建立微分方程数值解法,首先要将微分方程离散化首先要将微分方程离散化.一般采用以下几种方法一般采用以下几种方法:(1)(1)用差商近似导数用差商近似导数 建立数值解法的常用方法建立数值解法的常用方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法8(2)(2)用数值积分近似积分用数值积分近似积分实际上是矩形法实际上是矩形法宽宽高高数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法9(3)(3)用用TaylorTaylor多项式近似并可估计误差多项式近似并可估计误差数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法102 简单的数值方法简单的数值方法 Euler方法的基本原理方法的基本原理将将 在点在点 处进行处进行Taylor展开展开略去略去 项：项：然后用然后用 代替代替 ，即得，即得称上述公式为向前称上述公式为向前Euler 公式。公式。一、一、Euler方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法11若将若将 在点在点 处进行处进行Taylor展开展开略去略去 项：项：然后用然后用 代替代替 ，即得，即得称上述公式为向后称上述公式为向后Euler 公式。公式。向后向后Euler 公式为公式为隐式隐式格式，需要利用格式，需要利用迭代法迭代法求解求解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法12 Euler方法的方法的几何意义几何意义数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法13Y=y(x)ab数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法14解：解：向前向前Euler公式：公式：例例1：分别利用向前和向后分别利用向前和向后Euler方法方法求解初值问题求解初值问题的的数值数值解解（取步长为（取步长为 ）向后向后Euler公式：公式：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法15具体计算结果具体计算结果:数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法16数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法17利用数值积分将微分方程离散化利用数值积分将微分方程离散化得梯形公式得梯形公式:解决方法：有的可化为显格式，但有的不行解决方法：有的可化为显格式，但有的不行梯形方法为隐式算法梯形方法为隐式算法二、二、改进的改进的Euler方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法18梯形公式比梯形公式比EulerEuler法精度高一些法精度高一些,但计算量较大但计算量较大 实际计算中只迭代一次，这样建立的预测实际计算中只迭代一次，这样建立的预测校正系统称作改进的校正系统称作改进的Euler公式。公式。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法19数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法20例例解解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法21Euler近似解近似解精确解精确解0 1.0.1 1.10.2 1.191820.3 1.277440.4 1.358210.5 1.435130.6 1.50897y0-1y0.1-1.09545y0.2-1.18322y0.3-1.26491y0.4-1.34164y0.5-1.41421y0.6-1.483240 1.0.1 1.097740.2 1.187570.3 1.271290.4 1.350130.5 1.424990.6 1.49657改进改进Euler近似解近似解结果比较结果比较数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法22三、三、常微分方程数值解法的常微分方程数值解法的稳定性稳定性设一个数值方法以定步长设一个数值方法以定步长 求解实验方程求解实验方程得到线性差分方程的解得到线性差分方程的解 。当时。当时 ，若，若 ，则称该方法对步长为则称该方法对步长为绝对稳定绝对稳定的；否则称为不稳定的。的；否则称为不稳定的。将数值方法应用于实验方程，若对一切将数值方法应用于实验方程，若对一切都是绝对稳定的，则称区域都是绝对稳定的，则称区域 为该方法的为该方法的绝对稳定域绝对稳定域。上述定义表明，若数值方法可使任何一步产生的误差在后面的上述定义表明，若数值方法可使任何一步产生的误差在后面的计算中都能逐步削弱，则该方法为绝对稳定。计算中都能逐步削弱，则该方法为绝对稳定。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法23例如例如，对于向前，对于向前Euler法：法：将其应用于将其应用于实验实验方程方程当当 时，误差将逐步减弱，故此时方法稳定。时，误差将逐步减弱，故此时方法稳定。向前向前Euler法法绝对稳定域绝对稳定域：当当 因有误差变为因有误差变为 时，则有时，则有数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法24四、四、单步方法单步方法的局部误差的局部误差和和阶阶单步法单步法的一般形式的一般形式隐式隐式单步法单步法通常称通常称 为为增量增量函数函数显式显式单步法单步法称称 为某方法在点为某方法在点 的的整体截断整体截断误差误差设设 是是准确准确的，用某种方法计算的，用某种方法计算 时产生的截时产生的截断误差，称为该方法的断误差，称为该方法的局部局部截断截断误差，即误差，即(单步法：在计算单步法：在计算yn+1 时只利用时只利用yn)数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法25其中其中 为自然数，则称该方法是为自然数，则称该方法是 阶的或具有阶的或具有 阶精度。阶精度。如果给定方法的如果给定方法的局部截断局部截断误差为误差为如果一个如果一个 阶单步方法的阶单步方法的局部局部截断截断误差为误差为则称则称 为该方法的局部截断误差的为该方法的局部截断误差的主项主项。如向前如向前Euler方法的方法的局部局部截断截断误差误差一阶一阶方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法26 Euler方法的误差分析方法的误差分析对初值问题中的微分方程两端在区间对初值问题中的微分方程两端在区间 上积分上积分如果用如果用左矩形公式左矩形公式计算右端积分，并令计算右端积分，并令其中其中上述等式中如果用上述等式中如果用 代替代替 ，即得向前，即得向前Euler格式。格式。其其局部截断局部截断误差为误差为数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法27设设 关于关于 和和 均满足均满足Lipschitz条件，即条件，即和和数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法28其中其中而而整体截断整体截断误差为误差为数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法29 注意到注意到数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法30对于初值问题对于初值问题 ，如果，如果 关于关于 满足满足（向前向前Euler方法的方法的整体截断整体截断误差误差）Lipschitz条件，条件，为对应的为对应的Lipschitz常数，当常数，当时，向前时，向前Euler方法的数值解方法的数值解 一致收敛一致收敛于初值问题于初值问题 的精确解，且的精确解，且整体截断整体截断误差满足估计式误差满足估计式如果如果 ，Euler方法的方法的整体整体截断截断误差为误差为 数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法31一、一、Runge-Kutta方法的基本思想方法的基本思想3 龙格龙格-库塔（库塔（Runge-Kutta）方法）方法显式显式单步法的一般形式：单步法的一般形式：R-K方法是利用一些点的线性组合方法是利用一些点的线性组合构造构造增量函数增量函数，使得相应方法的使得相应方法的局部局部截断误差的截断误差的阶数阶数尽可能尽可能高高。二阶二阶Runge-Kutta方法方法确定参数确定参数 ，使得，使得与与 在点在点 的的Taylor展开式有尽可能多的展开式有尽可能多的相同项相同项。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法32比较两式的比较两式的相同项相同项得得方程组有无穷多解方程组有无穷多解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法33若取其一组解若取其一组解则得到则得到改进改进的的Euler公式（公式（二阶二阶方法）方法）若取其另一组解若取其另一组解则得到则得到二阶二阶的的Heun（休恩）公式。（休恩）公式。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法34二、显式二、显式Runge-Kutta方法及其稳定性方法及其稳定性和和设设 是一个正整数，代表使用函数值是一个正整数，代表使用函数值 的个数，的个数，是一些特定是一些特定的权因子（均为的权因子（均为实数），则称下列方法（公式）实数），则称下列方法（公式）为初值问题为初值问题 的的m级级显式显式RungeKutta公式，公式，其中其中数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法35类似前面的处理方法，可以得到类似前面的处理方法，可以得到四级四级方法：方法：m=4局部截断局部截断误差误差最常用的一种最常用的一种四阶四阶方法：方法：经典显式经典显式Runge-Kutta公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法36解：解：例例2：用用经典的四阶经典的四阶Runge-Kutta方法求解下列初值问题方法求解下列初值问题 。经典的四阶经典的四阶Runge-Kutta公式：公式：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法37数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法38数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法39注注：对于显式对于显式N级级R-K方法，最多只能得到方法，最多只能得到N阶方法。阶方法。上述方法的缺陷：计算非常复杂。上述方法的缺陷：计算非常复杂。可通过积分方法确定参数。可通过积分方法确定参数。例例2：确定如下三级三阶显式确定如下三级三阶显式Runge-Kutta公式中的参数：公式中的参数：解：解：对微分方程对微分方程 两边积分两边积分得得数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法40采用采用Simpson公式计算上式右端积分项公式计算上式右端积分项可设参数可设参数则有则有选择选择剩余参数剩余参数，使得，使得数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法41取取数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法42数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法43取取利用利用Taylor展开式展开式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法44代入代入当当 时，时，数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法45例例3：求求经典四阶的经典四阶的R-K方法的方法的绝对稳定域。绝对稳定域。解：解：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法46其其绝对稳定域绝对稳定域为为三、隐式三、隐式Runge-Kutta方法方法m级级隐式隐式RK方法的一般形式方法的一般形式其中系数的确定方法同其中系数的确定方法同显式显式RK方法完全类似方法完全类似数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法47(1)(1)一一级级二阶二阶的隐式的隐式中点中点方法：方法：(2)(2)二二级级四阶四阶的隐式的隐式R-K方法：方法：N级隐式级隐式R-K法法可以达到可以达到2N阶阶缺陷：需要求解缺陷：需要求解非线性方程（组）非线性方程（组）数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法48一、一、k步线性步线性多步法多步法4 线性线性多步法多步法/*Linear Mutistep Method and Predictor-Corrector Format*/所谓的线性所谓的线性多步法多步法，指的是某一步解的公式不仅，指的是某一步解的公式不仅与前一步的值有关，而且与前面若干步解的值有关的与前一步的值有关，而且与前面若干步解的值有关的方法。方法。对初值问题对初值问题 两边两边积分积分得得数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法49将将 换为节换为节点点取节点取节点 ，构造，构造 的的k+1个点的个点的Lagrange插值多项式：插值多项式：多步多步显式显式公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法50其中其中记记若函数值若函数值 已知，则得已知，则得r+1步步显式显式方法方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法51如如 时，可得二步时，可得二步显式显式阿达姆斯（阿达姆斯（Adams）格式）格式其中其中数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法52 Adams显式公式的显式公式的局部局部截断误差：截断误差：由由Lagrange插值余项知插值余项知其中其中（第（第二二积分积分中值中值定理）定理）k阶方法阶方法数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法53取节点取节点 ，构造，构造 的的k+1个点的个点的Lagrange插值多项式：插值多项式：多步多步隐式隐式公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法54其中其中记记则得到则得到r+1步步q+1阶的阶的隐式隐式方法方法如如 时，可得二步时，可得二步隐式隐式阿达姆斯（阿达姆斯（Adams）格式）格式梯形梯形公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法55常用的一种常用的一种预测预测-校正校正公式：公式：四阶四阶Adams预测预测-校正公式：校正公式：（显式显式）（隐式隐式）初始迭代值由初始迭代值由4阶阶R-K方法计算方法计算数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法56例例4：用用Adams预测预测-校正公式校正公式求解下列初值问题求解下列初值问题 。解：解：Adams预测预测-校正公式：校正公式：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法57R-K方法方法Adams预预-校法校法 精确精确解解 0 11.0000000000 0.11.0954461.0954451153 0.21.1832171.1832159566 0.31.2649121.2649110640 0.41.34164135711.3416407864 0.51.41421383341.4142135623 0.61.48323982421.4832396974 0.71.54919338041.5491933384 0.81.61245153641.6124515496 0.91.67331999931.6733200530 1.01.73205071981.7320508075数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法58 5 一阶方程组一阶方程组与与高阶高阶方程的数值解法方程的数值解法一、一、一阶微分方程组初值问题的一阶微分方程组初值问题的一般形式一般形式初始条件初始条件:数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法59写成写成向量向量的的形式形式:数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法60n=2对应的对应的Runge-Kutta公式公式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法61数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法62例例 考虑考虑LorenzLorenz模型：模型：其中参数b=8/3，c=10，r=28解：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法63图中，图中，x x1 1的图形为实线的图形为实线(蓝），蓝），x x2 2的图形为的图形为“*”“*”线（绿），线（绿），x x3 3的图形为的图形为“+”“+”线（红）。线（红）。MatlabMatlab软件计算数值解，计算结果如下图软件计算数值解，计算结果如下图若自变量区间取若自变量区间取0，10，计算结果如下：，计算结果如下：右图的坐标为右图的坐标为x x1 1，x x2 2，x x3 3。左图的坐标为左图的坐标为t t，x x1 1或或x x2 2或或x x3 3。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法64曲线呈震荡发散状三维图形的混沌状若自变量区间取若自变量区间取0，20、0，40，计算结果如下：，计算结果如下：数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法65观察结果：观察结果：1、该曲线包含两个“圆盘”，每一个都是由螺线形轨道构成。某些轨道几乎是垂直地离开圆盘中一个而进入另一个。2、随着t的增加，x(t)先绕一个圆盘几圈，然后“跳”到另一个圆盘中。绕第二个圆盘几圈，又跳回原来的圆盘。并以这样的方式继续下去，在每个圆盘上绕的圈数是随机的。数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法66作下列作下列变量代换变量代换可将其化为可将其化为一阶方程组一阶方程组的初值问题的初值问题:二、二、高阶高阶微分方程初值问题的微分方程初值问题的一般形式一般形式数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法67可类似可类似前述方法前述方法构造构造相应的数值方法求解相应的数值方法求解数值分析数值分析第第9 9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法68希望大家考出好成绩!预祝大家研究生涯愉快！谢谢