基于局部偏转吻切方法的多级压缩乘波体设计1)

郑晓刚 朱呈祥 尤延铖

(厦门大学航空航天学院,福建厦门 361005)

引言

高超声速飞行器具备全球实时侦察、快速部署和远程精确打击能力,是临近空间力量部署的主要载体,对国家战略安全产生重要影响,已然成为各国竞相争夺的空间技术焦点之一[1-5].常规气动布局的高超声速飞行器在飞行过程中通常会产生脱体激波,整体构型升阻比较低.随着飞行马赫数的增加,这一现象愈发明显,飞行器的综合性能急剧下降,出现“升阻比屏障”[6].乘波体构型则利用前缘附体激波近乎完美地将高压气流限制在飞行器下表面,获得比常规外形高得多的升阻比,被认为是突破常规“升阻比屏障”的有效途径,几乎已经成为高超声速飞行器设计的首选方案.

自1959 年Nonweiler[7]提出“Λ”型乘波体概念以来,有关乘波体设计方法的研究已经得到了长足的发展[8].为解决“Λ”型楔导乘波体下反角过大、容积率不足等问题,Jones 等[9]和Rasmussen 等[10]将基准流场由二维平面流动推广至轴对称锥形流动,提出了锥导乘波理论.而为了进一步提高对激波形状的控制,Sobieczky 等[11-12]和Rodi[13]在锥导乘波理论基础上,考虑激波展向曲率的可控设计,进一步发展了吻切锥、吻切轴对称以及吻切流场三类吻切乘波理论.随后,国内外一些学者以吻切理论为基础,获得了构型更为复杂、具有特定目的的一系列新型吻切类乘波构型[14],如变马赫数乘波体[15-16]、双后掠乘波体[17-18]以及双乘波构型[19-20]等.这些设计方法仅能实现二维、轴对称或准三维激波的逆向乘波设计,在附体激波的选择上仍存在一定局限.为此,Zheng 等[21-22]在吻切锥理论的基础上进一步发展了一种局部偏转吻切(local-turning osculating cones,LTOCs)方法,实现了全三维非轴对称流动的快速逆向求解,极大拓宽了乘波体的设计空间.

可以看到,上述方法多是基于一道激波设计的单级压缩乘波体,压缩效率不够.对于吸气式飞行器而言,需要进气道进行进一步的压缩,这无疑增加了整体的设计难度,也使得进气道内流动变得更为复杂.针对这一问题,国内外学者也已经进行了一定的探索.例如,美国X-51A 吸气式高超声速验证机便采用了两级平面压缩的乘波前体方案[23].刘嘉和王发民[24]提出了利用相交楔/锥基准流场构建两级压缩乘波前体的方法.吕侦军等[25-27]基于吻切锥乘波理论和Taylor-Maccoll 流动方程,首次提出了倾转圆锥绕流法的多级压缩乘波体设计方法,并通过试验和数值模拟的手段验证了该方法的可靠性.在此基础之上,同一课题组的Wang 等[28]进一步提出了多级锥乘波前体与Busemann 进气道相结合的一体化设计方法.该方法在提升一体化构型升阻比的同时,改善了进气道入口处的均匀性.贺旭照和倪鸿礼[29]和卫锋等[30]发展了一种吻切曲面锥乘波体设计方法.该方法采用具有直线激波和等熵压缩波系的曲面锥作为基本流场,利用第一道激波和马赫波对来流进行压缩,一定程度上克服了传统乘波体压缩量不足和容积率偏小的问题.

然而,现阶段绝大多数多级压缩乘波体的设计方法都是基于二维、轴对称或者吻切理论开展的,其激波形状较为单一,难以同日益复杂的进气道唇口相匹配,无法实现一体化构型的综合性能最优.针对这一问题,本文将多级压缩乘波体的概念进一步应用于LTOCs 方法,实现了多道复杂三维激波的逆向乘波设计.在接下去的内容中,本文以两级椭圆锥压缩乘波体为例,详细介绍了基于LTOCs 方法的多级压缩乘波体设计方法.随后,基于此法设计了3 种不同长短轴比的两级椭圆锥压缩乘波体.通过数值仿真验证了该方法的准确性,并对比分析了不同长短轴比对此类乘波体的性能影响.

