张 进， 姜志恒， 李 瑜， 温铁钝

(中国航发控制系统研究所 传感系统工程部， 江苏 无锡 214063)

0 引 言

研究表明[1]，光冰对于飞机表面的影响较霜冰更为显著。飞机表面结冰大多会带来很多不良后果，本文着重关注的温度传感器表面结冰，尤其是感温元件表面结冰，会导致感温元件与环境换热热阻值变大，温度响应迟钝，受冰层自身温度的影响，传感器测量温度会与待测环境温度出现较大偏差，严重影响温度传感器的正常工作。因此，开展发动机飞行环境下温度传感器防冰结构设计及防冰效果研究具有必要性。

研究飞机表面的结冰过程常用方法中[2]对于较为复杂的结霜、结冰问题，建立数理模型，采用数值模拟仿真方法，编程或者使用商业软件求解过冷水滴运动、撞击壁面及结冰的热物理现象是较为可行的研究途径[3]。本文基于计算流体动力学(computational fluid dynamics,CFD)计算仿真软件STAR-CCM+11.06，对某型温度传感器的防冰结构进行了模拟分析，获得防冰结构绕流流场，过冷水滴运动，防冰结构表面过冷水滴收集情况、冰层增长情况及冰层形状。

1 研究对象与求解设置

1.1 结构形式与网格划分

本文选取的研究对象为某型号航空发动机T1总温传感器，布置位置为发动机进气道入口附近，该设置具有较高的结冰风险，因此，传感器具备防冰结构的设计。传感器结构及外形尺寸如图1所示。

图1 传感器的结构形式及外形尺寸

可以看出，传感器的防冰结构为一段弯管，端部具有特别切口，其防水原理为：气流携带过冷水滴流过弯管结构，利用转弯过程的离心作用将气流中的水滴分离出去，使过冷水滴无法到达感温元件附近，避免了感温元件表面结冰，主要的结冰区域位于防冰结构管腔表面的离心侧。因此，建立相应的计算域模拟防冰结构的实际工作状态。传感器防冰结构CFD计算域如图2所示。

图2 传感器防冰结构CFD计算域

可知CFD计算域为一个143 mm×186 mm×417 mm的长方体区域，一端为气流入口，另一端为气流出口，传感器垂直安装于安装壁面上，计算域其余表面设置为对称边界。运用多面体网格对计算域进行网格划分，传感器防冰结构表面为重点网格细化区域，最小网格尺寸为0.5 mm，传感器防冰结构附近空间区域也需要网格细化，以便提高多相流分离区域的网格解析度。

1.2 仿真边界条件及求解设置

本文基于star-ccm+11.06模拟不同马赫数Ma(0.1,0.4)具有一定过冷水滴含量的气流流过传感器防冰结构时，传感器防冰结构的防冰效果。本文选取的计算边界条件如下：入口气流马赫数取0.4和0.1；入口气流温度取-20 ℃；入口过冷水滴温度为-20 ℃；过冷水滴粒径取10 μm；入口过冷水滴体积分数为3.24×10-7。

求解器设置为隐式非稳态求解，求解物理时间为25 s，平均时间步长为0.01 s，时间步长内迭代步数为10～100。

2 数理模型

运用star-ccm+11.06求解流场，非稳态、可压、粘性流动，湍流模型选取双方程模型realizable k-epsilon model，传感器表面设定为无滑移壁面，选用标准壁面函数。过冷水滴输运计算引入如下假设[4]：

1)过冷水滴尺寸很小，水滴运动不会影响绕流流场，流场计算为单向耦合；

2)过冷水滴的密度和体积在运动过程中不会发生变化，水滴不会发生蒸发及其他形式的质量迁移；

3)将过冷水滴等效为当量小球，且当量球径在运动中不发生变化；

4)过冷水滴的初始速度与来流速度一致。

计算中考虑过冷水滴的重力,浮升力及气流阻力，根据牛顿第二定律，过冷水滴的运动方程可写为

(1)

式中md为水滴质量，xd为水滴位移，ρd和ρa分别为水滴及空气密度，Vd为水滴体积，gn为当地重力加速度，Ad为水滴的迎风面积，Cd为水滴的阻力系数，u为当地空气流速，ud为水滴的运动速度。由于过冷水滴当量直径很小，气流流过的雷诺数(Re)较小，因此阻力系数Cd可取24/Re，再结合结冰工况阻力修正[5]，阻力系数Cd可写为

(2)

关于结冰量及冰层厚度的计算，过冷水滴在撞击壁面的位置凝结成冰，剩余部分随气流流向后方，在此需要考虑壁面控制体内的质量守恒及能量守恒。本文假设所有的过冷水滴在撞击壁面时完全凝结成冰且不再融化，另外假设冰层只沿撞击表面法线方向生长，冰层厚度是可以计算，可以建立单位结冰时间步长内壁面控制体积结冰量为

mice=mw+mrbi-mrbi-1-mes

(3)

式中mice为壁面结冰量，mw为撞击壁面的过冷水滴质量，mrbi为从上游控制体流入未凝水质量，mrbi-1为由控制体流出的未凝水质量，mes为蒸发或者升华的质量。

