现在让我们考虑第二个问题。“你的求解法对不对？”这也是个相当麻烦的问题。不管怎样，关于求解法的可靠性，曾经有很多数学家（也许不是数学家，或者被某些自认为是数学家的人鄙视为不是真正的数学，但我们这里不关心那些自命不凡的家伙）进行研究过，以至于很多人干脆认为数值分析就是研究“计算结果的误差”的学科。

让我先把一些不太相关的东西清理掉。比如有些人认为舍入精度的选择会影响计算稳定性，比如单精度和双精度计算会引起计算稳定性计算收敛性上的差距。我不是说这个想法是错误的，它只是不太得当，实践中我们选择单精度和双精度几乎从来不是因为这个。（比如当年有个人叫卢刚或者叫卢瑟，他用单精度算他的课题，然后老板认为用单精度是错误的，要用双精度重算，当然卢瑟同学没办法就重算了一下然后耽误了找工作。。。这个故事的结局我们最后再说）

最常见的物理学相关的数值模拟归结为解方程。什么是物理理论？就是给个系统的初始状态，然后计算一段时间后的状态的技术。由于时间是连续的，或者至少我们近似把它当成连续的，于是这问题化归为一个演化方程，通常是积分微分方程。如果问题足够简单，或者我们用某些技巧把积分处理掉，就会得到一个需要求解的微分方程。

那么接下来的问题就是，计算机本身并不能理解微积分，它只懂加减乘除。所以我们不得不用某种近似手段把求导化成除法，把积分变成求和，这个称为“离散化”。这样我们就得到第一个问题，这种替换肯定是近似的，你计算的结果精度必然受到替换准确度的限制。当然，你肯定也可以意识到，这个精度可以靠不断缩减离散的“步长”来提高。

然而现在你就会意识到一个问题，拿最简单的质点运动方程为例，我们的离散化无非就是把df/dt用(f(t+s)-f(t))/s来取代，其中s称为推进步长，是一个很小但并不是任意小的数。现在你可以时间上缩短每一个推进步的长度，每一步的误差都在降低；但是同时，如果你需要研究的总时间T是固定的，那么缩短每一个推进步s，总的执行步数T/s就要增加。那么，算完整个问题后，你的真正误差f(T)-f0(t)是否降低了？这不是一个很容易判定的问题。

当然，对于f函数是光滑的情况，常微分方程问题终究比较简单，数学家已经给了很明确的结论，我们可以用改变计算方法的手段，让“预计误差”是s^2的高阶无穷小,于是随着步长降低，最终误差肯定会随着s逐步减小。

这些看起来都非常简单也非常自然，根本用不着啥创新研究（当然！）。但如果你对计算问题的基本状况有一点点了解，你就能猜到问题出在什么地方了。在数值问题上，理论上的误差分析完成是一个重大的成就，但同时也仅仅相当于火星登陆计划建立了第一个地面观测站。实际上，即使对于常微分方程问题，这种误差分析的用处也不是很大。第一个问题是，实际上我们根本不可能进行足够细致的时间分析。这种情况一般发生在系统中有速度相差若干数量级的几个过程的情况。比如说原子核衰变，可能有若干谱系，一般长寿命放射性元素的平均衰变周期可以长达一百亿年；而中间产物的寿命可以从几亿年到几秒钟。如果你要详细分析所有的衰变产物演化，那么你需要跑上10^16步才能进行准确的分析。当然，拜最新的计算机技术所赐，我们可以在一个廉价CPU上进行每秒钟上百亿次计算，问题是，即使这样，时间步长充分解析仍然是不可能的。这种问题称为“刚性问题”。

既然时间步长无法充分解析，那么很明显，你可以适当加大一点时间步。“我们可以容忍一定的误差！”这是真理。但是，事情并不总是这么简单，原子核衰变问题本质是个强阻尼方程，随便你乱搞一气都不会出问题。当你考虑一个弹簧振子的时候，问题变得严重得多。弹簧振子方程就是dx/dt=v, dv/dt=-x。简单的“离散化”无非就是把求导变成差分。现在你考虑如果时间步长太大会是什么情况：首先粒子在平衡点左侧的时候，你计算到了一个向右的推力，然后把粒子推向平衡点，过了平衡点之后，再计算推力，开始降低粒子的速度。现在如果你的步长太大，粒子就可能在你的推动作用下一下子跑到超过右侧最远位置；这完全是你计算精度不够，和物理无关；然后到了右侧，你计算的推力还是太大，于是把粒子推回到左侧但推得更远，如此循环一会，粒子的振幅会越来越大，于是得到了一个完全反物理的结果。这一切毫不涉及物理，只是因为你的算法不能容忍大时间步长。这种情况，我们说结果发散了。这时候你也不用想什么精度问题了，因为结果已经彻底错误了。

