数值分析第八章常微分方程数值解法解析课件

数值分析

NumericalAnalysis第八章常微分方程数值解法郑州大学研究生课程（2014-2015学年第一学期）

数值分析

NumericalAnalysis第八章郑州大学2/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis第八章常微分方程数值解法

§8.1引言§8.2欧拉(Euler)法§8.3改进欧拉(Euler)方法§8.4单步法的稳定性2/66郑州大学研究生2011-2012学年课程3/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言问题提出

倒葫芦形状容器壁上的刻度问题.对于圆柱形状容器壁上的容积刻度,可以利用圆柱体体积公式其中直径D为常数.由于体积V与相对于容器底部的任意高度H的函数关系明确,因此在容器上可以方便地标出容器刻度。对于几何形状不是规则的容器,比如倒葫芦形状容器壁上如何标出刻度呢?3/66郑州大学研究生2011-2012学年课程数值分4/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言下表是经过测量得到部分容器高度与直径的关系.H00.20.40.60.81.0D00.110.260.561.041.17根据上表的数据,可以拟合出倒葫芦形状容器的图,建立如图所示的坐标轴后,问题即为:如何根据任意高度x标出容器体积V的刻度,由微元思想分析可知4/66郑州大学研究生2011-2012学年课程数值分5/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言其中x表示高度,直径D是高度x的函数,记为D(x),因此得到如下微分方程初值问题只要求解上述方程,就可求出体积V与高度x之间的函数关系,从而可标出容器壁上容积的刻度,但问题是函数D(x)无解析表达式,我们无法求出其解析解.5/66郑州大学研究生2011-2012学年课程数值分6/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

包含自变量x、未知函数y及未知函数的导数或微分的方程称为微分方程。在微分方程中,自变量的个数只有一个,称为常微分方程。自变量的个数为两个或两个以上的微分方程叫偏微分方程。微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。如果未知函数y及其各阶导数都是一次的,则称它是线性的,否则称为非线性的。6/66郑州大学研究生2011-2012学年课程数值7/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis常微分方程(ODEs未知函数是一元函数)

偏微分方程(PDEs未知函数是多元函数)

7/66郑州大学研究生2011-2012学年课程数值分8/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis同一个微分方程,具有不同的初始条件微分方程的定解条件：8/66郑州大学研究生2011-2012学年课程数值分9/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当x=0时,y=1,可得c=1时特解当x=0时,y=1,可得c=-1时特解两边积分通解本章重点讨论一阶常微分方程的初值问题，9/66郑州大学研究生2011-2012学年课程数值分10/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

在高等数学中，对于常微分方程的求解，给出了一些典型方程求解析解的基本方法，如可分离变量法、常系数齐次线性方程的解法、常系数非齐次线性方程的解法等。但能求解的常微分方程仍然是有限的，大多数的常微分方程是不可能给出解析解。10/66郑州大学研究生2011-2012学年课程11/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

待求解的问题：一阶常微分方程的初值问题/*Initial-ValueProblem*/:

解的存在唯一性（“常微分方程”理论）：只要f(x,y)在[a,b]R1上连续，且关于y

满足Lipschitz

条件，即存在与x,y无关的常数L

使则上述IVP存在唯一解。11/66郑州大学研究生2011-2012学年课程数值12/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis解析解法：（常微分方程理论）只能求解极少一类常微分方程；实际中给定的问题不一定是解析表达式，而是函数表，无法用解析解法。数值解法：递推法如何求解12/66郑州大学研究生2011-2012学年课程数值13/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis13/66郑州大学研究生2011-2012学年课程数值14/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis记号：14/66郑州大学研究生2011-2012学年课程数值15/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法推导Euler格式：★Taylor展开法★数值微分★数值积分法对微分方程的离散，可以有多种思路，但最基本的想法是“以直代曲”15/66郑州大学研究生2011-2012学年课程数值16/66

郑州大学研究生2014-2015学年课程数值分析

NumericalAnalysis16/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法(1)Taylor展开法方程初值问题Euler公式设给定等距剖分当步长h充分小时，略去h2项,得16/66郑州大学研究生2014-2015学年课17/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法(2)用差商近似导数差分方程初值问题向前Euler方法17/66郑州大学研究生2011-2012学年课程18/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法若用向后差商近似导数，即向后Euler方法18/66郑州大学研究生2011-2012学年课程数值19/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法（3）用数值积分方法19/66郑州大学研究生2011-2012学年课程数值20/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法若对积分用梯形公式，则得梯形欧拉公式20/66郑州大学研究生2011-2012学年课程21/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法

