童福林， 唐志共,*， 李新亮， 吴晓军， 朱兴坤

1.中国空气动力研究与发展中心 计算空气动力研究所， 绵阳 621000 2.中国科学院 力学研究所 高温气体动力学国家重点实验室， 北京 100190

压缩拐角激波与旁路转捩边界层干扰数值研究

童福林1， 唐志共1,*， 李新亮2， 吴晓军1， 朱兴坤2

1.中国空气动力研究与发展中心 计算空气动力研究所， 绵阳 621000 2.中国科学院 力学研究所 高温气体动力学国家重点实验室， 北京 100190

为了研究激波与旁路转捩边界层的干扰机理，采用直接数值模拟(DNS)方法对来流马赫数Ma∞=2.9，24° 压缩拐角内激波与转捩边界层的相互作用进行了系统的研究。考察了旁路转捩干扰下压缩拐角内分离区形态和激波波系结构的典型特征。比较了转捩干扰与湍流干扰流动结构的差异，并分析了造成差异的原因。研究了拐角内转捩边界层的演化特性，探讨了转捩干扰下脉动峰值压力和峰值摩阻的分布规律及形成机制。研究结果表明：相较于湍流干扰，两侧发卡涡串的展向挤压使得分离区起始点以V字型分布，且分离激波沿展向以破碎状态为主，激波脚呈现多层结构；拐角内的干扰作用急剧加速了边界层的转捩过程；转捩干扰下的拐角内峰值脉动压力以单峰结构出现在分离区的下游，同时干扰区内的强湍动能和高雷诺剪切应力使得其局部峰值摩阻系数要高于湍流干扰。

压缩拐角； 激波/边界层干扰； 旁路转捩； 脉动压力； 摩阻； 直接数值模拟

激波/边界层干扰是高速飞行器上普遍存在的典型流动问题，它会导致飞行器物面边界层出现大尺度非定常分离、强压力脉动以及局部干扰峰值压力和热流，严重影响飞行器的飞行安全和气动性能[1-2]。

自20世纪50年代以来，大批学者对该问题进行了广泛的研究，在非定常脉动压力和热流的预测[3-4]、激波/边界层干扰的控制[5]以及分离激波的低频运动特性[6-7]等方面都进行了详细的分析研究。但是，以往大部分的工作都是集中在激波/层流或激波/湍流干扰问题上，而对激波/转捩边界层干扰的风洞试验研究和数值模拟都相对缺乏。造成这种现象主要有两方面的历史原因，一方面是激波/转捩干扰在试验或是数值上都较难实现，风洞背景噪声、模型粗糙度等不确定因素都有可能导致边界层转捩机理发生质的改变，而在数值模拟方面，如何构造接近真实状态的转捩入口来流条件是制约激波/转捩干扰研究的主要困难；另一方面，在工程应用上曾一度认为激波/转捩干扰是激波/层流干扰和湍流干扰两类流动的中间状态，其流动结构、压力和热流分布规律等都应在两者之间，因而对激波/转捩边界层干扰的独特性和重要性缺乏足够认识[8]。实际上，转捩干扰问题要比湍流干扰或是层流干扰复杂得多，如果转捩出现在边界层再附点附近，转捩与边界层分离/再附同时发生，相互影响，可能会导致转捩干扰峰值热流要远高于湍流干扰[9]。因此，深入研究激波/转捩边界层的相互作用，无论是在工程应用还是理论研究中都具有十分重要的意义。

目前，激波与转捩边界层干扰的研究多以风洞试验手段为主，在数值模拟方面开展的工作还相对较少。风洞试验模型以入射激波与平板边界层(或轴对称旋成体)和压缩拐角两类典型构型为主。按照转捩干扰的发生方式[10]，可以分为触发转捩干扰(Triggered Transition Interaction)和自然转捩干扰(Natural Transition Interaction)。对于触发转捩干扰，转捩是由于激波与边界层干扰引起的；而自然转捩干扰则是通过控制激波与自然发生的转捩边界层相互干扰位置而形成。对于压缩拐角，由于激波是因流场中逆压梯度的存在而诱导产生的，试验较难对其干扰位置做到精确控制，只能通过改变来流参数，使得转捩在拐角内发生，并与诱导激波发生干扰，因此压缩拐角试验多以触发转捩干扰为主。对于入射激波干扰，触发转捩干扰和自然转捩干扰则都相对容易实现。Erdem等[10]的试验研究进一步比较了触发转捩干扰与自然转捩干扰的差别，研究发现，触发转捩干扰引起的局部压力和热流峰值要高于自然转捩干扰。Sandham和Schulein[9]试验测量了马赫数Ma=6时入射激波与转捩边界层的自然转捩干扰作用。结果表明，在相同雷诺数下转捩干扰的峰值热流要高于湍流干扰，尤其是在干扰的相对位置位于转捩的初期状态下。Schrijer[11]的高超声速钝锥-裙试验表明，如果转捩出现在边界层再附点附近，转捩干扰峰值热流也要高于湍流干扰。Xie等[12]利用风洞试验对压缩拐角激波/转捩边界层干扰的物理现象也进行了研究，分析了转捩干扰下的一些动力学特性。

