跌宕起伏的CFD收敛之旅

以下文章来源于LBM与流体力学 ，作者卢比与钢蛋

对于我们CFDer来说，计算的发散可以被称为“CFDer心中永远的痛”。为了更好的理解计算的发散问题，首先从剁手的起源——微分方程的数值解法说起…

现实中的流体无时无刻不在发生着变化，因此流动问题原则上都是非定常的。不过具体到某一工程问题，如果流场的变化幅度很小，或者流场的细节变化并非我们感兴趣的内容，则通常可简化为定常。比如对于管路流动，人们更想知道流体通过管道的压降及平均流场，而不太关心瞬态的流动细节，这样可使问题大大简化。

因此对于某一具体的工程问题，我们首先要界定是用定常还是非定常的算法。简而言之，定常算法只求解流场的某一稳定状态，而非定常计算则需包含跟时间相关的物理过程的历史细节。

对于喜欢记账的朋友，大概可以把定常计算理解为只关心每个月在不同领域的消费额度，而非定常计算则帮助消费者准确记录每时每刻的消费。

由于认知水平所限，人们感兴趣的大多数工程问题，都忽略了瞬态效应而近似为定常问题。比如飞机在平流层稳定的巡航，或者汽车在高速公路上定速行驶。因此定常问题的求解是人们在过去多年以来研究的聚焦点。

由于流体力学控制方程奇葩的数学性质，即便是定常计算，其数值求解的实现也不简单。亚声速不可压缩问题的控制方程为椭圆型，任意一点的扰动都会影响全场，而超声速流动的控制方程为双曲型，即流场中的扰动只会向下游传播。对于不同的流动问题，人们需要构造不同的数值格式——这导致求解方法丧失了普适性。

于是人们尝试在定常的方程中人为加入时间项，使得不同工况下控制方程的数学性质更为相似，便于构造统一的解法，这便是用非定常的思路来求解定常问题的由来。

在计算机条件欠佳的年代，人们首先发明了 “时间推进法”。具体思路是，通过一个假定的初始流场不断推进，最终逼近到流动稳定的时刻。这就是我们通常所说的显式算法。不过这里的时间步并非真实的物理时刻，因此被称为“伪时间步”。关于“显式算法”更粗暴的理解则是，此刻的结果可以直接通过上一步的结果推出来。如果用剁手的观点来看的话，则可以理解为老板上个月发多少工资，我这个月就买多少钱的东西。

显式算法由于构造简单，尤其是对计算机内存的需求小而广受欢迎。但显式算法由于需要兼顾稳定性，时间步不能太大，这大幅增加了计算时长，难以满足工程需求。就好比对于大多数人来说，如果只依靠工资消费的话，总是满足不了购物的需求和欲望，于是信用卡就应运而生了。

如果说“挣多少，花多少”是显式的剁手方式，那么利用信用卡借钱消费，再还钱的循环迭代模式则可以形容为隐式的剁手方式。隐式算法虽然在一定程度上解决了时间步的限制，但是其复杂的循环迭代模式，一不小心可能就引起信用卡逾期未还的发散。

为了更好的理解CFD计算中隐式的概念，我们一起来看看最常见的SIMPLE方法是如何让CFDer甘心先消费后付款的。

看着Spalding和Patankar这对师徒大神的迷之微笑，笔者不由的回忆起当年学习Patankar的《数值传热和流动》的痛苦经历。尽管教材一直保留着，但往事还是不要再提。

言归正传，作为一种压力修正解法，SIMPLE算法通过“先猜（消）想（费）后修（付）正（钱）”的方法得到压力场，并求解离散化的动量方程（N-S方程）。其基本思路如下：

• 基于给定的压力场，求解离散形式的动量方程，得出速度场。但压力场是假定的或不够精确的，由此得到的速度场一般不满足连续方程，因此，必须对给定的压力场加以修正。

• 修正压力使用的便是压力泊松方程，又称压力修正方程，此方程由离散形式的动量方程和连续性方程通过代数运算得出。

• 接着，使用修正后的压力场修正速度，使之满足连续性方程，但此时新的速度场又不一定满足动量方程。

• 然后不断迭代计算，并检查压力和速度的修正值是否收敛。若不收敛，则用修正后的压力开始新一轮的计算。如此反复，直到获得收敛的解。

  
  

通过上述的介绍可以看到SIMPLE方法巧妙的解决了N-S方程中的一个大问题：速度项通过动量方程和连续方程错综复杂的耦合在一起。而更繁琐的压力项，仅出现在动量方程中，不能得到直接求解。

在SIMPLE算法提出之后，一些改进算法也随之提出，其中比较知名的包括SIMPLER算法、SIMPLEC算法、PISO算法等，此处不再赘述。

通过上述对于SIMPLE方法的介绍，可以看到迭代法在工程上求解复杂方程具有很好的实用性。而紧紧扼住CFD计算收敛咽喉的“残差”，便是迭代法的产物。