1 多级压缩乘波体设计方法

1.1 局部偏转吻切方法

针对任意的三维激波曲面,一些学者尝试使用三维特征线法进行逆向求解[31-33].但是,由于特征线法需要不断迭代以及三维求解过程中特征线容易相交等问题,直接逆向求解是较为复杂且耗时的.一种合理的替代方式是同吻切理论那般对整个求解过程进行降维处理,将三维的求解过程简化为一系列二维过程.

基于这一思路,LTOCs 方法采用了同吻切锥理论一样的降维方式,将三维流动沿横向离散成若干切片,仅需求解各切面内的二维流场便可组合求得整个三维流动.其不同之处在于,吻切锥理论的切片是单一的平面,忽略了切面间的横向流动,难以实现三维复杂激波的精确求解.而LTOCs 方法中的切片则完整包含当地流线,这张切面是一种曲面,其本质上是一张流面(见图1).因此,LTOCs 方法中各切面内的求解方式将明显有别于常规吻切锥理论.在LTOCs 方法中,各流面内流动的求解过程大致分为3 步:首先,运用微元吻切面的概念离散给定的三维激波曲面,求得各流面内的激波曲线;随后,将各离散激波线绕X轴旋转至同一虚拟子午面内,利用文献[21]提出的非共轴二维特征线线法(non-coaxial method of characteristics,NMOC)进行求解;最后,根据“cross-stream marching (CSM)”方法将该虚拟子午面内求得的二维流场参数转换至三维笛卡尔直角坐标系内.激波后完整的三维流动则由这些流面内的流动组合而成.

图1 局部偏转吻切方法示意图[22]Fig.1 Schemaitc of LTOCs method[22]

1.2 多级压缩乘波体设计步骤

多级压缩乘波体的设计预期是该乘波体在设计条件下产生的多道激波相交于乘波体的底部截面内(即设计截面).设计激波的数目主要对乘波体的压缩性能有所影响,而设计的基本原理是一致的.本文在接下去的内容中以两级压缩乘波体为例介绍基于LTOCs 方法的多级压缩乘波体设计方法的概念.利用这一方法,只需激波参数合理,原则上可以设计出任意压缩级数的乘波体.基于LTOCs 方法的两级压缩乘波体设计方法大致可归纳为4 步.

(1)构造两道相交的非轴对称激波.椭圆锥激波是一种最为常见且具有标准数学表达式的非轴对称激波.在后续工作中,文章主要以椭圆锥激波作为非轴对称激波进行两级压缩乘波体设计.如图2 所示,两椭圆锥激波面相交于设计平面内的同一个长短轴比为a/b的椭圆上.在图示以O1点为坐标原点的坐标系内,第一道椭圆锥激波的表达式为

图2 两级椭圆锥压缩乘波体设计示意图Fig.2 Schemaitc for the design of dual elliptic cone compression waveriders

而第二道椭圆锥激波的表达式为

上述方程中,L1和L2分别表示两个椭圆锥的长度,即L1为线段O1O的长度而L2为线段O2O的长度.L1和L2必须满足任意激波点处的激波角度要大于当地马赫角,以确保所构造激波在物理上是可行的.

(2)指定设计平面内FCT (flow capture tube)曲线,并将其水平向前投影至第一道激波面上确定两级压缩乘波体的前缘型线.由于FCT 曲线对乘波体性能的影响并不是本文的研究重点,指定两级压缩乘波体的前缘线方程为二次曲线

其中A和R为两个无量纲常数.沿流向将FCT 曲线水平向前投影即可得到两级压缩乘波体的上表面.

(3)从前缘型线出发,流线追踪获得第一级压缩型面.将前缘型线离散成N点,利用LTOCs 方法追踪从这些离散点出发的流线,直至第二级压缩激波面.求出各流线与第二道激波面的交点,便可确定第一级压缩型面与第二道激波的交线,如图2 所示.交线上各离散点处的流动参数可在第一级压缩激波后的三维流场中插值得到.需要强调的是,第一级压缩型面是在第一道激波后流场内进行流线追踪得到的,其来流条件是均匀来流.