在转捩干扰的数值研究方面，Teramoto[13]采用大涡模拟计算了转捩边界层与激波的相互作用。在验证了计算结果的前提下,发现在转捩区内存在大尺度的纵向涡对、低速条带和发卡涡结构。Krishnan和Sandham[14]采用大涡模拟方法研究了飞行马赫数Ma=8的高超声速飞行器进气道内的激波和转捩边界层的相互作用。结果表明，处于转捩状态的边界层在与激波相互作用后会加快边界层的转捩。近年来，直接数值模拟(Direct Numerical Simulation，DNS)由于不引入任何湍流模型或亚格子应力模型，同时能够直接求出所有尺度的湍流脉动量时空演化信息，具有较高的可靠性，因此，DNS方法逐渐在激波/转捩边界层干扰机理的数值研究中得到应用。Tokura和Maekawa[15]对入射激波与空间发展的转捩边界层干扰问题进行了DNS计算，重点关注了转捩干扰下的低频特性以及其物理机制。此外，Ludeke和Sandham[16]还采用DNS方法研究了压缩拐角内流动在不同雷诺数和压缩角下转捩的发生过程。

本文采用DNS方法对来流马赫数Ma∞=2.9，24°压缩拐角内激波与旁路转捩边界层的相互作用进行系统研究。在笔者前期的论文[17]中，针对激波/转捩边界层干扰对分离泡结构的影响展开了研究，给出了转捩对分离泡尺度及形态影响的初步规律。与以往压缩拐角内激波/转捩边界层干扰研究的不同之处在于，本文在压缩拐角入口上游的平板采用添加吹吸扰动的方法来激发流动发生转捩，通过控制平板长度，使得进入拐角角部干扰区的转捩边界层处于旁路转捩过程的初期状态。因此，此时拐角内的转捩干扰类型属于自然转捩干扰，而非以往的触发转捩干扰。此外，为了进一步比较分析转捩干扰和湍流干扰下流动结构的差异性，还对相同来流条件下的湍流干扰进行DNS。为了便于比较和验证结果，选取的计算参数尽量接近Bookey等[18-19]的试验以及Wu和Martin[20]的DNS。

1 数值方法

控制方程采用一般曲线坐标系下的三维无量纲化后的可压缩Navier-Stokes方程组。流场变量采用无穷远来流参数进行无量纲化，长度变量采用单位mm进行无量纲化，方程具体形式如下：

(1)

式中：ξ、η和ζ分别为曲线坐标系流向、法向和展向；t为时间；守恒变量U和通量F的表达式分别为

其中：ρ为密度;u、v和w分别为流向、法向和展向速度;p为压力;e为总内能;其余参数的表达式分别为

J-1为坐标系的Jacobi变换系数，黏性应力和热流通量的计算式为

(2)

(3)

方程组中其他变量的具体计算公式参见文献[20]。同时，通量G和H的计算过程与F类似，这里不再赘述。

方程组中无黏项采用Martin等[21]优化构造的WENO(Weighted Essetntially Non-Oscillatory)格式以及Steger-Warming流通量分裂方法求解。该格式是在8阶中心格式网格基架点上优化得到的，由于采用了对称网格基架点，格式具有很高的波数分辨率和较低的数值耗散。黏性项采用8阶中心差分格式进行离散，时间推进采用3阶Runge-Kutta方法计算。

2 DNS设置

计算模型为含吹吸扰动带的上游平板与24° 压缩拐角，如图1所示。图中坐标系原点取为压缩拐角的拐点，其中x、y和z分别对应流向、法向和展向方向。

图1 压缩拐角计算模型示意图 Fig.1 Illustration of compression ramp computation model

计算域的流向跨度Lx由上游平板的流向跨度(Lx1+Lx2+Lx3)和压缩拐角流向跨度(Lx4+Lx5)两部分组成，其中压缩拐角的流向跨度包含角部区域(Lx4=35 mm)和斜面区域(Lx5=51.5 mm)。对于转捩干扰和湍流干扰工况，两者流向计算域的差别在于上游平板流向跨度的不同。上游平板流向计算域由Lx1、Lx2和Lx3三部分组成，其中Lx1为层流入口剖面距吹吸扰动带起始点的距离，Lx2为平板吹吸扰动带的长度，Lx3为吹吸扰动带终点位置距压缩拐角入口(x=-35 mm)的距离。计算中Lx1均为30 mm，Lx2均为20 mm。对于转捩干扰工况，Lx3=65 mm，对于湍流干扰工况，Lx3=250 mm。此外，对于转捩干扰和湍流干扰，计算域法向高度均为Ly=35 mm，展向宽度均为Lz=14 mm。

来流马赫数Ma∞=2.9，基于单位长度的来流雷诺数Re∞=5 581.4 mm-1，来流静温T∞=108.1 K，壁温Tw=307 K。计算域入口为层流边界层，实际计算过程中首先对相同来流条件下的二维平板层流边界层进行数值模拟，取距平板前缘200 mm处的层流解作为计算域层流入口条件。出口边界统一使用超声速出口无反射边界条件。物面边界为无滑移条件和等温壁。上边界取为简单无反射边界条件，展向为周期性条件。上游平板吹吸扰动带内壁面法向扰动速度分量与文献[22]相同，扰动形式为多频正弦波的壁面吹吸扰动。依据Gao[23]和Li[24]等的分析，在该多频正弦波扰动形式下触发的转捩类型为旁路(Bypass)型转捩，该转捩发展过程跳过了自然转捩中的不稳定扰动波的线性增长阶段而直接进入非线性增长阶段，相关详细的机理可进一步参阅文献[23-24]。