欧拉（Euler）方法是解初值问题的最简单的数值方法。初值问题的解y=y(x)代表通过点的一条称之为微分方程的积分曲线。积分曲线上每一点的切线的斜率等于函数在这点的值。

21/66郑州大学研究生2011-2012学年课程22/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysisEuler法的求解过程是:从初始点P0(即点(x0,y0))出发,作积分曲线y=y(x)在P0点上切线(其斜率为

),与x=x1直线相交于P1点(即点(x1,y1),得到y1作为y(x1)的近似值,如上图所示。过点(x0,y0),以f(x0,y0)为斜率的切线方程为

当时,得这样就获得了P1点的坐标。

22/66郑州大学研究生2011-2012学年课程23/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis同样,过点P1(x1,y1),作积分曲线y=y(x)的切线交直线x=x2于P2点,切线的斜率直线方程为当时,得23/66郑州大学研究生2011-2012学年课程24/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当时,得由此获得了P2的坐标。重复以上过程,就可获得一系列的点:P1,P1,…,Pn。对已求得点以为斜率作直线

取24/66郑州大学研究生2011-2012学年课程数值25/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

从图形上看,就获得了一条近似于曲线y=y(x)

的折线。这样,从x0逐个算出对应的数值解25/66郑州大学研究生2011-2012学年课程数值26/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法x0x1x2x3y0hhh欧拉折线法26/66郑州大学研究生2011-2012学年课程数值27/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法27/66郑州大学研究生2011-2012学年课程数值28/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例8.2.1

用欧拉法解初值问题

取步长h=0.2,计算过程保留4位小数

解:h=0.2,欧拉迭代格式

当n=0时，已知x0=0，y0=1，有

y(x1)=y(0.2)

y1=0.2×1(4－0×1)＝0.8当n=1时，已知x1=0.2，y1=0.8，有

y(0.4)y2

=0.2×0.8×(4－0.2×0.8)＝0.6144当n=2,时，已知x2=0.4，y2=0.6144，有

y(0.6)y3=0.2×0.6144×(4-0.4×0.6144)=0.461328/66郑州大学研究生2011-2012学年课程数值29/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis解：Euler公式为当h=0.5时29/66郑州大学研究生2011-2012学年课程数值30/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当h=0.25时30/66郑州大学研究生2011-2012学年课程数值31/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis00.50.751.010.25h=0.5h=0.2531/66郑州大学研究生2011-2012学年课程数值32/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性假设第n步是准确的32/66郑州大学研究生2011-2012学年课程数值33/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法局部截断误差称为局部截断误差33/66郑州大学研究生2011-2012学年课程数值34/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性定义若给定方法的局部截断误差满足则称该方法是P阶的，或称为具有P阶精度。34/66郑州大学研究生2011-2012学年课程数值35/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法整体截断误差35/66郑州大学研究生2011-2012学年课程数值36/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性36/66郑州大学研究生2011-2012学年课程数值37/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis由此知，当

37/66郑州大学研究生2011-2012学年课程数值38/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法

注

整体截断误差与局部截断误差的关系：

38/66郑州大学研究生2011-2012学年课程数值39/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式隐式欧拉法或向后欧拉法

/*implicitEulermethodorbackwardEulermethod*/xn+1点向后差商近似导数隐式或后退欧拉公式39/66郑州大学研究生2011-2012学年课程数值40/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式由于未知数yn+1

同时出现在等式的两边，故称为隐式/*implicit*/

欧拉公式，而前者称为显式/*explicit*/欧拉公式。隐式公式不能直接求解，一般需要用Euler显式公式得到初值，然后用Euler隐式公式迭代求解。因此隐式公式较显式公式计算复杂，但稳定性好（后面分析）。

隐式欧拉公式中的未知数yn+1

可通过以下迭代法求解：40/66郑州大学研究生2011-2012学年课程数值41/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

见上图，显然，这种近似也有一定误差，如何估计这种误差y(xn+1)

yn+1

？方法同上，基于Taylor展开估计局部截断误差。但是注意，隐式公式中右边含有f(xn+1

,yn+1),由于yn+1不准确，所以不能直接用y'(xn+1)代替f(xn+1

,yn+1)设已知曲线上一点Pn(xn,yn),过该点作弦线，斜率为(xn+1

,yn+1)点的方向场f(x,y)。若步长h充分小，可用弦线和垂线x=xn+1的交点近似曲线与垂线的交点。几何意义xnxn+1PnPn+1xyy(x)41/66郑州大学研究生2011-2012学年课程数值42/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