(4)从交线出发,在二级压缩激波后流线追踪生成第二级压缩型面.从第3 步得到的各离散交点出发,在第二道激波后同样利用LTOCs 方法追踪流线,直至设计平面为止,便可获得第二级压缩型面.不同于第一级压缩型面,在流线追踪生成第二级压缩面时,O2O范围内各离散激波点对应的来流条件是非均匀的,并且应当通过在第一道波后流场中三维插值得到.在生成多级压缩乘波体时,之后压缩面的处理方式与二级完全相同.

至此,两级压缩乘波体的设计已基本完成.根据这一方法,合理构建多级压缩激波系,理论上可设计出任意压缩级数的乘波体.

1.3 流线与第二道激波面相交判据

在两级压缩乘波体的设计过程中,需要计算各流面内从前缘型线出发的流线与第二道激波的交点位置,并利用空间插值计算求得对应点处的流动参数作为下一级型面计算时的来流参数.LTOCs 方法在计算各流面内的流动参数时是通过将各点旋转至同一虚拟子午面内进行计算的.因此,借助虚拟子午面可较为清晰的解释求交点过程.图3 给出了任意虚拟子午面内流线同两道压缩激波面间的相对位置关系.图中,A点为前缘离散点,坐标轴x为两道椭圆锥激波的轴线,即三维笛卡尔坐标系中的X.由于图示激波由流面内激波线旋转而成,故而其呈现为曲线状.

图3 虚拟子午面内流线与第二道激波相交示意图Fig.3 Schemaitc for the intersection of the streamline with the second shock in the fictional meridian plane

可以看到,判断某一微元流线段是否与第二道激波相交,可通过对比这一微元段流线上两个端点与第二道激波面间的相对位置关系来实现.若流线的前端点Pn-1位于第二道激波面外,而后端点Pn位于第二道激波上或激波内则表示此段流线与第二道激波面相交.由于Pn-1点在前一微元段内已经经过判断,因此,只需对当前点Pn的位置进行判断激波.对于任意一点,其空间坐标(Xn,Yn,Zn)必须满足如下判据,方能位于第二道椭圆锥激波上或激波内

此外,Pn的坐标还应满足如下位置关系

由式(4)和式(5)组合的判据便可判断第一级压缩型线是否与第二道激波面相交.利用点Pn处的微元吻切方程以及第二道激波面方程,便可计算得到此微元吻切面内的激波曲线,进而求得流线与第二道激波面的交点Pt.

交点Pt处的气动参数利用反距离权重(inverse distance weighted,IDW)插值法在第一道激波后流场内插值得到.IDW 插值法假定周边离散点对插值点都有着局部的影响,而且这种影响随着距离的增加而减弱.该方法的核心在于周边离散点对应权重函数的求解,各离散点的权重可表示为

式中,下标0 代表插值点,下标i代表周边离散点,而n则表示距离插值点最近的离散点数目.根据各离散点的权重即可求得插值点处的任意参数

在本文中,为同时兼顾计算精度与效率,取n为8.

随后,从交点Pt开始,利用LTOCs 方法对第二道椭圆锥激波面进行离散求解.离散时,激波点处的微元吻切面应当由该点合速度方向及曲率方向共同确定,而非水平来流方向与曲率方向所确定.各离散激波点的参数利用IDW 插值法求得,用作第二道激波求解的来流条件.

1.4 局部网格自适应加密方法

两级压缩乘波体设计方法的本质是在各流面内给定两道激波曲线条件下反设计流线.即在任意虚拟子午面内,给定激波曲线AB和PtB以及来流条件,反求这两道激波对应的物面形状与依赖域(见图4).激波AB上的参数由来流条件和Rankine-Hugoniot (R-H)方程求得,随后利用非共轴二维特征线法获得相应型面AD以及左行特征线边界BD.由激波角与马赫角之间的关系可知,左行特征线边界BD必定位于第二道激波PtB的上游.仅靠激波AB是无法确定物面DPt与区域BDPt的.本文参照文献[34]的做法,将激波AB延长用以确定上述区域的参数.

图4 第二级压缩面特性线网格跨度过大示意图Fig.4 Schemaitc for the large characteristic grids of the second compression surface

