《最新单步法的收敛与稳定PPT课件.ppt》由会员分享，可在线阅读，更多相关《最新单步法的收敛与稳定PPT课件.ppt（38页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、2 定义定义3 3 若一种数值方法（如单步法（4.1）对于固定的 , 当 时有 ，其中 是（1.1），（1.2）的准确解，则称该方法是收敛收敛的. nhxxn00h)(nnxyy)(xy 显然数值方法收敛是指 ，对单步法（4.1）有下述收敛性定理： 0)(nnnyxye 定理定理1 1 假设单步法（4.1）具有 阶精度，且增量函数 关于 满足利普希茨条件 p),(hyxy,),(),(yyLhyxhyx（4.2）又设初值 是准确的，即 ，则其整体截断误差整体截断误差 0y)(00 xyy).()(pnnhOyxy（4.3）9 9.4.2 9.4.2 绝对稳定性与绝对稳定域绝对稳定性与绝对稳定域
2、 定义定义5 5 若一种数值方法在节点值 上大小为 的扰动，于以后各节点值 上产生的偏差均不超过 ，则称该方法是稳定稳定的. ny)(nmym 以欧拉法为例考察计算稳定性. 例例4 4 考察初值问题 .1)0(,100yyy其准确解 是一个按指数曲线衰减得很快的函数，如图9-3所示. xexy100)( 用欧拉法解方程 得 yy10010.)1001(1nnyhy若取 ，则欧拉公式的具体形式为 025.0h,5.11nnyy计算结果列于表9-4的第2列. 可以看到，欧拉方法的解 （图9-3中用号标出）在准确值 的上下波动，计算过程明显地不稳定. ny)(nxy 但若取 则计算过程稳定. nny
3、yh5.0,005.01图9-311再考察后退的欧拉方法，取 时计算公式为 025.0h.5.311nnyy计算结果列于表9-4的第3列（图9-3中标以号），这时计算过程是稳定的. 0067.00625.5100.00233.0375.3075.00816.025.2050.02857.05.1025.049后退欧拉方法欧拉方法节点计算结果对比表12 例题表明稳定性不但与方法有关，也与步长 的大小有关，当然也与方程中的 有关. h),(yxf 为了只考察数值方法本身. 通常只检验将数值方法用于解模型方程的稳定性，模型方程为 ,yy（4.8）其中 为复数，对一般方程可以通过局部线性化化为这种形式
4、，例如在 的邻域，可展开为 ),(yx,)(,()(,(),(),(yyyxfxxyxfyxfyxfyyx略去高阶项，再做变换即可得到 的形式. uu 对于 个方程的方程组，可线性化为 ，这里mAyy 13 为 的雅可比矩阵 . Amm )(jiyf 若 有 个特征值 ，其中 可能是复数，所以，为了使模型方程结果能推广到方程组，方程（4.8）中 为复数. Amm,21i 为保证微分方程本身的稳定性，还应假定 . 0)Re( 先研究欧拉方法的稳定性. 模型方程 的欧拉公式为 yy.)1(1nnyhy（4.9）设在节点值 上有一扰动值 ，它的传播使节点值 产生大小为 的扰动值，假设用 按欧拉公ny
5、n1ny1nnnnyy*14式得出 的计算过程不再有新的误差，则扰动值满足 11*1nnnyy.)1(1nnh可见扰动值满足原来的差分方程（4.9）. 如果差分方程的解是不增长的，即有 ,1nnyy则它就是稳定的. 显然，为要保证差分方程（4.9）的解是不增长的，只要选取 充分小，使 h.11h（4.10）在 的复平面上，这是以 为圆心，1为半径的单位圆域. 称为欧拉法的绝对稳定域，一般情形可由下面定义. h)0, 1(15 定义定义6 6 单步法（4.1）用于解模型方程（4.8），若得到的解 ，满足 ，则称方法（4.1）是绝对稳定绝对稳定的. nnyhEy)(11)(hE 在 的平面上，使
6、的变量围成的区域，称为绝对稳定域绝对稳定域，它与实轴的交称为绝对稳定区间绝对稳定区间. h1)(hE 对欧拉法 ，其绝对稳定域已由（4.10）hhE 1)(.11h给出，绝对稳定区间为 .02h 在例5中 ，即 为绝对稳定区间. 01002,100h02.00 h 例4中取 故它是不稳定的，当取 时它是稳定的. 025.0h005.0h16 对二阶R-K方法，解模型方程（4.1）可得到 ,2)(121nnyhhy故 .2)(1)(2hhhE绝对稳定域由 得到，于是可得绝对稳定区间为 ，即 . 12)(12hh02h/20 h 类似可得三阶及四阶的R-K方法的 分别为 )(hE,!3)(!2)(
7、1)(32hhhhE17.!4)(!3)(!2)(1)(432hhhhhE由 可得到相应的绝对稳定域. 1)(hE 当 为实数时则得绝对稳定区间. 分别为 三阶显式R-K方法： 即 051.2h./51.20 h 四阶显式R-K方法： 即 078.2h./78.20 h 从以上讨论可知显式的R-K方法的绝对稳定域均为有限域，都对步长 有限制. 如果 不在所给的绝对稳定区间内，方法就不稳定. hh 例例5 5 分别取 及 用经典的四阶R-K方法（3.13）计算. , 1)0(),10(20yxyy1 .0h2 .0h18 解解 本例 分别为 及 ，前者在绝对稳定区间内，后者则不在，用四阶R-K方
8、法计算其误差见下表： h,20240 .31250 .6250 .1250 .2598. 42 . 01017. 01015. 01014. 01012. 01093. 01 . 00 . 18 . 06 . 04 . 02 . 043211hhxn 以上结果看到，如果步长 不满足绝对稳定条件，误差增长很快. h 对隐式单步法，可以同样讨论方法的绝对稳定性，例如对后退欧拉法，用它解模型方程可得 ,111nnyhy19故 .11)(hhE由 可得绝对稳定域为 ，它是以 为圆心，1为半径的单位圆外部，故绝对稳定区间为 . hhE11)(11h)0, 1(0h 当 时，则 ，即对任何步长均为稳定的.
9、 0 h0 对隐式梯形法，它用于解模型方程（4.8）得 ,21211nnyhhy20故 .2/12/1)(hhhE对 有 ，故绝对稳定域为 的左半平面，绝对稳定区间为 ，即 时梯形法均是稳定的. 0)Re(12/12/1)(hhhEh0h h021 9.5 9.5 线性多步法线性多步法 在逐步推进的求解过程中，计算 之前事实上已经求出了一系列的近似值 ，如果充分利用前面多步的信息来预测 ，则可以期望会获得较高的精度. 1nynyyy,101ny 这就是构造所谓线性多步法的基本思想线性多步法的基本思想. 构造多步法的主要途径多步法的主要途径是基于数值积分方法和基于泰勒展开方法，前者可直接由方程（
10、1.1）两端积分后利用插值求积公式得到. 本节主要介绍基于泰勒展开的构造方法. 22 9.5.1 9.5.1 线性多步法的一般公式线性多步法的一般公式 如果计算 时，除用 的值，还用到 的值，则称此方法为线性多步法线性多步法. kny1 knyiyin()2, 1 ,0k 一般的线性多步法公式可表示为 ,010kiinikiiniknfhyy（5.1）其中 为 的近似， 为常数, 及 不全为零，则称（5.1）为线性线性 步法步法. iny)(inxy,),(0ihxxyxffininininii,00k 计算时需先给出前面 个近似值 , 再由（5.1）逐次求出 . k110,kyyy,1kky
11、y23 如果 ，称（5.1）为显式显式 步法步法，这时 可直接由（5.1）算出；如果 ，则（5.1）称为隐隐式式 步法步法，求解时与梯形法（2.7）相同，要用迭代法方可算出 . kny0kk0kkkny （5.1）中系数 及 可根据方法的局部截断误差及阶确定，其定义为： ii 定义定义7 7 设 是初值问题(1.1)，(1.2)的准确解，线性多步法（5.1）在 上的局部截断误差为 )(xyknx. )()()();(1010kiinikiiniknnknxyhxyxyhxyLT（5.2）若 ，则称方法（5.1）是 阶的阶的， 则)(1pknhOT1pp24称方法（5.1）与方程（1.1）是相容
12、的相容的. 由定义7，对 在 处做泰勒展开，由于 knTnx,)(!3)()(!2)()()()(32 nnnnnxyihxyihxyihxyihxy.)(!2)()()()(2 nnnnxyihxyihxyihxy代入（5.2）得 ,)()()()()(2210 npppnnnknxyhcxyhcxyhcxycT（5.3）25其中 ),(11100kc),()1(2101211kkkkc2)!1(1)1(2(!11211121kqqkqqqpkqkkqc.,3,2q（5.4） 若在公式（5.1）中选择系数 及 ，使它满足 ii.0,0110ppcccc26由定义可知此时所构造的多步法是 阶的
13、，且 p).()(2)1(11pnpppknhOxyhcT（5.5）称右端第一项为局部截断误差主项局部截断误差主项， 称为误差常数误差常数. 1pc 根据相容性定义， , 即 ，由（5.4）得 1p010 cc., 1011110kikiikiik（5.6）故方法（5.1）与微分方程（1.1）相容的充分必要条件是（5.6）成立. 27 当 时，若 ，则由（5.6）可求得 1k01.1, 100此时公式（5.1）为 ,1nnnfhyy即为欧拉法. 从（5.4）可求得 ，故方法为1阶精度，且局部截断误差为 02/12c),()(21321hOxyhTnn 这和第2节给出的定义及结果是一致的. 28
14、 对 , 若 ，此时方法为隐式公式，为了确定系数 ，可由 解得 1k01100,0210ccc.2/1, 1100于是得到公式 ),(211nnnffhyy即为梯形法. 由（5.4）可求得 , 故 ，所以梯形法是二阶方法，其局部截断误差主项是 , 这与第2节中的讨论也是一致的. 12/13c2p12/)(3nxyh 对 的多步法公式都可利用（5.4）确定系数 ，并由（5.5）给出局部截断误差. 2kii,29 9.5.2 9.5.2 阿当姆斯显式与隐式公式阿当姆斯显式与隐式公式 考虑形如 kiiniknknfhyy01（5.7）的 步法，称为阿当姆斯阿当姆斯（Adams)方法方法. k 为显式
15、方法， 为隐式方法，通常称为阿当姆斯显式与隐式公式，也称Adams-Bashforth公式与Adams-Monlton公式. 0k0k 这类公式可直接由方程（1.1）两端积分（从 到 积分）求得. 1 knxknx 下面可利用(5.4)由 推出，对比(5.7) 01pcc30kiiniknknfhyy01与（5.1） ,010kiinikiiniknfhyy可知此时系数 . 1,01210kk 显然 成立，下面只需确定系数 ，故可令 ，则可求得 .00ck,10011kcck,10（若 ，则令 来求得 ）. 0k01kcc110,k 以 为例，由 ，根据（5.4）得 3k04321cccc31
16、.65)278(4,19)94(3,5)32(2, 13213213213210 若 ，则由前三个方程解得03,1223,1216,125210得到 的阿当姆斯显式公式是3k).51623(121223nnnnnfffhyy（5.8）由（5.4）求得 ，所以（5.8）是三阶方法，局部8/34c32截断误差是).()(835)4(43hOxyhTnn 若 ，则可解得 03.83,2419,245,2413210于是得 的阿当姆斯隐式公式为 3k),5199(2412323nnnnnnffffhyy（5.9）它是四阶方法，局部截断误差是 ).()(720196)5(53hOxyhTnn（5.10）
17、33 用类似的方法可求得阿当姆斯显式方法和隐式方法的公式，表9-5及表9-6分别列出了 时的阿当姆斯显式公式与阿当姆斯隐式公式，其中 为步数， 为方法的阶， 为误差常数. 4,3,2, 1kkp1pc720251)9375955(244483)51623(1233125)3(22221115912334122311211nnnnnnnnnnnnnnnnnnpffffhyyfffhyyffhyyfhyycpk公式阿当姆斯显式公式表341603)19106264646251(7205472019)5199(2443241)85(1232121)(22169123434123231212111nnn
18、nnnnnnnnnnnnnnnnnnnpfffffhyyffffhyyfffhyyffhyycpk公式阿当姆斯隐式公式表35 例例6 6 用四阶阿当姆斯显式和隐式方法解初值问题 .1)0(, 1yxyy取步长 . 1 .0h 解解 本题 . 从四阶阿当姆斯显式公式得到 nnhxxyfnnnn1.0, 1).24. 324. 09 . 07 . 39 . 55 .18(241)9375955(2412312334nyyyyffffhyynnnnnnnnnn对于四阶阿当姆斯隐式公式得到 36).324. 01 . 05 . 01 .229 . 0(241)5199(2412312323nyyyyffffhyynnnnnnnnnn由此可直接解出 而不用迭代，得到 3ny).324. 01 . 05 . 01 .22(9 .241123nyyyynnnn计算结果见表9-7，其中显式方法中的 及隐式方法中的 均用准确解 计算得到，对一般方程，可用四阶R-K方法计算初始近似. 3210,yyyy210,yyyxexyx)(37 从以上例子看到同阶的阿当姆斯方法，隐式方法要比显式方法误差小，这可以从两种方法的局部截断误差主项 的系数大小得到解释，这里 分别为 及 . )()1(11npppxyhc1pc720/251720/19