韩汉桥,张陈安,王发民
(中国科学院力学研究所 高温气体动力学国家重点实验室,北京 100190)
钱学森于1948年在美国火箭学会年会上提出了著名的钱学森弹道,它是一种飞行器再入后进行高超声速滑翔的弹道[1]。自钱学森弹道提出以来,高超声速滑翔的概念受到世界上几个航空航天大国的重视。近年来,随着临近空间高超声速滑翔飞行器(如CAV、HTV-1、HTV-2等)概念的提出,这类采用高超声速滑翔来实现远程机动、快速响应及到达、有效突防等目标的飞行器受到了更加广泛的关注。
高超声速滑翔飞行器在高空(H=45km~75km)以高马赫数(M∞=15~25)滑翔时,这一飞行区间的来流密度较低,具有较高的马赫数和总焓,飞行器所面临的高空多物理效应(主要是粘性干扰效应和真实气体效应)对其气动性能会产生较大的影响。国内外学者针对平板等简单外形和轨道飞行器的多物理效应问题进行了大量的研究。由于实验条件的限制,目前风洞实验还无法进行大规模的高空高马赫数下粘性干扰效应和真实气体效应方面的研究,对其研究主要采用 CFD方法。Maus[2-3]等人在1983年采用CFD方法研究了粘性干扰效应、真实气体效应等对航天飞机气动性能的影响,并结合地面风洞实验的结果外推出飞行条件的结果。Anderson[4]等人1991年研究了平板的非平衡化学反应效应对粘性干扰的影响,指出有限速率化学反应能够减小高超声速粘性干扰效应。庄逢甘等[5-6]上世纪80年代后期对美国航天飞机面临的空气动力学问题进行了一系列初步研究,包括粘性干扰效应和真实气体效应等。近年来程晓丽[7]、龚安龙[8]、叶友达[9-10]等针对高超 声速轨道飞行器也进行了真实气体效应和粘性干扰效应对飞行器气动特性影响的研究。其结果表明,对于航天飞机这种轨道飞行器,粘性干扰效应使轴向力系数增加,升阻比减小,并导致俯仰力矩系数减小(低头力矩);真实气体效应也会使轴向力系数增加,升阻比减小,但使俯仰力矩系数增大(抬头力矩),压心位置前移。以上这些研究的结论是针对简单模型或美国航天飞机得到的,而对于其它类型的高超声速飞行器,由于构型和飞行条件的差别,这些物理效应产生的影响不完全相同,前面的研究结果只能提供部分参考。
首先采用文献[11]中的方法生成一个锥导乘波体,然后结合热防护、有效体积和飞行器的气动操纵舵面这些因素进行工程化设计,最终设计出来的高超声速滑翔乘波飞行器如图1所示,飞行器的总长为5m,以该长度作为参考长度。
图1 高超声速滑翔乘波飞行器外形Fig.1 Model view of hypersonic gliding waverider
针对该飞行器首先选取高度H=40km,马赫数M∞=20的状态进行数值模拟,分别计算了迎角0°、2°、3°、4°、6°等工况。图2给出了计算所得的升阻比随迎角变化的曲线,可以看出在迎角2°~3°之间升阻比达到最大,并且计算所得的升阻比均在3.5以上,可以满足高超声速远程滑翔的需要。在后面对多物理效应的研究中,均取迎角3°作为典型工况进行数值计算。
图2 不同迎角时的升阻比曲线Fig.2 Lift-to-drag ratio at different angles of attack
在笛卡尔坐标下三维无源非定常N-S方程组的守恒积分形式为:
式中,Q=(ρ,ρu,ρv,ρw,ρe)T,F=FC-FV,FC表示对流矢通量,FV为粘性矢通量。对于完全气体,需补充的状态方程为:
对于平衡气体,状态方程以及输运特性的计算公式不同于完全气体的状态方程,其中热力学属性由Tannehill[12]提供的平衡气体热力学属性的拟合曲线求得。计算过程中空间离散采用 Harten-Yee二阶TVD格式[13],时间离散采用隐式 LU-SGS方法[14]。
在飞行高度较高(如H=70km)时,飞行器头部和前缘附近局部存在一定的稀薄气体效应,但其对整体气动力影响较小,故依旧采用N-S方程计算,壁面边界条件仍采用无滑移条件做近似处理。
高超声速飞行器在高空飞行时,将面临高马赫数、低雷诺数的飞行环境,此时激波非常贴体,而边界层比较厚,相应的位移厚度也较大,会使边界层外的无粘流动发生很大改变,这一改变反过来又会影响边界层的增长,这种边界层与外部无粘流动之间的相互作用就称之为粘性干扰[15]。
为了更直观地表示粘性干扰效应对高超声速气动力的影响,可以选择一个粘性干扰参数作为粘性干扰效应强弱的度量。本文采用了在航天飞机上应用比较成功的第三粘性干扰参数¯V′∞。
式中T′为参考温度,μ′为参考温度T′时的空气粘性系数。T′/T∞为边界层内的参考温度与来流温度之比:
步态选择方面,选取稳定性高的三角步态作为四足机器人在斜面上的运动步态,并采用LF-RH-RF-LH的腿部运动顺序以获得最优的稳定裕度[6]。本文的运动步态时序图如图3所示。其中,黑色矩形表示腿部处于支撑相,白色矩形表示腿部处于摆动相。
根据Sutherland粘性系数公式:
针对所设计的乘波布局滑翔飞行器,选取五个状态点进行计算,飞行高度分别为:H=32km、40km、50km、60km、70km,其中马赫数均为20,迎角均为3°粘性干扰参数随着高度增加而增大,粘性干扰效应随之增强。图3给出了所设计的高超声速滑翔飞行器的气动参数随粘性干扰参数的变化,图3(a)和(b)分别是轴向力系数CA和俯仰力矩系数Cm随粘性干扰参数的变化,图3(c)和(d)分别是轴向力系数的粘性干扰增量ΔCA和俯仰力矩系数的粘性干扰增量ΔCm随粘性干扰参数的变化,计算时把粘性效应分成了诱导压力和剪切效应两项,粘性干扰增量由¯V′∞时的总气动系数减去¯V′∞=0时的无粘值得到。从图中可以看出:粘性干扰效应的增强使该飞行器的轴向力系数增加,其中诱导压力产生的轴向力受粘性干扰效应的影响很小,轴向力的增加主要是因为剪切效应增强,摩阻增加;俯仰力矩系数也随着粘性干扰效应的增强而增大,俯仰力矩系数的粘性干扰增量为正值,即产生抬头力矩,使静稳定性减弱。其中对俯仰力矩系数影响较大的是诱导压力,剪切效应对俯仰力矩系数的影响尽管也随着粘性干扰效应的增强而增强,但相对诱导压力的影响来说比较小。而对于航天飞机这种轨道飞行器,俯仰力矩系数的粘性干扰增量为负值[2],这是由于对不同类型的构型,粘性干扰效应作用的区域和强度不同。
图3 气动参数随粘性干扰参数的变化Fig.3 Aerodynamic parameters at different viscous interaction parameters
图4 不同高度时的马赫数云图Fig.4 Mach number distributions at different attitudes
图4给出了高度32km、50km和70km(单位长度雷诺数分别为:5.53×10m 3.98×10m 3.42×104m-1)时计算得到的马赫数云图。在马赫数不变时,高度越高,雷诺数越低,边界层的厚度越大。从图中可以看出,随着高度的增加,边界层的厚度快速增加,与边界层外的无粘流场形成强烈的相互干扰,激波形状发生了改变,变得越来越脱体,粘性干扰效应的作用增强。图5为不同高度时迎风面和背风面中轴线上的无量纲压力分布。随着高度的增加,物面的无量纲压力升高,由于迎风面距头部较近的区域(0<x/l<0.3)处于较强的粘性干扰区,压力增加较多,尽管背风面也处于较强的粘性干扰区,但在距头部较近的区域粘性干扰引起的压力增量比迎风面小,这是造成俯仰力矩系数粘性干扰增量为正值,产生抬头力矩的原因。
图5 不同高度时迎风面和背风面中轴线上的无量纲压力分布Fig.5 Windward and leeward centerline pressure distributions at different attitudes
真实气体效应指的是高超声速气流经过激波压缩或者粘性阻滞而减速时,部分有向运动的动能转化为分子随机运动的动能,引起气体温度升高,随着飞行速度的增加,这种温升使得气体逐渐偏离完全气体模式,相继出现分子的振动能激发、离解、原子的电离及电子激发和光辐射等一系列复杂物理化学现象,以及由此产生的对流场结构和飞行器性能的影响[16]。
高超声速乘波布局滑翔飞行器,相比航天飞机这种轨道飞行器来说属于细长体,并且飞行时迎角较小,空气因激波强烈压缩而引起较强化学反应的区域仅仅集中在头部附近和前缘的局部。边界层内因粘性引起温度上升而产生的化学反应也会改变飞行器表面的压力分布,进而影响飞行器的气动特性。
高空高超声速飞行器周围流场中往往既有粘性干扰效应,也有真实气体效应,这种耦合作用给飞行器气动性能带来的影响值得研究。
本文针对所设计的高超声速滑翔飞行器分别采用平衡气体模型和完全气体模型进行了计算分析。图6给出了高度H=60km,迎角3°时不同马赫数下采用完全气体模型和平衡气体模型计算所得的x=l/2截面处的压力云图比较,通过压力分布可以看出真实气体效应使得激波后的压力降低,激波位置更靠近物面,并且随着马赫数的增加,真实气体效应增强,这种差别更加明显。
图6 高度60km时不同马赫数下完全气体和平衡气体的压力云图对比Fig.6 Real-gas effects on pressure distribution at different Mach numbers
真实气体效应引起的流场变化必然会对飞行器气动性能产生影响。图7给出了采用完全气体模型和平衡气体模型计算所得的气动参数比较。图7(a)是轴向力系数随马赫数的变化,采用平衡气体模型计算所得的轴向力较大,主要是因为边界层中化学反应吸收了大量的热量,使边界层内温度降低,减小了边界层的厚度,从而使摩阻增加。
图7(b)是法向力系数随马赫数的变化,可以看出随着马赫数的增大,采用平衡气体模型计算得到的法向力系数相比完全气体得到的值越来越小。为了进一步了解真实气体效应对飞行器气动参数的影响规律,图8给出了物面中心线上的无量纲压力分布,真实气体效应使迎风面和背风面的整体压力都有所降低,并且随着马赫数增大,真实气体效应增强,这种压力减小越明显,这与粘性干扰效应的作用相反,印证了真实气体效应会减弱粘性干扰效应的结论。对于背风面,几乎整个区域都处在较强的粘性干扰区内,采用平衡气体模型计算得到的背风面压力在整个区域内都低于完全气体的结果;对于迎风面,随着马赫数的增大,一方面边界层增厚,强粘性干扰区域范围增大,另一方面激波压缩增强,且边界层中的温度上升,化学反应增强,使下表面压力受化学反应影响而降低的区域有所增大。在迎风面和背风面真实气体效应及粘性干扰效应的共同作用下,随着马赫数的增大,迎风面压力的降低比背风面压力的降低要多,从而采用平衡气体模型得到的法向力系数相比完全气体得到的值越来越小。
图7 真实气体效应对气动参数的影响Fig.7 Real-gas effects on aerodynamic parameters at different Mach numbers
图8 物面中心线上的无量纲压力分布Fig.8 Pressure distributions on surface centerline
图7(c)是俯仰力矩系数随马赫数的变化,计算结果表明真实气体效应降低了俯仰力矩,即产生低头力矩,这与航天飞机的结论[2]相反。结合图8分析可知:边界层中因粘性耗散而产生的真实气体效应使其厚度减小,边界层诱导压力减小,并且迎风面靠近头部的强粘性干扰区的压力比后部弱粘性干扰区的压力降低得多,这种压力分布的改变正是产生低头力矩的原因。而航天飞机的构型与本文设计的高超声速乘波布局滑翔飞行器差别较大,由于航天飞机头部钝度较大,且大迎角飞行,迎风面因强激波压缩空气而产生的真实气体效应占很大一部分,真实气体效应对压力分布的影响规律将有所不同,这就造成俯仰力矩系数变化规律的不同。
本文针对所设计的高超声速乘波布局滑翔飞行器,在高度H=60km,来流马赫数M∞=15、20、25的情况下,分别采用平衡气体模型N-S方程、完全气体欧拉方程、完全气体N-S方程进行了数值模拟,通过三种方程计算得到的气动参数的比较分析来揭示高超声速多物理效应对乘波布局高超声速滑翔飞行器气动性能的影响。
图9给出了数值计算得到的多物理效应对轴向力系数、升阻比、俯仰力矩系数以及纵向压心位置的影响,以M∞=15的完全气体模型计算结果作为比较的参考基准值,其中纵向压心位置为纵向压心与质心之间的距离,压心在质心之后为正。这里多物理效应主要指马赫数效应、粘性效应和真实气体效应,其中粘性效应包含了粘性干扰效应。从图9(a)中可以发现影响轴向力系数的三种物理效应中,马赫数效应使飞行器轴向力系数减小,粘性效应使得轴向力系数增大,真实气体效应也使得轴向力系数有所增大,且粘性效应较其它效应的影响要大。从图9(b)中可以看出,马赫数增加使升阻比略有降低,真实气体效应也使升阻比有所降低,粘性效应使升阻比下降较大。图9(c)给出了俯仰力矩系数的变化,完全气体模型的计算结果中俯仰力矩系数随着马赫数的增加而增加,另外两种气体模型计算得到的俯仰力矩系数受马赫数的影响较小。如前面所分析的,真实气体效应使得俯仰力矩系数减小,而粘性效应使得俯仰力矩系数增加较多。图9(d)给出了纵向压心位置的变化规律,可以发现完全气体模型的计算结果中纵向压心随着马赫数的增加而前移,而另外两种气体模型计算得到的纵向压心位置受马赫数的影响较小。真实气体效应使纵向压心后移,粘性效应使纵向压心后移较多。
由此可见,不同的气动力系数对不同物理效应的敏感程度是不一样的,在进行飞行器设计时,可以根据物理效应对设计的目标气动力系数的影响程度,重点关注某一个或某两个物理效应的影响来设计飞行器。如在生成乘波飞行器时,升阻比是一个很重要的指标,从图9中可以看出对升阻比影响最大的是粘性效应,而马赫数效应和真实气体效应对升阻比的影响不大,因此可以重点关注粘性效应对生成乘波飞行器的影响。
图9 多物理效应对气动参数的影响Fig.9 Multi-physical effects on aerodynamic parameters
本文设计了一种高超声速远程滑翔乘波飞行器,利用数值方法研究了高空多物理效应对其气动性能的影响。结果表明:
(1)该飞行器的轴向力和俯仰力矩系数都随着粘性干扰效应的增强而增大,俯仰力矩系数随粘性干扰效应的变化规律与航天飞机这类轨道飞行器不同,主要是由于对不同类型的构型,粘性干扰效应作用的区域和强度不同。
(2)本文设计的这类飞行器的真实气体效应主要集中在边界层内,它使边界层变薄,物面压力减小,与粘性干扰效应的影响相反,这说明真实气体效应在一定程度上减弱了粘性干扰效应。对轴向力的影响与粘性干扰效应相同,都使其增大;对俯仰力矩系数的影响是使其减小,产生低头力矩,这与航天飞机的结果相反,主要是构型不同和产生真实气体效应的原因不同所致。
(3)综合比较了马赫数效应、粘性效应和真实气体效应对气动参数的影响,分析了不同的气动力系数对不同物理效应的敏感程度,可以为高超声速滑翔飞行器的设计及气动评估提供参考。
[1]钱学森.星际航行概论[M].北京:科学出版社,1963.
[2]MAUS J R,GRIFFITH B J,TOLBERT D G,et al.Understanding space shuttle flight data by use of wind tunnel and CFD results[R].AIAA 83-2745,1983.
[3]MAUS J R,GRIFFITH B J,Szema K Y,et al.Hypersonic Mach number and real gas effects on space shuttle orbiter aerodynamics[J].Journal of Spacecraft,1984,21(2):136-141.
[4]HALLGREN W F,ANDERSON J D.The effects of nonequilibrium chemistry on hypersonic viscous interaction[R].AIAA 91-3323,1991.
[5]庄逢甘,赵梦熊.航天飞机的空气动力学问题[J].流体力学实验与测量,1987141-8.
[6]庄逢甘,赵梦熊.航天飞机的粘性干扰效应--航天飞机的空气动力学问题之二[J].流体力学实验与测量,1988,2(1):1-11.
[7]程晓丽,苗文博,周伟江.真实气体效应对高超声速轨道器气动特性的影响[J].宇航学报,2007,28(2):259-264.
[8]龚安龙,周伟江,纪楚群等.高超声速粘性干扰效应相关性研究[J].宇航学报,2008,29(6):1706-1710.
[9]叶友达.高空高速飞行器气动特性研究[J].力学进展,2009,39(4):387-397.
[10]叶友达.空间高速飞行器气动特性研究与布局设计优化[J].力学进展,2009,39(6):683-694.
[11]耿永兵.高超声速乘波飞行器优化设计[D].[博士学位论文].北京:中国科学院研究生院,2006.
[12]TANNEHILL J C,MUGGE P H.Improved curve fits for the thermodynamic properties of equilibrium air suitable for numerical computation using time-dependent or shock-capturing methods[R].NASA CR-2470,1974.
[13]YEE H C,HARTEN A.Implicit TVD schemes for hyperbolic conservation laws in curvilinear coordinates[J].AIAA Journal,1987,25(2):266-274.
[14]YOON S,JAMESON A.Lower-upper symmetric-Gauss-Sediel method for the Euler and Navier-Stokes equations[J].AIAA Journal,1988,26(9):1025-1026.
[15]卞荫贵,徐立功.气动热力学[M].合肥:中国科学技术大学出版社,1997.
[16]HOLDEN M S,WADHAMS T P,CANDLER G V.Experiment studies in the LENS shock tunnel and expansion tunnel to examine real gas effects in hypervelocity flows[R].AIAA 2004-916,2004.