转捩干扰和湍流干扰工况分别对应拐角入口为Bypass转捩过程的非线性发展阶段(也是该转捩类型的初始阶段)和完全湍流阶段。表1给出了计算中的自由来流参数及不同工况下拐角入口(x=-35 mm)处的边界层参数，包括了边界层动量厚度θ、名义厚度δ、摩阻系数Cf及动量厚度雷诺数Reθ。为了便于比较，表1中还给出相同来流条件下Wu和Martin[20]的湍流干扰DNS以及Bookey等[19]的风洞试验拐角入口处边界层参数。要特别说明的是，表1中Tran表示拐角入口边界层处于转捩状态，对应为转捩干扰工况，Turb表示拐角入口边界层处于完全湍流状态，对应为湍流干扰工况，下文类似。由于风洞试验没有详细给出拐角上游的转捩过程以及转捩的触发方式，湍流干扰工况下拐角入口流动参数只能尽量接近试验，无法做到完全一致。另外，Wu和Martin[20]的数值模拟采用的是循环重构方法生成完全湍流入口条件，因此入口参数也存在一定的差异。从表1中可看出，湍流干扰工况入口处的边界层名义厚度和动量厚度与风洞试验及前人DNS基本接近，但摩阻系数要高些。在笔者的前期研究中，对湍流干扰工况下干扰区内的速度型、压力和摩阻分布进行了结果验证。研究表明，压力分布与Wu和Martin[20]的DNS数据基本重合，两者均在Bookey等[19]试验数据误差带(5%)的范围内。同时，由于角部入口处来流摩阻偏高，抑制了角部分离泡的发展，计算得到的分离区起始点相比靠后，但再附点位置与试验值及DNS相比基本重合，而且拐角斜面上摩阻系数无论是在分布规律还是量级上都与Wu和Martin[20]的DNS是一致的，具体结果可进一步参考文献[17]。对于转捩干扰工况，拐角入口处的边界层名义厚度约为湍流干扰的1/3，摩阻系数约为湍流干扰的2/3。

表1 自由来流参数及拐角入口(x=-35 mm)的边界层参数Table 1 Free-stream parameters and boundary layer parameters at ramp inlet (x=-35 mm)

表2分别给出了转捩干扰和湍流干扰工况的计算网格参数,表中Nx、Ny和Nz分别为计算域的流向、法向和展向网格点数。图2为转捩干扰工况计算网格的示意图，其中流向网格每隔10个网格点显示。网格采用与文献[19]相同的代数方法生成，流向网格在拐角角部区域(-35 mm≤x≤35 mm)内密集分布，法向网格往壁面附近进行了加密处理，展向网格均匀分布。以x=-35 mm 处的壁面量为度量(见表2)，各工况的流向网格尺度Δx+约为4.0，壁面法向网格尺度Δy+均小于0.5，展向网格尺度Δz+约为5.0。由于角部分离区对数值黏性的敏感性，计算采用的网格尺度远小于平板边界层湍流的直接数值模拟[23-24]。

表2 计算网格Table 2 Computational grids

图2 转捩干扰工况计算网格Fig.2 Computational grid for transitional case

3 结果分析与讨论

3.1 流场结构

压缩拐角流场内蕴涵着诸多复杂的流动现象，如边界层分离和再附、激波/激波干扰、激波与湍流的相互作用以及分离激波的低频振荡运动等。对于转捩干扰工况，由于转捩边界层的演化与激波干扰、边界层分离/再附等现象同时发生，相互作用及影响，使得转捩干扰下的流场结构与以往湍流干扰存在着明显的差异。

3.1.1 分离区形态

图3给出了时间平均后的物面极限流线分布情况。为了便于比较说明，图中还给出了湍流干扰下的结果。图中红色实线代表分离区起始点沿展向的分布，绿色虚线代表拐角拐点。从极限流线的整体分布趋势上来看，转捩干扰下分离泡起始点在展向上呈现“V”字型结构，两侧的分离区长度要比中间位置小得多，前者出现在x=-5 mm 附近，而中间位置则出现在x=-15 mm附近。湍流干扰下分离泡起始点沿展向的分布则较为均匀，约在x=-16 mm附近。另外，从转捩干扰的结果中还可以清楚看到，此时分离区内存在着两个尺度和方向各不相同的小分离泡，该流动结构则与湍流干扰完全不同。

图3 时间平均物面极限流线分布情况Fig.3 Distributions of time-averaged surface limiting streamlines

研究表明，造成上述分离区形态的差异，主要是由于转捩干扰下拐角入口处两侧发卡涡串在拐点附近的展向“挤压”作用，进一步加速了分离泡沿展向中部层流区往上游的发展，使得分离泡在展向上出现“V”字型结构，同时诱导了分离区内小尺度分离泡的出现,具体的流动机理分析可进一步参阅笔者的前期研究[17]。

3.1.2 激波结构