隐式欧拉法的局部截断误差：即隐式欧拉公式具有1阶精度。42/66郑州大学研究生2011-2012学年课程数值43/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式比较欧拉显式公式和隐式公式及其局部截断误差显式公式隐式公式43/66郑州大学研究生2011-2012学年课程数值44/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

若将这两种方法进行算术平均，即可消除误差的主要部分/*leadingterm*/而获得更高的精度,称为梯形法梯形公式/*trapezoidformula*/—显、隐式两种算法的平均注：的确有局部截断误差，即梯形公式具有2

阶精度，比欧拉方法有了进步。但注意到该公式是隐式公式，计算时不得不用到迭代法，其迭代收敛性与欧拉公式相似。44/66郑州大学研究生2011-2012学年课程数值45/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例8.2.3

对初值问题

证明用梯形公式求得的近似解为

并证明当步长h0时,yn收敛于精确解证明:解初值问题的梯形公式为∵

∴

整理成显式

反复迭代,得到∵

∴

45/66郑州大学研究生2011-2012学年课程数值46/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis公式局部截断误差精度显隐稳定性步数欧拉显式公式1阶显差单步欧拉隐式公式1阶隐好单步梯形公式2阶隐好单步欧拉法小结46/66郑州大学研究生2011-2012学年课47/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法47/66郑州大学研究生2011-2012学年课48/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法48/66郑州大学研究生2011-2012学年课49/66郑州大学研究生课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法

显式欧拉公式计算工作量小,但精度低。梯形公式虽提高了精度,但为隐式公式,需用迭代法求解,计算工作量大。综合欧拉公式和梯形公式便可得到改进的欧拉公式。

结合已有格式的优点，以得到计算方便、计算量减少且精度保持的数值格式49/66郑州大学研究生课程数值分析Numeri50/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法

先用欧拉公式(8.2)求出一个初步的近似值,称为预测值,它的精度不高,再用梯形公式对它校正一次,即迭代一次,求得yn+1,称为校正值,这种预测-校正方法称为改进的欧拉公式：称为Euler公式与梯形公式的预测—校正系统。50/66郑州大学研究生2011-2012学年课51/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法实际计算时，常改写成以下形式几何解释xnxn+1ABPn+1=(A+B)/2欧拉法梯形法改进欧拉法51/66郑州大学研究生2011-2012学年课52/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysispredictorcorrector52/66郑州大学研究生2011-2012学年课53/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法

可以证明,改进的欧拉公式的精度为二阶。这是一种一步显式格式,它可以表示为嵌套形式。53/66郑州大学研究生2011-2012学年课54/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例8.3.154/66郑州大学研究生2011-2012学年课55/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis55/66郑州大学研究生2011-2012学年课56/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法56/66郑州大学研究生2011-2012学年课57/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis改进欧拉法的算法57/66郑州大学研究生2011-2012学年课58/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.4单步法的稳定性稳定性:误差在以后各步的计算中不会无限制扩大.稳定性在微分方程的数值解法中是一个非常重要的问题。因为微分方程初值问题的数值方法是用差分格式进行计算的，而在差分方程的求解过程中，存在着各种计算误差，这些计算误差如舍入误差等引起的扰动，在传播过程中，可能会大量积累，对计算结果的准确性将产生影响。这就涉及到算法稳定性问题。

58/66郑州大学研究生2011-2012学年课59/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例：考察初值问题在区间[0,0.5]上的解。分别用欧拉显、隐式格式和改进的欧拉格式计算数值解。0.00.10.20.30.40.5精确解改进欧拉法

欧拉隐式欧拉显式

节点xi

1.00002.00004.00008.00001.6000101

3.2000101

1.00002.5000101