在第二道激波后同样存在这样的问题,激波PtB仅能确定物面PtP1以及区域BPtP1.P1点与设计平面间的物面仍需通过延长第二道激波求得.不同于第一道激波,延长后第二道激波的来流边界条件被分为两个部分.交点B之前的来流参数是经过第一道激波压缩后的流动参数,而交点B之后的来流参数则是自由来流参数.因此,在B点后的第一个微元吻切面内,激波点B与B1的参数存在着较大的梯度.这两点发出的左行特征线的间距逐渐增加,并最终导致物面上相邻P1和P2两点的步长跨度过大,壁面过渡不够光顺,如图4 所示.针对这一问题,本文提出了一种基于流向压力梯度的特征线网格自适应加密技术,其基本原理如图5 所示.这一自适应加密技术的核心思想是:在利用非共轴特征线法求解虚拟子午面内参数时,在每一个小三角形区域内对流线上两点的压力梯度进行判断.若两者间的压力梯度大于给定值,则在流线两点间线性插值引入新的离散点作为下一个单元过程的输入,以此达到加密流线网格点数的目的.以图示区域为例,流线点D1原本是在三角形C1DD1内,利用非共轴特征线理论求得,在流线DD1上的步长较大.而引入自适应加密技术后,对流线点B1与C1间的压力进行判断,并在两者间线性插值得到点.在三角形中求得点参数,接着在三角形中求得D1点.借此,流线DD1上的网格点数得到了增加.通过不断重复上述过程,最终对物面P1P2上的特征线网格点数进行加密处理,使得两点间的型面得以光顺过渡.

图5 基于压力梯度的特征线网格自适应加密技术Fig.5 Adaptive refinement of characteristic grids based on the streamwise pressure gradient

2 多级压缩乘波体设计方案

根据上述多级压缩乘波体设计方法,本文在来流马赫6、飞行高度27 km 的条件下,设计了3 个不同的两级椭圆锥压缩乘波体,分别命名为Case A,Case B和Case C.3 个乘波体对应的两道椭圆锥激波的顶点位置完全相同.如图6 所示,3 个乘波体在设计平面内的激波曲线(shock wave profile curve,SWPC) 为不同长短轴比(a/b) 的椭圆.Z半轴上,3 个椭圆激波的的半轴长度b是相同的.这也意味着:3 个乘波体在对称平面内的两道设计激波是完全一致的.

图6 各乘波体设计平面激波形状对比Fig.6 Comparison of SWPCs in the design plane

可以看到,Case A 对应椭圆激波的长短轴比为8/7,椭圆焦点位于Y轴之上.Case C 对应椭圆激波的长短轴比为7/8,椭圆焦点位于Z之上.Case B 对应激波的长短轴比为1,是标准圆形.因此,Case B 实际上就是标准的两级圆锥压缩乘波体.3 个乘波体对应的FCT 曲线方程相同,均为

由上述条件可知,3 个乘波体具有相同的长度,但宽度并不相同.这是由三者在设计平面内的激波形状差异造成的.此外,从图6 中可以看到,虽然3 个乘波体在对称面内的激波形状完全相同,但三者的高度存在一定差异.这一结果表明椭圆锥激波后,横向流动对乘波体型面有一定影响.

图7 给出了3 个乘波体下表面形状对比结果,图中浅绿色部分为第一级压缩面,浅蓝色部分为第二级压缩压缩面.按照设计预期,俯视图中下表面的前缘处将产生第一道激波,一级压缩面与二级压缩面相交处将会产生第二道激波,且这两道激波相交于乘波体底部截面内.从图中可以看到,3 个乘波体一级、二级压缩面的交线形状各不相同.其中,Case A 的交线在对称面附近呈上凸,Case B 的交线在对称面附近基本水平,而Case C 的交线在对称面附近则呈下凹状.经统计,Case A 一级、二级压缩面的面积S1和S2分别为0.455 m2和0.472 m2.Case B 两个压缩面的面积分别为0.401 m2和0.439 m2.而Case C 两个压缩面的面积分别为0.351 m2和0.415 m2.3 个乘波体的压缩面积差别较大,因此三者的性能参数差别也较为明显.3 个乘波体压缩面积的大小关系为:Case A＞Case B＞Case C.可以预见,这3 个乘波体中Case A 的压缩性能最好,而Case C 的压缩性能最差.

