前后缘同时可控的乘波体气动修型设计与分析

李永洲，李光熙，张堃元，马 元

(1. 西安航天动力研究所，陕西 西安 710100；2. 南京航空航天大学 能源与动力学院，江苏 南京 210016)

0 引言

乘波体是通过流线追踪和吻切轴对称等方法在给定的基准流场中来截取满足设计要求的流面作为压缩面，设计点时其激波完全附在前缘，可以保持较高的升阻比。另外，乘波体具有反设计的特点，方法不但灵活而且适应性好。因此，作为高超声速飞行器的理想构型，乘波体越来越受到关注[1]。

从乘波体的设计流程可以看出，基准流场的设计十分关键。早期的研究主要是一些相对简单、容易求解的流场，比如平面楔形[2]、圆锥[3]和旋成体流场[4]等。随着数值计算方法的发展，一些三维流场也被用于乘波体设计，如楔锥组合体[5]、倾斜圆锥和椭圆锥[6]、任意截面锥形体流场[7]和钝头椭锥体和组合尖锥体等[8]等。另外，反设计的基准流场也被研究，如给定激波形状的曲面外锥[9]和气动参数分布可控的曲面外锥[10]。此外，为了进一步提高乘波体设计的灵活性，早期的锥导和吻切锥[4]方法发展为和吻切内锥和吻切流场[11]等方法，对应的设计参数由前缘或后缘控制发展到可以同时控制出口截面的激波形状[12]。

针对腹部进气布局的高超声速飞行器前体与进气道的一体化设计要求，期望乘波体出口附近的型面变化较缓，流场比较均匀，此时需要同时指定后缘和前缘型线。本文作者先前提出了一种前后缘同时可控的乘波体设计方法[13]。虽然该方法可以同时满足乘波体前缘和后缘要求，但是混合函数使得下表面两侧过渡处产生了诱导激波和较大范围的高压区，出口流场均匀性变差，而且后部的激波变形不再附着在前缘，升阻比下降。此外，乘波体下表面是升/阻力的主要来源，其外形的变化会造成气动性能的显著变化[14]。因此，寻求一种既可以实现型面光滑过渡，又不破坏原有乘波特性且气动性能变化较小的新方法就显得尤为重要。本文在文献[13]基础上，根据超声速气动原理，结合流线追踪和几何重构技术提出了一种前后缘同时可控的乘波体气动修型设计方法，且通过数值仿真手段验证设计方法并分析气动修型的乘波体设计点和非设计点的气动特性。

1 马赫数可控基准流场特征

根据指定的压缩面马赫数分布，采用有旋特征线法反设计出外锥形基准流场[13]，具体设计参数是设计马赫数Ma=6.0、捕获半径Ri=250 mm、前缘半锥角δc=7°及反正切马赫数分布。基准流场的结构如图1所示，其分为两个区域，前缘激波DA的依赖区域DAB和后段等熵压缩区域ACB。AC为出口，AB为结尾特征线，DB段为激波DA所决定的壁面，约为总长度的81.2%，这与内锥形基准流场不同，后者前缘激波所决定的壁面比例约为25%～35%。根据超声速扰动沿特征线传播的特点，只要区域DAB内的流线不做修改，就不会改变乘波体设计状态的三维激波形状，即不会改变其乘波特性。

图1 外锥形基准流场的结构Fig.1 Structure of external conical basic flow field

图2给出了基准流场的马赫数分布，前缘激波不断弯曲，出口马赫数分布较均匀。基准流场总收缩比Rct=3.43，增压比为7.4时总压恢复系数为0.816，长度Lb/Ri=3.70，出口平均马赫数Mae约为4.10。

图2 外锥形基准流场的马赫数分布Fig.2 Mach number distribution of external conical basic flow field

2 乘波体气动修型设计

