轴对称后向台阶不稳定流动及压强振荡数值研究①

张翔宇，何国强，刘佩进

(西北工业大学燃烧、热结构与内流场重点实验室，西安 710072)

0 引言

后向台阶流动具有几何结构简单、分离点固定和来流条件容易控制等特点，是研究湍流边界层分离和再附流动的重要形式之一。目前，国内外对层流至湍流条件下的流动进行了大量的试验和数值模拟［1－3］，试验大多在风洞和水洞中进行，测量手段有热线、热膜探针和PIV等。目前，为提高战术固体火箭发动机性能，大多采用高能复合推进剂及复杂装药构型，尤其是翼柱形装药的固体火箭发动机在翼槽烧完时，也会形成具有自身特点的后向台阶流动。这种台阶结构是轴对称的，具有一定倾斜角度，并且由于侧向加质的作用，近壁面速度垂直于轴线，流场中主流速度较大，一般为完全发展的湍流流动。前柱后翼形装药形成的后向台阶会使流动在台阶附近形成剪切层，进而产生转角涡脱落现象，这是固体火箭发动机中3种典型的旋涡脱落现象之一［4－5］。当涡脱落频率与声场某阶固有频率接近时，就可能形成共振，这种声涡耦合产生的共振是大长径比固体火箭发动机燃烧室压强振荡的源头。Flandro和Jacobs［6］最早指出，流动经过障碍物时，形成的旋涡脱落会将部分平均流能量转化为声能，从而产生声能增益。Mason 和 Folkman［7］以及 Brown［8］等得出结论，旋涡脱落可诱发燃烧室的轴向基频。Culick 和 Magiawala［9］、Dunlap［10］通过实验研究表明，当旋涡脱落频率与燃烧室声场固有频率接近时，会引起压强振荡现象。Vuillot［11］、Lupoglazoff［12］设计了一种转角涡脱落致压强振荡的发动机，并通过试验及数值模拟证明了转角涡脱落产生的声涡作用，可引起与声场固有频率相耦合的压强振荡现象。Shanbhogue［13］等通过冷流试验，研究了转角涡脱落引起的压强振荡问题，得到了St数随流场结构参数变化的规律。

本文对战术发动机内形成的轴对称后向台阶流动进行模型化处理，设计了一种无需外加扰动的旋涡脱落致压强振荡的流场结构。运用稳定性理论得到声涡耦合模态预测方法，并采用大涡模拟(LES)方法，针对宽来流条件进行数值模拟。通过对比计算结果，验证了模态预测方法的准确性，分析了转角涡脱落产生的一般过程及旋涡多尺度运动规律。通过分析流场关键位置的压强振荡，研究了气动声学特性及旋涡运动与压强振荡的关系。

1 计算构型设计

1.1 声涡耦合模态预测方法

对于某一特定构型的发动机燃烧室，如果其平均流边界层的不稳定频率处于轴向声模态范围内，就有可能形成旋涡脱落引起的压强振荡现象。Lupoglazoff等通过求解不可压Orr－Sommerfeld方程，得到了后向台阶结构燃烧室的径向速度分布曲线。分析表明，具有拐点的双曲正切形速度分布是剪切层发展为离散旋涡的重要条件。应用空间增长率理论，分析流体动力学稳定性问题，旋涡波动量为

其中，α =αr+iαi∈C，ω∈R，B 是积分常数，f(y)、φ(y)、π(y)是涡量波动分布。任意给定ω，计算求解本征问题的解为α=αr+iαi，并引入以下无量纲角频r=2πfδ/ΔUr+iai。

通过处理可得到cr与实际频率的关系曲线。分析表明，在cr的临界值(cr)cr以上，边界层是稳定的。此外，计算中还需定义以下参数:无量纲旋涡位移速度将边界层厚度表示成燃烧室直径D的函数δ=KDD;将剪切层速度差ΔUe表示成平均流速度U的函数ΔUe=KUU=KU(a0M)。最终形成后向台阶不稳定流动需同时满足以下2个条件:

(1)流场第p阶固有模态fp位于剪切层不稳定范围内，即无量纲频率(cr)p小于临界频率(cr)cr;

本文采用半封闭声腔模型fp=(2p－1)a0/(4L)。

(2)旋涡运动距离l应至少包含一个涡波长度λv，即 λv＜l。

代入KV、KU表达式，对于给定构型，可认为二者为常数。联立以上二式，整理得:

式中的关键参数参考文献［11］中的结果，代入数值计算得到:

式(4)就是后向台阶流动中，给定长径比及台阶位置条件下主流速度与转角涡脱落可能诱发的声模态之间的关系。需指出的是，上述结论仅是发生不稳定的必要条件。尽管如此，经验表明，轻微的不稳定趋势，就有可能引发声涡耦合现象。这也就是为什么当满足必要条件时，就很容易发生不稳定，而制造完全稳定的发动机却十分困难。

