诱导轮超同步旋转空化传播机理

项 乐，李春乐，许开富，陈 晖

(西安航天动力研究所，陕西 西安 710100)

0 引言

诱导轮是安装在涡轮泵主泵上游的一种轴流泵，通过适当增压提升来流压力，从而提高主泵的空化性能。由于其安装角小、稠度大等结构特点，能够在一定程度的空化条件下正常工作，而不发生明显的扬程性能下降。然而即便在设计工况下，诱导轮内也可能发生空化诱发的流动不稳定现象，例如我国新一代高压补燃循环液氧煤油液体火箭发动机研制过程中，超同步旋转空化现象长期以来是导致涡轮泵振动量级过高的重要激振源之一，我国在研的500 tf级液氧煤油发动机涡轮氧泵水试和热试车中均出现了明显的超同步旋转空化现象，其特征频率幅值显著大于叶轮转频幅值，导致振动幅值偏高，目前尚无有效技术手段彻底抑制，因此亟需对诱导轮超同步旋转空化现象发生机理展开深入研究，在涡轮泵设计阶段对该现象进行抑制，从而提升发动机的可靠性。

超同步旋转空化是旋转空化现象的一种，是相对于同步旋转空化和次同步旋转空化而言的。旋转空化现象是指诱导轮3个叶片上空化区波动沿周向传播的一种流动不稳定现象，按传播的频率不同可以分为3类：空化区波动传播频率快于诱导轮转速(或者传播方向与诱导轮旋转方向相同)则是超同步旋转空化，对应的频率一般大于1；空化区波动传播频率与诱导轮转速一致，则是同步旋转空化，空化区波动相对于叶片是固定的，对应的频率为1；空化区波动传播频率小于诱导轮转速(或者传播方向与诱导轮旋转方向相反)则是次同步旋转空化，对应的频率小于1，目前比较常见的是超同步旋转空化现象。Kamijo等第一次在LE-7氧泵中发现超同步的轴振动，并通过一系列实验研究证实该振动是由超同步旋转空化引起。1998年日本HII火箭第8次发射失利，经过分析确认氢泵诱导轮旋转空化诱发的压力脉动频率与叶片的固有频率发生了共振，使叶片发生了疲劳断裂，导致发动机停车，引起这次发射事故。该事件也使得超同步旋转空化现象逐渐引起航天界的重视，Tsujimoto等对诱导轮内的空化不稳定现象开展了全面的研究，对旋转空化的特征进行详细描述，并提出了完整的空化不稳定发生判据。Choi等通过对入口壳体的结构修正实现了旋转空化较好的抑制效果，但是其机理并不清楚。Angelo等研究了多个诱导轮内发生的空化不稳定现象。陈晖等首次在我国某型号发动机试车数据中识别出旋转空化特征频率。Li X等研究了叶顶间隙大小对同步旋转空化的影响。Lettieri等基于可视化实验发现叶尖泄漏涡与旋转空化的发生有密切联系。Kim等利用PIV测量了诱导轮三叶片旋转空化下叶尖的流场，为揭示旋转空化的机理提供了重要的实验数据支撑。Xiang L等结合高速摄像和动态参数采集技术研究了某三叶片诱导轮中发生的空化不稳定现象，清晰揭示了超同步旋转空化工况下空化形态演变过程，填补了国内该领域的可视化试验研究空白。

上述实验研究为揭示旋转空化的发生规律提供重要的支撑，关于其发生机理许多学者提出了不同的观点。Iga等基于二维叶栅的仿真认为旋转空化是由潜在的“旋转失速”元素在特定的工况下被激发引起，但目前为止该解释尚未得到试验结果证实。Kimura等基于数值模拟发现泄漏涡与相邻叶片的前缘相互作用与旋转空化的发生有密切联系。而Tani等数值研究工作证实泄漏涡与旋转空化没有直接关系，叶尖空化溃灭诱发的速度散度是引发旋转空化的原因。

综上可知，目前关于旋转空化的发生机理仍存在较大争议，本文以某二维平板叶栅为研究对象，利用数值模拟研究了超同步旋转空化现象，对超同步旋转空化沿周向传播的机理做出了明确解释。

1 数值计算方法

1.1 控制方程

采用均相流模型模拟空化流动，其假设液相和气相为均匀混合的介质，两相处于局部平衡状态，具有相同的速度、压力等流场信息。利用气相体积分数或质量分数对两相进行区分，同时引入气相体积分数的输运方程来描述气液两相之间的质量交换，方程为