对于前后缘同时可控的乘波体，控制后缘的主要目的是希望出口附近的型面和流场满足与进气道一体化设计要求，而非针对整个下表面，因此只要对基准乘波体(即指定前缘乘波体)靠近出口附近的型面进行修型使之达到匹配要求即可。由第1节可知，只要区域DAB内的流线不做修改，就不会改变乘波体设计状态的三维激波形状，即可以保持原有的乘波特性，所以在指定后缘型线的条件下，只需对区域ABC内对应的流面进行修型以满足设计要求。这里有两种修型方法，一种是继续采用混合函数，与文献[13]的区别在于此时型面过渡从结尾特征线AB处开始，过渡区域只有原来的四分之一左右，型面变化更加剧烈，产生的诱导激波将更强，不但前缘激波变形更加严重，而且出口流场的均匀性更差；另一种是几何重构方法，在修型区域ABC内增加多条B样条控制曲线，然后通过三维曲面造型软件来生成满足要求的自由曲面。这种方法控制变量多，构造和分析灵活方便，不但可以避免型面剧烈变化产生的诱导激波破坏原有流场，而且比较容易满足总体的结构要求。因此，本文采用第二种方法进行气动修型，下面以超椭圆前缘转圆弧后缘乘波体为例介绍完整的设计流程，具体步骤如下：

步骤1 在外锥形基准流场中，从指定的前缘水平投影出发，向下游追踪流线生成基准乘波体Waverider-F的下表面，上表面取自由流面。该方法与传统方法的不同在于前缘不是指定轴向投影，这样不但可以灵活控制前缘后掠角的变化，而且可以更好控制机体面积和布置飞行器组件。

根据飞行器总体要求，前缘FCT水平投影采用类三角翼形的超椭圆曲线(图3)，长半轴为490 mm，短半轴为200 mm，指数m=1.1，超椭圆中心与锥尖的距离DF=925 mm。Waverider-F下表面(橙色区域)的前缘与指定前缘FCT完全重合，与设计预期相符。

图3 前缘水平投影在基准流场中的相对位置Fig.3 Relative location of horizontal projection of leading edge in the basic flow field

图4给出了Waverider-F的三维视图，前缘对应的轴向投影近似为余弦曲线，虽然容积率较大，但是两侧空间很小且出口附近的型面明显下凹。

图4 指定前缘的基准乘波体三维构型(Waverider-F)Fig.4 Three-dimensional configuration of original waverider with assigned leading edge

步骤2 在基准乘波体Waverider-F的基础上，根据超声速气动原理，结合几何重构技术实现前后缘同时可控乘波体的气动修型设计。

为了与进气道一体化设计，希望Waverider-F出口对称面附近的型面更加平缓，同时为了提高两侧的空间利用率，因此后缘TCT取为圆弧，见图5，半径300 mm，圆心与基准流场轴心的距离DT=85 mm。根据超声速流场的扰动传播特点，以基准流场中的结尾特征线AB(图1)为界，从前缘出发的流线追踪到此处截止，生成气动修型乘波体Waverider-M的前部型面OGEH(图6中的蓝色区域)，该型面决定三维入射激波形状。对应的后部型面EGFH(图6中红色区域)采用几何重构方法进行设计，由于对称性，选取一半即可：先确定两条B样条控制曲线，本文取对称面内的EF和1/4后缘展向长度处的E1F1，然后根据这些控制曲线和轮廓线GE、GF，使用曲面造型软件完成曲面重构和光顺性分析。最终生成的下表面见图6，流线追踪区域OGEH和修型区域EGFH过渡光滑，出口对称面附近型面平缓，两侧型面饱满，空间利用率提高。

图5 后缘轴向投影在基准流场中的位置Fig.5 Relative location of axial projection of trailing edge in the basic flow field

图6 气动修型的乘波体三维构型(Waverider-M)Fig.6 Three-dimensional configuration of aerodynamic modification waverider

步骤3 对乘波体Waverider-F和Waverider-M的无粘构型进行粘性修正，流线追踪区域采用文献[15]的附面层计算方法沿每条流线进行修正，见式(1)，式中a和b为系数，x为轴向坐标。对于后部修型区域，先按照公式(1)计算两条B样条控制曲线的修正量进行修正，然后采用几何重构方法再次生成三维曲面，这样便完成了整个型面的粘性修正。

δ(x)=ax+bxe-x

(1)

3 设计结果的数值分析

3.1 数值计算方法