1.2 研究对象的几何结构

为得到后向台阶流动涡脱落及气动声学的一般规律，需对计算模型进行适当简化。忽略侧向加质的影响，将后向台阶扩张角θ确定为90°。长径比参考战术发动机燃烧室特征，并进行适当缩比，整个流动区域为轴对称结构，总长1 200 mm;将后向台阶布置在流场中央，台阶前后段L1、L2等长，且均为600 mm;单侧后向台阶高度H为25 mm;台阶前后管路直径D1、D2分别为50 mm和100 mm。流动区域示意图见图1。

2 数值研究方法

在后向台阶流动中，台阶后缘会产生类似平板剪切层的大涡序列，整个流场中的旋涡结构十分复杂。工程应用中，普遍采用的雷诺时均(RANS)方法难以模拟这种复杂的流动特性和旋涡结构。因此，本文采用大涡模拟(LES)方法，对多尺度旋涡的运动过程进行数值模拟。

2.1 控制方程离散及求解

LES在湍流数值模拟中只计算大尺度脉动，将小尺度脉动对大尺度运动的作用建立亚格子模型，过滤小尺度脉动，并得到亚格子应力项，亚格子模型的选择至关重要。后向台阶流动中，剪切层与近壁面流动关系密切。因此，选用对近壁面区域进行修正的WALE(Wall－Adatping Local Eddy－viscosity)亚格子模型，从而避免了Smagorinsky模型难以准确预测近壁区流动及湍流转捩现象的不足［14］。计算中，不考虑化学反应以及凝相粒子的影响，控制方程采用纯气相可压N－S方程。空间参数离散采用具有二阶精度的BCD(Bounded Central Differencing)格式，时间上采用二阶隐式格式，时间步长为10－5s。为提高计算的收敛性，先采用空间二阶离散格式的RANS方法得到流场的稳态解，然后采用上述LES方法求解非定常流动。

2.2 计算网格与边界条件

由于LES受网格密度及人为划分的影响显著，须将可能出现大尺度回流区与涡脱落区的计算域进行高精度的空间网格剖分。在台阶前后、转角及壁面附近进行加密。计算采用二维轴对称方法，通过网格收敛验证，总网格数约为55 000(图2)。介质为常温常压可压气体，粘性模型采用Sutherland定律，入口边界给出空气质量流率、来流总压和初始温度，出口边界给出静压和回流温度，壁面采用无滑移绝热壁面条件。该计算模型是一个进口为声学反射边界的半封闭声腔。

2.3 计算方法验证

后向台阶流动数值模拟的关键是对旋涡运动及压强振荡的准确模拟。陈晓龙等［15］通过搭建冷流试验系统，在高压空气中加入示踪粒子，对转角涡脱落进行了可视化研究。

图3是应用上述数值方法对流场进行的二维大涡模拟结果与冷流试验中的高速摄影结果。在激光照射下，混合了粒子的流向气体呈亮白色，而照片中的“黑洞”现象是由于旋涡中高速运动气体将粒子甩出造成的。可看出，本文采用的数值方法得到的旋涡位置和结构与试验结果较为吻合。

此外，还对VKI实验室针对Ariane 5 P230助推器缩比模型的冷流试验进行了数值模拟校验。该实验器为轴对称结构［16］(图4)，通过高频压强传感器测量了流场不同位置的压强振荡。图5是计算域网格划分，与后向台阶流动数值模拟采用的网格剖分尺寸一致。

发动机头部压强振荡频率的试验结果和数值模拟结果见表1。由表1可看出，采用的二维轴对称LES数值方法得到的压强振荡频率与试验结果较为接近，前四阶频率与实验结果间的最大误差小于6%，可较好地对湍流场气动声学问题进行研究。

表1 LES结果与试验结果对比Table 1 Results of LES and cold experiments Hz

3 后向台阶流动的旋涡结构

3.1 转角涡脱落的形成过程

在后向台阶流动中，一般以雷诺数作为特征参数进行研究，定义Re=UH/ν。其中，U为来流速度;H为台阶高度;ν为介质的运动粘度。为方便分析流场特征，需将流动区域进行划分。

以Re=11 288(Ma=0.02)的算例为例，处理稳态计算结果绘制流函数分布图(图6)，在紧靠台阶后缘突扩区域形成一个低速角涡区，在其下游是一个狭长的近壁面回流区，接着流动再附着并继续发展。通过上述分析，可得如图7所示的后向台阶湍流分离流动流场的基本划分方法。

连续介质经过后向台阶时，会在台阶下游形成典型的流动剪切层，随着剪切层向下游发展，形成流向涡，这种大尺度旋涡运动对流动稳定性及气动声学问题都具有重要影响。

