ch常微分方程初值问题的数值解法PPT课件

1、第六章 常微分方程初值问题的数值解法6.1 欧拉方法 6.2 计算公式的误差分析 6.3 龙格库塔方法 6.4 向一阶方程组与高阶方程的推广第1页/共20页问题：) 15(,)(),(000bxxyxyyxfy求解数值求解方法：的近似值的函数值然后逐个求解出节点上生成节点作等距分割对区间)(,),(),(,), 1 , 0(,21000NixyxyxyNxbhniihxxbx第2页/共20页 欧拉公式与改进欧拉公式欧拉公式与改进欧拉公式)(,()()(:)(,()()()(11iiiiiiiiixyxfhxyxyxyxfxyhxyxy得由)26()1, 1 ,0()(),(001Nixyyyx

2、fhyyiiii由此可建立计算公式6.1 欧拉方法欧拉方法这称为欧拉公式欧拉公式例例6.1 以 h=0.1为步长，用欧拉法求常微分方程初值问题第3页/共20页)(,()()(:)(,()()()(,1111111iiiiiiiiixyxfhxyxyxyxfxyhxyxy得由类似地)36()1, 1 ,0()(),(00111Nixyyyxfhyyiiii由此可建立另一公式后退欧拉公式是一个隐式公式，通常采用迭代法求解。这称为后退欧拉公式后退欧拉公式第4页/共20页 梯形公式与改进欧拉公式梯形公式与改进欧拉公式欧拉公式与后退欧拉公式也可采用积分近似的方法推出)46()(,()()()(,()()

3、()(),(11111iiiiiixxiixxxxiidxxyxfxyxydxxyxfdxxyxyxyyxfy化为积分方程问题从而把微分方程问题转得对微分方程积分从而可导出欧拉公式得由)(,()()()(,()(,()(,(111iiiiiixxiixxxyxfhxyxyxyxfhdxxyxfdxxyxfiiii第5页/共20页式从而可导出后退欧拉公得由)(,()()()(,()(,()(,(111111111iiiiiixxiixxxyxfhxyxyxyxfhdxxyxfdxxyxfiiii)56()(),(),(2)(,()(,(2)()()(,()(,(2)(,(00111111111x

4、yyyxfyxfhyyxyxfxyxfhxyxyxyxfxyxfhdxxyxfiiiiiiiiiiiiiiiixxii从而导出梯形公式得由梯形公式梯形公式也是隐式单步法公式第6页/共20页11)(1)1(1)(1)(,|,iikikikiyxyyyyk的近似值值作为函数取时使得迭代计算到某个用梯形公式计算时，通常取欧拉公式的解作为迭代初值进行迭代计算，即采用下式), 1 ,0(),(),(2),()(11)1(1)0(1kyxfyxfhyyyxfhyykiiiiikiiiii) 1, 1 , 0(),(),(2),(,11111)1(1Niyxfyxfhyyyxfhyyyyiiiiiiiiii

5、ii则得到显式计算公式代替取直接简单地这称为改进欧拉公式改进欧拉公式第7页/共20页例例 仍取步长h ，采用改进欧拉法重新计算例 6.1 的常微分方程初值问题。这时改进欧拉公式为计算结果见表6-2（书125页）解解第8页/共20页6.2 计算公式的误差分析计算公式的误差分析 定义定义 若 yi+1 是 yi=y(xi) 从计算得到的近似解，则称y(xi+1) yi+1为所用公式的局部截断误差局部截断误差，简称为截断误差截断误差。 定理定理6.1 若单步法 yi+1 = yi+h (xi , yi , h) 的局部截断误差为 O (h p+1) ，且增量函数 (x , y , h) 关于 y 满

6、足李普希兹条件，即存在常数 L0，使对 成立不等式yy,| ),(),(|yyLhyxhyx则其整体截断误差 y(xi) yi=O(hp) 第9页/共20页截断误差的估计(基本假设： yi = y( xi ) )设 y(x)C 3 x0 , b , 则 )76()()(2)()()()(321 hOxyhxyhxyhxyxyiiiii（1）对欧拉公式，有)86()()()(2)()()()(,()(232111 hOhOxyhyxyxyhxyxyxfhxyyiiiiiiiii故因此，欧拉公式的局部截断误差为 O (h2)（2）对后退欧拉公式，有)96()()()(2)()(,()(232111

7、11 hOhOxyhyxyxyxfhxyyiiiiiii故因此，后退欧拉公式的局部截断误差为 O (h2)第10页/共20页（3）对梯形公式，注意到其公式可改写为),(),(21111iiiiiiiyxfhyyxhfyy故由式（6-9）和（6-9）得)106()()()(2)()(221),()(),()(21)(33232111111 hOhOxyhhOxyhyxfhyxyyxfhyxyyxyiiiiiiiiiiii因此，梯形公式的局部截断误差为 O ( h3 )第11页/共20页（4）对改进欧拉公式，有),(),(2)()(1111iiiiiiiiyxfyxfhyxyhxyy而由 ，故有)

8、,(),(),(yxfyyxfyyxfyyx 得)()(2)()()()()()()(,()()(,(),(32122)(,(11hOxyhxyhxyyhOxyhxyhOfyhfhxyxfxyhxyhxfyxfiiiiiixyxyxiiiiiiiii 所以与式（6-7）比较得 y(xi+1) yi+1 = O ( h3 ) 因此，改进欧拉公式的局部截断误差为 O ( h3 ) 第12页/共20页 定义定义 若一种求解常微分方程初值问题的数值计算方法的局部截断误差为 O ( hp+1 ) ，则称该方法为 p阶精度阶精度，或称该方法为 p阶方法阶方法。 由此定义知，欧拉方法与后退欧拉方法为一阶精度

9、，梯形法与改进欧拉方法为二阶精度。第13页/共20页6.3 龙格龙格-库塔方法库塔方法由中值定理，有),(,)(,()()()()(111iiiiiixxyfhyxxxyxy 因此，以上介绍的各种单步法本质上都是对平均斜率 f( , y( ) 进行近似，龙格-库塔据之提出了适当选取若干点上的斜率值作近似以构造高精度计算公式的方法，其基本思想是基于泰勒展式的待定系数法。第14页/共20页 二阶二阶R-KR-K公式公式问题问题：建立二阶精度的计算格式形为在 y(xi) = yi 的假设下，有)()()()(,()(,()()(,(2)(,(2)(,(1121hOfybfahxyhOfbhKfahx

10、yxfbhKxyahxfKxyxyxfKiiiixyxyxixyxyxiiiiiii)()()()(3)(,(222211hOfybfahxyhxyyiixyxyxiii故解解第15页/共20页)()()()(3)(,(222211hOfybfahxyhxyyiixyxyxiii而)(2121)()()()()2/()()()(3)(,(2321hOfyfhxyhxyhOxyhxyhxyxyiixyxyxiiiiii 根据格式为二阶精度，即 y(xi+1) yi+1 = O(h3) 比较两式系数得)136(2/12/112221ba 系数满足(6-13)的形为(6-12)计算格式统称为二阶R-

11、K公式。当令1=1/2时，解得 2=1/2 ，a=b=1，即为改进欧拉公式。若令 1=0，解得 2=1，a=b=1/2，则得另一计算公式)146()2/, 2/(),(12121hKyhxfKyxfKhKyyiiiiii变形欧拉公式变形欧拉公式第16页/共20页 四阶四阶 R-K R-K 公式公式每一步需计算的 f 值的个数1234567n8精度阶1234456n-2 1965年，Butcher研究发现显式R-K公式的精度与需要组合的斜率值的个数具有如下关系 可见，超过四阶精度的R-K公式效率并不高，实际计算通常选用如下四阶格式)156(),()2/, 2/()2/, 2/(),()22(63

12、42312143211hKyhxfKhKyhxfKhKyhxfKyxfKKKKKhyyiiiiiiiiii经典经典R-KR-K公式公式第17页/共20页这时经典R-K公式为 例例 取步长h ，采用经典R-K法计算例 6.1 的常微分方程初值问题。)()()2()2()2()2()22(63)(42)2/(31)2/(2143211hKyehxKKhyehxKKhyehxKyexKKKKKhyyihxiihxiihxiixiiiiiii 取 h=0.2 计算得到表6-4(书133页）。 与例和例比较可见，用经典R-K法计算得到的解比用欧拉法和改进欧拉法所得到的解精确得多。解解第18页/共20页 步长的自动选择步长的自动选择 对于 p 阶精度的计算格式，当取步长为 h 时，记 为从 y(xi) 计算得到的 y (