查看“︁冯诺依曼稳定性分析”︁的源代码

[[数值分析]]中, '''冯诺依曼稳定性分析''' (亦作'''傅立叶稳定性分析''') 用于验证计算线性偏微分方程时使用特定[[有限差分法]]的[[数值稳定性]]<ref>{{Cite web |url=http://books.google.co.in/books?id=y77n2ySMJHUC&pg=PA523&dq=von+Neumann+stability+analysis#PPA523,M1 |title=Analysis of Numerical Methods by E. Isaacson, H. B. Keller |accessdate=2011-05-20 |archive-date=2011-05-21 |archive-url=https://web.archive.org/web/20110521003030/http://books.google.co.in/books?id=y77n2ySMJHUC&pg=PA523&dq=von+Neumann+stability+analysis#PPA523,M1 |dead-url=no }}</ref>，该分析方法基于对数值误差的傅立叶分解。1947年英国研究人员[[約翰·克蘭克]]和[[菲利斯·尼科爾森]]在文章中对该方法进行了简要介绍<ref>
{{Citation 
 | surname1 = Crank
 | given1 = J.
 | surname2 = Nicolson
 | given2 = P.
 | title = A Practical Method for Numerical Evaluation of Solutions of Partial Differential Equations of Heat Conduction Type
 | journal = Proc. Camb. Phil. Soc.
 | volume = 43
 | year = 1947
 | pages = 50–67
 | id = {{doi|10.1007/BF02127704}}
}}</ref>，
尔后又出现在冯诺依曼合作的文章中
<ref>
{{Citation 
 | surname1 = Charney
 | given1 = J. G.
 | surname2 = Fjørtoft
 | given2 = R.
 | surname3 = von Neumann
 | given3 = J.
 | title = [http://atoc.colorado.edu/~dcn/ATOC7500/members/Reading/Charney_ndbvm_1950.pdf Numerical Integration of the Barotropic Vorticity Equation]
 | journal = Tellus
 | volume = 2
 | year = 1950
 | pages = 237–254
}}</ref> 。
[[洛斯阿拉莫斯国家实验室]]对该方法进行了进一步发展。

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

冯诺依曼稳定性方法只适用于满足 Lax–Richtmyer 条件 （[[Lax 等价定理]]） 的某些特殊差分法： 偏微分方程系统须线性，常系数，满足周期性边界条件，只有两个独立变量，差分法中最多含两层时间步<ref>
{{Citation 
 | surname1 = Smith
 | given1 = G. D.
 | title = Numerical Solution of Partial Differential Equations: Finite Difference Methods, 3rd ed.
 | year = 1985
 | pages = 67–68
}}</ref>。 由于相对简单，人们常使用冯诺依曼稳定性分析代替其他更为详细的稳定性分析，用以估计差分方法中对容许步长的限制。

== 方法描述 ==
冯诺依曼误差分析将误差分解为[[傅立叶级数]]。为了描述此过程，考虑一维[[热传导方程]]
:<math>
  \frac{\partial u}{\partial t} = \alpha \frac{\partial^2 u}{\partial x^2}
</math>
空间网格间隔为 <math>L</math>, 对网格作 [[FTCS]] （Forward-Time Central-Space，时间步前向欧拉法，空间步三点中心差分) 离散处理，
:<math>
  \quad (1) \qquad u_j^{n + 1} = u_j^{n} + r \left(u_{j + 1}^n - 2 u_j^n + u_{j - 1}^n \right)
</math>
其中 <math>r = \frac{\alpha\, \Delta t}{\Delta x^2}</math> 。<math>u_j^{n}</math> 为离散网格上的数值解，用于近似此偏微分方程的精确解 <math>u(x,t)</math> 。

定义[[舍入误差]] <math>\epsilon_j^n = N_j^n - u_j^n</math> 。
其中 <math>u_j^n</math> 是离散方程 (1) 式的精确解，<math>N_j^n</math> 为包含有限浮点精度的数值解。 因为精确解 <math>u_j^n</math> 满足离散方程, 误差 <math>\epsilon_j^n</math> 亦满足离散方程 <ref>{{cite book | title = Computational Fluid Dynamics: The Basics with Applications | author = Anderson, J. D., Jr. | authorlink=John D. Anderson|edition = | publisher = [[McGraw Hill]] | year = 1994 | isbn = }}</ref>：
:<math>
  \quad (2) \qquad \epsilon_j^{n + 1} = \epsilon_j^n + r \left(\epsilon_{j + 1}^n - 2 \epsilon_j^n + \epsilon_{j - 1}^n \right)
</math>
此式将确定误差的递推关系。方程 (1) 和 (2) 中，误差和数值解随时间具有一致的变化趋势。对于含周期性边界条件的线性微分方程，间隔 <math>L</math> 上的空间部分误差可展开为傅立叶级数
:<math>
  \quad (3) \qquad \epsilon(x) = \sum_{m=1}^{M} A_m e^{ik_m x}
</math>
其中[[波数]] <math>k_m = \frac{\pi m}{L}</math>，<math>m = 1,2,\ldots,M</math>， <math>M = L/\Delta x</math>。 通过假设误差幅度 <math>A_m</math> 是时间的函数，可以给出误差和时间的关系。 不难知单步中，误差随时间指数增长，因此 (3) 式可以写作
:<math>
  \quad (4) \qquad \epsilon(x,t) = \sum_{m=1}^{M} e^{at} e^{ik_m x}
</math>
其中 <math>a</math> 为常量。

由于误差所满足的差分方程是线性的（级数每一项的行为与整个级数一致），只估计一项的误差变化便足以估计整体趋势：
:<math>
  \quad (5) \qquad \epsilon_m(x,t) = e^{at} e^{ik_m x}.
</math>
为找出误差随时间步的变化, 将方程 (5) 式应用于离散后的误差表达式上
*:<math>
\begin{align}
 \epsilon_j^n & = e^{at} e^{ik_m x} \\
 \epsilon_j^{n+1} & = e^{a(t+\Delta t)} e^{ik_m x} \\
 \epsilon_{j+1}^n & = e^{at} e^{ik_m (x+\Delta x)} \\
 \epsilon_{j-1}^n & = e^{at} e^{ik_m (x-\Delta x)},
\end{align}
</math>
再代入到 (2) 式中，求解方程后可得
:<math>
  \quad (6) \qquad e^{a\Delta t} = 1 + \frac{\alpha \Delta t}{\Delta x^2} \left(e^{ik_m \Delta x} + e^{-ik_m \Delta x} -  2\right).
</math>

使用已知的指数三角关系式
:<math> \qquad \cos(k_m \Delta x) = \frac{e^{ik_m \Delta x} + e^{-ik_m \Delta x}}{2} </math> 和 <math> \sin^2\frac{k_m \Delta x}{2} = \frac{1 - \cos(k_m \Delta x)}{2} </math>
可以将方程 (6) 变作
:<math>
  \quad (7) \qquad e^{a\Delta t} = 1 - \frac{4\alpha \Delta t}{\Delta x^2} \sin^2 (k_m \Delta x/2).
</math>

定义涨幅因子
:<math>
  G \equiv \frac{\epsilon_j^{n+1}}{\epsilon_j^n},
</math>
则误差有限的充要条件为 <math>\vert G \vert \leq 1</math> 。 已知
:<math>
  \quad (8) \qquad G = \frac{e^{a(t+\Delta t)} e^{ik_m x}}{e^{at} e^{ik_m x}} = e^{a\Delta t}, 
</math>
联立 (7) 和 (8) 两式，易得稳定性条件为
:<math>
  \quad (9) \qquad \left\vert 1 - \frac{4\alpha \Delta t}{\Delta x^2} \sin^2 (k_m \Delta x/2) \right\vert \leq 1
</math>
即
:<math>
  \quad (10) \qquad \frac{\alpha \Delta t}{\Delta x^2} \leq \frac{1}{2}.
</math>
(10) 即为该算法的稳定性条件。 对于 FTCS 求解一维热传导方程，给定 <math>\Delta x</math> ， 所允许的 <math>\Delta t</math> 取值需要足够小以满足 (10) ，才能保证计算的数值稳定。

==參考資料==
{{reflist}}

[[Category:数值分析]]