LES可获得旋涡结构的运动演化过程，通过整个计算区域内的等涡量线图可观察剪切层的运动特征，包括旋涡的产生、脱落和合并过程。

图8描述了流动稳定后转角涡脱落形成的典型过程。图8(a)是一个旋涡脱落后新旋涡产生的初始过程，随着剪切层在后向台阶下游持续发展，由于流动自身的不稳定性，剪切层沿径向摆动，同时受到近壁面回流区的卷吸作用，沿流向运动了一定距离后，剪切层向壁面方向发展(X=0.67 m)。

远离台阶后，流动速度逐渐降低，向壁面移动的部分剪切层又受到后续发展的剪切层的挤压作用，此时新的涡核基本形成(X=0.68 m)，如图8(b)所示。

此后，随着旋涡向下游继续发展，新的涡量注入涡核，旋涡体积增加。当旋涡增长到一定体积时，完整的涡结构基本形成(X=0.69 m)，其受剪切层的影响也逐渐减小(图8(c))。

此后，逐渐形成了稳定的涡脱落(X=0.70 m)，涡量呈离散式分布(图8(d))。下文分析中，采用的旋涡脱落频率按上述过程进行判定。

3.2 转角涡脱落的形成过程

后向台阶流动的轴向速度分布特点是台阶下游流速低，且距离台阶越远，流速越低。由于流道结构突扩及壁面的粘性作用，管流中心流速大于近壁面处流速。由于湍流的多尺度性和随机性，旋涡运动也呈现类似特点。流场中旋涡运动的典型过程如图9所示。

台阶后缘由剪切层发展出离散的小尺寸旋涡(X=0.64 m)，其在向下游发展过程中，运动速度逐渐降低，被后续形成的旋涡追赶，并发生初始合并过程(X=0.74 m)，随着壁面涡量及离散涡量的注入(X=0.9 m)，进而形成大尺度旋涡向下游运动，直至离开流场(X=1.0 m)。

通过分析涡量云图，将不同来流条件下的旋涡多尺度运动频率列于表2。由于旋涡运动的随机性，取旋涡运动3～5个周期的平均值。可看出，随着来流速度增加，各尺度的旋涡运动频率均相应增加。

表2 多尺度旋涡运动频率Table 2 Frequency of the multi-scale vortex motion

4 压强振荡及气动声学特性

后向台阶流动中的旋涡脱落和旋涡合并过程，会引发流动的不稳定性。Powell方程指出［17］，在流速很低的条件下，仅存在涡的地方才会产生声。当旋涡运动与声场固有频率接近时，就可能发生声涡耦合现象，造成流场中持续的压强振荡。LES方法不但可模拟不同来流条件下的旋涡发展过程，还可得到流场各处压强随时间的变化规律。通过FFT分析，获得压强波动的特征频率和振幅。计算中，记录流场中关键部位的压强值，压强监测点分布如图10所示。

图10中，1位于流场进口壁面上，试验中高频振荡压强传感器一般布置在燃烧室头部，是研究燃烧不稳定问题的重要布测点;2位于台阶转角处;3表示布置在角涡区和回流区的监测点;4表示布置在再发展区的监测点。

对于不规则的声腔结构，ANSYS软件的模态分析，可确定固有频率及振型。试验证实，该方法对频率预测准确性很高。ANSYS分析声场前5阶轴向模态分别为 99.5、182.6、381.5、464.7、663.6 Hz。对各工况下数值模拟记录的监测点1压强值进行FFT分析，得到结果如图11所示，横坐标表示压强振荡频率，纵坐标表示振幅(波动量与平均压强之比)。将流场的结构及各工况来流参数代入式(4)进行计算，得到诱发声场前5阶轴向振型所需的来流Ma范围(表3)。并将上述结果列于表4。

表3 声涡耦合模态预测Table 3 Prediction of acoustic-flow coupling mode

表4 数值模拟与声涡耦合模态预测结果对比Table 4 Results of the mode prediction and LES

对比ANSYS模态分析结果与数值模拟结果可知，本文设计的后向台阶流动出现了较为明显的压强振荡现象。在设计来流范围下，流动最多可诱发声场的前5阶轴向频率。随着来流速度增加，各阶轴向振型相继被激发，同时压强振荡振幅也逐渐增大，在流场头部前两阶轴向频率的振幅最大，最大振幅接近平均压强的5‰。声涡耦合模态预测方法对不同速度下可能诱发的轴向声模态预测结果与数值模拟结果十分吻合，当来流速度接近激发某阶固有频率所需的速度下限时，该频率逐渐被激发。分析可知，激发轴向基频所需的来流速度最低。工况1对应的来流Ma=0.01，此时没有发生声涡耦合现象。因此，通过该方法，也可对后向台阶稳定流动的速度上限进行预估。