为了研究转捩干扰下的激波结构，图4分别给出了某无量纲时刻下z=2 mm和z=7 mm截面内角部区域的瞬时纹影图。图4中z=2 mm对应为左侧包含发卡涡串的xOy截面，而z=7 mm 对应为展向的中截面，位于两侧发卡涡串中间的层流区。同时，图5给出了相同时刻湍流干扰下各展向截面瞬时纹影图。

图4 转捩干扰下z=2,7 mm的xOy截面内瞬时纹影图Fig.4 Instantaneous numerical schlieren plots for z=2,7 mm plane in transitional case

图5 湍流干扰下z=2,7 mm的xOy截面内瞬时纹影图Fig.5 Instantaneous numerical schlieren plots for z=2,7 mm plane in turbulent case

为了增强流动的显示效果，构造了变量Ns[20]，即

Ns=C1exp(-C2(Φ-Φmin)/(Φmax-Φmin))

(4)

转捩干扰下的激波结构则与湍流干扰有着明显的差异。首先，由于转捩干扰下的分离区尺度要小于湍流干扰，这使得转捩干扰下的分离激波更为靠近折角物面。同时，尽管在转捩干扰分离区下游边界层外缘也延伸出了多道激波束，但从主激波的大幅度变形来看(如图4中粉色箭头所示)，此时激波束的强度要明显强于湍流干扰情况。另外，由于转捩干扰下入口的来流边界层为2.3 mm，约湍流干扰下的1/3，此时转捩边界层内的强脉动对分离激波脚产生了大尺度的弯曲和变形，激波脚呈现出与湍流干扰完全不同的多层分布结构。

对于转捩干扰工况，如图6(a)所示，分离激波在脚部附近出现了大范围破裂，激波阵面在展向上破碎成了很多小尺度的激波等值面，而且计算域左右两侧的破裂形式和破碎程度也有较大的差异。拐角上游入口处发卡涡串的尺度与激波阵面的破碎程度及形式有较强的相关性，即大尺度发卡涡串对应激波阵面大范围的强破碎(见图6(b)右侧)。同时，由于转捩干扰下分离区形状沿展向呈”V”型结构，而湍流干扰下的分离泡结构沿展向的变化则要平缓得多，因此分离区结构上的差异进一步加剧了转捩干扰下激波阵面的弯曲和破碎程度。另外，在干扰区下游，由于转捩干扰下激波阵面更为靠近斜面，此时激波面形状受再附激波的影响较大，与湍流干扰相比，展向形状有更大程度的变形和褶皱。

图6 转捩干扰下压力梯度瞬时等值面Fig.6 Instantaneous contours of magnitude of pressure gradient in transitional case

图7 湍流干扰下压力梯度瞬时等值面Fig.7 Instantaneous contours of magnitude of pressure gradient in turbulent case

3.2 转捩边界层特性

3.2.1 平均速度剖面

图8给出了转捩干扰下不同流向位置处Van-Direst变换后的平均速度沿物面法向的分布。为了便于比较，图中还给出不可缩平板湍流边界层的结果(黑色点划线)。从图中可以看到，在不同流向位置处，速度剖面在黏性底层(y+<5)内均满足线性关系式。然而在对数层内，尽管各流向位置处速度剖面与对数律存在较大的偏差，但是随着压缩拐角内转捩边界层的演化发展，在激波干扰区的下游，如x=40 mm处的测点，此时速度剖面在对数层内已逐步接近对数律，这表明转捩边界层在与激波相互作用后会加快边界层的转捩。

图8 转捩干扰下不同流向位置处平均流向速度沿物面法向的分布Fig.8 Van-Driest transformed mean streamwise velocity at various streamwise locations in transitional case

3.2.2 脉动强度