由于模型和流场的对称性，选取一半模型进行数值计算，采用ICEM软件进行网格划分，壁面附近的网格进行加密，总网格单元约为140万。采用压力远场、对称面和压力出口边界条件，壁面取绝热无滑移和固体边界。Fluent软件求解，无粘时采用二阶迎风格式求解欧拉方程。有粘时，通量差分采用Advection Upstream Splitting Method (AUSM)格式，湍流模型为两方程的Re-Normalization Group (RNG)k-ε模型，近壁处采用非平衡壁面函数。多项式拟合比热，Sutherland公式计算分子粘性系数，各残差至少下降4个数量级并且流量沿程守恒视为收敛。文献[16]校核了计算方法的可信度，表明其能较准确地模拟高超声速复杂流场。下文数值计算的来流条件：设计点马赫数6.0，静压2 549.22 Pa，静温221.55 K；接力点马赫数4.0，静压5 529.3 Pa，静温216.65 K。

3.2 设计点时乘波体的流场与性能

3.2.1 乘波体的无粘流场特征

对Waverider-F和Waverider-M无粘构型进行无粘计算，图7给出了其对称面和出口截面的马赫数分布，Waverider-F的前缘弯曲激波与基准流场的激波形状相同(图2)，出口处的圆弧激波紧贴前缘，出口马赫数分布均匀，约为4.07。Waverider-M虽然采用了几何重构方法，但是对称面内的激波形状和马赫数分布都与Waverider-F相同，出口激波也保持了圆弧形状而且基本贴紧前缘。另外，Waverider-M出口截面没有诱导激波，两侧的马赫数不断减小，而对称面附近的马赫数分布均匀，这有利于与进气道的匹配。

图8给出了Waverider-F和Waverider-M沿程横截面的马赫数分布，Waverider-F具有良好的乘波特性。横截面内的激波均为同心圆弧，体现了外锥形基准流场的波系特征。Waverider-M的波系特点与Waverider-F相同，下表面激波紧贴前体前缘，只是在出口截面两侧有少量偏移，横截面激波也为同心圆弧。上述流场特点表明本文的气动修型方法可行，在不改变乘波体乘波特性的条件下实现了前后缘同时可控。

图7 无粘设计点时乘波体对称面和出口马赫数分布Fig.7 Mach number distribution of waverider symmetry plane and exit at inviscid design point

图8 无粘设计点时乘波体沿程横截面的马赫数分布Fig.8 Mach number distribution of cross sections along the flow direction at inviscid design point

为了进一步分析流场，图9给出了二者下表面的静压分布和流线，Waverider-F的流线向外偏转，静压分布中间高两侧低，这体现了曲面基准流场的压缩特点。Waverider-M的静压分布在流线追踪区域与Waverider-F相同，绝大部分的修型区域静压分布连续，只有两侧由于过渡距离过小造成很小范围的高压区，这与文献[13]中存在大范围高压区不同，因此这里未产生诱导激波破坏原有流场的特征，也说明本文的几何重构方法有效灵活，可以保持原流场的波系特点。

图9 无粘设计点时乘波体下表面的静压分布与流线Fig.9 Static pressure distribution and streamlines of lower surface of waverider at inviscid design point

3.2.2 乘波体的有粘流场特征

对Waverider-F和Waverider-M粘性修正构型进行有粘计算。图10可以看出，二者对称面的激波形状和马赫数分布基本一致，出口的圆弧激波贴着前缘。二者与无粘时的波系和马赫数分布(图7)基本一致，Waverider-M出口对称面附近的马赫数分布仍然较均匀。

有粘时，Waverider-F和Waverider-M沿程横截面的马赫数分布如图11所示，虽然附面层不断发展，但是二者的激波形状与无粘时(图8)相同，横截面内激波呈同心圆分布且整体附着在前缘，Waverider-M出口两侧有少许偏移。

图10 有粘设计点时乘波体对称面和出口的马赫数分布Fig.10 Mach number distribution of symmetry plane and exit at viscous design point

图11 有粘设计点时乘波体沿程横截面的马赫数分布Fig.11 Mach number distribution of cross sections along flow direction at viscous design point