所以你会看到，当碰到这类问题的时候，你的计算结果可能会完全和物理没有关系。当然，对于弹簧振子，你总可以缩小时间步长来解决这个问题。可是如果弹簧振子耦合上一个慢速转动的车轮呢？这就又回到前面说的刚性问题。刚性问题并不是不能容忍精度误差，它真正的问题是快过程往往会引起发散，发散了的计算就没有意义，或者你两次计算得到的结果完全不同，或者是计算出的全是无穷大。

当然这问题总是可以解决的，这一般是采用所谓“隐格式”的技巧来完成。在下一部分我会回到这个事情，现在你只要意识到，即使对于谐振子，数值计算也不是那么简单，并且你需要首先在理论上校核你的算法本身有没有问题。

熟悉计算的应该马上意识到其实我没有对你说出大部分真相。在很多情况下，对于谐振子，即使没有发生上面出现的发散，结果在物理上也是本质上（！）错误的。没错，这就是我要说的第二个问题。拿几种最简单的数值差分方法去算一个简谐振荡，然后他们看上去都能算出差不离的结果。。。但是算足够多的步数之后，你会发现有些差分方法（比如你通常用的那种），计算出的振子振幅正在缓慢的减小！

为什么？因为我们并没有用无限短的步长，导致计算仍有很大的误差，而在很多情况下，这个误差会体现在振幅上。当然你马上就会意识到，这种计算有点扯。因为谐振子的基本特性就是能量守恒，振幅的减小在本质上就是破坏能量守恒，然后这个计算在物理上的意义可能就会大打折扣。更重要的是，在实践中我们的问题比谐振子要复杂很多，振荡可能仅仅是物理问题中的一环，那么能量守恒的破坏就可能导致系统长期行为出现严重的偏离甚至本质错误。

解决这个问题的方法很简单，那就是在设计计算方法的时候就人为强制守恒定律成立。这里所说的守恒定律不见得是能量守恒，也许是其他更复杂的东西。。。。在哈密顿系统的情况，有一系列技巧用于处理这个问题，统称为“辛算法”。（应该一提的是这问题在60年代NASA的研究人员就指出过，并也给出了解决方案的雏形，然而没有人重视他的工作，所以要到80年代之后才有人重新捡起这个问题，并且在本世纪由一群在微分方程数值问题上碰得头破血流的家伙发扬光大。Anyway，如果他们不是非要相信自己的“实际工程经验”而是记得看看过去那些“毫无用处的论文”的话，本来不至于如此。）

现在你看到，即使对于最简单的问题，数值模拟的结果可信度也是必须小心的。最简单的数值算法中都充斥着陷阱，稍微不留神就会掉进某个坑洞。更重要的是，和一般的工程开发不同，大部分数值模拟人员并不可能对整个问题有较为完整的理解：一个稍微复杂一点正规一点的数值模拟问题，如同我们前面强调的那样，就不会有很多人真的了解其中全部的陷阱。问题可能出在物理，数学，计算机，以及其他任何一部分，而有多少人真的可以有勇气说他了解了自己研究的问题从物理建模到计算机上实际运行的整个过程呢？即使你真的有这样的能力，你如何让别人相信你不是在吹牛呢？记住，在得出一个计算结果之前，一切都可能带来错误。数值的失败有一千个父亲，成功却没有。

我们已经看到问题所在，数值本身总是在某些限定条件下进行，于是必然有一定的误差，而这种误差的有限性（或者叫模拟的收敛性）以及守恒律问题是我们用来考核模拟的基本手段。一旦你的算法在这两个方面有所缺陷，那么带来的就是严重的误差。而且在很多问题中，这些限制并不像谐振子问题那么容易推断，这时候，你就会碰到几乎无法确定模拟的适用性的困境，甚至你根本无法判断模拟结果是不是有意义的。

（让我们回到一开始的那个可悲的例子。卢刚认为老板在人为刁难他，然后枪毙了他的老板以及同学。下一部分我会提到他那个问题的一些侧面，你会看到他发疯并不能完全怪他，当然老板刁难他也不能说是完全无理取闹，总之这是个杯具。）