(1)

(2)

式中:混合密度=+(1-)，、分别为液相、汽相密度;有效黏度=+，为湍流黏度，为混合黏度，=+(1-)，、分别为液相、汽相黏度。由于本文研究的介质为常温水，空化过程热效应较弱，可以作等温处理，因此这里不考虑能量方程。

1.2 空化模型

本文采用Singhal模型进行空化流动的仿真，其具体表达式为

(3)

(4)

(5)

1.3 数值计算与设置

图1 计算域Fig. 1 Computation domain

数值计算是基于ANSYS CFX平台，利用CEL语言将上述Singhal模型嵌入软件。常规二方程湍流模型由于高估了空化区尾缘的黏性，无法准确预测空化脱落等非定常现象。目前有研究者通过修正湍流模型中的涡黏系数表达式，实现对空化脱落过程的更准确预测。根据Byungjin等的研究，诱导轮空化不稳定是一种系统不稳定现象，与空化脱落这种局部不稳定现象没有直接关系。这里为了简化流场结构，提炼经典流动模型，不考虑空化区的脱落，选择SST-进行仿真计算，其具体输运方程和相关表达式见文献[20]。计算域网格数量约为2.2×10，具体网格细节如图2所示，非定常计算时间步长为0.000 1 s，根据文献[21-23]计算的经验，这里采用的网格数和计算时间步长满足了无关性要求。计算结果显示壁面值约20，满足壁面函数对网格尺寸的要求。

图2 网格细节Fig.2 Mesh details

2 结果与分析

2.1 定常计算结果

图3为仿真预测的空化性能曲线与实验结果的对比，可以看到，仿真预测的断裂点空化数与实验结果比较接近，但是扬程变化趋势有一定区别，这是由于简化的二维叶栅无法复现三维诱导轮内的复杂流动特性。

图3 仿真预测的空化性能与实验结果对比Fig.3 Comparison of predicted cavitation performance with experimental results

图中标注了旋转空化的发生范围，仿真预测的旋转空化范围大于实验结果，但是二维叶栅的仿真结果可为揭示超同步旋转空化的发生机理提供重要的流场信息。

从气相体积分数分布来看，空化数比较高时，空化区范围比较小，而且3个叶片上空化区均匀对称分布；随空化数降低，空化区范围逐渐扩大；在旋转空化发生范围内，3个叶片上空化区分布明显不均匀，下文将详细分析；直至3个叶片上空化区均发展至喉部位置，造成流道较严重阻塞，开始发生空化性能断裂，整体空化区发展规律与实验结果比较一致，因此整体来看，基于二维叶栅的仿真一定程度能够反映典型的诱导轮内空化流动发展规律。

2.2 非定常计算结果

对叶片吸力面上气相体积分数分布进行面积分可以得到空化区面积，图4给出了=0.09时3个叶片上升力和空化区面积随时间变化关系。可以看到，随时间步推进，升力和空化区面积均由较小的稳定值逐渐扩大，最终形成稳定的周期性波动，而且空化区波动明显沿叶片3-2-1方向传播，而这正是诱导轮旋转方向。

图4 升力和空化区面积随时间步变化关系(σ=0.09)Fig.4 Lifting force and cavitation area versus timesteps(σ=0.09)

图5为实验结果捕捉的超同步旋转空化现象，可以看到空化区波动以一定的速度沿周向传播，如图5(a)所示，假如关注每一圈最短的空化区，明显能看到最短空化区沿着叶片3-2-1的方向传播，大约在第6圈时，最短空化区又回到叶片3，即对于单个叶片而言，空化区波动的传播周期约为5圈。图5(b)更清晰地反映了该趋势，可以看到单个叶片表面空化长度的波动频率为0.18倍转频，而且空化区波动的传播方向与叶轮旋转方向一致，因此从绝对坐标系来看，其特征频率为1.18倍转频，详细的分析见文献[13]。

图5 超同步旋转空化Fig.5 Super-synchronous rotating cavitation

对比实验结果，可以得出仿真预测的正是超同步旋转空化现象。同时注意到升力和空化区波动呈明显的相关性，即空化区面积越大，相应的升力也越大。升力由叶片表面静压值积分所得，因此其绝对值取决于压力面和吸力面的压差，由于空化主要发生在叶片的吸力面，对于空化区面积较大的叶片，其吸力面低压区范围较大，压差也更大，因此叶片升力较大，下文将重点分析叶片表面空化区的变化。