能量守恒方程可写为

Qs=Qw+Qrbi+Qhot-Qrbi-1-Qes-Qh

(4)

式中Qs为过冷水滴结冰释放的凝结潜热；Qw为进入控制体的水滴带入的能量；Qrbi，Qrbi-1分别为从上游控制体流入及流向下游控制体的未凝水携带的能量；Qhot为由固体壁面传导入壁面控制体的热量；Qes为水蒸发或者水滴升华带走的能量。由式(3)、式(4)，可得壁面任意位置的结冰量。

关于冰层增长引起壁面偏移的计算，上文中模型能够完成结冰量及冰层厚度的计算，随着计算时间的推移，冰层厚度增加会影响结冰区域的壁面形态。本文中的温度传感器防冰结构是一个近似于弯管的结构，弯管内壁面会由于结冰发生管径变化，且这一变化会进一步影响管内流动和冰层积累，因此，必须引入网格模型来处理由于冰层增长导致的壁面偏移，使用网格变形模型morphing在每个时间步长内根据计算的冰层厚度增量去调整壁面的网格节点，使相应的网格表面能够向流体域内部推移。

3 计算结果与分析

3.1 不同来流马赫数对防冰效果影响分析

由于本文中传感器防冰结构主要利用弯管处的惯性离心效应来分离空气中的过冷水滴，那么这种分离效果将受到来流速度的影响。来流马赫数分别为0.4，0.1的绕流流场云图如图3所示。

图3 2 s时不同来流马赫数形成绕流流场剖面

由图3可以看出，随着来流马赫数的减小，弯管内的流速显著降低。管内气流中过冷水滴的分离效果会产生显著差别，来流马赫数分别为0.4,0.1的过冷水滴体积分数分布云图如图4所示。

图4 2 s时不同来流马赫数形成过冷水滴浓度剖面

可以看出，不同来流马赫数(Ma≥0.1)情况下，弯管防冰结构均能实现较好的分离效果，过冷水滴均无法到达感温元件表面，防冰效果较好；另外，随着马赫数的降低，过冷水滴的分离效果变差，过冷水滴撞击弯管内壁面的区域逐渐增大，过冷水滴浓度较高区域逐渐接近感温元件表面。

2 s时来流马赫数分别为0.4，0.1情况下，传感器及防冰结构表面水滴收集率分布如图5所示。

图5 2 s时不同来流马赫数下防冰结构表面水滴收集率分布

可以看出，不同马赫数下，防冰结构外表面的水滴收集分布特点及位置比较接近，数值上的差异是由于不同气流速度带来的收集流量不同。较低马赫数下，防冰结构弯管内壁面的水滴收集面积较大。

防冰结构弯管内壁面纵轴线上沿传感器纵轴方向(以安装平面为零点)的水滴收集率分布如图6所示。

图6 2 s时不同来流马赫数下弯管内壁面纵轴向水滴收集率

可以看出，来流马赫数为0.4时，弯管内壁面水滴收集率主要集中于防冰结构的前端，而来流马赫数为0.1时，弯管内壁面水滴收集率分布较为均匀，但是范围较大，从防冰结构前端延伸至气流泻出孔上沿。由于气流速度的不同，马赫数为0.4时的内壁面水滴收集率均值远大于马赫数为0.1时的情况。

3.2 定来流马赫数冰层堆积形态随时间变化分析

本文选取0.4来流马赫数作为条件，研究防冰结构弯管内表面冰层堆积形态随时间的变化情况。考虑壁面冰层堆积引起的壁面推移，2 s与25 s时，防冰结构弯管内壁面的变化情况(未考虑弯管外壁面冰层增长)如图7所示。

图7 不同时刻防冰结构弯管内壁面的变化情况

可以看出，0.4马赫数下，防冰结构弯管内壁面的结冰区域集中在弯管前端，且管道内壁迎风面与背风面的差异较大。不同时刻弯管内部纵剖面的变化对水滴收集效果的影响如图8所示。

图8 不同时刻弯管纵剖面轮廓及过冷水滴浓度分布

可以看出，随着时间推移，弯管前端内壁面冰层逐渐增厚，导致壁面向流体域内部推移，从而前端内壁面等效迎风面逐渐增大，使局部水滴收集率逐渐提高，因此，单位时间冰层增长厚度逐渐增大，冰层增长情况呈现非线性。另外，过冷水滴撞击区域逐渐向弯管前端偏移，更加远离感温元件区域。

4 结束语

1)不同来流马赫数(Ma≥0.1)情况下，弯管防冰结构均能达到较好的分离，使过冷水滴均无法到达感温元件表面，防冰效果较好；但随着来流马赫数的减小，水滴分离效果变差，马赫数低于0.1时，感温元件表面有结冰风险；

2)来流马赫数为0.4时，冰层生长集中在弯管前端，由于冰层生长引起壁面偏移，局部水滴收集率逐渐变大，从而加速冰层生长，结冰区域逐渐向前端偏移。依据计算结果推算，来流马赫数0.4持续作用1 min之内，弯管封冻点将出现在弯管前端，从而影响到传感器正常工作。该防冰结构不具备自我除冰能力，需要改进设计，增加主动除冰功能。