1、数值分析数值分析数值分析数值分析第九章第九章 常微分方程数值解常微分方程数值解第一节第一节 求解初值问题数值方法的基本原理求解初值问题数值方法的基本原理第二节第二节 高精度的单步法高精度的单步法 第三节第三节 线性多步法线性多步法第四节第四节 一阶微分方程组的解法一阶微分方程组的解法第五节第五节 边值问题的打靶法和差分法边值问题的打靶法和差分法数值分析数值分析数值分析数值分析考虑一阶常微分方程的初值问题考虑一阶常微分方程的初值问题/*Initial-Value Problem*/:0)(,),(yaybaxyxfdxdy只要只要 f(x,y)在在a,b R1 上连续，且关于上连续，且关于 y
2、满足满足 Lipschitz 条条件件，即存在与，即存在与 x,y 无关的常数无关的常数 L 使使对任意定义在对任意定义在 a,b 上的上的 y1(x)和和 y2(x)都成立，则上述都成立，则上述IVP存存在唯一解。在唯一解。|),(),(|2121yyLyxfyxf 要计算出解函数要计算出解函数 y(x)在一系列节点在一系列节点 a=x0 x10,使得使得),(hyxhQyynnnn 111()()pO hp (,)nnQ xyh(,)(,)Q x y hQ x y hL yy对一切对一切 成立成立,则该方法收敛则该方法收敛,且有且有 yy和和)(pnhOe 由该定理可知整体截断误差总比局部
3、截断误差低一阶由该定理可知整体截断误差总比局部截断误差低一阶 对 Euler 方法，(,)(,)Q x y hf x y,当(,)f x y 关于y满足 Lipschitz 条件时是收敛的。数值分析数值分析数值分析数值分析对改进的对改进的Euler法法,),(,(),(),(yxhfyhxfyxfhyxQ 21于是有于是有 1(,)(,)(,)(,)21(,(,)(,(,)2Q x y hQ x y hf x yf x yf xh yhf x yf xh yhf x y 设设L为为f关于关于y的的Lipschitz常数常数,则由上式可得则由上式可得(,)(,)(1)2hQ x y hQ x y
4、 hLL yy限定限定h即可知即可知Q满足满足Lipschitz条件条件,故而改进的故而改进的Euler法法收敛收敛.数值分析数值分析数值分析数值分析3.稳定性稳定性设在计算ny时有误差n，即实际计算得第n步的 解为nnnyy，nnnyy。则第n+1 步的 解为111nnnyy，111nnnyy。如果有1/1nn，说明计算中的舍入误差可 以得到控制，数值方法是稳定的，否则是不稳定的。一般分析时为简单起见，只考虑一般分析时为简单起见，只考虑模型方程模型方程yy 常数，可以常数，可以是复数是复数数值分析数值分析数值分析数值分析一般分析时为简单起见，只考虑一般分析时为简单起见，只考虑模型模型方程方程
5、 yy 当步长取为当步长取为 h 时，将某算法应用于上式，并假时，将某算法应用于上式，并假设只在初值产生误差设只在初值产生误差 ，则若此误差以后，则若此误差以后逐步衰减，就称该算法相对于逐步衰减，就称该算法相对于 绝对稳定绝对稳定，的全体构成的全体构成绝对稳定区域绝对稳定区域。我们称。我们称方法方法A 比方法比方法B 稳定稳定，就是指，就是指 A 的绝对稳定区域比的绝对稳定区域比 B 的的大大。000yy h h h数值分析数值分析数值分析数值分析例：例：考察显式欧拉法的稳定性考察显式欧拉法的稳定性 1(1)nnnnyyh yh y 000yy 11110(1)()(1)(1)nnnnnnny
6、yhyyhh|1|1h 0-1-2ReImg例：例：考察梯形的稳定性考察梯形的稳定性 可见绝对稳定条件是：可见绝对稳定条件是：1(1)nnyh y 显式欧拉法的稳定性条件是显式欧拉法的稳定性条件是11()2nnnnhyyyy112112nnhh11212nnhyyh数值分析数值分析数值分析数值分析111nnyyh 111nnh 可见绝对稳定区域为：可见绝对稳定区域为：|1|1h210ReImg注：注：一般来说，隐式欧拉法的绝对稳定性比同阶一般来说，隐式欧拉法的绝对稳定性比同阶的显式法好。的显式法好。解：解：111(,9 4)(0,1,)nnnnyyhf xyn 讨讨论论隐隐式式欧欧拉拉法法 的
7、的绝绝对对稳稳定定区区和和局局例例部部截截断断误误差差。用于模型 yy，得到11nnnyyhy 数值分析数值分析数值分析数值分析考察其局部截断误差。在点1nx以h为增量作 Taylor 展开 211()()()()2nnnhy xy xhy xy 11121111()(,()(,)()2nnnnnnnEy xyhh f xy xf xyy111(,)(0,1,)nnnnyyhf xyn得到 2111()()(,()()2nnnnhy xy xhf xy xy 两式相减，得两式相减，得数值分析数值分析数值分析数值分析利用(,)f x y关于y满足李普希茨条件 111111(,()(,)()nnn
8、nnnf xy xf xyL y xy 1121111211()(,()(,)()2()()2nnnnnnnny xyhh f xy xf xyyhhL y xyy211()()2(1)nnhy xyyhL211()()nny xyO h隐式隐式欧拉公式是一阶方法欧拉公式是一阶方法数值分析数值分析数值分析数值分析例：例：对于常微分方程初值问题对于常微分方程初值问题证明证明隐式隐式欧拉公式欧拉公式是一阶方法。是一阶方法。00(,)()yf x yy xy 1112311()()()()()(,)2 nnnnnnnnnEy xyhy xhy xyxO hyhf xy解：解：2111(,)()()(
9、)()nnnnnf xyy xy xhyxO h111(,)0,1,2,.nnnnyyhf xyn1112322()()()()()()()()2 nnnnnnnnnEy xyhy xhy xyxOhyh y xh y xOh2()Oh隐式隐式欧拉公式是一阶方法欧拉公式是一阶方法数值分析数值分析数值分析数值分析第二节第二节 高精度的单步法高精度的单步法 在高精度的单步法中在高精度的单步法中,应用最广泛的是应用最广泛的是Runge-Runge-Kutta(Kutta(龙格龙格-库塔库塔)方法方法一一、基本原理基本原理由微分中值定理，在区间1,nnx x上存在01使 1()()()nnny xy
10、xy xhh其中，1,()(,)nnhxxy xf x y 于是，1()()()(,()nnnnny xy xhy xhhf xh y xh 令(,()nnkf xh y xh表示1,nnx x上函数()y x的 平均斜率。在1,nnx x上用不同的方法计算出平均斜率，就可得到求解常微初值问题（9-1）式的不同数值方法。数值分析数值分析数值分析数值分析111112121 (,)11()22(,)(,)nnnnnnnnnnEulerEuleryyhkEulerkfxyyyhkkEulerkfxykfxh yhk 如如果果将将公公式式与与改改进进公公式式写写成成下下列列形形式式：公公式式 改改进进
11、公公式式11(,)()(,)(,)nnf x yy xyf x yf x y以上两组公式都使用函数在某些点上的以上两组公式都使用函数在某些点上的值的线性组合来计算的近似值。值的线性组合来计算的近似值。EulerEuler公式：每步计算一次的值，为一阶方法。公式：每步计算一次的值，为一阶方法。改进Euler公式：需计算两次的值，二阶方法。改进Euler公式：需计算两次的值，二阶方法。数值分析数值分析数值分析数值分析Runge-Kutta法的一般形式法的一般形式11111(,)(,)(2,3,)pnniiinniininijjjyyhc kkfxykfxa h yhb kip 其中，ic，ja，j
12、tb都是待定的参数，选取这些参数的 原则是使 R-K 方法式有p阶精度。数值分析数值分析数值分析数值分析二、二阶龙格库塔方法二、二阶龙格库塔方法11111R-K(,)(,)(2,3,),(,)()pnniiinniininijjjiijinnnyyhc kkfxykfxa h yhb kipabcxyTaylory xxTaylor 方方法法的的一一般般形形式式为为其其中中，都都是是参参数数，确确定定它它们们的的原原则则是是使使近近似似公公式式在在处处的的展展开开式式与与在在处处的的展展开开式式的的前前面面项项尽尽可可能能多多地地重重合合。数值分析数值分析数值分析数值分析11122122211
13、2R-K()(,)(,)nnnnnnyyh c kc kkf xykf xa h yhb k 当当p p时时，方方法法的的一一般般形形式式为为 112221122321(,)(,)(,(,)(,)(,)(,)(,)(,)()nnnnnnnnnnnnnnnxnnynnnnxyTayloryyh c f xyc f xa h yhb f xyyh c f xycf xya hfxyhb fxyf xyO h 上式在处的展开式为上式在处的展开式为12232221()(,)(,)(,)(,)()nnnxnnynnnnyccf xyhc a fxybfxyf xyhO h 数值分析数值分析数值分析数值分
14、析123123()()()()()()2(,)(,)(,)(,)()2nnnnnnnnnxnnynnnny xxTaylorhy xy xhyxyxO hyf xyhhfxyfxyf xyO h 在在处处的的展展开开式式为为112232221()(,)(,)(,)(,)()nnnnxnnynnnnyyccf xyhc a fxybfxyf xyhO h 12222211 1/2 1/2 ),ccc ac bO 3 3有有 无无 穷穷 多多 组组 解解，每每 一一 组组 解解 得得 一一近近 似似 公公 式式，局局 部部 截截 断断 误误 差差 均均 为为(h h这这 些些 方方 法法 统统 称
15、称 二二 阶阶 方方 法法。数值分析数值分析数值分析数值分析12222211 1/2 1/2 ),ccc ac bO 3 3有有 无无 穷穷 多多 组组 解解，每每 一一 组组 解解 得得 一一近近 似似 公公 式式，局局 部部 截截 断断 误误 差差 均均 为为(h h这这 些些 方方 法法 统统 称称 二二 阶阶 方方 法法。122211121211,1,2()/2(,)(,)nnnnnnccabE uleryyh kkkfxykfxhyhk 取取此此 为为 改改 进进公公 式式。近近 似似 公公 式式 为为 122211212110,1,2(,)(2,2)nnnnnnccabyyhkkf
16、xykfxhyhk 取取此此 为为 常常 用用 的的 二二 阶阶 公公 式式，称称 为为 中中 点点 公公 式式。数值分析数值分析数值分析数值分析(0)(0)(0)(1)(1)(1)(2)(2)(2)()(1);2,3,jjjjyffffyffxyffyffxyffyffjxy 一般地有一般地有数值分析数值分析数值分析数值分析三、三阶龙格库塔方法三、三阶龙格库塔方法11231213123 (4)6(,)(,)22(,2)nnnnnnnnpRKRKhyykkkkf xyhhkf xykkf xh yhkhk 类类似似地地，对对，即即三三个个点点，通通过过更更复复杂杂的的计计算算，可可导导出出三三
17、阶阶公公式式。常常用用的的三三阶阶公公式式为为：数值分析数值分析数值分析数值分析四、四阶龙格库塔方法四、四阶龙格库塔方法1123412132434 (22)6(,)(,)2 2(,)22(,)nnnnnnnnnnpRKRKhyykkkkkf xyhhkf xykhhkf xykkf xh yhk 对对，即即四四个个点点，可可导导出出四四阶阶公公式式。常常用用的的四四 阶阶公公式式为为：数值分析数值分析数值分析数值分析 h=0.2,x=0 x=12(01);(0)1.xyyxyy 设设取取步步长长从从直直到到用用四四阶阶龙龙格格库库塔塔方方法法求求解解初初值值问问题题例例：数值分析数值分析数值分
18、析数值分析112341211322433(22);62;2;222;222().nnnnnnnnnnnnnnhyykkkkxkyyxhhkykhykxhhkykhykxhkyh kyh k 由由 经经 典典 的的 四四 阶阶 龙龙 格格 库库 塔塔解解公公 式式 得得：数值分析数值分析数值分析数值分析2RK RKRK RK）方方法法的的导导出出基基于于T Ta ay yl lo or r展展开开，故故要要求求所所求求问问题题的的解解具具有有较较高高的的光光滑滑度度。当当解解充充分分光光滑滑时时，四四阶阶方方法法确确实实优优于于改改进进E Eu ul le er r法法。对对一一般般实实际际问问
19、题题，四四阶阶方方法法一一般般可可达达到到精精度度要要求求。如如果果解解的的光光滑滑性性差差，则则用用四四阶阶方方法法解解的的效效果果不不如如改改进进E Eu ul le er r法法。两点说明两点说明:1RKRK46RK5）当当p=1,2,3,4p=1,2,3,4时时，公公式式的的最最高高阶阶数数恰恰好好是是p,p,当当p4p4时时，公公式式的的最最高高阶阶数数不不是是p p，如如p=5p=5时时仍仍为为，p=p=时时公公式式的的最最高高阶阶数数为为。数值分析数值分析数值分析数值分析R-K方法的绝对稳定区域/2211233223443 11 ,22111 )22411 )24 (,)()k(
20、)()()()()nnnnnnnRKhhhhf x yyhhyyyyhkkyyhhkkhyyhhhkk 代代入入公公式式：将将234112341 2)6111 12624(2()()()nnnhyykkkkyhhh 数值分析数值分析数值分析数值分析234111112624()()()nnhhhh 则则234111 112624()()()hhhh 绝绝对对稳稳定定区区域域：2 1-3-2 -1 0 -1 -2 数值分析数值分析数值分析数值分析五、变步长的龙格五、变步长的龙格库塔方法库塔方法()1()51115(2)15(2)11(2)11()11,(),2,2()2,2()1.()16hnhn
21、nnnhnhnnhnnhnnhyy xychhxxhychy xycy xyy xy n n以以经经典典四四阶阶龙龙格格库库塔塔公公式式为为例例。从从节节点点x x 出出发发，以以 为为步步长长求求一一近近似似值值将将步步长长折折半半，即即取取为为步步长长从从跨跨两两步步到到，求求一一近近似似值值每每跨跨一一步步的的截截断断误误差差是是因因此此有有由由上上两两式式 (2)(2)()11111().15hhhnnnny xyyy 数值分析数值分析数值分析数值分析记1()(1)211115nnyy，表示步长分半前后两次计算结果的 绝对偏差。下面讨论用来判断步长的选择是否合适。事先给定误差限要求，通过步长分半后计算出。以下分两种情况来选择步长：若，再将步长分半 一次计算下去。若，反复将步长分半，直到为止。以最后一 次的步长为合适的步长计算下去。数值分析数值分析数值分析数值分析习题九：习题九：P350-5，8