图6为不同空化数下3个叶片无量纲空化区面积(以叶片吸力面面积为基准)随时间变化关系，可以看到，当空化数较高或较低(=0.18、=0.02)时，空化区均未出现波动；当空化数较高时(=0.18)，3个叶片上空化区较小，其发展都不受相邻叶片的限制，呈现自由发展的特征；而当空化数降低至一定程度时(=0.02)，3个叶片上空化区范围都较大，受相邻叶片的限制，难以呈现出沿周向传播的特征。而当空化数处于一定范围内，空化区面积呈明显的周期性波动，对于单个叶片，空化区波动频率分别为28 Hz(=0.15)、26 Hz(=0.14)、22 Hz(=0.09)、9 Hz(=0.05)，该结果也充分表明了空化区是否与相邻叶片存在相互作用是旋转空化现象发生与否的先决条件。

图6 不同工况的空化区波动(仿真结果)Fig.6 Cavitation fluctuation atvarious conditions (simulated results)

根据Iga等的研究，空化区沿周向传播的频率比可定义为PVR=(+)，其中为空化区沿叶片旋转方向的传播速度。上述4个空化数下空化区的周向传播频率比PVR分别为1.336、1.312、1.264、1.108，这与旋转空化的特征频率范围相符，而且PVR随着空化数减小而不断降低，这也是旋转空化的典型特征，进一步证实了数值仿真较好地预测了二维叶栅内发生的超同步旋转空化现象。

进一步以=0.09计算结果为例来揭示空化区变化规律，图7(a)给出了连续一段时间内气相体积分数和液流角(flow angle)的云图，其中间隔时间Δ=0.048 s。首先从图中能直观地看出空化区波动沿叶片3-2-1方向(即叶片旋转方向)传播。叶片前缘的液流角与空化区发展呈现出较好的一致性，即液流角越大，相应的空化区也越大，表明空化区的增长或缩短是液流角的变化引起的。图7 (b)进一步给出了空化区与液流角相互作用关系，其中液流角取自叶尖前缘的监控点，相同颜色曲线代表相同叶片，F_B1和V_B1分别表示叶片1前缘液流角和叶片1表面空化区面积，其他类同。

图7 空化区和液流角的变化(σ=0.09)Fig.7 Variation of cavitation area and flow angle(σ=0.09)

可以看到，二者均呈现出非常规律的周期性，但是波动曲线形状有所不同；对于同一个叶片(如图中蓝色虚线框)，空化区的发展滞后于液流角的变化，即只有当液流角在最大值维持一定时间，空化区才会显著增长；而液流角一旦开始减小，空化区面积随即减小；如果关注两个相邻叶片(如图中黑色虚线框)，可以看到空化区和液流角的变化呈负相关，叶片3上空化区面积从最大值开始减小时，叶片1前缘液流角开始增加；而叶片3空化区从最小值开始增大时，导致叶片1前缘液流角减小，同时注意到，只有当叶片3空化区增大至一定程度时才会导致叶片1前缘液流角大幅减小，因此空化区发展不仅受液流角的影响，而且会反过来显著影响相邻叶片的液流角。

总结空化区面积和液流角相互作用关系可以得出3个重要特征：

1)叶片前缘液流角变化直接影响叶片表面空化区的发展，且两者成正相关；

2)空化区变化会影响相邻叶片的前缘液流角，且两者呈负相关；

3)只有当空化区增长至一定程度才会影响到相邻叶片的前缘液流角。

2.3 超同步旋转空化传播机理

为了更清晰揭示空化区变化对流场的扰动，图8给出不同时刻瞬态计算结果的差值。

图8 不同时刻空化区差值和速度矢量差值分布(σ=0.09)Fig.8 Cavitation area difference and velocity vector difference at various moments(σ=0.09)

例如图8(a)为图7中和两个时刻瞬态计算结果的差值，其中负值代表空化区减小，正值代表空化区扩大。从速度矢量差值来看，大部分流场速度几乎不受影响，但是在空化区尾部，空化区的变化会对流场产生显著的影响。