为了更好的理解迭代和残差，我们一起看看下面这个简单的例子：

下表更进一步的给出了15次迭代的过程和计算结果。可见，随着迭代的不断推进，x和y的值越来越逼近其真实值1和2。而残差的绝对值也随着迭代不断减小。

如果把变量的残差值画成随迭代步的变化曲线，即残差曲线。当然上述的方程非常简单，所以其收敛曲线也很完美，计算的残差随着迭代的进行而不断减小。本节提到的残差为每个变量的绝对残差，实际工程计算中还可能会使用相对残差、最大残差、平均残差、均方根残差等参数。

虽然很多人使用信用卡进行隐式消费，但是只要能够规划好收入和消费的关系，按时还款，便可以保证消费的收敛。那么对于我们经常使用的CFD软件，其收敛标准除了迭代的残差之外，还有什么额外的参数需要考量呢？

• 迭代残差判据，刚刚已经举例说明了，主要是指计算迭代过程，当前计算步的物理量和前一步的差异，需要达到某一设定的标准。比如小伙伴们通常将此值设置为10-4，10-5等，残差越小意味着结果越逼近收敛；

• 通量守恒判据，意味着计算域进出口的物理量要达到平衡。最简单的便是质量守恒，即流入多少，也流出多少。通常情况下，在计算起始阶段，流场距离最终的稳定状态较远，进出口可能出现不平衡的回流，导致进出口的流量不一致；

• 监控参量判据，指的是在计算中实时检测某些敏感位置的物理量如压力或速度，如果随着计算的进行，这些物理量不再变化，基本可认为计算趋于收敛。

上述三大判据便是我们判断计算收敛或者发散的依据了，而小伙伴们也可以根据实际情况联立使用或者单独使用。

无论是消费还是数值计算，大家都希望结果是收敛的。不过，出于各种各样的原因，总有很多人忘记或者无法按时支付信用卡的欠款，如果长时间的逾期未还，可能就要出现发散的问题了。

对于计算过程中的发散，小伙伴们都不陌生，其表现为残差的急剧增大并导致计算崩溃。而在流场中则表现为某些区域内变量的大幅波动并超出了计算机可容许的数位。

计算的发散通常与以下这些因素相关：

• 边界条件

边界条件是首先要检查的问题，比如开口系统中只设置了进口，而没有出口；或者进出口都使用了静压边界条件，而计算域中又有剧烈的压力变化等，均容易导致数值溢出而发散。

• 初始流场

对于确定的物理问题，给定了初始流场，也就意味着我们要通过计算弥补初始流场和真实物理场之间的鸿沟。因此给定一个好的初场极大程度上影响着收敛性与收敛速度。毕竟，历史上能够从开局一个碗到打下江山的只有一个朱元璋。有的CFDer为了获取更接近真实的初始流场，通常会单独进行此预计算，比如人为增加粘度尽快获得稳定流场，再将粘度恢复为正常值进行计算。

• 网格问题

如果计算没算几步就开始发散，则很可能是网格质量太差，甚至出现了奇点。如果把网格形容为道路，数值求解所构造的矩阵方程组就是道路上驰骋的车辆。老司机常说的“宁走十里远，不走一里险”便是其中的道理。

• CFL数

无论是显式还是隐式，都会涉及到调节CFL数，CFL数过大则会导致计算不稳定。CFL数的意义在于调节时间步长，使得求解数值解的依赖域完全包含真实解的依赖域，也可理解为，计算不能漏掉流场中扰动的传播。

• 松弛因子

不同于CFL数，松弛因子（Relaxation Factors）影响的并非时间步长，而是将该时间步的计算结果与上一步结果的差值作适当缩减，以避免由于差值过大而引起非线性迭代过程的发散。因此调节松弛因子也可以调节计算的收敛性。

• 另一种不收敛

除了由于残差剧增或者物理量大幅变化引起的发散之外，还有一种不收敛却是由于所求解物理场本身的属性引起的。比如机翼的大攻角工况下，吸力面有大面积的流动分离，这种流动现象本身就有很强的瞬态特征，强制使用定常算法自然容易引起发散。

传统CFD方法使用的隐式迭代不可避免的会碰到收敛和发散的问题，尽管各路大神发明了各种“小额贷”或者“还款小助手”，不断改善计算的稳定性，但本质上无法完全消除计算迭代所导致的发散的风险。

反观LBM方法，其天然瞬态与接近物理本原的特性，决定了它是一种稳定的显式方法。而且，LBM中的时间推进就是真正的物理时刻，而非前文提到的“伪时间步”。因此，LBM计算中，下一时刻的状态完全由上一时刻推算得出，每个物理时间步内并无额外的迭代。

如果把传统的CFD方法形容为信用卡消费，那么LBM方法则是典型的储蓄卡消费。因此，优秀的LBM方法通常具有非常好的鲁棒性，几乎不会出现计算发散的问题。