图7 各乘波体设计下表面形状对比Fig.7 Comparison of lower surfaces

表1 进一步给出了3 个乘波体的几何参数对比.可以看到,3 个乘波体的长Lw,宽W和高H与前文一致,符合设计需求.在体积Vol,水平投影面积Sp和垂直投影面积Sb方面,3 个乘波体间的关系与压缩面积关系相同,为Case A＞Case B＞Case C.而在容积率方面,三者的关系恰好相反,为Case C＞Case B＞Case A.相较于两级圆锥压缩乘波体Case B 而言,Case A 的容积率下降约3.29%,Case C 的容积率增加约2.35%.

表1 各乘波体几何参数对比Table 1 Comparison of geometric parameters

3 多级压缩乘波体设计方法验证

3.1 数值方法与算例验证

本文使用商业软件Fluent 对前文3 种乘波体模型进行数值模拟.采用有限体积法对控制方程进行空间离散,求解器为“density-based”,差分格式选用二阶AUSM (advection upstream splitting method)的迎风格式.流体假设为定比热比γ等于1.4 的理想气体.有黏计算时,湍流模型选用k-ωSST 模型,分子黏性采用Sutherland 公式计算,以保证计算精度和强间断计算的稳定性.计算过程中,监测残差变化、计算域进出口质量流量差以及乘波体表面的平均压力值.当残差曲线至少下降3 个数量级、进出口流量差值与进口流量之比小于0.01 且乘波体表面的平均压力保持不变时即判定计算收敛.

为了验证上述数值方法的可靠性,本文选用文献[35]中的双椭球模型进行三维高超声速流动的数值模拟,进而与文献中的实验数据进行对比.双椭球模型的具体尺寸如图8 所示.本次算例验证对比了不同攻角条件下模型对称面上下壁面的沿程压力分布,计算条件为来流马赫数4.94,单位长度雷诺数5.26 × 107m-1,计算攻角则为-5°和0°.

图8 双椭球几何模型示意图Fig.8 Schematic of double ellipsoid model

考虑到模型的对称性,为了节约计算资源和计算时间,本文仅对半模进行结构性网格划分,如图9所示.网格划分时在壁面附近以及其他预计流动参数梯度较大的区域,如两椭球交界处进行网格加密处理.第一层网格间距为1 μm,以保证y+小于2.5,壁面法向网格长度增长率为1.1,网格总数在180 万左右.

图9 双椭球模型结构网格示意图Fig.9 Schematic of structured-grids for the double ellipsoid model

图10 给出了不同攻角条件下,双椭球模型对称面上下壁面的无量纲压力对比结果.其中,横坐标X代表以双椭球模型前缘点为起始点的流向坐标,P为当地静压,而Pt2则表示双椭球模型的前缘驻点压力.从图中可以看到,CFD 计算得到的对称面内上下壁面的压力分布与实验结果基本一致.除第二个椭球体前缘形成的激波附近外,上下壁面的最大压力误差均控制在10%以内.可见,本文所采用的数值计算方法的误差较小,可作为后续工作的研究手段.

图10 双椭球模型对称面压力分布对比Fig.10 Comparison of pressure distribution on the symmetric plane

3.2 设计方法验证

为了验证本文基于LTOCs 方法发展的多级压缩乘波体设计方法的可靠性,本节以Case C 为例详细对比了该模型在无黏、设计条件下的激波形状与壁面压力分布.在构建结构网格时,乘波体表面布置C 型网格,并在乘波体前缘及其他流动梯度较大的区域进行加密处理,总网格数量为500 万左右.

图11 给出了Case C 在对称面内的压力云图,以展示无黏CFD 计算得到的激波形状同设计激波之间的差异.图中,红色虚线代表对称面内第一道设计激波,而白色虚线则表示第二道设计激波.由前文设计参数可知,对称面内两道设计激波均为直线,且激波角度分别为14.04°和21.80°.由图可知,无黏CFD 计算求得的两道激波在乘波体底部截面相交,且与预设计激波基本重合,验证了本文提出的多级压缩乘波体设计方法的正确性.

