计算方法习题PPT课件

1、上页上页下页下页第第9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法 9.1 引言引言 9.2 简单的数值方法简单的数值方法 9.3 龙格龙格- -库塔方法库塔方法上页上页下页下页其中f (x,y)是已知函数，(1.2)是定解条件也称为初值条件。初值条件。本章主要讨论一阶常微方程的初值问题本章主要讨论一阶常微方程的初值问题00( , )(1.1)( ) (1.2)yf x yaxby xy 9.1 引引 言言上页上页下页下页则称则称f(x, y) 关于关于y满足满足利普希茨利普希茨(Lipschitz)条件，条件，L称为称为f(x, y) 的的利普希茨常数利普希茨常数。121212

2、( ,)( ,),f x yf x yL yyyyR 定义：定义： 若存在实数若存在实数L0, 使得使得上页上页下页下页定理定理00( , )|, , ,(1.1)(1.2) , ( )fDx yaxb yRyxa byRxa by x 设设 在在区区域域上上连连续续，关关于于 满满足足利利普普希希茨茨条条件件，则则对对任任意意，常常微微分分方方程程初初值值问问题题式式和和式式当当时时存存在在唯唯一一的的连连续续可可微微解解. .上页上页下页下页 所谓所谓数值解法数值解法, 就是寻求解就是寻求解y(x)在一系列离散节点在一系列离散节点121iixxxx 上的近似值上的近似值 y1, y2, ,

3、 yi , yi+1 ,. 相邻两个节点的间距相邻两个节点的间距hi=xi+1- -xi 称为称为步长步长. 今后如不特别说明，总是假定今后如不特别说明，总是假定 hi=h(i=1,2,)为为常数常数, 这时节点为这时节点为xi=x0+ih(i=0,1,2,) (等距节点等距节点).上页上页下页下页 初值问题的初值问题的数值解法一般数值解法一般采取采取“步进式步进式”，即，即求解过程顺着节点排列的次序一步一步地向前推进求解过程顺着节点排列的次序一步一步地向前推进. 描述这类算法，只要给出用已知信息描述这类算法，只要给出用已知信息yi, yi- -1, yi- -2,计计算算yi+1的递推公式的