图9和图10分别给出了无量纲化流向脉动速度均方根Urms和脉动压力均方根prms的分布云图。从图9中可以看到，转捩干扰和湍流干扰下，流向脉动速度均方根的分布规律和峰值大小均差别明显。对于湍流干扰工况，流向脉动速度主要出现在分离泡上方的剪切层内，无量纲峰值大小约为0.2，而对于转捩干扰，分离泡尺度明显小于湍流干扰，此时流向脉动速度则主要出现角部靠近拐点的区域内(-10 mm

图9 流向脉动速度均方根分布云图 Fig.9 Contours of root mean square of streamwise fluctuation velocity

图10 脉动压力均方根分布云图 Fig.10 Contours of root mean square of fluctuation pressure

3.2.3 速度条带结构

高低速条带结构是壁湍流中的重要拟序结构[27]。所谓速度条带结构指的是在壁湍流的近壁区，高速区和低速区呈条带状交替分布。图11给出了转捩干扰下t=2 752时刻距物面y+=5的xOz截面内瞬时流向速度Us分布。图12对应为图11中拐角角部区域的局部放大图。从图中可以清楚看到，在角部干扰区的上游(-40 mm25 mm)，条带结构很快又重新恢复，而且沿展向呈均匀分布(如图12 中右侧箭头所示)。这也表明压缩拐角内的干扰作用将进一步加剧转捩边界层的转捩过程。

3.2.4 湍动能及其输运方程

可压缩湍动能的输运方程为[25,28]

(5)

图11 转捩干扰下距物面y+=5的xOz截面内瞬时流向速度分布Fig.11 Instantaneous streamwise velocity contour in xOz plane at y+=5 in transitional case

图12 转捩干扰下距物面y+=5的xOz截面内瞬时流向速度分布(图11的局部放大图)Fig.12 Zoomed visualization of Fig.11, i.e., instantaneous streamwise velocity contour in xOz plane at y+=5 in transitional case

图13 不同流向位置处时空平均湍动能沿物面法向分布Fig.13 Distributions of mean turbulent kinetic energy at different streamwise locations

图14分别给出了转捩干扰下不同流向位置处湍动能输运方程中各项(TKE budget)沿物面法向的分布情况。图14(a)～图14(f)分别对应x=-35.0，-21.0，-7.0，6.5，20.0，45.0 mm，图14(a)和图14(b)对应拐角内干扰区的上游(定义为1区)，图14(c)和图14(d)对应拐角内干扰区(定义为2区)，图14(e)和图14(f)对应拐角内干扰区的下游(定义为3区)。

湍动能生成项P表征了Reynolds应力通过平均运动的变形率向湍流脉动运动输入的能量。如图14所示，随着转捩边界层往下游的发展，1区内湍动能的生成项较小，峰值区域只出现在距壁面0.1～1.0 mm之间。而2区内的湍动能生成项峰值则急剧升高，峰值大小约为0.006，出现在边界层外缘2 mm附近。在3区内，边界层内靠近物面的位置(<0.1 mm)出现了一个局部峰值，约为0.004。此外，在转捩边界层外的流场中，由于激波的剪切作用，也产生了一个局部的湍动能生成峰值，如图14中(e)和图14(f)中黑色箭头所示。

黏性耗散项ε代表湍流脉动运动引起的黏性耗散。从图14中可以看到，该项以近壁区内为主，不同流向位置处的峰值大小差别较大，其中2区内黏性耗散项峰值达到0.015，约为拐角入口处的10倍。黏性扩散项D的分布规律与耗散项ε的分布类似，该项也以近壁区为主，峰值大小各有不同，其中以2区内的近壁值为最大，约为拐角入口处的6倍。

图14 转捩干扰下不同流向位置处湍动能方程各项的分布Fig.14 Turbulent kinetic energy budget at different streamwise locations in transitional case

脉动密度和脉动速度的四阶相关项Ts表征由脉动运动的相关所产生的湍动能的扩散。压力输运项S是脉动压力和脉动速度相关的梯度项。在1区内，Ts和S的值均非常小。但在2区，由于湍动能的生成急剧增加，而湍能的耗散仍以近壁区耗散为主，此时四阶相关项Ts和压力输运项S发挥了主要的平衡机制，把远离壁面的湍动能输运到近壁区耗散。另外，在3区内，由于激波的压缩作用，压力输运项在转捩边界层外还会出现一个局部的峰值。

从图14中还可以看到，压力膨胀项П在近壁区内其值非常小，仅在激波区域发挥作用，如图14(e) 和图14(f)所示。压力膨胀项表征了可压缩效应所产生的脉动密度对湍动能增长率的影响。这表明，除了主激波区域，其他区域内的本质压缩性效应是非常弱的，研究结果与Li等[25]的湍流干扰下的结论是一致的。

3.3 物面脉动压力分布

此外，从图16中还可以明显看到，对于转捩干扰，脉动压力均方根的分布以单峰值结构为主，峰值位置出现在0 mm

图15 物面瞬时脉动压力分布云图 Fig.15 Contours of instantaneous wall fluctuation pressure

图16 物面脉动压力均方根沿流向的分布 Fig.16 Distributions of root mean square of wall fluctuation pressure

3.4 物面摩阻系数分布

图17 时间平均物面摩阻系数分布云图 Fig.17 Contours of time-averaged wall skin-friction coefficient

图17分别给出了转捩干扰和湍流干扰下的时间平均物面摩阻系数〈Cf〉分布云图。图中In为拐角入口x=-35 mm。在分离区起始点Sp之前，转捩干扰下物面摩阻系数沿展向呈现中间低，两边高的非均匀分布，该分布规律与转捩干扰入口处的拟序涡结构是相对应的(参见文献[17])。在分离区再附点Rp之后，转捩干扰下摩阻系数在展向上呈现”W”型结构，而湍流干扰的结果在展向上呈现”V”型分布结构。

以往压缩拐角的研究表明[30-31]，拐角内的Görtler流向涡结构与摩阻系数展向分布规律密切相关。为了进一步研究转捩干扰下的Görtler流向涡结构对摩阻系数展向分布的影响，图18给出了压缩拐角内x=21.9 mm处垂直于斜面的截面(η，z)内的时间平均流场。图中从上往下依次为截面内的流线分布、流向和法向脉动速度以及物面摩阻的展向分布。左图对应为转捩干扰结果，右图对应为湍流干扰结果，其中截面内流线分布用脉动流向速度进行了渲染。图19还分别给出了拐角斜面上的脉动流向速度等值面，红色对应脉动的正值，蓝色对应脉动的负值。

该截面内流向和法向脉动速度的计算方式如下：首先将截面内笛卡儿坐标系下时间平均后的流向和法向速度转换为截面内坐标系(ξ，η)，其中ξ平行于拐角斜面，而η垂直于拐角斜面。然后将转换后的流向和法向速度减去其对应的展向平均值，即可得到截面内转换后的流向和法向脉动速度，如图18所示。

可以看到，对于转捩干扰，截面内存在着两对Görtler流向涡，而且左侧涡对明显强于右侧涡对。在Görtler涡的上抛和下洗作用下，截面内流向脉动速度沿展向呈现正负交替的分布结构。总体上来看，在转捩干扰下，沿展向出现了3个流向脉动正值区和2个流向脉动负值区(见图19(a))，流向脉动速度的正值与负值区展向位置与摩阻沿展向的W型分布是一一对应的，即脉动正值区域对应于摩阻展向分布的波峰，而脉动负值区域对应于摩阻展向分布的波谷。对于湍流干扰工况，截面内只存在一对Görtler流向涡，同样在Görtler涡的上抛和下洗作用下，流向脉动速度沿展向也呈现正负交替分布结构，但此时沿展向只出现了2个流向脉动正值和1个流向脉动速度负值，而且正负值的展向位置与湍流干扰下摩阻展向分布的V字型结构也是相对应。上述研究结果也进一步表明，转捩干扰下压缩拐角内Görtler流向涡数是造成拐角斜面摩阻系数沿展向W型分布的主要因素。

图20分别给出了转捩干扰和湍流干扰下压缩拐角内时空平均摩阻系数{Cf}沿流向的分布情况。为了比较，图中还给出了Wu和Martin[20]的湍流干扰DNS计算结果。从湍流干扰计算结果与Wu和Martin的结果比较来看，两者的差别主要体现在入口处摩阻值略有不同，这主要是由于入口湍流生成方法的差异造成，但其他区域(如拐角干扰区内以及干扰区下游等)的摩阻分布，两者吻合得较好，这也验证了计算程序的准确性和可靠性。

对于转捩干扰工况，在干扰区下游，摩阻分布沿拐角斜面急剧升高，随后缓慢降低，在10 mm

从物面摩阻系数的定义来看：

(6)

式中:下标w代表壁面处的值，∞代表来流值；μ为黏性系数。在一阶近似的假设下，式(6)还可以进一步写成

(7)

图18 x=21.9 mm处(η,z)截面内的流场结果(左图为转捩干扰，右图为湍流干扰) Fig.18 Flowfields of plane (η，z) at x=21.9 mm (Left and right subfigures denote transitional and turbulent cases, respectively)

对于转捩干扰和湍流干扰工况，来流条件、壁面温度Tw和网格间距Δyw均相同，因此，影响两种工况下摩阻分布差异的主要因素只是物面边界层内近壁区的流向速度大小。

依据FIK理论[32]给出的摩阻与雷诺应力的关系，摩阻还可以进一步分为4部分组成，包括层流、湍流、各向异性和瞬态，其中以湍流部分占主导作用，而该项又主要是通过对雷诺剪切应力〈-ρu″v″〉的加权积分来体现。这是因为雷诺剪切应力是流场中流向动量沿物面法向输运能力的重要表征，其正值对应将高动量的流体输运到物面附近，相反，负值对应将低动量的流体从物面输运到边界层外缘。

图21分别给出了转捩干扰和湍流干扰下拐角内雷诺剪切应力的分布云图。可以清楚看到，转捩干扰和湍流干扰下在激波区，由于激波的强剪切作用，流场中均出现了局部峰值。但在转捩干扰工况下，在拐角斜面10 mm

图19 压缩拐角斜面上流向脉动速度等值面Fig.19 Isosurface of streamwise fluctuation velocity at ramp surface

图20 时空平均摩阻系数沿流向分布Fig.20 Distributions of mean skin-friction coefficient

图21 雷诺剪切应力分布云图Fig.21 Contours of Reynolds shear stress

4 结 论

采用直接数值模拟(DNS)方法对来流马赫数Ma∞=2.9，24°压缩拐角内激波与转捩边界层的相互作用进行了研究。通过改变拐角上游含吹吸扰动带平板的长度，使得拐角内诱导激波与Bypass型转捩边界层的非线性阶段产生干扰。系统地考察了转捩干扰下压缩拐角内分离区形态和激波波系等典型结构的流动特征。通过与相同来流条件下湍流干扰工况的对比，研究了转捩干扰与湍流干扰的流动差异性以及相应物理机理。分析了激波干扰对转捩边界层演化特性的影响作用。探讨了转捩干扰下局部脉动峰值压力和摩阻系数的分布规律及形成机制。研究结果表明：

1) 与湍流干扰相比，由于入口处发卡涡串的展向挤压作用，分离区起始点沿展向呈现“V”型分布。同时，相较于湍流干扰的褶皱状态，分离激波在转捩边界层外缘以破碎状态为主，此时激波脚表现为多层结构。

