Spring-ICE 结冰算法述评-（2）水滴轨迹计算

320科技工作室 2020年6月17日浏览：2614 评论：8 收藏：1

Spring-ICE系列前述两篇文章：

今天要聊的内容，既常规又不常规。

说其常规，是因为水滴轨迹的计算是所有结冰计算的基础算法，不可或缺，是必备技能。说其非常规是因为，尽管目前文献对这个计算过程描述的已经足够详细，但是还有很多的技术细节，非实际操作而不可揣摩，正是这些具体的技术细节决定着计算的准确性、效率，甚至是整个结冰算法的设计逻辑。

我这段时间感悟，任何工程的东西，越涉及到基础的和实践的，自己必须要亲自去试一试，做一做，且不能是简单的感受下，而要做出效果，做出感悟。要一竿子插到底，摸到边。

如果我做复材，我就要亲自去做一做材料成型，研究试件加工，然后试验测试，最后是结果处理。反复的想，整个过程和我仿真或者理论，哪些是对应的，哪些是不同的，这些不同会有什么影响，有些试验现象的机理是什么，能不能从理论层面找到依据。

如果我们现在的身份不是一个学生，不是面向自己的课题。而是一个工程师面对一个项目，整个项目不是一个人能够完成的，那么作为组织者，是不是给别人分个工，等着大家的结果就完事了？我想不是，即便没有精力面面俱到，那么至少应该精通某一方面，然后瞄准项目的目标，搞清楚大家的工作和目标的联系，然后对重点问题深入研究，要搞懂要把控。总而言之，必须实践，必须讲求细节。

1 水滴轨迹计算到底要干啥？

从结冰计算整个层面来说，所谓结冰就是水滴撞在机翼表面然后冻结的现象。水滴轨迹的计算，就是要搞清楚水滴会撞击到机翼表面的哪些区域，每个区域水滴撞击的量（即水滴收集量）是多少，为下一步的结冰物理计算提供输入。

从算法的层面，水滴轨迹计算需要做以下工作：

（1）控制方程建立。说白了，水滴在流场中的运动，满足什么方程。

（2）求解方程，计算轨迹。水滴在流场中运动，肯定要有流场的计算，然后是水滴运动方程的计算。

（3）水滴运动状态的判断。这个意思是说，要判断水滴在运动中是否已经撞击到了翼面，是否已经冲出了我们预设的计算边界，如果是，这个时候需要结束计算，然后存储这个结果，为水滴收集量的计算作输入。

2 水滴轨迹计算控制方程

首先是水滴运动控制方程，这个文献里都有，通常的拉格朗日法基本是下面这种形式：

这个方程，是考虑了水滴的浮力、重力、阻力、升力，从而根据牛顿第二定律建立的。说白了，就是高中的“受力分析”，然后列方程。

说着简单，要求解出来却要费点劲。好在文献里都说了，这个方程进一步化简为一阶常微分方程组，根据初值，采用龙格库塔法进行时间推进求解。这里涉及的是本科数值计算的内容。

转眼间，我们就用到了牛顿第二定律和数值计算方法。别急，还有重头戏。上面的方程中，需要用到流场速度。流场计算用什么算法？

流场算法的选择直接影响到需不需要做网格，结冰后网格需不需要写更新算法，这是结冰算法大框架的问题，十分重要。我就在这个问题上徘徊很久。最早我是准备写NS方程求解，为此我还摸索一套自己想的网格自动生成和更新方法。最后卡在NS方程上。

回过头看，当时NS方程应该是写出来了，但是网格域的问题，还有别的基本认知问题导致我认为自己走不下去了。被迫放弃了NS方程，在师兄的指导下，选择了面元法。

从现在看，面元法的选择无疑是对的，首先它足够简单，计算效率高。早期的空气动力学还管它加鳞片布源法，就是根据空气动力里“源”和“汇”的概念，在翼面上布置“源”，每个源对其他位置的流场都有一个影响，通过叠加这些影响再综合其在翼面上法向速度为0这样的边界条件，运用线性代数方法，求解这组方程，得到每个源的强度。然后再将每个源对流场中特定点影响求出来，叠加就是这个点的流场速度。

最简单的面元法，不考虑涡，两组线性代数方程，简单易上手。

除了简单，通过我上面的描述可知，面元法不需要划分流场网格，那么冰形长出来后，只需要把边界往外推进即可，不用重构网格，无疑解决了冰形更新带来的很多问题。

有人立马就要问了，那你这么算出来的结果准不准？

要回答这个问题，首先我们要搞清楚水滴的计算是两相流的问题，水滴自身的惯性使得它不会完全按照流线运动，也就是说，水滴的运动对流场是敏感的，但是不是我们想象中那么敏感，这就给我们的流场计算留出了裕度。另外一个方面，冰形的冰角位置和水滴收集量大的地方不是绝对的对应，因为水滴在翼面上会流动，这个流动对冰角位置的影响是很大的，而这个流动由结冰热力学方程确定的。也就说，即便前面流场的计算存在一定的误差，导致水滴的运动出现一定的偏差，但是最后的结冰热力学计算，可以把冰角的位置修正过来。

这个对我很有启发，作为一个涉及多个学科的平台开发，在能不能算的准的问题上，一定要首先考虑好，最后的结果有没有兜底的修正项，找到这个兜底的，前面遇到一些影响效率的东西或者很复杂的东西，可以适当舍弃。