6.25001021.56251023.90631039.76561041.00002.50006.25001.56261013.90631019.76561011.00004.97871022.47881031.23411046.14421063.0590107§8.4单步法的稳定性59/66郑州大学研究生2011-2012学年课60/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.4单步法的稳定性60/66郑州大学研究生2011-2012学年课61/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis定义若某算法在计算过程中任一步产生的误差在以后的计算中都逐步衰减，则称该算法是绝对稳定的/*absolutelystable*/。一般分析某算法的稳定性时，为简单起见，只考虑模型方程或试验方程/*testequation*/§8.4单步法的稳定性61/66郑州大学研究生2011-2012学年课62/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis记数值误差为：引进试验方程：

62/66郑州大学研究生2011-2012学年课63/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis向前欧拉公式的稳定性试验方程的欧拉公式若每步计算有舍入误差，则

63/66郑州大学研究生2011-2012学年课64/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当λ为复数时当λ为实数时64/66郑州大学研究生2011-2012学年课65/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis本章小结

本章介绍了常微分方程初值问题的基本数值解法。包括单步法和多步法。单步法主要有欧拉法、改进欧拉法和龙格—库塔方法。多步法是亚当姆斯法。它们都是基于把一个连续的定解问题离散化为一个差分方程来求解，是一种步进式的方法。用多步法求常微分方程的数值解可获得较高的精度。实际应用时，选择合适的算法有一定的难度，既要考虑算法的简易性和计算量，又要考虑截断误差和收敛性、稳定性。65/66郑州大学研究生2011-2012学年课数值分析

NumericalAnalysis第八章常微分方程数值解法郑州大学研究生课程（2014-2015学年第一学期）

数值分析

NumericalAnalysis第八章郑州大学67/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis第八章常微分方程数值解法

§8.1引言§8.2欧拉(Euler)法§8.3改进欧拉(Euler)方法§8.4单步法的稳定性2/66郑州大学研究生2011-2012学年课程68/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言问题提出

倒葫芦形状容器壁上的刻度问题.对于圆柱形状容器壁上的容积刻度,可以利用圆柱体体积公式其中直径D为常数.由于体积V与相对于容器底部的任意高度H的函数关系明确,因此在容器上可以方便地标出容器刻度。对于几何形状不是规则的容器,比如倒葫芦形状容器壁上如何标出刻度呢?3/66郑州大学研究生2011-2012学年课程数值分69/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言下表是经过测量得到部分容器高度与直径的关系.H00.20.40.60.81.0D00.110.260.561.041.17根据上表的数据,可以拟合出倒葫芦形状容器的图,建立如图所示的坐标轴后,问题即为:如何根据任意高度x标出容器体积V的刻度,由微元思想分析可知4/66郑州大学研究生2011-2012学年课程数值分70/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言其中x表示高度,直径D是高度x的函数,记为D(x),因此得到如下微分方程初值问题只要求解上述方程,就可求出体积V与高度x之间的函数关系,从而可标出容器壁上容积的刻度,但问题是函数D(x)无解析表达式,我们无法求出其解析解.5/66郑州大学研究生2011-2012学年课程数值分71/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

包含自变量x、未知函数y及未知函数的导数或微分的方程称为微分方程。在微分方程中,自变量的个数只有一个,称为常微分方程。自变量的个数为两个或两个以上的微分方程叫偏微分方程。微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。如果未知函数y及其各阶导数都是一次的,则称它是线性的,否则称为非线性的。6/66郑州大学研究生2011-2012学年课程数值72/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis常微分方程(ODEs未知函数是一元函数)

偏微分方程(PDEs未知函数是多元函数)

7/66郑州大学研究生2011-2012学年课程数值分73/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis同一个微分方程,具有不同的初始条件微分方程的定解条件：8/66郑州大学研究生2011-2012学年课程数值分74/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当x=0时,y=1,可得c=1时特解当x=0时,y=1,可得c=-1时特解两边积分通解本章重点讨论一阶常微分方程的初值问题，9/66郑州大学研究生2011-2012学年课程数值分75/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

在高等数学中，对于常微分方程的求解，给出了一些典型方程求解析解的基本方法，如可分离变量法、常系数齐次线性方程的解法、常系数非齐次线性方程的解法等。但能求解的常微分方程仍然是有限的，大多数的常微分方程是不可能给出解析解。10/66郑州大学研究生2011-2012学年课程76/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.1引言

待求解的问题：一阶常微分方程的初值问题/*Initial-ValueProblem*/:

解的存在唯一性（“常微分方程”理论）：只要f(x,y)在[a,b]R1上连续，且关于y

满足Lipschitz

条件，即存在与x,y无关的常数L

使则上述IVP存在唯一解。11/66郑州大学研究生2011-2012学年课程数值77/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis解析解法：（常微分方程理论）只能求解极少一类常微分方程；实际中给定的问题不一定是解析表达式，而是函数表，无法用解析解法。数值解法：递推法如何求解12/66郑州大学研究生2011-2012学年课程数值78/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis13/66郑州大学研究生2011-2012学年课程数值79/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis记号：14/66郑州大学研究生2011-2012学年课程数值80/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法推导Euler格式：★Taylor展开法★数值微分★数值积分法对微分方程的离散，可以有多种思路，但最基本的想法是“以直代曲”15/66郑州大学研究生2011-2012学年课程数值81/66

郑州大学研究生2014-2015学年课程数值分析

NumericalAnalysis81/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法(1)Taylor展开法方程初值问题Euler公式设给定等距剖分当步长h充分小时，略去h2项,得16/66郑州大学研究生2014-2015学年课82/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法(2)用差商近似导数差分方程初值问题向前Euler方法17/66郑州大学研究生2011-2012学年课程83/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法若用向后差商近似导数，即向后Euler方法18/66郑州大学研究生2011-2012学年课程数值84/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法（3）用数值积分方法19/66郑州大学研究生2011-2012学年课程数值85/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法若对积分用梯形公式，则得梯形欧拉公式20/66郑州大学研究生2011-2012学年课程86/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法

欧拉（Euler）方法是解初值问题的最简单的数值方法。初值问题的解y=y(x)代表通过点的一条称之为微分方程的积分曲线。积分曲线上每一点的切线的斜率等于函数在这点的值。

21/66郑州大学研究生2011-2012学年课程87/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysisEuler法的求解过程是:从初始点P0(即点(x0,y0))出发,作积分曲线y=y(x)在P0点上切线(其斜率为

),与x=x1直线相交于P1点(即点(x1,y1),得到y1作为y(x1)的近似值,如上图所示。过点(x0,y0),以f(x0,y0)为斜率的切线方程为

当时,得这样就获得了P1点的坐标。

22/66郑州大学研究生2011-2012学年课程88/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis同样,过点P1(x1,y1),作积分曲线y=y(x)的切线交直线x=x2于P2点,切线的斜率直线方程为当时,得23/66郑州大学研究生2011-2012学年课程89/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当时,得由此获得了P2的坐标。重复以上过程,就可获得一系列的点:P1,P1,…,Pn。对已求得点以为斜率作直线

取24/66郑州大学研究生2011-2012学年课程数值90/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

从图形上看,就获得了一条近似于曲线y=y(x)

的折线。这样,从x0逐个算出对应的数值解25/66郑州大学研究生2011-2012学年课程数值91/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法x0x1x2x3y0hhh欧拉折线法26/66郑州大学研究生2011-2012学年课程数值92/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法27/66郑州大学研究生2011-2012学年课程数值93/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例8.2.1

用欧拉法解初值问题

取步长h=0.2,计算过程保留4位小数

解:h=0.2,欧拉迭代格式

当n=0时，已知x0=0，y0=1，有

y(x1)=y(0.2)

y1=0.2×1(4－0×1)＝0.8当n=1时，已知x1=0.2，y1=0.8，有

y(0.4)y2

=0.2×0.8×(4－0.2×0.8)＝0.6144当n=2,时，已知x2=0.4，y2=0.6144，有

y(0.6)y3=0.2×0.6144×(4-0.4×0.6144)=0.461328/66郑州大学研究生2011-2012学年课程数值94/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis解：Euler公式为当h=0.5时29/66郑州大学研究生2011-2012学年课程数值95/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis当h=0.25时30/66郑州大学研究生2011-2012学年课程数值96/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis00.50.751.010.25h=0.5h=0.2531/66郑州大学研究生2011-2012学年课程数值97/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性假设第n步是准确的32/66郑州大学研究生2011-2012学年课程数值98/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法局部截断误差称为局部截断误差33/66郑州大学研究生2011-2012学年课程数值99/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性定义若给定方法的局部截断误差满足则称该方法是P阶的，或称为具有P阶精度。34/66郑州大学研究生2011-2012学年课程数值100/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法整体截断误差35/66郑州大学研究生2011-2012学年课程数值101/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法欧拉方法的收敛性36/66郑州大学研究生2011-2012学年课程数值102/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis由此知，当

37/66郑州大学研究生2011-2012学年课程数值103/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法

注

整体截断误差与局部截断误差的关系：

38/66郑州大学研究生2011-2012学年课程数值104/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式隐式欧拉法或向后欧拉法

/*implicitEulermethodorbackwardEulermethod*/xn+1点向后差商近似导数隐式或后退欧拉公式39/66郑州大学研究生2011-2012学年课程数值105/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式由于未知数yn+1

同时出现在等式的两边，故称为隐式/*implicit*/

欧拉公式，而前者称为显式/*explicit*/欧拉公式。隐式公式不能直接求解，一般需要用Euler显式公式得到初值，然后用Euler隐式公式迭代求解。因此隐式公式较显式公式计算复杂，但稳定性好（后面分析）。

隐式欧拉公式中的未知数yn+1

可通过以下迭代法求解：40/66郑州大学研究生2011-2012学年课程数值106/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

见上图，显然，这种近似也有一定误差，如何估计这种误差y(xn+1)

yn+1

？方法同上，基于Taylor展开估计局部截断误差。但是注意，隐式公式中右边含有f(xn+1

,yn+1),由于yn+1不准确，所以不能直接用y'(xn+1)代替f(xn+1

,yn+1)设已知曲线上一点Pn(xn,yn),过该点作弦线，斜率为(xn+1

,yn+1)点的方向场f(x,y)。若步长h充分小，可用弦线和垂线x=xn+1的交点近似曲线与垂线的交点。几何意义xnxn+1PnPn+1xyy(x)41/66郑州大学研究生2011-2012学年课程数值107/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