～时刻，叶片3上空化区扩大，由于空化区内气泡的生长，在空化区尾缘诱发较大的顺时针涡扰动(如图中白色箭头)，当涡扰动靠近喉部附近时，引起正的流向速度分量，导致叶片1前缘液流角大幅减小，叶片1上空化区随即减小，与此同时大量气泡溃灭导致空化区尾缘形成较大的逆时针涡扰动，引起负的流向速度增量，导致叶片2前缘冲角显著增大，因此叶片2上空化区开始增长；～时刻，由于前缘液流角较大，叶片2上空化区不断扩大，但是诱发的顺时针涡扰动尚未达到喉部位置，即不足以影响叶片3前缘的冲角，此时叶片3上空化区仍然在增长直至最大，同时导致叶片1前缘冲角继续减小，相应的空化区也减小，进一步导致叶片2前缘冲角增大，促进了叶片上空化区的发展；～时刻，叶片2上空化区进一步增长，诱发的顺时针涡扰动达到喉部位置，导致叶片3前缘冲角大幅减小，因此叶片3上空化区也显著减小，并导致叶片1前缘冲角增大，由此叶片1上空化区开启了新一轮生长周期。可以看到，上述过程完成了最大空化区从叶片3向叶片2的传播，进一步重复该过程将会形成沿着叶片旋转方向的周期性传播。

进一步总结上述分析，可以得到旋转空化的传播机理，如图9所示。

图9 超同步旋转空化传播机理Fig.9 Propagation mechanism of super-synchronous rotating cavitation

假设在时刻，3个叶片表面空化区呈非对称分布，其中叶片3上空化区最长，叶片1上空化区最短，叶片2上空化区刚好发展至足以影响叶片3的程度，结合文献[6]的研究，可认为当空化区长度达到65%的叶片间距，就会对相邻叶片产生影响，开始出现各类空化不稳定现象。叶片2上空化区持续增长，其空化区尾缘诱发的顺时针涡扰动到达喉部位置，使叶片3前缘冲角大幅减小，相应空化区也大幅减小，同时在空化区尾缘诱发逆时针涡扰动，导致叶片1上空化区开始增大；在时刻，由于冲角增大，叶片1上空化区显著扩大，但是其长度尚未达到65%的叶片间距，叶片2上空化区继续增长，直至达到最大；在此过程中，受叶片2上空化区生长的影响，叶片3上空化区持续减小至最小；在时刻，叶片1上空化区开始影响叶片2以后，继续生长至最大；受其影响叶片2上空化区减小至最小，并进一步促进了叶片3上空化区的发展，由此完成了空化区波动沿叶片旋转方向的传播。

上述讨论清晰揭示了空化区与叶片冲角相互作用过程，解释了空化区波动沿周向传播的机制。图10为旋转空化工况下两个相邻叶片上的空化区结构，可以看到，叶片3几乎不受叶片2上空化区的影响，冲角较大(如图中虚箭头所示)，空化区范围也较大，几乎发展至叶片1前缘，从叶片1上空化区的形态能看出，受上游叶片3空化区的影响，此时叶片1前缘的冲角显著减小，因此实验结果一定程度证实了图9中传播机理的合理性。文献[12]利用PIV测量了发生旋转空化时诱导轮叶尖的液流角分布，首次通过实验证实了较大的空化区会导致其下游相邻叶片冲角为负值，从而抑制泄漏空化的发展；而较小的空化区会增大下游相邻叶片的冲角，促进叶尖泄漏空化的发展，进一步证实了图9中的传播机理。

图10 超同步旋转工况下空化结构Fig.10 Cavitation structure at SSRC condition

3 结论

本文利用数值模拟研究了二维叶栅中的超同步旋转空化现象，其中叶栅取自真实诱导轮，根据实验结果中旋转空化发生工况的转速、流量换算成来流速度和液流角作为边界条件，进行非定常计算，得出结论如下：

1)仿真预测的空化断裂点与实验结果比较接近，且叶栅内空化发展规律与实验结果一致，表明建立的二维叶栅计算模型能够较可靠地反映真实诱导轮内空化流动特征；

2)在一定空化数范围内3个叶片上空化区呈现规律的周期性波动，波动频率在1.1～1.4之间，且随空化数减小而降低，空化区波动传播方向与诱导轮旋转方向一致，证实仿真比较准确地捕捉到超同步旋转空化现象；

3)计算了不同时刻瞬态流场差值，分析了空化区变化对流场的扰动。发现空化区与叶片前缘冲角成正相关，同时空化区扩大会在尾缘诱发顺时针涡扰动，空化区减小会在尾缘诱发逆时针涡扰动，当扰动达到叶片喉部位置时，会对相邻叶片产生影响，因此空化区扩张会导致相邻叶片冲角减小，相应空化区也减小，空化区缩减会导致相邻叶片前缘冲角增大，相应空化区也扩大，在空化区和叶片前缘冲角这种相互作用机制下，空化区波动以一定速度沿周向传播。