常微分方程数值解法简介.ppt





《常微分方程数值解法简介.ppt》由会员分享,可在线阅读,更多相关《常微分方程数值解法简介.ppt(93页珍藏版)》请在淘文阁 - 分享文档赚钱的网站上搜索。
1、科学计算与数学建模科学计算与数学建模中南大学数学科学与计算技术学院第第77章章 常微分方常微分方程数值解法简介程数值解法简介 简单的数值方法与基本概念 2数值线性多步法3非线性高阶单步法Runge-Kutte法45一阶方程组和高阶方程的初值问题6常微分方程边值问题的数值解法1 实际问题的微分方程模型 第七章常微分方程数值解法简介 微分方程在科学和工程技术中有很广泛的应用。许多实际问题的数学模型都可以用微分方程来描述,归结为常微分方程的定解问题。很多偏微分方程问题,也可以化为常微分方程问题来近似求解,但是求出所需的解绝非易事。实际上,除了极特殊情形外,人们不可能求出微分方程的解析解,只能用各种近
2、似方法得到满足一定精度的近似解。在常微分方程中已经熟悉了级数解法和Picard逐步逼近法,这些方法可以给出解的近似表达式,称为近似解析方法。另一类方法只给出解在一些离散点上的值,称为数值方法。数值方法应用范围更广,特别适合用计算机计算,本章主要介绍常用的常微分方程数值解法。函数是事物的内部联系在数量方面的反映,如何寻找变量之间的函数关系,在实际应用中具有重要意义.在许多实际问题中,往往不能直接找出变量之间的函数关系,但是有时却容易找出变量的改变量之间的关系,从而建立描述问题的微分方程模型.1 实际问题的微分方程模型 例7.1 将初始温度 的一碗汤放置于环境温度 保持在的桌上,10 分钟后测得汤
3、的温度为1000C。如果汤的温度低于550C才可以喝,问再过20分钟后这碗汤能喝了吗?设物体在 时刻的温度为,从,温度从,注意到热量总是从温度高的物体向温度低的物体传导,因而,所以温度差 恒正,又因物体将随时间而逐渐冷却,则温度的改变量为 两边除以,并令 得温度变化速度为解:为了解决这一问题,需要了解有关热力学的一些基本规律.热量总是从温度高的物体向温度低的物体传导的;在一定的温度范围内,一个物体的温度变化速度与这个物体的温度和其所在介质温度的差值成正比。其中 是比例常数.从而得出描述物体冷却过程的微分方程模型为 容易求出这个一阶微分方程初值问题的解为 根据问题所给的条件知,当时,得到将,代入
4、,得(7.1.1)(7.1.2)从而得到这碗汤的温度随时间变化的函数关系为 于是,将 代入计算得到再过20 min汤的温度,这说明再过20 min后这碗汤能喝了.不过,并不是所有的微分方程模型都可求出解析解。例如,看似简单的微分方程,自德国数学家Wilhelmvon Leibniz提出100多年后才被法国数学家Joseph Liouville证明它没有解析解,只能借助于数值的方法求数值解.(7.1.3)例7.2 某地区发现一种有免疫性的传染病,为了控制疫情扩散对该地人 群进行隔离处理.为了分析受感染人数的变化规律,需要建立描述传染 病传播过程的数学模型.解 设该地区的总人数为常数,任意时刻病人
5、、健康人和病人治愈后 移出感染系统的移出者的比例分别为,病人的日接触率,日治愈率,则容易得出从 时刻,病人和健康人的改变量为 每个方程两边除以,并令,化简后得(7.1.4)其中(对任意的 t).式()就是描述病人和健康人的比例 和 随时间变化的微分方程模型,这是一个微分方程组的初值问题.但是,这一初值问题的解析解是无法求出的,因此不能直接利用 和 的解析式来分析和解决问题。在数学建模课程中学到的大量数学模型都是用微分方程形式给出的,各类微分方程本身和它们的解所具有的特性在常微分方程及数学物理方程中有所解释.虽然,求解微分方程有许多解析方法,但解析方法只能够求解一些特殊类型的方程,在实际应用中人
6、们更关心的是某些特定的自变量在某一个定义范围内的一系列离散点上的近似值.这样一组近似解称为微分方程在该范围内的数值解,寻找微分方程数值解的过程称为微分方程的数值解法。2 简单的数值方法与基本概念 设 在区域 上连续,求 满足 其中 是已知常数,这就是一阶常微分方程的初值问题.为使问题()的解存在、唯一且连续依赖初值,即初值问题()适定,还必须对右端项 加以适当限制,通常要求 关于 是已知函数,且满足Lipschitz条件,即存在常数L,使7.2.1 常微分方程初值问题(7.2.1)(7.2.2)对所有 及 成立。本章总假定满足条件()。1.Euler方法的导出与几何意义 最简单的数值解法是Eu
7、ler法。将区间 作N等分,小区间的长度 称为步长,点列 称为节点,。由已知初值,可算出 在 的导数。7.2.2 Euler法及改进的Euler法 其中,并略去二阶小量,得 下面用3种方法导出Euler法.本章用 表示函数 在 点 的精确值,表示 的近似值。就是 的近似值。利用 可算出,如此下去可算 出 在所有节点上的值,它的一般递推公式为()1)幂级数展开法 利用Taylor展式()这就是Euler法。实际上,初值问题()的解是xy平面上过点 的一条积分曲线。按Euler法,过初始点 作经过此点的积分曲线的切线(斜率为),沿切线取点(按式()计算)作为 的近似.然后,过 作经过此点的积分曲线
8、的切线(斜率为),沿切线取点(按式()计算)作为 的近似.如此下去,即得一条以 为顶点的折线,这就是用Euler法得到的近似积分曲线,如图7-1所示.从几何图形上看,越小,此折线逼近积分曲线越好,因此也称Euler法为Euler折线法。Euler法有明显的几何意义:2)数值微分法利用向前差商近似导数从而得出Euler法的一般递推公式3)数值积分法将初值问题()写成等价的积分形式:取,得用左矩形公式作为右端积分的近似,并用 替代,即得从而也可得出Euler法的一般递推公式为(7.2.5)2.改进的Euler方法 由Euler方法的数值积分导出法可知只要给出式()右端定积分的一种近似计算方法,就可
9、得出初值问题()的一种数值求解方法。如果用右矩形公式作为式()右端积分的近似,则可得从而也可得出一般递推公式为 称式()为后退Euler法。(7.2.6)显然,改进的Euler法比Euler法精度更高。后退Euler法和改进的Euler法,由于未知数 同时出现在等式的两边,故称为隐式算法;如果未知数 由已知量直接计算(即不出现在等式右端),则称为显式算法。对于隐式算法,每步计算需要解关于 的方程,而这样的方程往往是非线性的,通常将初值取为,用迭代法求解,一般只需迭代几步即可收敛。一般先用显式公式计算一个初值,再用隐式公式迭代求解。如果用梯形公式近似替代式()右端的定积分,则可得 从而得出一般递
10、推公式为 称式()为改进的Euler法。(7.2.7)如果先用显式Euler公式作预测,算出 再将 代入隐式梯形公式的右边作校正,得到从而可得 这种方法称为预估校正法。可以看到它是显示格式,迭代求解过程比隐式公式的简单。后面将看到,它的稳定性高于显式Euler法。如果在区间 上对初值问题()的方程两边积分,则有 并用中矩形求积公式近似替代右端的定积分,则得出一般递推公式为 称式()为Euler中点公式。(7.2.8)和,这样的方法称为双步法(或如果计算 的近似值 时只用到前一节点的值,则从初值样的方法称为单步法;而Euler中点公式计算到前两个节点的值多步法.多步法附加初值才能逐一计算出以后各
11、节点的值。出发可逐一计算出以后各节点的值,这时需要用二步法);计算时需要用到前面多个节点值的方法称为7.2.3 截断误差与算法精度的阶 从Euler法的几何意义得知,由Euler法所得的折线明显偏离了积分曲线,可见此方法非常粗糙(即误差太大)。现在分析一下求解初值问题()的数值方法误差的来源。为使问题简化,不考虑因计算机字长限制引起的摄入误差。在假设第 步计算是精确的前提下(即),第 步计算 的截断误差 称为局部截断误差。若某算法的局部截断误差为,即为 的同阶无穷小,则称该算法有 阶精度。若局部截断误差则称 为误差主项,为误差主项系数。1.Euler法的局部截断误差 由Euler法的一般递推公
12、式和 的Taylor展式,得 所以,Euler法的局部截断误差为,即Euler法为1阶精度算法,其误差主项为。2.后退Euler法的局部截断误差 同理,由后退Euler法的一般递推公式,和 的Taylor展式,得 所以,后退Euler法的局部截断误差也是,即后退Euler法也是1阶精度算法,其误差主项为。3.改进Euler法的局部截断误差 由改进Euler法的一般递推公式可得其局部截断误差为 所以,改进Euler法的局部截断误差为,即改进Euler法是2阶精度算法,其误差主项为。4.整体截断误差 当然人们更关心的是近似解 的误差,即,称为整体截断误差.由 的Taylor展式与Euler法的一般
13、递推公式相减,得 记,因 关于 满足Lipschitz条件,所以存在常数L,使得 以此递推,得 注意到,于是 上式右端依赖于初始误差 和局部截断误差的上界R。对于Euler法,可取(C 是与n 无关的常数),若,即,则 所以,比局部截断误差低一阶。用同样的方法可以证明改进Euler法的整体截断误差的阶为,也比局部截断误差低一阶。5.Euler算法的稳定性 在实际计算中,由于测量误差、舍入误差等因素的影响,初值 往往不能精确给出,其误差将依次传递下去。如果传递误差能够被控制,精确地说,传递误差连续依赖于初始误差,则称算法稳定;否则就说算法不稳定。显然,不稳定的算法是不能用的。下面仅考察Euler
14、法的稳定性。设从初值 和 算出的节点值分别为 和,则满足 两式相减,并令,则从而(因)。这说明 连续依赖初始误差,即Euler法稳定。同样可证改进的Euler法也稳定。例7.3 取步长,分别用Euler法、后退Euler法和中点法 求解初值问题:解 步长,所以各节点(1)因为 利用Euler法的计算公式 可得(2)利用后退Euler法的计算公式可解得的显示表达 于是,可得(3)由中点法的计算公式 可知需要两个初值。在此,我们利用后退欧拉法计算的结果 再依次计算得 而该初值问题的解析解是 用它计算各节点的函数值,可得 将上述3种方法计算的结果同精确值对照,可以看出它们确实都是精确值的近似值,只是
15、误差不一样。Euler法的误差较小,后退Euler法误差偏大,中点法误差最小。3 线性多步法 用Euler法计算节点 的近似值 时只用到前一节点的值,是线性的单步法.为了提高解的精度,需要构造线性多步法,其一般形式为(7.3.1)其中,和 是常数,且,和 不同时为0。按公式(7.3.1)计算 时,要用到前面 个节点的值,因此式(7.3.1)称为多步法(或 步法)。又因为方程(7.3.1)关于 是线性的,所以称为线性多步法。显然,若,则线性多步法(7.3.1)是显式的;若,则线性多步法(7.3.1)是隐式的.用线性多步法进行计算 时,除需要给定 外,还要附加初值,这可以用其他方法计算。由于多步法
16、每计算一步用到的信息更多,因此希望由此构造出精度更高的算法。7.3.1 数值积分法 将微分方程 在 上积分,得(7.3.2)适当选取 个节点,作被积函数 的 次Lagrange插值多项式 并用 近似代替,就可得到形如式(7.3.1)的线性多步法。插值节点的不同取法就可得出不同的线性多步法。1.Adams外插法 取 为节点,构造 的Lagrange插值多项式,则其中是插值余项.将式()代入式(),得 舍去余项并用 代替,即得 由此可见,Adams法的局部截断误差为(7.3.3)(7.3.4)(7.3.5)(7.3.6)下面给出Adams法()的具体形式。假设前 个节点处的函数值已知,即 的近似值
17、 已算出,从而函数值 也已知。这样,就可以利用 个数据点,,构造被积函数 的 次插值多项式其中 是Lagrange插值基函数.从而可得(7.3.6)记,则有 这就是求解初值问题的Adams显式公式。是关于 和 的线性表达式,所以它是线性 步法。在上述Adams显式公式的推导中,选用了,作为 插值节点,但构造的 次插值多项式 是代替区 间 上的未知函数,因此属于“外插”,称为 Adams外插法,也称为Adams-Bashorth法.显然,当 时,Adams外插法就是Euler法。(7.3.6)2.Adams内插法 如果将Adams外插法推导过程中的节点改为,公式(7.3.7)就相应地变为 由于式
18、(7.3.8)右端含有(可能是非线性表达式),所以式(7.3.8)属于隐式公式,称为Adams隐式公式,且插值区间包含了积分区间,因此属于“内插”,称为Adams内插法,也称为Adams-Moulton法。显然,当 时,Adams内插法就是改进的Euler法。表7-1和表7-2分别给出了当 时Adams显式公式和隐式公式的系数值。(7.3.8)表7-1Adams外插法系数值 0 1 2 3 4 513-123-16 555-59 37-91901-2774 2616-1274 2514277-7923 9982-7298 2877-475表7-2Adams内插法系数值 0 1 2 3 4 51
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 微分方程 数值 解法 简介

限制150内