2) 从平均速度剖面、脉动强度和速度条带结构的演化特性来看，压缩拐角内激波与转捩边界层的相互作用急剧加速了边界层的转捩过程。同时，转捩干扰下湍动能的输运机制与湍流干扰类似，但干扰区内的湍动能强度明显高于湍流干扰工况，而且更为靠近物面近壁区。

3) 转捩干扰下拐角干扰区内物面脉动压力呈现强烈的三维随机特征，相较于湍流干扰，脉动压力峰值出现在干扰区下游，且以单峰结构为主。转捩干扰下物面摩阻沿展向呈”W”型分布，这是由拐角内Görtler流向涡对的分布结构所决定的。此外，转捩干扰下干扰区内的强湍动能和高雷诺剪切应力，使得转捩干扰下的局部峰值摩阻系数要高于湍流干扰。

致 谢

感谢国家超级计算天津中心，中国科学院网络中心超级计算中心以及山西吕梁超算中心提供计算机时。

[1] HOLDEN M S. Reviews of aerothermal problems associated with hypersonic flight: AIAA-1986-0267[R]. Reston: AIAA, 1986.

[2] ANDREOPOULOS J, MUCH K. Some new aspects of the shock wave/boundary layer interaction in compression ramp flows[J]. Journal of Fluid Mechanics, 1987, 180: 405-428.

