冯诺依曼稳定性分析

数值分析中, 冯诺依曼稳定性分析 (亦作傅立叶稳定性分析) 用于验证计算线性偏微分方程时使用特定有限差分法的数值稳定性^[1]，该分析方法基于对数值误差的傅立叶分解。1947年英国研究人员約翰·克蘭克和菲利斯·尼科爾森在文章中对该方法进行了简要介绍^[2]，尔后又出现在冯诺依曼合作的文章中 ^[3] 。洛斯阿拉莫斯国家实验室对该方法进行了进一步发展。

数值稳定性

数值稳定性与数值误差密切相关。使用有限差分方法进行计算时，若任意时间步的误差不会导致其后计算结果的发散，则可称该有限差分法是数值稳定的。如果误差随着进一步计算降低最终消失，该算法被认为稳定；若误差在进计算中保持为常量，则认为该算法“中性稳定”。但如果误差随着进一步计算增长，结果发散，则数值方法不稳定。数值方法的稳定性可以通过冯诺依曼稳定性分析得到验证。稳定性一般不易分析，特别是针对非线性偏微分方程。

冯诺依曼稳定性方法只适用于满足 Lax–Richtmyer 条件（Lax 等价定理）的某些特殊差分法：偏微分方程系统须线性，常系数，满足周期性边界条件，只有两个独立变量，差分法中最多含两层时间步^[4]。由于相对简单，人们常使用冯诺依曼稳定性分析代替其他更为详细的稳定性分析，用以估计差分方法中对容许步长的限制。

方法描述

冯诺依曼误差分析将误差分解为傅立叶级数。为了描述此过程，考虑一维热传导方程

\frac{\partial u}{\partial t} = α \frac{\partial^{2} u}{\partial x^{2}}

空间网格间隔为 $L$ , 对网格作 FTCS （Forward-Time Central-Space，时间步前向欧拉法，空间步三点中心差分) 离散处理，

(1) u_{j}^{n + 1} = u_{j}^{n} + r (u_{j + 1}^{n} - 2 u_{j}^{n} + u_{j - 1}^{n})

其中 $r = \frac{α Δ t}{Δ x^{2}}$ 。 $u_{j}^{n}$ 为离散网格上的数值解，用于近似此偏微分方程的精确解 $u (x, t)$ 。

定义舍入误差 $ϵ_{j}^{n} = N_{j}^{n} - u_{j}^{n}$ 。其中 $u_{j}^{n}$ 是离散方程 (1) 式的精确解， $N_{j}^{n}$ 为包含有限浮点精度的数值解。因为精确解 $u_{j}^{n}$ 满足离散方程, 误差 $ϵ_{j}^{n}$ 亦满足离散方程 ^[5]：

(2) ϵ_{j}^{n + 1} = ϵ_{j}^{n} + r (ϵ_{j + 1}^{n} - 2 ϵ_{j}^{n} + ϵ_{j - 1}^{n})

此式将确定误差的递推关系。方程 (1) 和 (2) 中，误差和数值解随时间具有一致的变化趋势。对于含周期性边界条件的线性微分方程，间隔 $L$ 上的空间部分误差可展开为傅立叶级数

(3) ϵ (x) = \sum_{m = 1}^{M} A_{m} e^{i k_{m} x}

其中波数 $k_{m} = \frac{π m}{L}$ ， $m = 1, 2, \dots, M$ ， $M = L / Δ x$ 。通过假设误差幅度 $A_{m}$ 是时间的函数，可以给出误差和时间的关系。不难知单步中，误差随时间指数增长，因此 (3) 式可以写作

(4) ϵ (x, t) = \sum_{m = 1}^{M} e^{a t} e^{i k_{m} x}

其中 $a$ 为常量。

由于误差所满足的差分方程是线性的（级数每一项的行为与整个级数一致），只估计一项的误差变化便足以估计整体趋势：

(5) ϵ_{m} (x, t) = e^{a t} e^{i k_{m} x} .

为找出误差随时间步的变化, 将方程 (5) 式应用于离散后的误差表达式上

$\begin{matrix} ϵ_{j}^{n} & = e^{a t} e^{i k_{m} x} \\ ϵ_{j}^{n + 1} & = e^{a (t + Δ t)} e^{i k_{m} x} \\ ϵ_{j + 1}^{n} & = e^{a t} e^{i k_{m} (x + Δ x)} \\ ϵ_{j - 1}^{n} & = e^{a t} e^{i k_{m} (x - Δ x)}, \end{matrix}$

再代入到 (2) 式中，求解方程后可得

(6) e^{a Δ t} = 1 + \frac{α Δ t}{Δ x^{2}} (e^{i k_{m} Δ x} + e^{- i k_{m} Δ x} - 2) .

使用已知的指数三角关系式

\cos (k_{m} Δ x) = \frac{e^{i k_{m} Δ x} + e^{- i k_{m} Δ x}}{2}

和

\sin^{2} \frac{k_{m} Δ x}{2} = \frac{1 - \cos (k_{m} Δ x)}{2}

可以将方程 (6) 变作

(7) e^{a Δ t} = 1 - \frac{4 α Δ t}{Δ x^{2}} \sin^{2} (k_{m} Δ x / 2) .

定义涨幅因子

G \equiv \frac{ϵ_{j}^{n + 1}}{ϵ_{j}^{n}},

则误差有限的充要条件为 $| G | \leq 1$ 。已知

(8) G = \frac{e^{a (t + Δ t)} e^{i k_{m} x}}{e^{a t} e^{i k_{m} x}} = e^{a Δ t},

联立 (7) 和 (8) 两式，易得稳定性条件为

(9) | 1 - \frac{4 α Δ t}{Δ x^{2}} \sin^{2} (k_{m} Δ x / 2) | \leq 1

即

(10) \frac{α Δ t}{Δ x^{2}} \leq \frac{1}{2} .

(10) 即为该算法的稳定性条件。对于 FTCS 求解一维热传导方程，给定 $Δ x$ ，所允许的 $Δ t$ 取值需要足够小以满足 (10) ，才能保证计算的数值稳定。

參考資料

Template:Reflist

[1] Template:Cite web

[2] Template:Citation

[3] Template:Citation

[4] Template:Citation

[5] Template:Cite book

[1]

[2]

[3]

[4]

[5]

冯诺依曼稳定性分析

数值稳定性

方法描述

參考資料

导航菜单

搜索