数值模拟导论 ——第十四讲 多步法 Ⅱ 雅可比 ·怀特 合作者 Deepak Ramaswamy, MichalRewienski, and Karen Veroy 概要 多步法的小时间步问题 回顾局部切断误差最小化 一个不收敛的实例 稳定 +连续导出收敛 考察大时间步问题 两个时间刻度实例的绝对稳定性 振荡器 SMA-HPC ?2003 MIT 一般表示法 非线性差分方程 : k步多步法 : () () () , d xt f xt ut dt = () () 00 ?? , kk lj lj jjlj ax t f x utβ ?? ? == =? ∑∑ 多步系数 基本方程 多步法 离散点的解 时间间断 SMA-HPC ?2003 MIT 多步法 简单问题分析 () () () 0 , 0 C d vt vt v v dt λλ= =∈ 00 ?? kk lj lj jj ax t vβλ ? ? == =? ∑∑ Cλ∈ 标量ODE: 多步法标量公式: 对于所有 都必须考虑λ的值 SMA-HPC ?2003 MIT 多步法 收敛分析 收敛定义 ( ) 0, ? 0 max 0 l T l t tvvt ?? ∈ ?? ? ?? ? →??→当时 定义：对于给定的任何初始条件，用多步法 求解 [0,T]上的初值问题是收敛的。 SMA-HPC ?2003 MIT 多步法 收敛分析 收敛的两个条件 1）局部条件： “一步法 ”误差较小 （连续性） 可以用泰勒展开式证明 2）全局条件：单步误差不会快速增加 （稳定性） 多步法（ K>1）需要详细分析 SMA-HPC ?2003 MIT 多步法 收敛分析 全局误差等式 00 ?? 0 kk lj lj jj av t vβλ ?? == ??= ∑∑ () () 00 ? kk l jlj j lj jj d avt t vt e dt β ?? == ??= ∑∑ ( ) ? ll l E vt v≡ ? ( ) ( ) ( ) 1 00 11 ll lkl kk atEatE atEeλβ λβ λβ ?? ? ? + ?? + + ?? =" 多步法公式： 精确解基本上和多步法一致： 全局误差： 差分方程将局部切断误差和全局误差联系起来 局部切断误差（ LTE） SMA-HPC ?2003 MIT 多步法 使得局部切断误差足够小 () () 1p p d vt t vt pt dt ? =? = () () d vt vt dt λ= () () () () 1 00 kj kj kk pp jj d vt vt dt akjt t pkjt eβ ? ? ? == ? ??? ?? = ∑∑   局部切断误差： 局部切断误差 若 注：不能根据 式推导 () () 00 kk l jlj j lj jj d avt t vt e dt β ?? == ??= ∑∑ SMA-HPC ?2003 MIT 多步法 使得局部切断误差足够小 精度约束 k=2实例 () () 1 00 0 kk pp jj ak j pk jβ ? == ?? ? ??= ?? ?? ∑∑ 0 1 3 4 p p p p = = = = 0 1 2 0 1 2 1 1 1 0 0 0 0 2 1 0 -1 -1 -1 0 4 1 0 -4 -2 0 0 8 1 0 -12 -3 0 0 16 1 0 -32 -4 0 0 a a a β β β ?? ???? ?? ???? ?? ?? = ?? ???? ?? ?? ?? ???? ?? ?? 0 i a = ∑ 精度拘束： 对于 k＝ 2，产生系数的 5x6维矩阵方程 注：总存在 SMA-HPC ?2003 MIT 多步法 使局部切断误差较小 精度约束k=2实例, 生成方法 1 2 0 1 2 1 1 0 0 0 1 1 0 -1 -1 -1 2 1 0 -4 -2 0 4 1 0 -12 -3 0 8 1 0 -32 -4 0 16 a a β β β ??? ? ??? ??? ??? ? ??? ??? ??? =? ?? ? ? ??? ??? ??? ? ? ????? 首先介绍一个规范化的问题，例如 0 1a = 01 2 0 1 2 14 1 1, 0, 1, , , 333 aaa βββ== =?= = = 5 ()Ct? 01 2 0 1 2 1, 4, 5, 0, 4, 2aaa β ββ== =?= = = 4 ()Ct? 用两步法求解最小的局部切断误差 满足所有五个精度约束 局部切断误差＝ 用两步直接法求解最小的局部切断误差 只能满足四个精度 局部切断误差＝ SMA-HPC ?2003 MIT 多步法 使局部切断误差较小 前欧拉法，后欧拉法，梯 形捕捉法和 “最 ”直接法的 局部切断误差图 SMA-HPC ?2003 MIT 多步法 使局部切断误差较小 前欧拉法，后欧拉法，梯形捕捉 法和 “最 ”直接法的全局误差图 SMA-HPC ?2003 MIT 多步法 使局部切断误差较小 前欧拉法，后欧拉法，梯形捕捉 法和 “最 ”直接法的全局误差图 SMA-HPC ?2003 MIT 多步法 多步法的稳定性 差分方程 为什么 2步直接法不能收敛？ 多步法差分方程 全局误差 我们在使局部切断误差变小的过程中，全局 误差怎么会变得这样大？ SMA-HPC ?2003 MIT 多步法 多步法的稳定性 稳定性定义 ( ) ( ) ( ) 1 00 11 ll lkl kk atEatE atEeλβ λβ λβ ?? ? ? + ?? + + ?? =" () N 0, 0, 0 max max ll TT tt T tEC e t ?? ?? ∈∈ ?? ?? ?? ?? ?? ?→ ≤ ? 间断独立 多步法差分方程 定义：当 全局误差受持续时间段内局 部切断误差之和的限制 . 时， 多步法是稳定的。 稳定性含义 : SMA-HPC ?2003 MIT 差分方程的补充 离散和 根的关系 1 01 ll lkl k ax ax ax u ?? +++ =" 0 0, , 0 k xx ? = =" 0 l llj j xhu ? = = ∑  离散和 () () 1 , 10 1 01 + + =0 q M Q l m l qm q q k m k k a h a l zaz γ? ? == ? ↑ = ∑∑ 0 " 多个根 方程 的根 给出一个 k阶差分方程初始条件为 0 , x通过式 和 u联系起来 SMA-HPC ?2003 MIT 差分方程的补充 离散和 边界条件 () () ( ) , 1 , 10 q qm M Q lj m lj qm q qm R xlj uγ? ? ? == =? ∑∑  1 q ? < , m ax j qm j l RC u≤ 3与独立 ( ) 1 q ? ε< + ,0 m ax l j qj e RC u ε ε ↑ ≤ 边界不等根 若 ,那么 若 ,那么 SMA-HPC ?2003 MIT 多步法 多步法的稳定性 稳定性定理 定理：当且仅当方程 1 01 0 kk k az az a ? + ++=" 1.数值小于 1 2.数值等于 1,但是不相等 的根为下列两者之一时 ,多步法是稳定的 : SMA-HPC ?2003 MIT 多步法 多步法的稳定性 ( ) ( ) ( ) 1 00 11 ll lkl kk atEatE atEeλβ λβ λβ ?? ? ?+?? ++?? =" 0t? → ( ) ( ) ( ) 1 0011 0 ll kk atzatz atλβ λβ λβ ? ? ? + ?? + + ?? =" 1 tκ+? N 0, 0, 0, () max max max lt T ll l TT T ll tt t CT eCe E Ce e tTt κκ? ?? ?? ?? ∈∈ ∈ ?? ?? ?? ?? ? ?? ?? ?? ≤≤ ?? 给定的多步法差分方程 若 ,当时,方程 z数值上小于 1或 z根不同并极限 在 再由差分方程的补充式 的根 (κ >0)内 SMA-HPC ?2003 MIT 多步法 多步法的稳定性 SMA-HPC ?2003 MIT 多步法 多步法的稳定性 最直接法 01 2 0 1 2 1, 4, 5, 0, 4, 2aaa β ββ== =?= = = 最直接 2步法 “最 ”直接法非常不稳定 SMA-HPC ?2003 MIT 多步法 多步法的稳定性 Dahlquist 首次稳定性界限 对于一个稳定的，直接 k步多法， 可以满足的精度约束的最大值小于或 等于 k（注意由 2k个系数）。对于非直 接法，精度约束在 k为偶数时为 k＋ 2， 或 k为奇数时精度约束为 k＋ 1。 SMA-HPC ?2003 MIT 多步法 收敛分析 收敛性、稳定性和连续 性的条件 0 tt? <? 0 p 0 p () 0 1 1 0, max p l T l t eCt + ?? ∈ ?? ? ?? ?≤? 2 0, 0, max max ll TT ll tt T EC e t ?? ?? ∈∈ ?? ?? ?? ?? ?? ?≤ ? () 0 0, max p l T l t E CT t ?? ∈ ?? ? ?? ≤? 1)局部条件：一步误差较小（连续性） 对于 ，精度约束达到 （ 必须大于 0） 1)全局条件：一步误差缓慢增加（稳定性） 方程的根在单位圆的内部或在圆上 收敛结论： SMA-HPC ?2003 MIT 大时间步问题 多步法 采用后欧拉法 ,很容易在动态变化剧烈段 采用 小时间步 ,在动态衰减缓慢时采用大时间步 两个时间常数的回路 问题 SMA-HPC ?2003 MIT 两个时间步常数回路 的前欧拉法 大时间步稳定性 多步法 前欧拉法对小时间步较精确 ,而当时间步增大 时变得不稳定 . SMA-HPC ?2003 MIT 大时间步稳定性 多步法 标量 ode问题的前欧拉法 , 后欧拉法和梯形捕捉法 () () 0 (), 0 C d vt vt v v dt λλ= =∈ 1 ? l v + = ? l v + ? l tvλ?= ( ) ? 1 l tvλ+? 11tλ+?> 1 ? l v + = ? l v + ? l tvλ? () 1 1 ?? 1 ll vv tλ + ?= ?? () 1 1 1 tλ < ?? 1 ? l v + = ? l v + () ( ) () 11 10.5 ?? ? ? 0.5 10.5 l ll l l t tv v v v t λ λ λ ++ +? ?+?= ?? 标量 ODE: 前欧拉法 : 若 后欧拉法 : 若 梯形捕捉法 : ,即使λ <0解会变大 ,即使λ >0解会减小 SMA-HPC ?2003 MIT 绝对稳定性的前欧拉 法大时间步区域 大时间步稳定性 多步法 (1 )z tλ= +? 前欧拉法 SMA-HPC ?2003 MIT 电路实例中的前欧拉 法大时间步的稳定性 大时间步稳定性 多步法 0.1, 2.1, 0.1t λ? ==??电路实例 SMA-HPC ?2003 MIT 电路实例中的前欧拉 法大时间步的稳定性 大时间步稳定性 多步法 1.0, 2.1, 0.1t λ? ==?? 电路实例 SMA-HPC ?2003 MIT 绝对稳定性的后欧拉 法大时间步区域 大时间步稳定性 多步法 1 (1 )ztλ ? =?? 后欧拉法 SMA-HPC ?2003 MIT 电路实例中的后欧拉 法大时间步的稳定性 大时间步稳定性 多步法 0.1, 2.1, 0.1t λ? ==?? 电路实例 SMA-HPC ?2003 MIT 电路实例中的后欧拉 法大时间步的稳定性 大时间步稳定性 多步法 1.0, 2.1, 0.1t λ? ==?? 电路实例 SMA-HPC ?2003 MIT 稳定性定义 大时间步的稳定性 多步法 () 0 0 k kj jj j atzλβ ? = ?? = ∑ tλ? 多步法的绝对稳定性区域 的根在单位圆内时 A-稳定 : 一种方法若它的绝对稳定性区域包含整 个复杂平面的左上角 ,则为 A-稳定 . Dahlquist的第二稳定性界限 : 收敛阶大于 2时不存在 A-稳定多步法 , 且 的值当方程 梯形法最精确 . SMA-HPC ?2003 MIT 振动杆和振动块问题 数值试验 多步法 为什么前欧拉法振幅增加 ,后欧拉法衰减 ,梯 形法保持振幅 ? SMA-HPC ?2003 MIT 前欧拉法大时间步振 动实例 大时间步稳定性 多步法 (1 )z tλ= +? 前欧拉法 SMA-HPC ?2003 MIT 后欧拉法大时间步实例 大时间步稳定性 多步法 1 (1 )ztλ ? =?? 后欧拉法 SMA-HPC ?2003 MIT 梯形法大时间步振动 实例 大时间步稳定性 多步法 (1 0.5 ) (1 0.5 ) t z t λ λ + ? = ?? 梯形法 SMA-HPC ?2003 MIT 大时间步问题 多步法 两个时间常数的稳定性问题 (电路中 ) 前欧拉法 :稳定 ,不精确 ,时间步大小有限制 . 振动问题 前欧拉法产生一个不稳定的差分方程而 不考虑时间步大小 . 梯形捕捉法时最精确的 A-稳定 m步法 . 后欧拉法时 A-稳定 ,任何时间步都适合 梯形法映射虚轴 . 后欧拉法产生一个稳定的 (衰减的 )差分方程 而不管时间步大小 . SMA-HPC ?2003 MIT 局部切断误和精度 差分方程稳定性 稳定 +连续导出收敛 小结 多步法的小时间步问题 研究大时间步问题 两个时间刻度实例的绝对稳定性 振荡器 没有谈及的问题 朗格 -库拉定理 , 高阶 A-稳定法 SMA-HPC ?2003 MIT