[3] EDWARDS J R. Numerical simulations of shock/boundary layer interactions using time dependent modeling techniques: A survey of recent results[J]. Progress in Aerospace Sciences, 2008, 44(6): 447-465.

[4] KNIGHT D, LONGO J, DRIKAKIS D, et al. Assessment of CFD capability for prediction of hypersonic shock interactions[J]. Progress in Aerospace Sciences, 2012, 48-49(2): 8-26.

[5] DELERY J M. Shock wave/turbulent boundary layer interaction and its control[J]. Progress in Aerospace Sciences, 1985, 22(4): 209-280.

[6] CLEMENS N T, NARAYANASWAMY V. Low frequency unsteadiness of shock wave turbulent boundary layer interactions[J]. Annual Review of Fluid Mechanics, 2014, 46(1): 469-492.

[7] GAITONDE D V. Progress in shock wave boundary layer interactions[J]. Progress in Aerospace Sciences, 2015, 72: 80-99.

[8] DOLLING D S. Fifty years of shock-wave/boundary-layer interaction research: what next?[J]. AIAA Journal, 2001, 39(8): 1517-1530.

[9] SANDHAM N D, SCHULEIN E. Transitional shock wave/boundary layer interactions in hypersonic flow[J]. Journal of Fluid Mechanics, 2014, 752: 349-382.

[10] ERDEM E, KONTIS K, JOHNSTONE E. Experiments on transitional shock wave boundary layer interactions at Mach 5[J]. Experiments in Fluids, 2013, 54(10): 1598-1620.

[11] SCHRIJER F J. Experiments on hypersonic boundary layer separation and reattachment on a blunted cone flare using quantitative infrared thermograph: AIAA-2003-6967[R]. Reston: AIAA, 2003.

[12] XIE S F, GONG J, JI F. Research of flow field characteristics of hypersonic shock wave/transitional boundary layer interaction: AIAA-2015-3569[R]. Reston: AIAA, 2015.

[13] TERAMOTO S. Large eddy simulation of transitional boundary layer with impinging shock wave[J]. AIAA Journal, 2005, 43(11): 2354-2363.

[14] KRISHNAN L, SANDHAM N D. Shock wave/boundary layer interactions in a model Scramjet intake[J]. AIAA Journal, 2009, 47(7):1680-1691.

[15] TOKURA Y, MAEKAWA H. DNS of a spatially evolving transitional turbulent boundary layer with impinging shock wave: AIAA-2011-0729[R]. Reston: AIAA, 2011.

[16] LUDEKE H, SANDHAM N. Direct numerical simulation of the transition process in a separated supersonic ramp flow: AIAA-2010-4470[R]. Reston: AIAA, 2010.

[17] 童福林， 李新亮， 唐志共, 等． 转捩对压缩拐角激波/边界层干扰分离泡的影响[J]． 航空学报, 2016, 37(10): 2909-2921.