图11 无黏条件下Case C 对称面激波形状对比Fig.11 Comparison of shock waves on the symmetric plane for Case C under inviscid conditions

基于LTOCs 方法和无黏CFD 计算的Case C 表面压力分布情况如图12 所示.图中,上半部分为基于LTOCs 方法得到的壁面压力分布,标记为“LTOCs”;下半部分为基于无黏CFD 计算得到的壁面压力分布,标记为“CFD”.可以发现,两种方法的结果无论是在数值还是分布规律上均表现出良好的一致性.两种结果的压力等值线基本重合,误差主要出现在一、二级压缩型面的交线附近.这主要是因为CFD 计算中,激波总是存在一定的厚度.

图12 Case C 无黏壁面压力分布对比Fig.12 Comparison of pressure distribution on the lower surface for Case C under inviscid conditions

表2 进一步给出了设计条件下,Case C 无黏气动参数的对比结果,表中的相对误差均是相对于CFD 结果的.可以看到,两种方法求得的无黏气动力基本一致,各参数的误差均要小于0.3%.由此可见,本文基于LTOCs 方法提出的多级压缩乘波体设计方法能够以较高的精度准确复现预设计多道激波及其波后流场.

表2 Case C 无黏气动力参数对比Table 2 Comparison of inviscid aerodynamic parameters for Case C

4 多级压缩乘波体性能分析

至此,基于LTOCs 方法的多级压缩乘波体设计概念及其精度均已得到验证.在接下去的内容中,将重点对比3 个两级椭圆锥压缩乘波体在无黏和有黏条件下的气动性能和压缩性能,探讨分析不同长短轴比对此类两级椭圆锥压缩乘波体性能的影响.

4.1 无黏条件

图13 给出了无黏、设计条件下,3 个乘波体在不同流向截面内的壁面压力分布.不难发现,此类乘波体在第二级压缩面上的壁面压力可分为两个部分,前半部分的压缩效应明显大于后半部分.这主要是因为前者对应的激波位于第一道激波内,来流先经过了第一道激波的压缩.而后者对应的激波位于第一道激波外,其来流为水平均匀气流.由于3 个乘波体在对称面内的激波形状完全相同,三者在此处的壁面压力分布差别不大.综合来看,Case A 的平均增压比最大,Case B 次之而Case C 的平均增压比最小.此外,从图中数据可以判断,在乘波体底部截面内Case B 的流场畸变度最小,Case A 次之,而Case C 的流场畸变最大.

图13 无黏条件下两级椭圆锥压缩乘波体不同流向截面内壁面压力分布对比Fig.13 Comparison of pressure distribution in different streamwise planes under inviscid conditions

总压恢复系数 σ 为气流压缩后总压与来流总压的比值,反映的是气流经过多道激波压缩后的能量损失.流场均匀度 ε 则是描述乘波体底部截面内流场均匀程度的参数

3 个两级椭圆锥压缩乘波体在无黏、设计条件下的压缩与气动性能参数如表3 所示.相比于标准两级圆锥压缩乘波体Case B,长短轴比a/b大于1 的两级压缩乘波体Case A 的流量系数增加约3.02%,增压比增加约5.48%,可见Case A 的压缩性能要优于Case B.压缩性能的增加相应也导致Case A 的压缩后平均马赫数下降约1.27%,总压恢复系数减少约1.74%.Case A 高压缩性能也意味着其壁面压力的增加,相应的升阻力系数也会增长,Case A 的升阻比较Case B 增大约2.39%.长短轴比a/b小于1 的两级压缩乘波体Case C 则恰巧相反.其流量系数比Case B 减小约5.22%,增压比下降约8.12%,压缩后平均马赫数增加约1.86%,总压恢复系数增加约1.52%.在升阻系数方面,Case C 的升阻力系数均要小于Case B,而Case C 的升阻比下降约2.58%.此外,在流场均匀度方面,Case A 与CaseC 的流场均匀程度均要差于Case B.这主要是因为Case B 对应的两道激波为标准圆锥,是等强度激波,而Case A和Case B 的激波为椭圆锥激波,是变强度激波.Case A 压缩后的流场均匀程度要强于Case C,流场均匀程度提升约25%.

表3 乘波体无黏性能参数对比Table 3 Comparison of invisicd performance parameters for three waveriders