4、递推公式. 若计算若计算yi+1时只用到前一点的值时只用到前一点的值yi，称为，称为单步法单步法. 若若计算计算yi+1时时用到用到yi+1前面前面 k 点的值点的值yi, yi-1, yi-k+1，称为，称为k步法步法. 上页上页下页下页9.2 简单的数值方法简单的数值方法上页上页下页下页)(,()()()()(1iiiiiixyxfhxyxyhxyxy 由向前差商公式由向前差商公式hxyxyxyiii)()()(1 一、欧拉法一、欧拉法11 ,(),(),iiiiyyy xy xEuler 用用近近 似似 代代 替替有有公公 式式),(iiiyxhfy 1iy, 2 , 1 , 0 i得得

5、上页上页下页下页每步计算只用到1ny100021111(,)(,) (,)nnnnyyhf xyyyhf xyyyhf xy依上述公式逐次计算可得：依上述公式逐次计算可得：ny上页上页下页下页 几何意义是用一条初始点重合的折线来近似表示曲线。几何意义是用一条初始点重合的折线来近似表示曲线。 y P1 (x1, y1) P2 (x2, y2) y=y(x) 0 x0 x1 x2 x3 x )(1iixxh欧拉法欧拉法的几何意义的几何意义000(,)P xy上页上页下页下页1+1+1+1()()()()(, ()iiiiiiy xy xhy xy xh f xy x 由向后差商公式由向后差商公式1

6、+1()()()iiiy xy xy xh Euler于于是是 有有后后退退公公式式+1+1(,)iiiyhf xy 1iy, 2 , 1 , 0 i得得二、后退欧拉法二、后退欧拉法上页上页下页下页 后退的欧拉公式与欧拉公式有着本质的区别后退的欧拉公式与欧拉公式有着本质的区别, 后后者是关于者是关于yi+1的一个直接计算公式，这类公式称作是的一个直接计算公式，这类公式称作是显式的显式的；前者公式的；前者公式的右端含有未知的右端含有未知的yi+1，它实际上，它实际上是关于是关于yi+1的一个函数方程的一个函数方程, ,这类方程称作是这类方程称作是隐式的隐式的. . 隐式方程通常用迭代法求解，而迭

7、代过程的实隐式方程通常用迭代法求解，而迭代过程的实质是质是逐步逐步显式化显式化. .上页上页下页下页 设用欧拉公式设用欧拉公式(0)1( , )iiiiyyhf x y 给出迭代初值给出迭代初值 ，用它代入，用它代入后退欧拉公式式的后退欧拉公式式的右右端，使之转化为显式，直接计算得端，使之转化为显式，直接计算得(0)1iy (1)(0)111(,),iiiiyyhf xy 然后再用然后再用 代入代入后退欧拉公式的右端，又有后退欧拉公式的右端，又有(1)1iy (2)(1)111(,).iiiiyyhf xy 如此反复进行，得如此反复进行，得(1)( )111(,) (0,1, ).kkiiii

8、yyhf xyk 上页上页下页下页由于由于f(x, y)对对y满足满足Lipschitz条件条件，则，则(1)( )111111(,)(,)kkiiiiiiyyh f xyf xy ( )11.kiihL yy 由此可知，只要由此可知，只要hL1，迭代法，迭代法就收敛到解就收敛到解(1)( )111(,) (0,1,).kkiiiiyyhf xyk +1iy上页上页下页下页三三 单步法的局部截断误差与阶单步法的局部截断误差与阶 初值问题初值问题(1.1),(1.2)的的单步法单步法可用可用一般形式一般形式表示为表示为11( , ),iiiiiyyhx y yh 其中多元函数其中多元函数 与与f

9、(x, y )有关，当有关，当 含有含有yi+1时，方法时，方法是是隐式隐式的，若不含的，若不含yi+1则为则为显式方法显式方法，所以显式单步，所以显式单步法可表示为法可表示为 (x, y, h)称为称为增量函数增量函数，例如对欧拉法有，例如对欧拉法有1( , ),iiiiyyhx y h ( , , )( , ).x y hf x y 上页上页下页下页 定义定义1 设设y(x)是初值问题是初值问题(1.1),(1.2)的的准确解准确解, 称称11()( )( , ( ), ).iiiiiTy xy xhx y xh 为显式单步法的为显式单步法的局部截断误差局部截断误差. . Ti+1之所以称

10、为之所以称为局部的局部的，是假设在，是假设在xi前各步没有前各步没有误差误差. .当当yi= =y(xi)时，计算一步，则有时，计算一步，则有11111()() ( , )()( )( , ( ), ).iiiiiiiiiiiy xyy xyhx y hy xy xhx y xhT 所以，局部截断误差可理解为用显示单步法所以，局部截断误差可理解为用显示单步法计计算一步的误差算一步的误差上页上页下页下页11222()()(, ()= ()+()( )()()2( )().2iiiiiiiiiTy xy xhf xy xhy xhy xyy xhy xhyO h 定义定义2 若显式单步法的若显式单

11、步法的局部截断误差局部截断误差满足满足则称显示单步法具有则称显示单步法具有p阶精度或是阶精度或是p阶方法阶方法. .欧拉法的局部截断误差欧拉法的局部截断误差为为11()piTO h 一一阶方法阶方法上页上页下页下页11112111122()()(, ()() ()()( )()2( )O().2iiiiiiiiiTy xy xhf xy xhy xy xhy xyhy xhyh 以上定义对隐式单步法也是适用的以上定义对隐式单步法也是适用的. .例如，对后例如，对后退欧拉法退欧拉法，其局部截断误差为其局部截断误差为一一阶方法阶方法上页上页下页下页四、四、 梯形公式梯形公式 将将 区间区间 积分积

12、分 1( , ),方程在上iiyf x yx x 11 ( , )(0,1 , )iiiixxxxy dxf x y dxi )(,()(,(2)()(111 iiiiiixyxfxyxfhxyxy1113(,)(,()(,()2()12iixiiiixihfx y dxfxy xfxy xhy 忽略高阶项忽略高阶项, ,有有对右端采用梯形公式对右端采用梯形公式, , 有有上页上页下页下页) 1,., 0(),(),(2111 niyxfyxfhyyiiiiii注：梯形公式注：梯形公式局部截断误差局部截断误差 ， 即梯形公式具有即梯形公式具有2 阶精度，比欧拉公式有了进步。但注阶精度，比欧拉公

13、式有了进步。但注意到该公式是意到该公式是隐式隐式公式公式.31()iTO h 11 ,(),(), iiiiyyy xy x 用用代代 替替有有 梯梯 形形 公公 式式 上页上页下页下页五、五、 改进欧拉公式改进欧拉公式 预测预测: 先用先用显式显式欧拉公式作欧拉公式作预测预测，算出，算出),(1iiiiyxfhyy 校正校正: 再将再将 代入代入隐式隐式梯形公式的右边作梯形公式的右边作校正校正，得到，得到1 iy),(),(2111 iiiiiiyxfyxfhyy改进欧拉公式改进欧拉公式 这是一种显式格式这是一种显式格式, ,它可以表示为嵌套形式它可以表示为嵌套形式。 )1,., 0(),(

14、,),(211 niyxfhyxfyxfhyyiiiiiiii改进欧拉公式改进欧拉公式上页上页下页下页 )(21),(),(11cpipiiciiipyyyyxhfyyyxhfyy注：注：此法亦称为此法亦称为预测预测- -校正法校正法。可以证明该算法具。可以证明该算法具有有 2 2 阶精度。阶精度。或者表示成下列平均化形式或者表示成下列平均化形式 1,.,2 , 1 , 0ni改进欧拉公式改进欧拉公式上页上页下页下页例：用改进欧拉公式求解初值问题例：用改进欧拉公式求解初值问题 要求取步长要求取步长h=0.2，计算，计算y(1.2)及及y(1.4)的近似的近似值，小数点后至少保留值，小数点后至少

15、保留5位位. 解解 设设f(x,y)=- -y- -y2sinx , x0=1, y0=1, xi=x0+ih=1+0.2i, 改进欧拉公式为欧拉公式为 2dsin0d(1)1yyyxxy1111( ,) ( ,)(,)2iiiiiiiiiiyyhf x yhyyf x yf xy上页上页下页下页 于是有于是有 由由y0=1计算得计算得110.631706(1.2)0.715 489yyy220.476965(1.4)0.526112yyy212211110.2(sin)0.1(sinsin)iiiiiiiiiiiiiyyyyxyyyyxyyx上页上页下页下页由微分中值定理，有由微分中值定理，

16、有*1*1()( )(), 01 ()( )iiiiiiiy xy xy xhf xh y xhKhy xy xhK K* *称为曲线称为曲线y( (x) )在区间在区间 xi, , xi+1+1 上的平均斜率，只要知道平上的平均斜率，只要知道平均斜率，就可计算均斜率，就可计算y(xi+1+1).).因此只要对平均斜率提供一种近似因此只要对平均斜率提供一种近似算法，则由算法，则由(1)(1)式可导出一种相应的求解公式。式可导出一种相应的求解公式。(1) 9.3 龙格龙格库塔方法库塔方法上页上页下页下页欧拉公式可改写成欧拉公式可改写成 ),(),(111iiiiiiyxfKyxfKhKyy 取取

17、局部截断误差为局部截断误差为 。 改进的欧拉公式又可改写成改进的欧拉公式又可改写成 )(2hO2111212112121),(),()2121(KKKhKyxfKyxfKKKhyyiiiiii 取取)(3hO111(), (),*,iiiiiiy xyy xy KKyyhK 设设则则局部截断误差为局部截断误差为 。 斜率一定取K1， K2 的平均值吗？步长一定是h 吗？即第二个节点一定是xi+1吗？上页上页下页下页将改进的欧拉公式推广为将改进的欧拉公式推广为),(),(12122111phKyphxfKyxfKKKhyyiiiiii l ll l其局部截断误差为其局部截断误差为 ，)()(31

18、11hOxyyRiii 21,1221 pl ll ll l这里有这里有 3 个未知个未知数，数， 2 个方程。个方程。存在存在无穷多个解无穷多个解。所有满足上式的格式统称为。所有满足上式的格式统称为2阶龙格阶龙格 - 库塔公式库塔公式。其中其中2阶龙格阶龙格 - 库塔格式库塔格式21, 121 l ll lp注意到，注意到， 就是改进的欧拉公式。就是改进的欧拉公式。 问题问题: 为获得更高的精度，应该如何进一步推广？为获得更高的精度，应该如何进一步推广？上页上页下页下页改进的Euler 公式推广为二阶Runge-Kutta公式带来这样的启示：若在若在xi, xi+1上多取几个点的斜率值，然后将它上多取几个点的斜率值，然后将它们的线性组合作为平均斜率的近似值，则有可能们的线性组合作为平均斜率的近似值，则有可能构造出具有更高精度的计算公式。构造出具有更高精度的计算公式。-Runge-Kutta方法的基本思想。上页上页下页下页显式显式N N级龙格级龙格- -库塔公式库塔公式其中其中l li ( i = 1, , N )， i ( i = 2, , N ) 和和 ij ( i = 2, , N; j = 1, , i 1 ) 均为待定系