图12为有粘时二者静压与极限流线分布，可以看出，与无粘时类似(图9)，极限流线整体向两侧偏转，静压分布相同，只是此时附面层的存在使得Waverider-M最外侧的流线向外偏转更加严重，这会对出口处激波产生少许影响。

上述有粘流场特征验证了粘性修正方法的有效性，核心区保持了无粘流场特点。

图12 有粘设计点时乘波体下表面的静压与极限流线Fig.12 Static pressure distribution and limit streamlines of lower surface of waverider at viscous design point

3.2.3 乘波体的总体性能

Waverider-F无粘构型的容积率η约为0.23，有粘构型的容积率约为0.22。Waverider-M的容积率更高，二者分别为0.24和0.23。

表1给出了设计点时乘波体的总体性能参数，包括无粘和有粘性能，其中：L为升力；CL为升力系数；D为总阻力；CD为总阻力系数；CDp为压差阻力系数；L/D为升阻比。俯仰力矩系数为CM式(2)的参考点选为前体上表面与对称面交线的终点。下标e表示出口截面，出口截面定义为出口处激波与后缘线围成的区域，如果激波没有紧贴后缘，过前后缘交点做水平线与激波相交。σe为总压恢复系数，pe/p0为增压比，各参数均按照质量加权平均。

CM=M/(0.5·ρ0v02A·Lw)

(2)

式中：M为俯仰力矩(抬头力矩为正)；ρ0和v0分别为来流密度和速度；A为水平投影面积；长度Lw/Ri=1.96。

无粘时，Waverider-F的乘波特性良好，升阻比达到3.46，俯仰力矩系数为0.040，出口增压比和马赫数与基准流场几乎相等。与Waverider-F相比，Waverider-M出口总压恢复系数近似相等，这是因为二者前缘激波形状相同，激波后都是等熵压缩的缘故，但是Waverider-M几何重构区域型面的变化造成压缩量明显增加。Waverider-M的俯仰力矩不变，升力和阻力均增加，而阻力增加的比例更高(46%)，从而升阻比降为2.90。

有粘时，相对无粘参数，Waverider-F和Waverider-M出口增压比变化很小，粘性造成了二者的总压恢复系数分别降低了6.2%和7.5%；在升力下降和总阻力增加情况下，Waverider-F和Waverider-M的升阻比分别下降15.9%和12.8%。相对Waverider-F，Waverider-M的性能变化趋势与无粘时一致，升阻比由2.91降为2.53，出口总压恢复系数降低了1.6%。

上述性能分析进一步验证了气动修型方法的可行性和粘性修正方法的有效性。

表1 设计点时乘波体的总体性能参数

3.3 接力点时乘波体的有粘流场与性能

在接力点对Waverider-F和Waverider-M粘性修正构型进行有粘计算，图13给出了此时对称面和出口的马赫数分布，二者对称面的激波形状和马赫数分布相同，出口都是圆弧激波且与前缘有间距。Waverider-F出口马赫数按照圆弧均匀分布，主流区变化范围为3.01～3.12，Waverider-M对称面附近的马赫数分布与Waverider-F近似相等。

图13 有粘接力时乘波体对称面和出口的马赫数分布Fig.13 Mach number distribution of symmetry plane and exit at viscous relay point

图14给出了有粘接力点时二者沿程横截面的马赫数分布，由图进一步看出，Waverider-F和Waverider-M前部的激波仍然可以完全乘波，这也体现了弯曲激波压缩的优势，后部的激波开始逐渐偏离前缘，造成两侧气流上溢。另外，Waverider-F和Waverider-M横截面的激波形状相同且均为同心圆弧。

图14 有粘接力点时乘波体沿程横截面的马赫数分布Fig.14 Mach number distribution of cross sections along the flow direction at viscous relay point

图15给出有粘接力点时下表面的静压与极限流线分布，由图可以看出，下表面的极限流线和静压分布与设计点时(图12)基本一致，Waverider-M的几何重构区域绝大部分静压分布连续，只在出口附近的两侧存在很小的高压区。

图15 有粘接力点时乘波体下表面的静压与极限流线Fig.15 Static pressure distribution and limit streamlines of lower surface of waverider at viscous relay point

