2. 四川大学制造科学与工程学院, 成都 610065
2. College of Manufacturing Science and Engineering, Sichuan University, Chengdu 610065, China
电弧加热器湍流平板试验技术是针对高超声速飞行器热防护材料考核的一种地面试验技术，主要用于表面为平板或以平板为基板的防热材料试验[1]。由电弧加热器产生的高温气流，经转接过渡段时形成层流边界层并在喷管喉道附近转捩为湍流边界层，在超声速矩形型面喷管出口处，与气流有一定迎角地放置平板试件，两者密接齐平无缝隙。试件上的边界层是喷管型面壁上边界层的自然延续，从而在试件上得到充分发展的湍流边界层流动，而试件迎角造成的前缘斜激波，可以提高试件上模拟的参数[2]。这种试验方法能够模拟试件飞行条件下的表面压力、热流和剪切力等气动参数，可以承担大面积防热材料烧蚀试验研究、天线窗口周围局部烧蚀试验研究、控制翼局部烧蚀试验研究以及槽孔和缝隙局部烧蚀试验研究等，是电弧加热设备上最常用的试验方法之一。
在正式试验之前，需要对电弧加热器和试件当地流场进行调试和校测，以确认流场特性和当地参数分布，这也是确保试验质量的关键环节。这个过程通常是由外形与试件相同的传感器校测模型来完成，需要对各项控制参数进行反复调试和匹配。近年来随着CFD技术的发展和其明显的优势，国外相关研究机构已越来越多地将CFD计算应用到电弧加热器流场预测中，它能够获得校测试验难以测量的局部区域，以及整个流场其他参数的分布，可以指导试验调试的实施，减少流场校测环节的探索过程[3]。NASA的Ames研究中心基于N-S方程开发了CFD程序DPLR，应用于大量的电弧加热设备试验中[3-11]，并将其作为地面试验的固有流程和标准配置[4]。应用分析显示，即便DPLR对试件表面热流的预测不确定度较大，仍不失为一种有效的辅助分析工具。文献[5]对改建的平板试验设备TPTF进行不同迎角的流场模拟，与校测试验数据误差小于15%，获得了多种工况下平板试验的流场特性。文献[6]利用ADSI程序对已有数据库进行插值，以预测数据库以外的流场状态，分析结论指出基于CFD数据库插值的结果与实际测量偏差小于10%，且优于以试验数据库进行插值的结果。国内CFD技术发展较快，但在电弧加热试验中并未大量广泛应用，仅限于单一状态或问题的模拟。张友华、涂建强等对湍流导管试验流场进行了工程估算和数值模拟[12-13]，以指导试验方法。李媛、于胜春等利用商业软件对超声速喷管自由射流场进行了分析[14-15]，杨鸿对湍流平板烧蚀试验流场的工程估算、数值模拟和测试数据进行了比较，证明了CFD计算对试验前的状态估算是完全可行的[16-17]，隆永胜利用不同气体状态方程计算了半椭圆喷管流场，得到了湍流和真实气体效应对流场结构的影响[18]。
随着高超声速飞行器的发展，迫切需要进行更精确的防热设计，防热材料的地面考核试验也要求尽可能模拟到更宽范围的弹道参数变化。因此，近年来在电弧加热器平板试验上发展了变迎角试验技术，即在同一车次中可连续变化试件迎角，以增大表面气动参数的匹配范围。这种方法提升了试验的模拟能力，但同时也增加了状态调试的复杂性，使得调试和校测环节占据整个试验的大部分周期和费用。本文基于FLUENT软件，采用CFD方法对电弧加热器湍流平板试验的不同压力和迎角流场进行模拟和分析，研究试件表面气动参数特性分布，以期为试验的实施提供参考，帮助前期的试验策划。
1 计算方法 1.1 控制方程与离散方法计算采用的控制方程为一般直角坐标系中无量纲可压缩N-S方程[19]：
(1)
对控制方程采用二阶迎风精度，基于微元中心有限体积法 (控制容积法) 进行空间离散求解，基于密度的隐式耦合求解方法，即同时求解连续方程、动量方程、能量方程及输运方程的耦合方程组，然后逐一地求解湍流标量方程。
1.2 物理模型及边界条件电弧湍流平板试验的原理如图 1所示，据此建立的计算域如图 2所示。喷管内流区和射流外流区组成整个流体计算域，喷管为出口尺寸100mm×32mm、名义马赫数2.3的二元矩形超声速型面喷管，平板模型长度为150mm。喷管左边入口的高温高压空气经内部型面加速后，作用到与喷管出口齐平安装的平板模型上，其后半自由地射流到大气环境中。为了减小喷管出口射流对边界的影响，将与喷管出口临近的边界左移20倍特征长度的距离，并设置为压力远场条件，上下空间远场边界距离射流区60倍特征长度，喷管入口采用压力入口条件，右边出口采用压力出口条件，模型表面采用固定无滑移等温壁条件。
|
| 图 1 电弧湍流平板试验原理 Fig.1 Scheme of arc-heated turbulent flat plate testing |
|
| 图 2 计算域及边界条件 Fig.2 Domain and boundary condition |
为了保证计算精度和收敛速度，整个计算域采用四边形结构网格，保证网格物面正交并向物面和变化梯度较大的区域按指数规律加密。由于涉及模型表面附面层加热，近壁面的网格划分尤为重要，除了确保足够的网格数量和密度之外，需要根据采用的非平衡壁面处理方法，将壁面第一层无量纲网格距离Yplus控制在30~100之间[20]。Yplus按照如下公式估算：
(2)
(3)
式中：τw为壁面剪切力，Δy在非平衡壁面处理方法中为四分之一近壁面第1和第2个网格点之间距离。划分的网格结构如图 3所示。
|
| 图 3 生成的网格 Fig.3 Grid generation |
湍流模型采用两方程的RNG k-ε湍流模型。RNG k-ε模型的输运方程采用与标准k-ε模型相同的形式, 然而它的导出依据了更为严格的统计技术--重整化群理论, 显著提高了解决快速变形流动问题的精度。在该模型中, 模式常数由重整化群理论给出: C1ε= 1.42, C2ε=1.68。
2 测试试验及结果对比为验证计算方法，进行了有限车次的湍流平板流场校测试验。试验在CARDC (中国空气动力研究与发展中心) 的20MW管式电弧加热器上进行。试验方法及装置布置如图 1及前文所述。
利用平板校测模型 (见图 4) 测量表面热流和压力分布。表面压力是在需要的位置上，法向方向开Φ1mm的测压孔，并连接微型压力传感器进行测量。测量热流的量热塞嵌入校测板中，并与表面平齐，采取措施隔离量热塞与校测板之间的传热。在量热塞瞬态感受流场热流的过程中，可以忽略量热塞侧向的传热，认为表面热流只向其法向方向进行一维热传导，根据一维热传导理论，表面热流可以近似为：
|
| 图 4 校测板 Fig.4 Calibration plate |
(4)
因此，通过测量量热塞的温升速率即可获得表面热流。校测板上测点位置分布如图 5所示，通过调整校测板的位置，也可改变测点在流场中的相对位置。
|
| 图 5 测点分布 Fig.5 Distributing of measurement spot |
电弧加热器点火前，校测板处于流场外部，待电弧加热器运行且状态稳定后，送进机构快速将校测板送入指定的流场位置进行测量，校测板通常在流场中停留1~2s，获得数据后再由送进机构快速送出流场完成测试，图 6为试验过程中的流场照片。
|
| 图 6 试验过程照片 Fig.6 Photo of testing |
验证试验选用了2个典型试验状态condA和condB进行测试并与数值模拟结果进行比较。condA的工况为：总温1700K，总压1.4MPa，迎角0°；condB的工况为：总温3200K，总压1.3MPa，迎角19°。图 7和8显示了2个状态下水冷校测板表面的热流、压力的计算结果和测试数据。从图中可以看出，测试和计算结果在分布趋势上有着较好的一致性，相对误差小于20%，说明文中的计算方法能够有效预测试验条件下校测板的表面气动热参数。
|
| 图 7 表面热流计算结果与测试数据对比 Fig.7 Comparison of computed surface heat flux of the calibration plate with the test data |
|
| 图 8 表面压力计算结果与测试数据对比 Fig.8 Comparison of computed surface pressure of the calibration plate with the test data |
流场模拟分为2个部分：第1部分给定入口气流状态，改变模型迎角进行模拟，为了凸显模型迎角的影响，设定状态使喷管出口静压约为0.1MPa即环境大气压的气流状态；第2部分给定模型迎角，改变入口气流状态进行模拟。
第1部分 (定入口气流状态) 计算：给定入口总压1.4MPa，总温1500K的入口条件，分别对模型迎角0°、3°、6°、9°、12°、15°和18°的情况进行了稳态计算。图 9和10为不同模型迎角情况下喷管出口附近气流马赫数和静压云图。
|
| 图 9 不同迎角气流马赫数云图 Fig.9 Computed flowfield Mach number contours of various angle |
|
| 图 10 不同迎角气流静压云图 Fig.10 Computed flowfield pressure contours of various angle |
从计算结果可以看出模型无迎角时，出射气流在模型上均匀发展，流场内部结构一致。在有迎角时，由于流动方向强迫改变，出射气流在模型头部形成斜激波，使得模型位置马赫数减小，压力增加，激波在气流边界上反射为膨胀波，触及模型壁面后又反射为膨胀波。随着模型迎角的增大，斜激波角度和强度随之增大，反射波系更为明显，并在模型后形成较强的激波串。可见计算方法较好地捕捉了超声速气流激波/激波干扰的现象，波系明显，且波系结构符合基本的气动力学规律和试验结果。
利用修正牛顿理论计算不同迎角校测板前部10mm处表面压力，并与数值计算对比 (见图 11)。2种方法对压力的计算一致性较好，相对误差小于5%。
|
| 图 11 不同迎角校测模型表面压力数值计算与工程计算 Fig.11 CFD computed and inferred surface pressure of the calibration plate with various angle |
为了了解流场对模型表面参数的影响，将不同迎角情况下模型表面冷壁热流和压力分布列于图 12和13。从图中可以看出，0°迎角时，模型表面热流和压力沿气流方向缓慢下降，下降幅度较小，基本可认为整个模型表面参数是均匀的。有迎角时，由于前部斜激波的作用，模型表面热流和压力升高，且受到反射波的影响，尾部壁面参数出现骤降。随着迎角的增大，模型表面参数随之升高，然而激波反射位置逐渐前移，使得模型表面的均匀区逐渐减小。由此可以看出，通过增大模型迎角，能够明显提升表面参数，从0°~18°，均匀区的冷壁热流大约增加至原来的260%。然而这是以牺牲均匀区范围为代价的，18°迎角的均匀区只有约70mm。需要注意的是，此次计算选用的是较早期的小口径喷管，对于后期配置的大口径喷管，这种情况会改善很多，甚至不会出现，这取决于模型长度和喷管出口高度。
|
| 图 12 不同迎角校测模型表面热流分布 Fig.12 Computed surface heat flux of the calibration plate with various angle |
|
| 图 13 不同迎角校测模型表面压力分布 Fig.13 Computed surface pressure of the calibration plate with various angle |
值得注意的是，从15°迎角开始，模型前部出现了一段低热流区域，在18°迎角时这种现象更为明显，区域长度约为25mm。分析认为这是由于大迎角情况下，气流在模型前部的拐角处出现了附面层分离，并且从图 14中的流线能够看出喷管与模型相连的拐角位置的确没有气流贴附壁面。在更大迎角情况下，由于附面层分离再附而导致的低参数区域可能更长，鉴于此，在进行电弧加热器湍流平板试验时，模型迎角不宜超过20°，这也给本文的试验方法提供了参考和指导。
|
| 图 14 拐角处流场流线图 Fig.14 Computed pathline of flowfield on corner |
第2部分 (定迎角) 计算：为了捕捉波系干扰，设定模型的迎角为9°，入口气流总温为1500K，分别对入口气流总压为0.8、1.0、1.4、1.6、1.8和2.0MPa的流场进行稳态数值计算。图 15和16为不同入口总压喷管出口附近气流马赫数和静压云图。从计算结果可以看出，入口压力1.4MPa时，由于喷管出口静压约为0.1MPa，与环境压力相近，喷管口只有由模型引起的斜激波，当入口压力大于1.4MPa时，喷管口上缘出现膨胀波，使气流加速、压力减小，气流边界膨胀。膨胀波向下与激波交汇，由于此区域气流受激波影响压力较大，膨胀波经交汇点后依然是膨胀波，并在模型壁面反射为新的膨胀波，此后又与气流边界上反射的膨胀波交汇。波系交汇、反射以及干扰形成菱形区，随着入口压力的增大，这种现象更加明显。当入口压力小于1.4MPa时，喷管口上缘出现激波，使气流压力增加，气流边界收缩。上下两条激波交汇，并分别在模型壁面和气流边界上反射为激波和膨胀波。由于气流的上边界收缩，这种情况的菱形区域更加复杂。可见计算方法有效捕捉了超声速气流波系交汇、干扰和反射现象。
|
| 图 15 不同总压流场马赫数云图 Fig.15 Computed flowfield Mach number contours of various inlet total pressure |
|
| 图 16 不同总压流场静压云图 Fig.16 Computed flowfield pressure of various inlet total pressure |
将不同入口压力情况下模型表面冷壁热流和压力分布列于图 17和18中。从图中可以看出，当入口压力大于1.4MPa时，模型表面参数能够得到提升，但由于膨胀波的反射干扰，反射点之后的表面参数骤降，均匀区范围大幅减小至50mm，而且若以此方法来提升表面参数，将牺牲大量的加热器功率消耗。随着入口压力的增加，均匀区范围变化很小，但反射点之后的表面参数下降幅度逐渐增大。当入口压力小于1.4MPa时，模型前部表面参数降低，但由于2道激波的反射干扰，表面参数在反射点之后大幅度跳跃变化，且均匀区范围只有不到40mm，无法满足试验要求。鉴于此，在进行电弧加热器湍流平板试验时，要求喷管出口静压大于等于环境压力。
|
| 图 17 不同弧室压力模型表面热流分布 Fig.17 Computed surface heat flux of the calibration plate with various total pressure |
|
| 图 18 不同弧室压力模型表面压力分布 Fig.18 Computed surface pressure of the calibration plate with various total pressure |
(1) 计算方法能够正确捕捉流场激波反射、干扰和附面层分离等现象，典型试验状态的计算结果与试验数据相符，可以有效预测试验条件下的流场结构和参数。
(2) 湍流平板试验中，改变入口气流状态和迎角均可以改变模型表面参数，2种方法各有其局限性。前者的局限性在于均匀区较小，且耗费大量的电弧功率，优势在于提升空间大，且随着参数的提高均匀区变化小；后者的局限性在于随着迎角的增加，均匀区也会逐渐减小 (比前者大)，且过大的迎角将导致附面层分离，优势在于变化效率高，控制灵活经济，而且对于大尺寸喷管情况，模型表面参数可以避免波系干扰，均匀性较好。当然，对于不同的具体情况，需要结合2种方法以满足试验的需求。
(3) 根据试件表面参数分布特性均匀性要求，喷管出口静压应约等于或大于环境压力，模型迎角大于15°时，应充分考虑附面层分离的影响。
(4) 本文仅针对电弧湍流平板试验的中低温度和二维流场条件，后续将进一步开展涉及真实气体效应的高温气流三维流动研究。
| [1] | 中国空气动力研究与发展中心. GJB7050-2010导弹防热材料电弧加热器湍流平板试验方法[S]. 北京: 总装备部军标出版发行部, 2010: 1. China Aerodynamics Research & Development Center. GJB 7050-2010 Plat test method in turbulent flow for thermal protection materials of missile at arc-heater[S]. Beijing:Military Standard Publication and Distribution Section of General Armament Department, 2010:1. |
| [2] | 张志成. 高超声速气动热和热防护[M]. 北京: 国防工业出版社, 2003: 278-279. Zhang Z C. Hypersonic aerothermodynamic and thermal protection[M]. Beijing: National Defence Industry Press, 2003: 278-279. |
| [3] | Tahir G, Kristina S. Computational analysis of arc-jet wedge calibration tests in IHF 6-inch conical nozzle[R]. AIAA-2009-1348, 2009. |
| [4] | Tahir G, Antonella I A, Kristina A S. Computational simulations of panel test facility flow:compression-pad arc-jet tests[R]. AIAA-2011-3635, 2011. |
| [5] | John A B, Tahir G, Frank C, et al. Calibration of the truncated panel test arc-jet facility[R]. AIAA-2009-4090, 2009. |
| [6] | Grant P, Dinesh P, Imelda T S. Arc jet facility test condition predictions using the ADSI code[R]. AIAA-2015-2663, 2015. |
| [7] | Balboni J A, Dinesh P, Imelda T S. Consolidating NASA's arc jets[R]. AIAA-2015-2667, 2015. |
| [8] | Balboni J A, Dinesh P, Imelda T S. On laminar to turbulent transition of arc-jet flow in the NASA Ames panel test facility[R]. AIAA-2012-3304, 2012. |
| [9] | Balboni J A, Dinesh P, Imelda T S. Computational simulation of high enthalpy arc heater flows[R]. AIAA-2006-1183, 2006. |
| [10] | Balboni J A, Dinesh P, Imelda T S. CFD analysis framework for arc-heated flowfields, Ⅱ:shear testing in arc-jets at NASA ARC[R]. AIAA-2009-4081, 2009. |
| [11] | Balboni J A, Dinesh P, Imelda T S. Computational simulations of the 10-MW TP3 arc-jet facility flow[R]. AIAA-2015-3103, 2015. |
| [12] | 张友华, 陈连忠. 超声速湍流导管烧蚀流场稳定性研究[J]. 宇航材料工艺, 2010(4): 64–67. Zhang Y H, Chen L Z. Stability of ablative flow field of supersonic turbulent duct[J]. Aerospace Materials & Technology, 2010(4): 64–67. |
| [13] | 涂建强, 陈连忠, 许考. 高温含水气流条件下燃烧室材料考核的电弧加热试验模拟方法[J]. 实验流体力学, 2015, 29(4): 81–85. Tu J Q, Chen L Z, Xu K. Testing of combustor chamber material in arc jet flow mixing with transverse injected water[J]. Journal of Experiments in Fluid Mechanics, 2015, 29(4): 81–85. |
| [14] | 李媛. 基于Fluent的传热风洞流场品质研究[J]. 机械制造与自动化, 2012, 41(5): 144–149. Li Y. Study of flow quality of heat transfer wind-tunnel base on fluent[J]. Machine Building Automation, 2012, 41(5): 144–149. |
| [15] | 于胜春, 汤龙生. 固体火箭发动机喷管及羽流流场的数值分析[J]. 固体火箭技术, 2004, 27(2): 95–97. Yu S C, Tang L S. Numerical analysis of nozzle and plume flow field of a solid rocket motor[J]. Journal of Solid Rocket Technology, 2004, 27(2): 95–97. |
| [16] | 杨鸿. 超声速湍流平板烧蚀试验技术研究[D]. 长沙: 国防科技大学, 2005. Yang H. The research on the supersonic turbulent flat ablationexperiment technique[D]. Changsha:National University of Defense Technology, 2005. |
| [17] | 杨鸿, 陈伟芳, 柳森. 电弧湍流平板烧蚀矩形喷管研制及应用[J]. 实验流体力学, 2006, 20(1): 27–30. Yang H, Chen W F, Liu S. The development and application on rectangular nozzle of arc turbulent flat plate ablation test[J]. Journal of Experiments in Fluid Mechanics, 2006, 20(1): 27–30. |
| [18] | 隆永胜, 胡振震, 袁竭, 等. 椭圆喷管设计与数值模拟[J]. 实验流体力学, 2015, 29(3): 80–86. Long Y S, Hu Z Z, Yuan J, et al. Design and numerical simulation of an elliptical nozzle[J]. Journal of Experiments in Fluid Mechanics, 2015, 29(3): 80–86. |
| [19] | 邹宁. 超声速喷管设计及其数值模拟和实验研究[D]. 南京: 南京航天航空大学, 2009. Zou N. Supersonic nozzle design and its investigation with numerical and experimental methods[D]. Nanjing:Nanjing University of Aeronautics and Astronautics, 2009. |
| [20] | 王福军. 计算流体动力学分析[M]. 北京: 清华大学出版社, 2004: 161-172. Wang F J. Computational fluid dynamic analysis[M]. Beijing: Tsinghua University Press, 2004: 161-172. |