由此可见,设计截面椭圆的长短轴比a/b对此类两级压缩乘波体的压缩和气动性能有着显著的影响,可通过调整这一参数设计出满足不同性能需求的多级压缩乘波构型.

4.2 黏性条件

以上关于3 个两级椭圆锥压缩乘波体的设计与分析都是在无黏条件下的,并未考虑黏性带来的影响.然而,在实际应用中,黏性问题是不可避免的.因此,有必要研究上述无黏构型在黏性条件下的综合性能,以此探讨本文基于LTOCs 方法提出的多级压缩乘波体设计方法的工程实用性.黏性计算同样采用结构网格.以Case C 为例,其表面布置有C 型网格,在壁面、前缘及其他流动梯度较大的区域进行加密处理.第一层网格间距为10 μm,以保证y+小于3,壁面法向网格长度增长率为1.1,网格总数在650 万左右.

图14 给出了Case C 在有黏、设计条件下波系结构.如图所示,在考虑空气黏性后,乘波体的前缘型线附近出现了一定程度的溢流,乘波体下表面的高压气流上洗,导致前缘型线附近存在脱体激波.黏性的影响在一定程度上相当于增加了乘波体下表面的厚度,一级和二级压缩面的厚度同时增加,导致两道压缩激波的角度均有所增大.与无黏条件下的激波系相比,有黏条件下整体的激波结构差别不大,两道椭圆锥激波在乘波体底部截面内大致相交,尤其是在对称面附近.整体而言,Case C 在黏性下仍具有较佳的乘波特性.其余两个乘波体在黏性条件下的波系结构与之类似,便不再赘述.

图14 黏性条件下Case C 激波系结构Fig.14 Configurations of shock waves for Case C under viscous conditions

表4 给出了3 个乘波体在黏性条件下的气动性能参数.从中不难看出,相比于无黏条件,三者的升力系数有所下降,阻力系数有所增加,致使三者的升阻比下降明显.以Case C 为例,其黏性条件下的升力系数下降约1.79%,阻力系数增加约10.53%,而升阻比下降约12.79%.可见,升阻比大幅下降的主要原因在于黏性的存在极大地增加飞行器的阻力.

表4 乘波体有黏气动性能参数对比Table 4 Comparison of viscous aerodynamic parameters for three waveriders

另外,在有黏情况下,3 个乘波体升阻比之间的关系依然与无黏时相同,按升阻比由大到小依次为:Case A＞Case B＞Case C.

5 结论

针对常规单级压缩乘波前体压缩量不足的问题,本文基于LTOCs 方法提出了多级压缩乘波体设计方法,可实现针对多道变强度激波的快速逆向乘波设计.文章以两级椭圆锥压缩乘波体为例详细介绍了这一方法,并以此为基础设计分析了3 种不同长短轴比的两级压缩乘波体,得到结论如下.

(1)无黏、设计条件下,两级椭圆锥压缩乘波体能够准确复现两道预设计椭圆锥激波形状.基于LTOCs 方法得到的乘波体壁面压力分布同无黏CFD 计算的结果基本一致,且两种方法求得无黏气动力参数也基本相同,误差仅为0.3%左右.这些结果均证明了基于LTOCs 方法的多级压缩乘波体设计方法的准确性.

(2)长短轴比对此类乘波体的性能影响较大.无黏设计条件下,相比于两级圆锥压缩乘波体Case B,长短轴比a/b大于1 的Case A 升阻比增加约2.39%,增压比增加约5.48%,流量系数增加约3.02%,总压恢复系数下降约1.74%,容积率下降约3.29%;而长短轴比a/b小于1 的Case C 升阻比下降约2.58%,增压比下降约8.12%,流量系数下降约5.22%,总压恢复系数增加了1.52%左右,容积率上升约2.35%.

(3)考虑黏性后,此类乘波体的激波系结构变化不大,两道椭圆锥激波在底部截面基本相交,且与设计激波也基本重合,误差主要出现在乘波体的两侧.黏性对乘波体阻力的影响远大于对升力的影响,从而导致了升阻比相比无黏状态有所下降.综合来看,此类乘波体在黏性条件下仍具备较好的乘波特性.