TONG F L, LI X L, TANG Z G, et al. Transition effect on separation bubble in a compression ramp[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(10): 2909-2921 (in Chinese).

[18] RINGUETTE M J, BOOKEY P, WYCKHAM C, et al. Experimental study of a Mach 3 compression ramp interaction atReθ=2 400[J]. AIAA Journal, 2009, 47(2): 373-385.

[19] BOOKEY P, WYCKHAM C. SMITS A J, et al. New experimental data of STBLI at DNS/LES accessible Reynolds numbers: AIAA-2005-0309[R]. Reston: AIAA, 2005.

[20] WU M, MARTIN M P. Direct numerical simulation of supersonic turbulent boundary layer over a compression ramp[J]. AIAA Journal, 2007, 45(4): 879-889.

[21] MARTIN M P, TAYLOR E M, WU M, et al. A bandwidth optimized WENO scheme for the effective direct numerical simulation of compressible turbulence[J]. Journal of Computational Physics, 2006, 220(1): 270-289.

[22] PIROZZOLI S, GRASSO F. GATSKI T B. Direct numerical simulation and analysis of a spatially evolving supersonic turbulent boundary layer atM=2.25[J]. Physics of Fluids, 2004, 16(3): 530-545.

[23] GAO H, FU D X, MA Y W, et al. Direct numerical simulation of supersonic turbulent boundary layer flow[J]. Chinese Physics Letters, 2005, 22(7): 1709-1712.

[24] LI X L, FU D X, MA Y W, et al, Acoustic calculation for supersonic turbulent boundary flow[J]. Chinese Physics Letters, 2009, 26(9): 094701.

[25] LI X L, FU D X, MA Y W, et al. Direct numerical simulation of shock wave turbulent boundary layer interaction in a supersonic compression ramp[J]. Science China: Physics, Mechanics & Astronomy, 2010, 53(9): 1651-1658.

[26] WU P P, MILES R B. Megahertz visualization of compression-corner shock structures[J]. AIAA Journal, 2001, 39(8): 1542-1546.

[27] LEE C B, WU J Z. Transition in wall-bounded flows[J]. Applied Mechanics Reviews, 2008, 61(3): 0802.

[28] PIROZZOLI S, BERNARDINI M. Direct numerical simulation database for Impinging shock wave/turbulent boundary layer interaction[J]. AIAA Journal, 2011, 49(6): 1307-1312.

[29] SKOTE M, HENNINGSON D S. Direct numerical simulation of a separated turbulent boundary layer[J]. Journal of Fluid Mechanics, 2002, 374(5): 379-405.

[30] LOGINOV M S, ADAMS N A, ZHELTOVODOV A A. Large eddy simulation of shock wave/turbulent boundary layer interaction[J]. Journal of Fluid Mechanics, 2006, 565(1): 135-169.

[31] DAWSON D M, LELE S K. Large eddy simulation of a three dimensional compression ramp shock turbulent boundary layer interaction: AIAA-2015-1518[R]. Reston: AIAA, 2015.

[32] FUKAGATA K, IWAMOTO K, KASAGI N. Contribution of Reynolds stress distribution to the skin friction in wall bounded flows[J]. Physics of Fluids, 2002, 14(11): L73-L76.

Numericalstudyofshockwaveandbypasstransitionalboundarylayerinteractioninasupersoniccompressionramp

TONGFulin1，TANGZhigong1,*,LIXinliang2，WUXiaojun1,ZHUXingkun2

1.ComputationalAerodynamicsInstituteofChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China2.StateKeyLaboratoryofHighTemperatureGasDynamics,InstituteofMechanics,ChineseAcademyofSciences,Beijing100190，China

Adirectnumericalsimulation(DNS)ofshockwaveandbypasstransitionalboundarylayerinteractionfora24°compressionrampatMachnumberMa∞=2.9isconducted.Theintricateflowphenomenaintheramp-corner,includingseparationbubblecharacteristicsandshockwavebehavior,havebeenstudiedsystematically.TheDNSresultsoftransitionalinteractionarecomparedwiththecorrespondingturbulentinteractionandthereasonsforthedifferencesareanalyzed.Theevolutionofthetransitionalboundarylayerintherampisresearched.Thefluctuationofwallpressureanddistributionofskinfrictioncoefficientintransitionalinteractionareinvestigatedindetail.Resultsindicatethatthedistributionofcoherentvortexstructuresisnon-uniforminthespanwisedirectionandtheseparationbubbleisreducedtoaV-shapebythemutualinteractionsofthehairpinvorticeschains.Theshockfrontsaredestroyedbadlyandevenbreakdownbytheinteraction.Themultiplelayerofshockfootsisobservedobviously.Theinteractionsrapidlyacceleratetheevolutionoftransitionandgreatlyamplifytheintensityoffluctuations.Thepeakofwallpressurefluctuationsappearswithsingle-peakstructureatthedownstreamofseparationregion.AndtheovershootofskinfrictioninducedbytransitionalinteractionisexplainedbythestrongReynoldsshearstressandhighturbulentkineticenergy.

compressionramp;shockwaveandboundarylayerinteraction;bypasstransition;fluctuationpressure;skinfriction;directnumericalsimulation

2016-01-13；Revised2016-03-14；Accepted2016-03-17；Publishedonline2016-03-301726

URL：www.cnki.net/KCMS/detail/11.1929.V.20160330.1726.006.html

s：NationalNaturalScienceFoundationofChina(91441103,11372330,11472278)

2016-01-13；退修日期2016-03-14；录用日期2016-03-17； < class="emphasis_bold">网络出版时间

时间：2016-03-301726

www.cnki.net/KCMS/detail/11.1929.V.20160330.1726.006.html

国家自然科学基金 (91441103,11372330,11472278)

*

.Tel.：0816-2463133E-mail515363491@qq.com

童福林， 唐志共， 李新亮， 等． 压缩拐角激波与旁路转捩边界层干扰数值研究J． 航空学报,2016,37(12):3588-3604.TONGFL,TANGZG,LIXL,etal.NumericalstudyofshockwaveandbypasstransitionalboundarylayerinteractioninasupersoniccompressionrampJ.ActaAeronauticaetAstronauticaSinica,2016,37(12):3588-3604.

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2016.0096

V211.3； O354.3

A

1000-6893(2016)12-3588-17

童福林男， 博士研究生， 助理研究员。主要研究方向： 可压缩湍流直接数值模拟， 高超声速气动热和热防护。Tel.: 0816-2463133E-mail: wowo2020@sohu.com

唐志共男， 博士， 研究员， 博士生导师。主要研究方向： 高超声速空气动力学。Tel.: 0816-2463133E-mail: 515363491@qq.com

*Correspondingauthor.Tel.：0816-2463133E-mail515363491@qq.com