表2给出了有粘接力点时总体性能参数，与设计点(表1)相比，Waverider-F和Waverider-M的升力有所增加，升阻比相对降低7.6%和8.3%。对于出口截面的性能，二者的增压比都下降了约一半，总压恢复系数高达0.92。此外，接力点和设计点时，相对Waverider-F，Waverider-M的变化趋势一致，只是接力点时出口参数变化幅度更小，总压恢复系数仅下降了0.4%。总体而言，二者在接力点也保持了较高的气动性能。

表2 有粘接力点时乘波体的总体性能参数

4 结 论

1)根据气动原理，确定了前缘激波的依赖区，采用流线追踪和几何重构技术实现了前后缘同时可控的乘波体气动修型设计，数值仿真结果验证了设计方法的可行性。

2)气动修型的乘波体型面过渡光滑，只在出口两侧有很小的高压区，可以保持基准乘波体的波系结构和乘波特性，出口对称面附近的流场均匀，接力点时前部也可以乘波，具有较高的升阻比。

3)气动修型的乘波体具有较高的容积率和预压缩效率。相对基准乘波体，其容积率、升力和压缩效率更高，而升阻比降低。无粘构型容积率提高了4.3%，设计点的出口总压恢复系数近似相等，增压比提高了45.8%。

参考文献：

[1] DUNCAN LUNAN M A. Waverider，a revised chronology： AIAA 2015-3529 [R]. Reston： AIAA，2015.

[2] NONWEILER T R F. Aerodynamic problems of manned space vehicles [J]. Journal of the royal aeronautical society，1959，63： 521-530.

[3] JONES J U，MOORE K C，PIKE J，et al. A method for designing lifting configurations for high supersonic speeds using axisymmetric flow field [J]. Archive of applied mechanics，1968，37(1)： 56-72.

[4] SOBIECZKY H，ZORES B，WANG Z. High speed flow design using the theory of osculating cones and axisymmetric flows [J]. Chinese journal of aeronautics，1999，12(1)： 1-8.

[5] TAKASHIMA N. LEWIS MJ. A cone-wedge waverider configuration for engine-airframe integration [J]. Journal of aircraft，1995，32(5)： 1142-1144

[6] RASMUSSEN M L. Waverider configurations derived from inclined circular and elliptic cones [J]. Journal of spacecraft and rockets，1980，17(5)： 537-545.

[7] LOBBIA M A，SUZUKI K. Experimental investigation of a Mach 3.5 waverider designed using computational fluid dynamics [J]. AIAA journal，2015，53(6)： 1590-1601.

[8] LIU C Z，BAI P，CHEN B Y，et al. Rapid design and optimization of waverider from 3D flow： AIAA 2016-3288 [R]. Reston： AIAA，2016.

[9] JONES K D，SOBIECZKY H，SEEBASSA R. et al. Waverider design for generalized shock geometries [J]. Journal of spacecraft rockets，1995，32(6)： 957-963.

[10] 李永洲，张堃元. 基于马赫数分布可控曲面外/内锥形基准流场的前体/进气道一体化设计[J].航空学报，2015，36(1)： 289-301.

[11] RODI P .E. The osculating flowfield method of waverider geometry generation, AIAA 2005-0511 [R]. Reston： AIAA，2005.

[12] LEWIS M J，CHAUFFOUR M L.Shock-based waverider design with pressure gradient corrections and computational simulations [J]. Journal of aircraft，2005，42(5)： 1350-1352.

[13] 李永洲，孙迪，张堃元. 前后缘型线同时可控的乘波体设计[J]. 航空学报，2017，38(1)： 120153.

[14] MAXWELL J R.Hypersonic waverider stream surface actuation for variable design point operation： AIAA 2016-4706 [R]. Reston： AIAA，2016.

[15] DRAYNA T W，NOMPELIS I，CANDLER G V. Hypersonic inward turning inlets： design and optimization： AIAA 2006-297 [R]. Reston： AIAA，2006.

[16] ZHANG K Y. Research progress of hypersonic inlet inverse design based on curved shock compression system:AIAA 2015-3647 [R]. Reston： AIAA，2015.