对比相同工况下流场中不同压强监测点的FFT分析结果发现，监测点2、3与1的压强振荡频率一致，但各阶频率的振幅更低，流场头部压强振荡最为强烈。处于再发展区的监测点4各阶频率的振幅更小，并出现了小于轴向基频的低振幅的宽峰振荡，此频率随来流速度增加而逐渐增加。

工况2～5监测点4的FFT分析结果如图12所示。结合表2给出的多尺度旋涡运动结果分析发现，低于轴向基频的频率与大尺度旋涡运动频率接近。在来流速度较低的情况下(工况2)，再发展区中大尺度旋涡运动频率对压强振动的影响较大。随来流速度增加，声涡耦合现象逐渐增强，当声场某阶固有频率被激发后，压强振荡频率与声场固有频率相耦合，流场中的压强振荡以声场轴向各阶频率为主。当大尺度旋涡运动频率与基频接近时，会造成基频的多峰现象(工况4、5)。表2的结果也表明，当近台阶小尺寸旋涡形成频率接近声场某阶固有频率时，该频率及更低阶频率就有可能被激发。这说明后向台阶流动中，转角处旋涡脱落及旋涡运动是诱发声涡耦合现象的源头。

5 结论

(1)后向台阶下游形成剪切层，由于流动不稳定性和近壁面回流区的作用，剪切层会发展为涡量离散分布的旋涡结构。随来流速度增加，旋涡脱落频率增大，脱落位置也更靠近后向台阶。小尺寸旋涡被后续形成的旋涡追赶，并发生初始合并过程，随壁面涡量及离散涡量的注入，最终形成大尺度旋涡，直至离开流场。

(2)从稳定性理论推导出后向台阶流动的声涡耦合模态预测方法，可对不同速度下可能激发的声场频率进行合理预测。当来流速度增加至激发某阶固有频率所需的速度下限时，该频率逐渐被激发。通过该方法，还可对后向台阶稳定流动的速度上限进行预估。

(3)流场头部压强波动振幅最大，频率与声场轴向模态相耦合。再发展区出现低于轴向基频的宽峰振荡受大尺度旋涡运动的影响，该频率与基频接近时，会造成基频的多峰现象。当近壁面小尺寸旋涡形成频率接近声场某阶固有频率时，该频率及更低阶频率就有可能被激发。

［1］Armaly B F，Durst F，Pereira J C F，et al.Experimental and theoretical investigation of backward－facing step flow［J］.Journal of Fluid Mechanics，1983，127:473－496.

［2］Lee T，Mateescu D.Experimental and numerical investigation of 2D backward facing step flow［J］.Journal of Fluids and Structures，1998，12:703－716.

［3］刘应征，朴英守，成亨镇.后台阶分离再附湍流流动的实验研究［J］.上海交通大学学报，2005，39(5):810－817.

［4］Fabignon Y，Dupays J.Instability and pressure oscillations in solid rocket motors［J］.Aerospace Science and Technology，2003，7:191－200.

［5］Telaro M，Avio Spa，Stella F，et al.Pressure oscillations in solid rocket motors［R］.AIAA 2006－4423.

［6］Flandro G A，Jacobs H R.Vortex－generated sound in cavities［R］.AIAA 1973－1014.

［7］Mason D R，Folkan S L，Behring M A.Thrust oscillations of the space shuttle solid rocket booster motor during static tests［R］.AIAA 79－1138.

［8］Brown R S，et al.Vortex shedding studies［R］.AIAA 80－1092.

［9］Culick F E C，Magiawala K.Excitation of acoustic modes in a chamber by vortex shedding［J］.Journal of Sound and Vibration，1979，64(3):455－475.

［10］Dunlap R，Brown R S.Exploratory experiment on acoustic oscillations driven by periodic vortex shedding［J］.AIAA Journal，1981，19(3):408－409.

［11］Vuillot F.Vortex－shedding phenomena in solid rocket motors［J］.Journal of Propulsion and Power，1995，11(4):626－639.

［12］Lupoglazoff N ，Vuillot F.Numerical simulation of vortex shedding phenomenon in 2D test case solid rocket motors［R］.AIAA 1993－3006.

［13］Shanbhogue S J，et al.Aeroacoustic instability in a solid rocket motor with finocyl grain［R］.AIAA 2003－4632－4636.

［14］张兆顺，崔桂香，许春晓.湍流理论与模拟［M］.北京:清华大学出版社，2005.

［15］Chen Xiao－long，He Guo－qiang，Liu Pei－jin.Experimental investigation on pressure oscillation induced by flow instability in SRM［C］//60th IAC，Daejeon，Korea，2009.

［16］Anthoine J.Experimental and numerical investigations of nozzle geometry effect on the instabilities in solid－propellant boosters［R］.AIAA 2000－3560.

［17］Powell A.Theory of vortex sound［J］.J.Acoust.Soc.Am.，1964，36:177－195.