隐式欧拉法的局部截断误差：即隐式欧拉公式具有1阶精度。42/66郑州大学研究生2011-2012学年课程数值108/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.2欧拉(Euler)法向后欧拉公式比较欧拉显式公式和隐式公式及其局部截断误差显式公式隐式公式43/66郑州大学研究生2011-2012学年课程数值109/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis

若将这两种方法进行算术平均，即可消除误差的主要部分/*leadingterm*/而获得更高的精度,称为梯形法梯形公式/*trapezoidformula*/—显、隐式两种算法的平均注：的确有局部截断误差，即梯形公式具有2

阶精度，比欧拉方法有了进步。但注意到该公式是隐式公式，计算时不得不用到迭代法，其迭代收敛性与欧拉公式相似。44/66郑州大学研究生2011-2012学年课程数值110/66郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis例8.2.3

对初值问题

证明用梯形公式求得的近似解为

并证明当步长h0时,yn收敛于精确解证明:解初值问题的梯形公式为∵

∴

整理成显式

反复迭代,得到∵

∴

45/66郑州大学研究生2011-2012学年课程数值111/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis公式局部截断误差精度显隐稳定性步数欧拉显式公式1阶显差单步欧拉隐式公式1阶隐好单步梯形公式2阶隐好单步欧拉法小结46/66郑州大学研究生2011-2012学年课112/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法47/66郑州大学研究生2011-2012学年课113/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法48/66郑州大学研究生2011-2012学年课114/66郑州大学研究生课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法

显式欧拉公式计算工作量小,但精度低。梯形公式虽提高了精度,但为隐式公式,需用迭代法求解,计算工作量大。综合欧拉公式和梯形公式便可得到改进的欧拉公式。

结合已有格式的优点，以得到计算方便、计算量减少且精度保持的数值格式49/66郑州大学研究生课程数值分析Numeri115/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法

先用欧拉公式(8.2)求出一个初步的近似值,称为预测值,它的精度不高,再用梯形公式对它校正一次,即迭代一次,求得yn+1,称为校正值,这种预测-校正方法称为改进的欧拉公式：称为Euler公式与梯形公式的预测—校正系统。50/66郑州大学研究生2011-2012学年课116/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysis§8.3改进欧拉(Euler)方法实际计算时，常改写成以下形式几何解释xnxn+1ABPn+1=(A+B)/2欧拉法梯形法改进欧拉法51/66郑州大学研究生2011-2012学年课117/66

郑州大学研究生2011-2012学年课程数值分析

NumericalAnalysispredictorcorrector52/66