0 引 言
1904年德国流体力学大师Prandtl提出了边界层的概念，为人们认识粘性流体及解决工程中的流体力学问题打开了一扇大门。边界层是存在于粘性流体和物面之间的一薄层流体，这一薄层流体对飞行器气动特性的影响一直是空气动力学研究的经典问题，特别是在高超声速流动中，由飞行器前缘弓形激波产生的熵层会凸显边界层的粘性效应。
高超声速二维钝楔是二元高超声速进气道的抽象模型，由于在进气道的流动中具有强激波、粘性效应显著及激波/边界层干扰的流动特征，边界层在高超声速进气道的研究设计中扮演着重要角色，对进气道性能(压缩性能、总压恢复性能、起动性能和耐反压性能等)有重要影响。因此，了解二维钝楔边界层参数的分布规律对进气道的型面设计及性能预测有重要意义。早在1966年，Benson 等[1] 就采用数值方法研究了一种马赫数Ma=8轴对称进气道的前缘钝度及壁面温度对边界层厚度的影响，得到了同实验一致的结果。John和Leonard[2] 针对四种高超声速轴对称等熵进气道开展了风洞试验，通过测量轴对称进气道的边界层位移厚度分布，详细研究了壁面温度、马赫数及雷诺数对边界层位移厚度的影响。Varner 和 [CM(22*2]Adams[3] 针对双钝锥外形开展了边界层参数的风洞 试验，详细研究了熵吞和低雷诺数效应对边界层参数的影响，并同理论结果进行了对比。Drayna等[4]为设计高超声速内转式进气道，针对轴对称Busemann流动采用数值方法研究了不同头部半径时湍流边界层位移厚度分布情况，为进气道边界层修正提供了准确数据。
针对高超声速方转圆进气道传统设计的缺陷，Gollan and Ferlemann[5]开发了一套高超声速方转圆进气道的设计分析工具，该工具可以快速计算进气道复杂三维粘性流场的边界层位移厚度，允许设计者实时考察进气道设计特征并进行型面修正，以提高进气道性能。同样，为设计高超声速方转椭圆进气道，Smart[6] 针对基于无粘理论设计的进气道开展粘性流场的计算，研究了压力、边界层位移厚度、摩擦力系数沿流向的变化规律，最终对理论型面进行边界层修正得到了进气道物理型面。Jenn[7] 针对一种轴对称尖锥外形计算了马赫数Ma=4.63，攻角α=12°时的边界层参数，包括边界层边缘的压力、密度、速度以及摩擦阻力系数、位移厚度，通过与试验数据对比验证了数值算法的精度。Martellucci[8] 针对半锥角为7.25°的尖锥在Ma=8时的湍流边界层参数进行了试验测量，并同理论数据进行了对比，在规律上获得了较好的一致性。国内专门针对边界层厚度的系统研究并不多，大多集中在对边界层转捩过程的研究[9,10,11]，程兴华等[12] 通过采用薄激波层理论计算高超声速二维钝楔激波形状，分析了边界层外缘熵的分布规律，研究了熵层对壁面主要气动参数的影响。
对边界层的研究主要有试验测量和数值计算两种方法。采用试验的方法，测量中会存在很多干扰因素，而且成本较高。当前CFD技术发展已经十分成熟，可以很好地解决复杂流动问题。因此，用数值方法研究高超声速二维钝楔边界层参数不失为一种高效、经济的方法。通过数值计算可以全面、准确地模拟边界层的流动结构和特征，进而可以对边界层参数进行精确预测，了解边界层参数的分布规律。
本文采用数值方法对四个二维钝楔的边界层参数进行了系统研究，得出了边界层位移厚度、边界层总压恢复、压力梯度沿流向的变化规律，对指导二元高超声速进气道型面设计具有一定参考价值。 1 控制方程及数值方法 1.1 控制方程
使用自主研发的非结构流场解算器——CACFD软件[13]——开展边界层参数计算，其控制方程是可压缩流粘性气体动力学方程组，积分型式如下：
其中：
1.2 湍流模型
在CACFD软件中，软件作者对几种湍流模型曾进行了比较研究，认为对于超声速流动，在计算边界层和压缩拐角流动方面，Menter k-ω SST湍流模型具有较好的模拟能力，故本文计算中选择Menter k-ω SST湍流模型开展研究，Menter k-ω SST湍流模型的守恒形式是：
P为湍动能生成项，湍流粘性定义：
N-S方程数值求解中，时间推进采用LU-SGS方法，在空间方向采用有限体积法离散，空间通量计算采用二阶精度的Roe-FDS格式。流动方程、湍流模型方程的离散均选择二阶精度。
计算域离散采用结构化网格，总数为8.5万，在需要捕捉激波和粘性的位置进行了加密。近壁面网格y+=0.1~0.4，边界层内布置有30个网格点。网格分布及计算域如图 1。
 |
| 图 1 计算域及网格Fig. 1 Computational domain and mesh |
计算中采用的边界条件有：
(1) 远场：给定来流的马赫数、 静温和单位雷诺数；
(2) 超声速压力出口：外插；
(3) 壁面边界：采用无滑移、绝热壁和等温壁条件；
(4) 对称面边界。
计算的收敛准则为：连续方程、动量方程、能量方程以及k-ω方程的残差下降5个数量级，且气动力保持恒定。 2 计算模型
研究的前体型面包括二维多级楔面及二维曲线型面，前缘倒圆半径有R=2mm及R=5mm两种(倒圆圆弧和Y轴及型面同时相切)，总计四种外形。型面外形及分析的边界层站位位置(以总长进行无量纲化)如图 2~图 5。计算状态包括不同的马赫数Ma、高度H、攻角α，如表 1。
 |
| 图 2 楔面外形(R=2mm)Fig. 2 Wedge shape(R=2mm) |
 |
| 图 3 楔面外形(R=5mm)Fig. 3 Wedge shape(R=5mm) |
 |
| 图 4 曲线外形(R=2mm)Fig. 4 Curved shape(R=2mm) |
 |
| 图 5 曲线外形(R=5mm)Fig. 5 Curved shape(R=5mm) |
| No. | Ma | H/km | α/(°) | Temperature of Wall |
| 1 | 4 | 17 | 0 | Adiabatic/Isothermal |
| 2 | 4.5 | 18.5 | 2 | Adiabatic/Isothermal |
| 3 | 5.5 | 23 | 3 | Adiabatic/Isothermal |
| 4 | 6.6 | 27 | 4 | Adiabatic/Isothermal |
位移厚度(以总长进行无量纲化):
总压恢复系数:
阻力系数:
沿程压力梯度:
3.2 流场云图
图 6至图 9给出了部分流场的马赫数云图和静压云图，从流场云图可见，前体第一道激波为脱体的弓形激波，楔面的第二道激波为平面激波，曲线型面前面为连续的弱压缩激波区，后面为膨胀区。
 |
| 图 6 压力云图(R=2mm,Ma=4,绝热壁)Fig. 6 Pressure contour(R=2mm,Ma=4,Adiabatic wall) |
 |
| 图 7 压力云图(R=2mm,Ma=6.5,绝热壁)Fig. 7 Pressure contour(R=2mm,Ma=6.5,Adiabatic wall) |
 |
| 图 8 压力云图(R=2mm, Ma=4,等温壁)Fig. 8 Pressure contour(R=2mm,Ma=4,Isothermal wall) |
 |
| 图 9 压力云图(R=2mm,Ma=6.5,等温壁)Fig. 9 Pressure contour(R=2mm,Ma=6.5,Isothermal wall) |
图 10给出了不同型面、不同钝化半径及壁温条件下的阻力系数对比曲线。本文计算得到的阻力系 数为水平线与钝头相交以上部分的阻力除以来流动压及各自参考长度。由图可见，整体情况是阻力系数 随马赫数增大而增大，在相同马赫数下，随钝化半径 增大而增大。在相同半径时，等温壁计算的阻力系数 与绝热壁基本相当，略大一些，这是因为等温壁壁温小于绝热壁，粘性系数小，边界层薄，但壁面法向速度梯度大，最终是等温壁的摩阻大一些，而二者压阻相同。
 |
| 图 10 阻力系数曲线Fig. 10 Curves of drag coefficient |
与常规的平板边界层不同，本文研究的两个型面流场相对复杂：前缘钝度引起熵增，压缩拐角及等熵压缩区存在激波，曲线型面后段流场中还有膨胀波，这些复杂流动使得速度及密度沿壁面法向没有明显的界限，这就给边界层的外边界的确定带来很大困难。通过调研，了解到国外的做法[14,15,16]是使用总温来确定边界层的外边界，因为总温只有在粘性层内变化较大。结合实际情况，对于等温壁(图 11(a))，横坐标r为垂直方向的距离，总温曲线没有明显峰值时，定义总温T0=0.995T0∞处为边界层外边界。对于绝热壁(图 11(b))，总温曲线在边界层内有明显峰值，定义峰值右侧的总温T0=0.999T0max处为边界层外边界。这样就解决了如何确定复杂流场边界层外边界的问题，为进一步的边界层参数计算提供了边界层外边界的理论判据。
 |
| 图 11 总温法向分布曲线Fig. 11 Vertical distributing curves of total temperature |
图 12~图 15给出了两级楔面前体外形在不同马赫数下边界层位移厚度、边界层总压恢复系数、壁 面压力梯度的沿程分布情况。从曲线来看，边界层位移厚度对前缘钝化半径及壁温较为敏感。R=5mm 的边界层位移厚度在第一楔面小于R=2mm的边界层位移厚度，而在第二楔则快速增长而超过R=2mm 的边界层位移厚度，这样的规律与文献[4]研究 结果的规律一致。压力梯度的影响是：R=2mm的边界 层位移厚度在逆压梯度作用下增速减小，而R=5mm 的边界层位移厚度对这样的逆压梯度不是很敏感。壁温的影响是：在相同钝化半径下，绝热壁的边界层位移厚度明显大于等温壁的边界层位移厚度。
 |
| 图 12 Ma=4时楔面边界层参数对比Fig. 12 Comparison of boundary layer parameters (Wedge, Ma=4) |
 |
| 图 13 Ma=4.5时楔面边界层参数对比Fig. 13 Comparison of boundary layer parameters (Wedge, Ma=4.5) |
 |
| 图 14 Ma=5.5时楔面边界层参数对比Fig. 14 Comparison of boundary layer parameters (Wedge, Ma=5.5) |
 |
| 图 15 Ma=6.5时楔面边界层参数对比Fig. 15 Comparison of boundary layer parameters (Wedge, Ma=6.5) |
图 16~图 19给出了曲线型面前体外形在不同 马赫数下边界层位移厚度、边界层总压恢复系数、壁 面压力梯度的沿程分布情况。与两级楔面前体类似，边界层位移厚度对前缘钝化半径及壁温较为敏感。R=5mm的边界层位移厚度在前面的顺压梯度和等 压力梯度时小于R=2mm的边界层厚度，而在后面 的膨胀区则快速增长而超过R=2mm的边界层位移 厚度。压力梯度的影响是：R=2mm的边界层位移 厚度在等压力梯度作用下增速不断放缓，在后面的膨 胀区则快速增长。而R=5mm的边界层位移厚度对 压力梯度不是很敏感。壁温的影响是：在相同钝化半径下，绝热壁的边界层位移厚度随着站位靠后明显大于等温壁的边界层位移厚度。边界层总压恢复对前缘钝化半径较为敏感，在相同站位，随前缘半径增加，边界层位移厚度增加，但总压恢复下降。
 |
| 图 16 Ma=4时曲面边界层参数对比Fig. 16 Comparison of boundary layer parameters (Curved,Ma=4) |
 |
| 图 17 Ma=4.5时曲面边界层参数对比Fig. 17 Comparison of boundary layer parameters (Curved,Ma=4.5) |
 |
| 图 18 Ma=5.5时曲面边界层参数对比Fig. 18 Comparison of boundary layer parameters (Curved,Ma=5.5) |
 |
| 图 19 Ma=6.5时曲面边界层参数对比Fig. 19 Comparison of boundary layer parameters (Curved,Ma=6.5) |
本文针对两种二维钝楔型面开展了边界层参数的数值研究，得到如下结论：
(1) 本文计算结果规律性较好，说明所采用的数值计算方法可以反映实际高空状态二维前体湍流边界层参数的物理规律；
(2) 针对边界层参数的三种影响因素：前缘钝化半径、壁面温度、壁面压力梯度，本文进行了较为详细的研究。研究发现，前缘半径及壁面温度对边界层参数影响较大，在两者一定的情况下，逆压梯度使边界层位移厚度增速放缓，顺压梯度会使边界层位移厚度快速增长；
(3) 边界层会对超燃冲压发动机进气道性能产生重要影响，本文研究结果可为进气道设计提供指导性意见，具体的量值预测精度还有赖于试验(特别是飞行试验)的进一步验证。
| [1] | Benson J L, Maslowe S A. Bluntness and boundary-layer displacement effects on hypersonic inlet flowfields[J]. Journal of Spacecraft and Rockets, 1966, 3(9):1394-1401. |
| [2] | Stroud J F, Mille L D. Hypersonic inlet boundary-Layer research[R]. AIAA-65-605, 1965. |
| [3] | Varner M O, Adams J C. Variable edge entropy and low Reynolds Number effects on hypersonic turbulent boundary layers[R]. AIAA-80-0131, 1980. |
| [4] | Drayna T W, Nompelis Ioannis, Candler G V. Hypersonic inward turning inlets:design and optimization[R]. AIAA 2006-297, 2006 |
| [5] | Gollan R J, Ferlemann P G. Investigation of REST-class hypersonic inlet designs[R]. AIAA 2011-2254, 2011. |
| [6] | Smart M K. Design of three-dimensional hypersonic inlets with rectangular to elliptical shape transition[R]. AIAA-98-0960, 1998. |
| [7] | Jenn A A. Preliminary design estimates of high-speed streams and boundary layer properties[R]. AIAA-92-0640, 1992. |
| [8] | Martellucci Anthony. Effects of mass transfer on hypersonic turbulent boundary layer properties[R]. AIAA-1972-184, 1972. |
| [9] | Cao Wei. A study of the transition prediction of hypersonic boundary layer on plane and wedge flow[J]. Acta Aerodynamic Sinica, 2009, 27(5):516-523. (in Chinese)曹伟. 高超声速边界层的转捩问题[J]. 空气动力学学报, 2009, 27(5):516-523. |
| [10] | Xu Ding, Ma Huiyang. Engineering transition models for hypersonic boundary layer[J]. Journal of the Graduate School of the Chinese Academy of Sciences, 2009, 26(1):43-48. (in Chinese)许丁, 马晖扬. 高超声速边界层工程转捩模式研究[J]. 中国科学院研究生院学报, 2009, 26(1):43-48. |
| [11] | Pan Honglu, Ma Handong, Wang Qiang. Large eddy simulation of transition in a hypersonic blunt wedge boundary layer[J]. Acta Aeronautica et Astronautica Sinica, 2007, 28(2):269-274.(in Chinese)潘宏禄, 马汉东, 王强. 高超声速钝楔边界层转捩大涡模拟[J]. 航空学报, 2007, 28(2):269-274. |
| [12] | Cheng Xinghua, Yang Tao, Chang Zhongdong. Influence of entropy layer on aerodynamics parameters for hypersonic 2-D blunt wedge[J]. Journal of Aerospace Power, 2012, 27(6):1362-1367.(in Chinese)程兴华, 杨涛, 常中东. 熵层对高超声速二维钝楔气动参数的影响[J]. 航空动力学报, 2012, 27(6):1362-1367. |
| [13] | 刘周, 白鹏, 周伟江. CACFD大型气动数值模拟软件平台[C].中国力学学会学术大会, 郑州:2009. |
| [14] | McInville, Moore R M, Frank G. Incorporation of boundary layer heating predictive methodology for NSWCDD application code[R]. AIAA-94-2001, 1994. |
| [15] | Liechty Derek S, Berry, et al. Comparison of methods for determining boundary layer edge conditions for transition correlations[R]. AIAA 2003-3590, 2003. |
| [16] | Albertson, Cindy W, Bauserman. Total temperature probes for high-temperaturehypersonic boundary-layer measurements[R]. NASA TM-4407, 1993. |