刚性方程

摘自Wikipedia——刚性方程

1. 定义

在数学领域中,刚性方程(stiffness equation)是指一个微分方程,其数值分析的解只有在时间间隔很小时才会稳定,只要时间间隔略大,其解就会不稳定。目前很难去精确地去定义哪些微分方程是刚性方程,然而粗略而言,若此方程式中包含使其快速变动的项,则其为刚性方程。
在积分微分方程时,若某一区域的解曲线的变化很大,会希望在这个区域的积分间隔密一些,若另一区域的曲线近似直线,且斜率接近零,会希望在这个区域的积分间隔松一些。不过针对一些问题,就算曲线近似直线,仍然需要用非常小的积分间隔来积分,这种现象称为“刚性”。有时可能会出现两个不同问题,一个有“刚性”,另一个没有,但两个问题却有同一个解的情形。因此“刚性”不是解本身的特性,而是微分方程的特性,也可以称为是刚性系统

2. 举例

考虑如下初值问题:

y=15y(t),t0,y(0)=1\begin{array}{c} y^{'} = -15y(t), t \geq 0, y(0) = 1 \end{array}

其精确解是:

y(t)=e15t\begin{array}{c} y(t) = e^{-15t} \end{array}

易知,当t{t \rightarrow \infty}时,y(t)0{y(t) \rightarrow 0}。当我们用数值分析的方法来解决该问题时,我们希望数值解也能满足该性质。

  • 若以欧拉方法来求数值解,则使用不同的步长将会得到不同的结果。第一种,步长h=1/4{h=1/4}的欧拉法强烈的震荡并且很快离开了图的边界。当将步长减半为 h=1/8{h=1/8}时,得到的结果在图的范围以内。但是它依然在0{0}附近震荡,并且不可能表示精确的解。

  • 而梯形法(即两阶段亚丹士-莫耳吞法)表达为

yn+1=yn+12h(f(tn,yn)+f(tn+1,yn+1))\begin{array}{c} y_{n+1} = y_n + {1 \over 2}h(f(t_n, y_n) + f(t_{n+1}, y_{n+1})) \end{array}

其求得的结果比欧拉法的结果要好很多。如上图所示,数值结果单调地减少到零,如同精确解一样。

3. 特性

刚性系统的特色是该系统所有特征值的实部均为负数,并且其中特征值实部绝对值中,最大和最小的比值远大于1。

4. 稳定区域

  • 将龙格-库塔法应用至测试方程 y=ky{y^{'} = ky} ,可以得到如 yn+1=ϕ(hk)yn{y_{n+1} = \phi(hk)y_n} 的形式,并可归纳出 yn=(ϕ(hk))ny0{y_n = (\phi(hk))^ny_0} ,其中 ϕ\phi 称为稳定性函数。因此

limnyn=0\begin{array}{c} \lim_{n \rightarrow \infty}y_n = 0 \end{array}

的条件等价于

ϕ(hk)<1\begin{array}{c} |\phi(hk)| \lt 1 \end{array}

这启发了绝对稳定区域(有时简称为稳定区域)的定义,亦即集合{zCϕ(z)<1}{ \{z \in \mathbb {C}|\,\,|\phi(z)| \lt 1\} }

  • 若一个方法的稳定区域包含 {zCRe(z)<0}{ \{z \in \mathbb {C} |\,\,\mathrm {Re} (z) \lt 0 \} }(即左半平面),则称该方法为A-稳定。
欧拉法 梯形法

本博客所有文章除特别声明外,均采用 CC BY-SA 4.0 协议 ,转载请注明出处!