考虑烧蚀表面形貌影响的气动特性精确预示方法_苏伟.pdf

2023年第1期 导 弹 与 航 天 运 载 技 术（中英文）No.1 2023 总第392期 MISSILES AND SPACE VEHICLES Sum No.392 收稿日期：2022-09-06；修回日期：2022-09-14 文章编号：2097-1974(2023)01-0021-05 DOI：10.7654/j.issn.2097-1974.20230105 考虑烧蚀表面形貌影响的气动特性精确预示方法 苏 伟，董 超，薛 普，刘文伶，杨天鹏（北京航天长征飞行器研究所，北京，100076）摘要：高速飞行器在飞行过程中防热层烧蚀形成的粗糙表面会影响飞行器的气动特性。针对该问题，对烧蚀表面的近壁流动数值模拟方法进行了研究。采用等效砂粒粗糙度和粗糙函数相结合的方法，对粗糙表面形貌以及对边界层湍流流动的影响进行表征。基于高速粗糙平板风洞试验，开展了3种粗糙壁面湍流模型的数值模拟和验证。最后，以10尖锥标模外形开展了不同壁面粗糙条件下的流场数值模拟，分析了壁面粗糙对气动特性的影响。研究结果表明，壁面粗糙会对法向力系数、轴向力系数和压心系数产生较大的影响。关键词：烧蚀表面；近壁流动；数值模拟 中图分类号：V211.3 文献标识码：A Study on Numerical Simulation Method of Near Wall Flow on Ablated Rough Surface Su Wei,Dong Chao,Xue Pu,Liu Wen-ling,Yang Tian-peng(Beijing Institute of Space Long March Vehicle,Beijing,100076)Abstract:During the flight of high-speed flight vehicle,the ablation of the heat shield will form rough surface,which will affect the aerodynamic characteristics of the vehicle.The numerical simulation method of near wall flow on ablated surface topography and boundary layer flow is characterized by the combination of equivalent sand roughness and roughness function.Based on the supersonic rough plate wind tunnel experiment,three numerical simulation methods of rough wall turbulence model are used to simulate flow field and be verified.Finally,taking the 10 degree conical standard model as the research object,the flow field simulation under different rough wall conditions is carried out,and the influence of wall roughness on aerodynamic characteristic parameters is analyzed.The results show that the wall roughness has a great impact on the normal force coefficient,the axial force coefficient and the pressure center coefficient.Key words:ablated rough surface;near wall flow;numerical simulation 0 引 言 高速飞行器具有强机动、超远程、强突防等特点，是当前世界航空航天强国研究的重点方向1。高速飞行器在大气层内飞行时，气动热环境严酷，防热材料表面受气动加热作用会发生烧蚀，使表面变得粗糙。这种粗糙表面会影响边界层的流动特性，产生气动特性预示偏差。对于长时间大气层内巡航飞行的飞行器，气动预示偏差中轴向力预示偏差的时间积分效应更加显著，易导致航程预示散差大的问题，严重制约飞行器性能的提升，因此须对烧蚀粗糙表面流动特征进行研究，发展高效的数值模拟方法，提高飞行器气动预示精度。Cebeci 等2通过修正内层的混合长度值引入粗糙壁面影响。Wilcox 等3在k-量方程模型的基础上修改了涡耗散方程的边界条件；Spalart 等4在原始 S-A 模型的基础上进行修正，提出 Boeing 和 ONERA 两种湍流模型粗糙修正方法，与试验结果对比吻合较好；Hellsten和Laine5对k-SST湍流模型的限制器进行了修改，取得了较好的模拟效果；Knopp 等6在k-SST模型的基础上也进行了粗糙修正方法研究；刘通7研究了粗糙度对表面摩阻系数、斯坦顿数和升力系数的影响。目前关于粗糙壁面流动的研究主要针对的是低速不可压缩流动和高速可压缩边界层转捩问题，对于高速可压缩流动湍流问题研究较少。因此本文将对高速粗糙表面流动模拟方法进行研究，并分析壁面粗糙对 导 弹 与 航 天 运 载 技 术（中英文）2023年 22 气动特性的影响。1 烧蚀表面形貌表征方法 1.1 等效砂粒粗糙高度 高速飞行器多采用非金属热防护材料，烧蚀形成粗糙的典型烧蚀形貌如图 1 所示。为了便于描述，本文采用 Nikuradse8提出的等效砂粒粗糙度方法。该方法建立了等效砂粒粗糙高度ks与粗糙单元形状参数的关系式，并作为不同粗糙度的统一度量，适用于大面积的粗糙表面流动模拟。图1 阿波罗11号返回舱防热层烧蚀形貌9 Fig.1 Ablation Morphology of the Thermal Protective Layer of Apollo 11 Re-entry Capsule9 基于ks的粗糙高度雷诺数+ss/ku kv=，其中 ww/u=（1）式中 w为壁面剪切应力。一般认为+sk5 时粗糙元完全位于粘性底层，被称为水力光滑,对边界层流动的影响可以忽略不计，粘性力占据主导作用；570时被称为完全粗糙，在摩擦阻力中砂粒受到的压力占主导作用，摩擦阻力显著增大。1.2 粗糙函数 Clauser 湍流速度亏损律表明，壁面上存在的小尺度连续粗糙将使得壁面律对数层平均流向速度剖面()Uy+向下偏移，如图 2 所示。图2 湍流速度亏损 Fig.2 Turbulent Velocity Deficiency 无量纲公式为 1lnUyAU+=+（）-（2）式中 ，A分别为 Von Karman 常数和积分常数，=0.41，A=5.5；U+为粗糙函数，表示粗糙表面产生的速度剖面偏移量；y+为摩擦雷诺数，/yu yv+=。U+为粗糙高度的函数：s1lnUkB+=+（）（3）其中常数B与粗糙元分布密度相关，结果表明，当+sk在一定范围内时，公式与试验结果吻合较好10。基于式（1）和式（2）可以得到粗糙壁面的壁面律为 s1lnyUABk+=+（）-（4）通常认为粗糙壁面形貌给定后，粗糙函数为常数。Nikuradse 给出的粗糙函数公式为 s1ln8.5UkA+=+（）-（5）Sigal等11对普通粗糙元形状参数与等效砂粒粗糙高度+sk之间的映射关系式进行了修正，给出的U+的值与粗糙元的形状和分布相关：17.35(1.634log1.01.44.892.24.8913.259.55 1.0-0.686log13.25100.0U+=）（）-（6）9.7 1.00.794log16.020.0U+=（）-（7）其中式（6）和式（7）分别适用于 2D 和 3D 粗糙元情况，3D 粗糙元对波阻的影响小于 2D 粗糙元，因此曲线相对 2D 粗糙元下移，且-1.6sffSASA=（）（8）式中 S为光滑壁面（无粗糙元）的面积；Sf为表面所有粗糙单元迎风面的总面积；Af为单个粗糙元的迎风投影面积；As为单个粗糙元的迎风浸润面积。在该式中，S/Sf为粗糙度密度参数，Af/As为粗糙度形状参数。2 粗糙壁面湍流模型 2.1 B-L 零方程湍流模型 B-L 模型由内外两层组成，内层遵循 Prandtl-Van Driest 模型，外层遵循 Clauser 模型：tinnercrossovertcrossovertouter()()yyyy=（9）内层涡粘系数t,inner的计算方法为 2t,innerelR=（10）式中 为涡量的模。对于光滑壁面，混合长度l：1exp(/)lkyyA+=-（11）对于粗糙壁面，Cebeci2给出的混合长度l修正为()1exp()/)lk yyyyA+=+-+（12）苏 伟等 考虑烧蚀表面形貌影响的气动特性精确预示方法 23第1期 式中 +ss0.9(/)/6yukkv=-。外层涡粘系数touter()的计算方法与原始模型相同。2.2 S-A 一方程湍流模型 S-A 一方程湍流模型公式为 jb1t2j2b1t2v2t2W1W2b2b22jj11(1)()1 1(1)uCfstxCfffCfdCCxxvvvvvvvvv+=-+-+-+（13）式中 3v133v1fC=+,vv=,17.1vC=，湍流粘性系数可以定义为1tvfv=。对于光滑壁面，边界层条件采用狄里克莱条件wall0v=。对于粗糙壁面，Spalart4采用壁面偏移的方法并且采用 Neumann 边界条件：sw,kcnydyvvvv=+?（）（14）式中 y=0.03 ks。2.3 k-SST 两方程湍流模型 k-SST 模型的无量纲形式为 jktkjjj(/)kkkuPktxxx *+=-+（15）2jtjjj12jj(/)12(1)uPtxxxkFxx+=-+-（16）式中 2ktP=，2P=，湍流粘性系数为 1t12max(,)a kaF=（17）式中 a1=0.31，为平均应变率的模，1/2ijij2W W=（）,jii,jji1()2uuWxx=-。对于光滑壁面，湍流能量方程采用狄里克莱边界条件，湍流扩散速率方程采用 Menter12提出的边界条件：wall0k=，1wall21160()d=（18）式中 下标 1 表示靠近壁面第 1 层网格格心处的值。对于粗糙壁面，湍流能量方程边界采用与光滑壁面相同的狄里克莱边界条件，Wilcox3湍流扩散速率方程中对涡耗散率进行了修正以提高近壁面的涡粘系数：2+wallrs1()uS k=（19）Sr的值根据 Schlichting 试验数据得到：2ssrss50/)25100/25kkSkk+=（20）为了防止公式中的限制器在接近壁面时被启用，Hellsten等5在以上公式基础上增加混合函数F3进行了如下修正：41t321231501tanhmax(,)akFaSF Fyv=-，（）（21）3 计算结果与分析 3.1 控制方程 数值模拟控制方程是三维 Navier-Stokes 方程。在直角坐标系下，无量纲化积分形式的 Navier-Stokes 方程为 a1dddvQFFtRe+=SS （22）式中 Q为流场张量；F为无粘张量；vF为粘性张量。采用有限体积法进行求解，空间离散方法采用二阶精度 Roe 格式，时间推进方法采用 LU-SGS 隐式方法。3.2 高速平板 采用Ekoto 等13高速平板风洞试验对粗糙壁面湍流模型的计算精度进行验证。该试验在德州农工大学高速风洞完成，自由来流马赫数为 2.86，总压和总温分别为 690 kPa 和 313 K，壁温Tw=292 K。模型由有机玻璃加工而成。粗糙壁面模型粗糙元为立方体，通过加工 0.79 mm 深的沟槽产生。图