跨声速翼型风洞开孔壁的简化模型数值模拟

RAE2822翼型跨声速绕流1.算例概述RAE2822翼型在跨声速流动中上翼面激波与边界层的分离现象一直是个难题，对于考核计算精度和湍流模型的综合作用是一个很好的算例。

本算例主要测试HyperFlow非结构解算器对2D的一般翼型计算能力，主要考核计算格式、湍流模型的计算能力和网格收敛性。

该算例一般采用湍流计算，测试人曾尝试用层流计算，但计算失败，理所当然的结果。

测试环境：Windows 7 PC，串行计算2.计算条件马赫数（Ma）单位长度雷诺数（Re）攻角（°）侧滑角（°）壁温0.73 6.5e6 2.79 ——参考长度（m）参考面积(m2) 参考点1 1 0，0，0对于RAE2822翼型，Cooke等进行了多种状态实验，本测试计算条件对应于Case 9，测量攻角为3.19○，由于风洞壁面干扰效应，修正后对应的计算攻角为2.79°[1]。

3.计算网格图 1 RAE2822初始计算网格（粗网格）图 2 RAE2822中等网格（上翼面和尾部加密网格）图3 RAE2822密网格（上翼面及激波附件的物面层继续加密）RAE2822初始计算网格四边形/三角形混合网格，二者之间用各向异性四边形均匀过渡，见图1，计算域为x方向100倍弦长半径、y方向50倍弦长半径。

在物面上分布了305个点，壁面第一层网格约为5.0E-06。

网格量为：26004三角形，10032四边形第二套网格是根据初始网格计算结果后进行加密的网格，主要在上上翼面和尾缘进行了网格加密，见图2。

物面分布412个点，第一层网格为1.0E-06。

第二套网格量为：25522三角形，18360四边形第三套网格（密网格）将上翼面加密、且在激波附件的边界层结构网格也进行了加密，见图3。

物面分布464个点，第一层高度1.0E-06。

网格量为：37099三角形，22032四边形4.参数设置参数值备注visname 1eq-sa 粘性类型uns_scheme_name roe 无粘项离散格式uns_limiter_name NONE 限制器eps_vencat 限制器参数ivencat vencat限制器类型uns_vis_name test 粘性项离散格式uns_gradient_name ggcell 梯度重构alf_l，alf_n 0.00001 熵修正系数tscheme 4（LUSGS）时间推进格式cfl_start，cfl_end 2 – 20cfl_nstep 1000sweeps 5 LUSGS扫描次数5.测试结果5.1.残差收敛SA SSTsa-ggnode-vencat sst-ggnodelaplacian图4 粗网格计算残差sstsa图5 中等网格的残差、气动力收敛sstsa图6 密网格的残差、气动力收敛上图可以看到，对于粗网格，采用sa-ggnode加限制器后，残差能很快收敛，但是气动力系数（见后文）却不见得改善。

哈尔滨空气动力研究所:中国航空工业空气动力研究院隶属于中国航空工业第一集团公司，于2000年7月由哈尔滨空气动力研究所（627所）和沈阳空气动力研究所（626所）合并组建，注册地为哈尔滨市，地址在原哈尔滨军事工程学院院内。

气动院现有职工713人，专业技术人员488人，其中，高级技术职务140人，中级技术职务204人，研究生38人，大学本科283人，享受国家政府特贴专家22人。

气动院是国家第一批授予流体力学硕士研究生招生培养权单位，2002年10月国家人事部和全国博士后管委会批准设立博士后科研工作站。

自改革开放以来，先后有118人次出国培训、技术合作和技术考察。

获得国家级科技奖17项，部级科技奖100项。

气动院拥有先进的科研设备，现有低速风洞两座，亚跨音速风洞三座。

经国防科工委批复，填补国内空白的低速增压风洞2002年开始建设。

气动院充分利用自己的技术实力，在非军品科研生产方面取得了长足发展。

以传感器技术、计算机应用和工业自动控制等为主，在油田、烟草、制药、制革、橡胶、铁路和煤炭等行业均取得较好经济效益和社会效益。

1. 风洞设备2.风洞试验技术3.气动力设计与理论研究风洞风洞（wind tunnel），是能人工产生和控制气流，以模拟飞行器或物体周围气体的流动，并可量度气流对物体的作用以及观察物理现象的一种管道状实验设备，它是进行空气动力实验最常用、最有效的工具。

简介风洞实验是飞行器研制工作中的一个不可缺少的组成部分。

它不仅在航空和航天工程的研究和发展中起NF-3低速风洞翼型实验着重要作用，随着工业空气动力学的发展，在交通运输、房屋建筑、风能利用和环境保护等部门中也得到越来越广泛的应用。

用风洞作实验的依据是运动的相对性原理。

实验时，常将模型或实物固定在风洞内，使气体流过模型。

这种方法，流动条件容易控制，可重复地、经济地取得实验数据。

为使实验结果准确,实验时的流动必须与实际流动状态相似,即必须满足相似律的要求。

第８卷㊀第６期２０２３年１１月气体物理ＰＨＹＳＩＣＳＯＦＧＡＳＥＳＶｏｌ.８㊀Ｎｏ.６Ｎｏｖ.２０２３㊀㊀ＤＯＩ:１０.１９５２７/ｊ.ｃｎｋｉ.２０９６￣１６４２.１０８８基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析李泳德ꎬ㊀郭㊀力ꎬ㊀季㊀辰(中国航天空气动力技术研究院ꎬ北京１０００７４)ＣｈａｒａｃｔｅｒｉｚａｔｉｏｎｏｆＴｒａｎｓｏｎｉｃＡｅｒｏｄｙｎａｍｉｃＤａｍｐｉｎｇｏｆＲｏｃｋｅｔｓＢａｓｅｄｏｎＣＦＤ/ＣＳＤＣｏｕｐｌｉｎｇＬＩＹｏｎｇ￣ｄｅꎬ㊀ＧＵＯＬｉꎬ㊀ＪＩＣｈｅｎ(ＣｈｉｎａＡｃａｄｅｍｙｏｆＡｅｒｏｓｐａｃｅＡｅｒｏｄｙｎａｍｉｃｓꎬＢｅｉｊｉｎｇ１０００７４ꎬＣｈｉｎａ)摘㊀要:随着新型大推力火箭的发展ꎬ弯曲模态频率的不断降低ꎬ以及流动分离和跨声速飞行时产生的激波震荡等因素ꎬ其在跨声速飞行过程中更容易出现非定常振动发散ꎮ文章以某带助推的运载火箭模型为研究对象ꎬ通过数值计算获取火箭强迫振动时的气动阻尼ꎬ并对影响火箭气动阻尼的因素进行了分析ꎮ包括结构节点位置㊁振动振幅大小㊁脉动压力等ꎮ研究表明:助推主要起到增大气动阻尼的作用ꎻ前节点主要影响收缩段的气动阻尼ꎻ振动振幅大小和脉动压力对气动阻尼的影响可忽略不计ꎮ关键词:气动阻尼ꎻ数值计算ꎻ跨声速ꎻ气动弹性ꎻ运载火箭㊀㊀㊀收稿日期:２０２３￣０９￣２５ꎻ修回日期:２０２３￣１０￣２３第一作者简介:李泳德(１９９５￣)㊀男ꎬ工学硕士ꎬ助理工程师ꎬ主要研究方向为气动弹性分析ꎮＥ￣ｍａｉｌ:５６２０６４１６９＠ｑｑ.ｃｏｍ通信作者简介:季辰(１９８２￣)㊀男ꎬ工学博士ꎬ研究员ꎬ主要研究方向为气动弹性力学ꎮＥ￣ｍａｉｌ:ｊｉｃｈｅｎ１６７＠ｈｏｔｍａｉｌ.ｃｏｍ中图分类号:Ｖ４７５.１㊀㊀文献标志码:ＡＡｂｓｔｒａｃｔ:Ｗｉｔｈｔｈｅｄｅｖｅｌｏｐｍｅｎｔｏｆｎｅｗｈｉｇｈ￣ｔｈｒｕｓｔｒｏｃｋｅｔｓꎬｔｈｅｄｅｃｒｅａｓｉｎｇｆｒｅｑｕｅｎｃｙｏｆｔｈｅｂｅｎｄｉｎｇｍｏｄｅｓｏｆｔｈｅｒｏｃｋｅｔꎬａｓｗｅｌｌａｓｔｈｅｆａｃｔｏｒｓｓｕｃｈａｓｆｌｏｗｓｅｐａｒａｔｉｏｎａｎｄｓｈｏｃｋｏｓｃｉｌｌａｔｉｏｎｓｇｅｎｅｒａｔｅｄｄｕｒｉｎｇｔｒａｎｓｏｎｉｃｆｌｉｇｈｔｍａｋｅｉｔｍｏｒｅｐｒｏｎｅｔｏｎｏｎ￣ｃｏｎｓｔａｎｔｖｉｂｒａｔｉｏｎ.Ｉｎｔｈｉｓｐａｐｅｒꎬａｌａｕｎｃｈｖｅｈｉｃｌｅｍｏｄｅｌｗｉｔｈｂｏｏｓｔｗａｓｔａｋｅｎａｓｔｈｅｒｅｓｅａｒｃｈｏｂｊｅｃｔꎬａｎｄｔｈｅａｅｒｏｄｙ￣ｎａｍｉｃｄａｍｐｉｎｇｏｆｔｈｅｒｏｃｋｅｔｄｕｒｉｎｇｆｏｒｃｅｄｖｉｂｒａｔｉｏｎｗａｓｏｂｔａｉｎｅｄｔｈｒｏｕｇｈｎｕｍｅｒｉｃａｌｃａｌｃｕｌａｔｉｏｎ.Ｔｈｅｆａｃｔｏｒｓａｆｆｅｃｔｉｎｇｔｈｅａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｆｔｈｅｒｏｃｋｅｔｗｅｒｅａｎａｌｙｚｅｄꎬｉｎｃｌｕｄｉｎｇｔｈｅｐｏｓｉｔｉｏｎｏｆｓｔｒｕｃｔｕｒａｌｎｏｄｅｓꎬｔｈｅｍａｇｎｉｔｕｄｅｏｆｖｉｂｒａｔｉｏｎａｍ￣ｐｌｉｔｕｄｅꎬｐｕｌｓａｔｉｎｇｐｒｅｓｓｕｒｅａｎｄｓｏｏｎ.Ｔｈｅｓｔｕｄｙｓｈｏｗｓｔｈａｔｔｈｅｂｏｏｓｔｍａｉｎｌｙｐｌａｙｓｔｈｅｒｏｌｅｏｆｉｎｃｒｅａｓｉｎｇａｅｒｏｄｙｎａｍｉｃｄａｍｐ￣ｉｎｇａｎｄｔｈｅｆｒｏｎｔｎｏｄｅｍａｉｎｌｙａｆｆｅｃｔｓｔｈｅａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｆｔｈｅｃｏｎｔｒａｃｔｉｏｎｓｅｃｔｉｏｎ.Ｔｈｅｖｉｂｒａｔｉｏｎａｍｐｌｉｔｕｄｅｓｉｚｅａｎｄｔｈｅｐｕｌｓａｔｉｎｇｐｒｅｓｓｕｒｅｈａｖｅａｎｅｇｌｉｇｉｂｌｅｅｆｆｅｃｔｏｎｔｈｅａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇ.Ｋｅｙｗｏｒｄｓ:ａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇꎻｎｕｍｅｒｉｃａｌｃａｌｃｕｌａｔｉｏｎꎻｔｒａｎｓｏｎｉｃꎻａｅｒｏｅｌａｓｔｉｃｉｔｙꎻｌａｕｎｃｈｖｅｈｉｃｌｅ引㊀言通常情况下人们认为气动力对火箭的振动起到阻尼作用ꎬ即气动阻尼为正值ꎮ然而随着大推力火箭发展ꎬ火箭的长细比逐渐加大ꎬ导致弯曲刚度越来越小ꎬ同时为了满足有效载荷的外形要求ꎬ火箭头部整流罩尺寸不断加大ꎬ后续箱体的直径却保持不变ꎬ形成了典型的锤头体外形ꎮ国内外大量的火箭研制经验表明[１￣９]ꎬ对于此类锤头体外形火箭的气动设计ꎬ必须要进行动态气动载荷与动态气弹稳定性分析ꎬ否则设计的疏忽可能会导致火箭结构出现毁灭性的破坏进而导致发射失败ꎮ目前常用的衡量气弹稳定性的方法是通过风洞试验来获取气动阻尼系数ꎮ早在１９６３年ꎬ美国国家航空航天局Ａｍｅｓ研究中心(ＮＡＳＡＡｍｅｓＲｅ￣ｓｅａｒｃｈＣｅｎｔｅｒ)采用半刚性模型开展试验研究[１０]ꎬ获取火箭头部的气动阻尼来评估其稳定性ꎬ但这只能用来模拟火箭弯曲振型前节点之前部分的结构动力学特性ꎮ直到兰利研究中心(ＮＡＳＡＬａｎｇｌｅｙＲｅｓｅａｒｃｈＣｅｎｔｅｒ)开发了全弹性模型气动阻尼试验气体物理２０２３年㊀第８卷技术ꎬ其可以模拟整体的结构动力学特性以及气动外形ꎬ并应用于多款运载火箭研制[１１￣１５]ꎮ国内ꎬ中国航天空气动力技术研究院对气动阻尼问题开展过较多的研究[１６￣２０]ꎬ从模型设计方法㊁模型制作工艺㊁试验机构设计和数据处理等诸多方面ꎬ逐步改进实现了从半刚性模型到全弹性模型的过渡ꎬ并在多个型号上得到验证ꎮ然而通过风洞试验研究气动弹性问题ꎬ技术难度大ꎬ试验成本高ꎬ同时几乎不可能开展全尺寸试验ꎮ因此通过数值计算的方法开展相关研究是另一种重要的手段ꎮ刘子强等[２１]实现了通过数值计算确定气动阻尼系数的技术和方法ꎬ并与试验结果进行对比ꎬ证实了该方法的可靠性ꎮ冉景洪等[２２]通过模态数据结合准定常理论的方法分析了减阻杆加后体这一弹性结构的气动阻尼ꎬ结果表明减阻杆造成的分离流会对后体的气动阻尼系数产生影响ꎮ朱剑等[２３]针对新一代捆绑式运载火箭发展了非结构网格下的气动阻尼计算方法ꎬ并分析了攻角㊁Ｍａｃｈ数等参数对气动阻尼的影响ꎮ本文在之前的计算方法[２３]的基础上采用ＩＤＤＥＳ模型ꎬ考虑脉动压力的影响ꎬ通过强迫振动的方式ꎬ针对捆绑式运载火箭的某一特定模态进行数值计算仿真ꎬ研究前节点位置ꎬ振动振幅ꎬ脉动压力等参数对气动阻尼的影响规律ꎮ１㊀计算方法图１为本文所用的捆绑式运载火箭的计算模型ꎬ是典型的锤头体结构ꎮ在跨声速阶段ꎬ其头部会产生激波造成激波边界层干扰ꎬ而在锤头体外形的过渡段会出现气流分离ꎮ为探究各部分气动阻尼的变化ꎬ将整个箭体分为头部㊁过渡段㊁弹身３个部分ꎮ图１㊀表面网格及区域划分Ｆｉｇ.１㊀Ｓｕｒｆａｃｅｇｒｉｄａｎｄｒｅｇｉｏｎｄｉｖｉｓｉｏｎ１.１㊀流场仿真模型本文分别用Ｒｅｙｎｏｌｄｓ平均法(Ｒｅｙｎｏｌｄｓ￣ａｖｅｒａｇｅｄＮａｖｉｅｒ￣ＳｔｏｋｅｓꎬＲＡＮＳ)和改进的延迟分离涡模拟(ｉｍｐｒｏｖｅｄｄｅｌａｙｅｄｄｅｔａｃｈｅｄ￣ｅｄｄｙｓｉｍｕｌａｔｉｏｎꎬＩＤ￣ＤＥＳ)[２４￣２５]进行计算ꎬ在ＲＡＮＳ方程中ꎬ将变量分为平均值和波动值两部分ꎬ对于速度分量有ｕｉ＝ｕｉ＋ｕᶄｉ其中ꎬｉ＝１ꎬ２ꎬ３ꎬｕｉ和ｕᶄｉ分别代表平均量和波动量ꎬ对于压强和其他标量也采用类似的形式ꎬ将这种形式代入连续性方程和动量方程中ꎬ并写成张量形式∂ρ∂ｔ＋∂∂ｘｉ(ρｕｉ)＝０(１)∂∂ｔ(ρｕｉ)＋∂∂ｘｊ(ρｕｉｕｊ)＝∂ｐ∂ｘｉ＋∂∂ｘｊμ∂ｕｉ∂ｘｊ＋∂ｕｊ∂ｘｉ－２３δｉｊ∂ｕｋ∂ｘｋæèçöø÷éëêêùûúú＋∂∂ｘｊ(－ρｕᶄｉｕᶄｊ)(２)其中ꎬｉꎬｊꎬｋ可分别取１ꎬ２ꎬ３ꎻρ是密度ꎻｔ是时间ꎻ当ｉ＝ｊ时δｉｊ取０ꎬ否则取１ꎮ式(１)㊁(２)是ＲＡＮＳ方程ꎬ由方程可知ＲＡＮＳ方法将湍流脉动对平均流动的作用模化为Ｒｅｙｎｏｌｄｓ应力项即－ρｕᶄｉｕᶄｊꎬ之后采用湍流模型进行封闭ꎬ本文采用的湍流模型为ＳＳＴｋ￣ω模型ꎬ其输运方程为∂∂ｔ(ρｋ)＋∂∂ｘｉ(ρｋｕｉ)＝∂∂ｘｊΓｋ∂ｋ∂ｘｊæèçöø÷＋Ｇｋ－Ｙｋ∂∂ｔ(ρω)＋∂∂ｘｉ(ρωｕｉ)＝∂∂ｘｊΓω∂ω∂ｘｊæèçöø÷＋Ｇω－Ｙω其中ꎬｋ和ω分别代表湍流动能和湍流耗散率ꎬΓｋ和Γω分别代表ｋ和ω的有效扩散系数ꎬＧｋ和Ｇω分别代表ｋ和ω的生成率ꎬＹｋ和Ｙω分别代表ｋ和ω的耗散率ꎮ因此ＲＡＮＳ方法只能计算大尺度的平均流动ꎬ本文采用ＩＤＤＥＳ方法计算脉动压力对气动阻尼的影响ꎮＩＤＤＥＳ方法是由分离涡模拟(ｄｅｔａｃｈｅｄ￣ｅｄｄｙｓｉｍｕｌａｔｉｏｎꎬＤＥＳ)方法改进而来ꎬ其本质思想与ＤＥＳ方法相同ꎬ是想以网格尺度和模型中的特征尺度隐式划分ＲＡＮＳ和大涡模拟(ｌａｒｇｅ￣ｅｄｄｙｓｉｍｕｌａｔｉｏｎꎬＬＥＳ)区域ꎬ使其既能处理ＲＡＮＳ方法无法得到的脉动场ꎬ也能降低ＬＥＳ方法在模拟高Ｒｅｙｎｏｌｄｓ数流动时所需的计算资源ꎮ区别在于当边界层较厚或者分离区域较窄时ꎬＤＥＳ方法会出现如模型应力损耗(ｍｏｄｅｌｅｄｓｔｒｅｓｓｄｅｐｌｅｔｉｏｎꎬＭＳＤ)ꎬ网格诱导分离(ｇｒｉｄ￣ｉｎｄｕｃｅｄｓｅｐａｒａｔｉｏｎꎬＧＩＳ)以及对数层不匹配(ｌｏｇａｒｉｔｈｍｉｃ￣ｌａｙｅｒｍｉｓｍａｔｃｈꎬＬＬＭ)问题[２４]ꎬ而ＩＤＤＥＳ模型通过改良计算区域划分ꎬ结合延迟分离涡模拟(ｄｅｌａｙｅｄｄｅｔａｃｈｅｄ￣ｅｄｄｙｓｉｍｕｌａｔｉｏｎꎬＤＤＥＳ)和０３第６期李泳德ꎬ等:基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析壁面模型大涡模拟(ｗａｌｌ￣ｍｏｄｅｌｅｄｌａｒｇｅ￣ｅｄｄｙｓｉｍｕｌａ￣ｔｉｏｎꎬＷＭＬＥＳ)ꎬ定义新的长度尺度解决了这些问题ꎬ具体公式详见文献[２５]ꎮ流场网格如图２㊁图３所示ꎬ边界层采用棱柱层结构ꎬ并调整第１层网格高度使得ｙ＋小于１ꎬ远场部分采用六面体结构网格ꎬ与边界层的过渡层采用非结构网格ꎮ整体网格单元数量为４.２ˑ１０６ꎮ图２㊀ｙ方向截面网格示意图Ｆｉｇ.２㊀Ｓｃｈｅｍａｔｉｃｄｉａｇｒａｍｏｆｃｒｏｓｓ￣ｓｅｃｔｉｏｎａｌｇｒｉｄｉｎｔｈｅｙ￣ｄｉｒｅｃｔｉｏｎ图３㊀ｘ方向截面网格示意图Ｆｉｇ.３㊀Ｓｃｈｅｍａｔｉｃｄｉａｇｒａｍｏｆｃｒｏｓｓ￣ｓｅｃｔｉｏｎａｌｇｒｉｄｉｎｔｈｅｘ￣ｄｉｒｅｃｔｉｏｎ物面边界条件为无滑移壁面条件ꎬ远场采用压力远场边界条件ꎬ湍流模型采用ＳＳＴｋ￣ω模型ꎬ采用密度基求解ꎬ气体黏性采用Ｓｕｔｈｅｒｌａｎｄ定律ꎬ空间离散采用２阶迎风格式ꎬ对流通量采用Ｒｏｅ格式ꎮ１.２　结构分析模型结构与流场耦合分析过程中ꎬ结构部分可以采用模态方法描述ꎮ结构模态可以通过有限元方法与结构模态试验方法获得ꎮ本文采用有限元分析结果获得的模态ꎬ图４所示为结构的前３阶模态ꎬ本文只分析计算结果中气动阻尼最小的第２阶模态ꎮ(ａ)ｆ＝１.２００Ｈｚ(ｂ)ｆ＝２.４６０Ｈｚ(ｃ)ｆ＝２.９５７Ｈｚ图４㊀结构的前３阶模态Ｆｉｇ.４㊀Ｆｉｒｓｔｔｈｒｅｅｍｏｄｅｓｏｆｔｈｅｓｔｒｕｃｔｕｒｅ由于火箭结构外形简单ꎬ一般不考虑其扭转影响ꎬ因此可以将其简化为简单的梁模型ꎬ这样就可以给出其模态振动方程ｑ㊆ｉ＋２ｂｉωｉｑ˙ｉ＋ω２ｉｑｉ＝ｆｉ(３)式中ꎬｑｉ为第ｉ阶模态的广义位移ꎬｂｉ为第ｉ阶模态的结构阻尼系数ꎬωｉ为第ｉ阶模态的固有频率ꎬ１３气体物理２０２３年㊀第８卷ｆｉ为第ｉ阶模态下质量归一化的广义气动力ꎮ若将ｆｉ按照Ｔａｙｌｏｒ展开并略去高阶项ꎬ可以将其转化为气动阻尼项与气动刚度项的形式ꎬ则式(３)可写为ｑ㊆ｉ＋２(ｂｉ＋Ｂｉ)ωｉｑ˙ｉ＋(Ｋｉ＋１)ω２ｉｑｉ＝０(４)式中ꎬＢｉ为气动阻尼系数ꎬＫｉ为气动刚度系数ꎬ研究表明[２６]ꎬ气动刚度相对于结构刚度为小量可以忽略不计ꎬ而在计算中结构阻尼往往设置为０ꎬ因此气动阻尼可以直接反映其气弹稳定性ꎮ１.３㊀气动阻尼分析原理气动阻尼的分析可以采用强迫振动或者自由振动的方式进行ꎬ这两种方法获得的时域数据不同ꎬ提取气动阻尼的方式也不同ꎮ强迫振动方法初始演化过程较短ꎬ因此计算量较小ꎬ同时能够分析某一种振动形式的气动阻尼ꎬ明确该振动形式是收敛还是发散ꎮ分析过程中能够获得不同部位与部件的气动阻尼ꎮ但是对于多模态相互作用引起的发散(例如颤振)较难预测ꎮ自由振动方法需要一定的自由演化时间才能够对时域数据进行分析ꎬ不过自由振动方法能够获得最能够吸收能量的模态及其振动频率ꎮ对于本研究所关注的问题ꎬ气动载荷对结构振动的过程中气动阻尼的影响较大ꎬ而对气动刚度与气动质量影响较小ꎬ即结构的固有振动频率受到来流的影响较小ꎬ其稳定性问题主要由气动阻尼的正㊁负引起ꎬ所以采用强迫振动方法分析ꎮ强迫振动下结构做简谐模态振动ｑｉ(ｔ)＝Ａｓｉｎ(ωｉｔ)式中ꎬＡ表示振动的振幅ꎬ将其代入计算气动力的公式中[２１]并做正交积分可得Ｂｉ＝ʏｌ０Ｂｘ(ｘ)ｄｘ＝－１ＭｉＡω２ｉＴʏｌ０ʏｔ０＋Ｔｔ０Ｇ(ｘꎬｔ)ｃｏｓ(ωｉｔ)ｄｔｄｘ(５)式中ꎬＭｉ为第ｉ阶模态的模态质量ꎬＴ为整数倍周期ꎬＧ为广义气动力ꎮ根据式(５)便可以得到局部或分区域的气动阻尼ꎮ１.４㊀耦合计算流程首先进行模态分析ꎬ以确定结构的模态频率与振型ꎬ用以设计强迫振动的频率和振幅ꎮ非定常流场计算前先进行定常流场计算ꎬ来加快非定常计算的演化速度并增强收敛性ꎬ结构节点位移通过径向基函数(ＲＢＦ)插值方法[２７]映射到气动网格节点上ꎬ来进行网格的变形ꎬ这里径向基函数选用ＷｅｎｄｌａｎｄＣ２ꎬ如下所示φ(ｘ)＝(１－ｘ)４(４ｘ＋１)最后将计算出来的广义力提取出来ꎬ截取演化完毕的整数倍周期ꎬ进行气动阻尼计算ꎮ耦合计算流程图如图５所示ꎮ图５㊀耦合计算流程图Ｆｉｇ.５㊀Ｆｌｏｗｃｈａｒｔｏｆｃｏｕｐｌｅｄｃａｌｃｕｌａｔｉｏｎ２㊀结果分析与讨论２.１㊀流场分析结果计算的来流Ｍａｃｈ数范围为０.７~１.２ꎮ其中中截面的压力分布如图６所示ꎮ可以看出在头部出现了膨胀波以及跨声速激波ꎬ在过渡段存在流动分离ꎬ随着Ｍａｃｈ数的增大ꎬ头部低压区域逐渐扩张ꎬ并且能明显看到ꎬ在流动再附的位置产生了再附激波ꎮ(ａ)Ｍａ＝０.７０２３第６期李泳德ꎬ等:基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析(ｂ)Ｍａ＝０.７５(ｃ)Ｍａ＝０.８０(ｄ)Ｍａ＝０.８５(ｅ)Ｍａ＝０.８８(ｆ)Ｍａ＝０.９０(ｇ)Ｍａ＝０.９２(ｈ)Ｍａ＝０.９６(ｉ)Ｍａ＝０.９８３３气体物理２０２３年㊀第８卷(ｊ)Ｍａ＝１.００(ｋ)Ｍａ＝１.０５(ｌ)Ｍａ＝１.１０图６㊀不同Ｍａｃｈ数下的中截面压力分布Ｆｉｇ.６㊀ＰｒｅｓｓｕｒｅｄｉｓｔｒｉｂｕｔｉｏｎｉｎｔｈｅｍｉｄｄｌｅｓｅｃｔｉｏｎａｔｄｉｆｆｅｒｅｎｔＭａｃｈｎｕｍｂｅｒｓ２.２　气动阻尼分布通过上述流场分析ꎬ可以看出火箭不同部位流动结构并不相同ꎬ在头部与箭身上ꎬ流动主要为附着流动ꎬ而在过渡段会出现较为复杂的波系结构以及流动分离ꎮ针对不同的流动结构随流向站位ｘ的变化ꎬ设该位置上广义力与广义位移的相位差为φ(ｘ)ꎬ并且简谐振动没有引入其他模态的广义力ꎬ则广义力的表达式为Ｇ(ｘꎬｔ)＝Ｆｇｅｎ ｓｉｎ[ωｔ＋φ(ｘ)]＋Ｆ０(６)其中ꎬＦｇｅｎ为广义力的振动幅度ꎬＦ０为广义力的常数偏移量ꎮ将式(６)代入到式(５)中得到Ｂ(ｘ)＝－ＦｇｅｎＭＡω２Ｔʏｔ０＋Ｔｔ０ｓｉｎ[ωｔ＋φ(ｘ)]ｃｏｓ(ωｔ)ｄｔ其中ꎬ广义力的常数偏移量Ｆ０的积分为０ꎬ因此省略ꎮ通过将等式中的正弦函数部分进行和差化积得到Ｂ(ｘ)＝－ＦｇｅｎＭＡω２Ｔʏｔ０＋Ｔｔ０ｓｉｎ(ωｔ)ｃｏｓ[φ(ｘ)]ｃｏｓ(ωｔ)ｄｔ＋[ʏｔ０＋Ｔｔ０ｓｉｎ[φ(ｘ)]ｃｏｓ(ωｔ)ｃｏｓ(ωｔ)ｄｔ](７)式(７)中第１部分在整个周期中的积分为０ꎬ只有第２部分保留ꎬ因此得到Ｂ(ｘ)＝－Ｆｇｅｎｓｉｎ[φ(ｘ)]ＭＡω２Ｔʏｔ０＋Ｔｔ０ｃｏｓ２(ωｔ)ｄｔ(８)式(８)中积分部分恒为正值ꎬ决定整个气动阻尼的部分只有相位角φ(ｘ)的正弦值ｓｉｎ[φ(ｘ)]ꎬ为了能够更加直观地获得相位角与气动阻尼Ｂ之间的关系ꎬ须将符号转化为对应的正弦函数转角ꎬ根据正弦关系ꎬ此转角为πꎬ因此得到Ｂ(ｘ)＝－Ｆｇｅｎ(ｘ)ｓｉｎ[φ(ｘ)＋π]ＭＡω２Ｔʏｔ０＋Ｔｔ０ｃｏｓ２(ωｔ)ｄｔ(９)图７为气动阻尼变化曲线ꎬ可以看出随着Ｍａｃｈ数的增大ꎬ整体气动阻尼先增大后减少ꎬ在Ｍａｃｈ数为０.９８时达到最大值ꎬ过渡段与箭体的气动阻尼变化趋势与整体基本相同ꎬ而头部区域则不同ꎬ是随着Ｍａｃｈ数的增大一直增大ꎬ只是增长速率变缓ꎮ图７㊀有助推时气动阻尼变化曲线Ｆｉｇ.７㊀Ａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｃｈａｎｇｅｃｕｒｖｅｗｉｔｈｂｏｏｓｔ根据式(９)ꎬ得到相位角与气动阻尼Ｂ之间的关系为:当φ(ｘ)ɪ(－πꎬ０)时ꎬ相位角滞后ꎬ气动阻尼Ｂ为负值ꎻ当φ(ｘ)ɪ(０ꎬπ)ꎬ相位角提前ꎬ４３第６期李泳德ꎬ等:基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析气动阻尼Ｂ为正值ꎻ为当φ(ｘ)＝０时ꎬ无相位角差别ꎬ气动阻尼Ｂ为０ꎮ在过渡段上ꎬ复杂的波系结构以及流动分离ꎬ使得气动力与结构位移之间会出现较为明显的迟滞现象ꎬ从而导致相位角φ(ｘ)ɪ(－πꎬ０)ꎬ由此在过渡段上产生了负的气动阻尼ꎮ计算过程中的广义力与广义位移随时间变化曲线如图８所示ꎬ可以看出所有工况计算结果都表现良好ꎬ需要注意的是在非定常计算初期ꎬ演化的不完全导致广义力存在一些突变异常的结果ꎬ计算气动阻尼时须剔除ꎬ选择后面演化完全的周期ꎮ本文计算了９个周期ꎬ剔除了第１个周期出现的错误结果ꎬ采用后８个周期进行气动阻尼分析ꎮ强迫运动振幅为芯级直径的０.５％ꎮ(ａ)Ｍａ＝０.７０㊀㊀㊀(ｂ)Ｍａ＝０.７５(ｃ)Ｍａ＝０.８０㊀㊀㊀(ｄ)Ｍａ＝０.８５(ｅ)Ｍａ＝０.８８㊀㊀㊀(ｆ)Ｍａ＝０.９０５３气体物理２０２３年㊀第８卷(ｇ)Ｍａ＝０.９２㊀㊀㊀(ｈ)Ｍａ＝０.９６(ｉ)Ｍａ＝０.９８㊀㊀㊀(ｊ)Ｍａ＝１.００(ｋ)Ｍａ＝１.０５㊀㊀㊀(ｌ)Ｍａ＝１.１０图８㊀不同工况下的广义力与广义位移随时间变化曲线Ｆｉｇ.８㊀Ｔｉｍｅｄｅｐｅｎｄｅｎｔｃｕｒｖｅｓｏｆｇｅｎｅｒａｌｉｚｅｄｆｏｒｃｅａｎｄｇｅｎｅｒａｌｉｚｅｄｄｉｓｐｌａｃｅｍｅｎｔｕｎｄｅｒｄｉｆｆｅｒｅｎｔｏｐｅｒａｔｉｎｇｃｏｎｄｉｔｉｏｎｓ２.３㊀气动阻尼影响因素２.３.１㊀有无助推对气动阻尼的影响捆绑式运载火箭相比于传统的运载火箭ꎬ最大的区别就是在尾部四周捆绑了助推器ꎬ使得其流场特性变得复杂ꎬ因此须分析其对气动阻尼的影响ꎮ图７㊁图９分别为有无助推时气动阻尼变化曲线ꎬ可以看出随着Ｍａｃｈ数的增大整体气动阻尼先增大后减少ꎬ在Ｍａｃｈ数为０.９８时达到最大值ꎬ过６３第６期李泳德ꎬ等:基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析渡段与箭体的气动阻尼变化趋势与整体基本相同ꎬ而头部区域则不同ꎬ是随着Ｍａｃｈ数的增大一直增大ꎬ只是增长速率变缓ꎮ对比两个图可知ꎬ助推主要起增大气动阻尼的作用ꎮ还可以看出有无助推情况下头部的气动阻尼变化很小ꎬ意味着在箭体尾部施加控制很难影响到头部的气动阻尼ꎬ特别是在超声速流场中ꎮ图９㊀无助推时气动阻尼变化曲线Ｆｉｇ.９㊀Ａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｃｈａｎｇｅｃｕｒｖｅｗｉｔｈｏｕｔｂｏｏｓｔ２.３.２㊀前节点位置影响为了考察前节点位置变化对气动阻尼的影响ꎬ在保持振动频率不变㊁头部最大振型位置与振幅不变的条件下移动前节点ꎬ变化后的振型如图１０所示ꎮ(ａ)Ｆｒｏｎｔｎｏｄｅａｆｔｅｒｔｈｅｔｒａｎｓｉｔｉｏｎｒｅｇｉｏｎ(ｂ)Ｆｒｏｎｔｎｏｄｅｉｎｔｈｅｔｒａｎｓｉｔｉｏｎｒｅｇｉｏｎ(ｃ)Ｆｒｏｎｔｎｏｄｅｂｅｆｏｒｅｔｈｅｔｒａｎｓｉｔｉｏｎｒｅｇｉｏｎ图１０㊀前节点变化后的振型Ｆｉｇ.１０㊀Ｖｉｂｒａｔｉｏｎｍｏｄｅａｆｔｅｒｔｈｅｃｈａｎｇｅｏｆｆｏｒｍｅｒｎｏｄｅ根据对计算结果的分析分别获得了不同前节点位置的整体气动阻尼对比与过渡段气动阻尼对比ꎬ如图１１㊁图１２所示ꎬ可以看出前节点位置的改变并没有影响整体气动阻尼随Ｍａｃｈ数增大而增大的趋势ꎬ且前节点在过渡段上与过渡段前的整体气动阻尼相差不大ꎬ而前节点在过渡段后的整体气动阻尼要高于另两种情况ꎬ因此过渡段与头部放在同一侧有助于提高气动阻尼ꎮ过渡段的气动阻尼会随着前节点的变化发生剧烈改变ꎬ前节点在过渡段前后随Ｍａｃｈ数增大的变化规律相反ꎬ节点前后的振动相位变化导致不同节点位置过渡段的振动相位不同ꎬ进而导致气动阻尼发生变化ꎮ图１１㊀不同节点位置的整体气动阻尼Ｆｉｇ.１１㊀Ｏｖｅｒａｌｌａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇａｔｄｉｆｆｅｒｅｎｔｎｏｄｅｐｏｓｉｔｉｏｎｓ图１２㊀不同节点位置的过渡段气动阻尼Ｆｉｇ.１２㊀Ａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｆｔｈｅｔｒａｎｓｉｔｉｏｎｒｅｇｉｏｎａｔｄｉｆｆｅｒｅｎｔｎｏｄｅｐｏｓｉｔｉｏｎｓ２.３.３㊀强迫振动振幅大小对气动阻尼的影响为了考察强迫振动振幅大小对气动阻尼的影响ꎬ在保证流场结构不发生改变的前提下ꎬ振动振幅分别为原来的一半和两倍ꎬ根据工程经验ꎬ如果振幅超过芯级直径的５％ꎬ则须考虑流场结构改变所造成的影响ꎮ图１３㊁图１４分别为不同振幅下的整体与头部气动阻尼ꎮ７３气体物理２０２３年㊀第８卷图１３㊀不同振幅下整体气动阻尼Ｆｉｇ.１３㊀Ｏｖｅｒａｌｌａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇａｔｄｉｆｆｅｒｅｎｔａｍｐｌｉｔｕｄｅｓ图１４㊀不同振幅下头部气动阻尼Ｆｉｇ.１４㊀Ａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｆｔｈｅｈｅａｄｒｅｇｉｏｎａｔｄｉｆｆｅｒｅｎｔａｍｐｌｉｔｕｄｅｓ可以发现改变振幅无论是对整体气动阻尼还是头部气动阻尼来说变化都很小ꎬ这意味着气动阻尼的大小主要取决于气动力与结构振动的相位差ꎬ不依赖于振动幅度的大小ꎮ２.３.４㊀脉动压力对气动阻尼的影响为了模拟出脉动压力的影响ꎬ采用ＩＤＤＥＳ方法对火箭气动阻尼进行计算ꎬ计算来流Ｍａｃｈ数为０.９２ꎬ计算过程中的广义力与广义位移如图１５所示ꎬ相较于图８可以看出广义力随时间变化曲线并不光滑ꎬ脉动压力的存在导致广义力由多个频率叠加而成ꎮ由于第２阶模态的频率为２.４６Ｈｚꎬ而由分离流㊁激波振荡等引起的脉动压力频率往往远大于此频率ꎬ因此这里选择３.５Ｈｚ为分界ꎬ将高于３.５Ｈｚ的部分视为由抖振脉动压力引起的广义力ꎬ低于３.５Ｈｚ的部分视为强迫振动引起的广义力ꎬ通过低通滤波把高于３.５Ｈｚ的广义力滤掉ꎬ可以获得由强迫振动引起的广义力与广义位移变化曲线ꎬ如图１６所示ꎬ通过此广义力计算的气动阻尼为２.０８ɢꎮ同样地ꎬ进行高通滤波将低于３.５Ｈｚ的广义力滤掉ꎬ可以获得由抖振脉动压力引起的气动阻尼为(２.９４ˑ１０－３)ɢꎬ由此得到脉动压力引起的气动阻尼变化为０.１４％ꎬ可以忽略不计ꎮ同时使用ＲＡＮＳ方法计算的气动阻尼为２.０７ɢꎬ与ＩＤＤＥＳ的计算结果相比误差约为(２.９４ˑ１０－３＋２.０８－２.０７)/２.０７ʈ０.４８％ꎬ这说明针对气动阻尼的模拟ꎬ抖振引起的脉动压力对气动阻尼的计算结果影响很小ꎬ起主要作用的还是广义力的变化ꎬ该变化由强迫振动引起的结构边界变化所导致ꎮ图１５㊀基于ＩＤＤＥＳ的广义力与广义位移变化曲线Ｆｉｇ.１５㊀ＶａｒｉａｔｉｏｎｃｕｒｅｓｏｆｇｅｎｅｒａｌｉｚｅｄｆｏｒｃｅａｎｄｇｅｎｅｒａｌｉｚｅｄｄｉｓｐｌａｃｅｍｅｎｔｂａｓｅｄｏｎＩＤＤＥＳ图１６㊀滤波后的广义力与广义位移变化曲线Ｆｉｇ.１６㊀Ｖａｒｉａｔｉｏｎｃｕｒｅｓｏｆｇｅｎｅｒａｌｉｚｅｄｆｏｒｃｅａｎｄｇｅｎｅｒａｌｉｚｅｄｄｉｓｐｌａｃｅｍｅｎｔｖａｒｉａｔｉｏｎｃｕｒｖｅａｆｔｅｒｆｉｌｔｅｒｉｎｇ３㊀结论本文通过数值计算方法研究了火箭的气动阻尼特性ꎮ根据流动特征分析与理论推导ꎬ发现火箭过渡段几何外形的收缩导致该区域出现复杂的分离与激波结构ꎬ从而造成了气动力相对于结构振动８３第６期李泳德ꎬ等:基于ＣＦＤ/ＣＳＤ耦合的火箭跨声速气动阻尼特性分析相位的滞后ꎬ导致了该区域为气动负阻尼ꎬ即气动不稳定性的主要来源ꎮ在此机理的基础上ꎬ分析了前节点位置㊁振动振幅㊁脉动压力等因素对气动阻尼的影响规律ꎮ可以得出以下结论:１)助推增加了正阻尼区域的面积ꎬ从而相对于没有助推的构型起到了增加气动阻尼的作用ꎮ２)前节点位置的改变对过渡段气动阻尼影响很大ꎬ节点前后的振动方向相反ꎬ导致节点在过渡段前后的气动阻尼变化规律也截然相反ꎬ将过渡段与头部区域放在节点的同一侧有助于增加气动阻尼ꎮ３)在不改变流场结构的前提下ꎬ改变振动的振幅ꎬ气动力也会产生相应幅度的变化ꎬ因此结构振幅对气动阻尼的影响可忽略不计ꎮ４)高频部分的广义力对气动阻尼的贡献很小ꎬ即结构振动引起的广义力变化对气动阻尼起主要作用ꎬ而脉动压力对计算气动阻尼影响不大ꎬ可忽略不计ꎮ参考文献(Ｒｅｆｅｒｅｎｃｅｓ)[１]㊀ＣｏｅＣＦ.Ｓｔｅａｄｙａｎｄｆｌｕｃｔｕａｔｉｎｇｐｒｅｓｓｕｒｅｓａｔｔｒａｎｓｏｎｉｃｓｐｅｅｄｓｏｎｔｗｏｓｐａｃｅ￣ｖｅｈｉｃｌｅｐａｙｌｏａｄｓｈａｐｅ[Ｒ].ＮＡＳＡＴＭＸ￣５０３ꎬ１９６１.[２]ＣｏｌｅＳＲＪｒꎬＨｅｎｎｉｎｇＴＬꎬＲａｉｎｅｙＡＧ.ＮＡＳＡｓｐａｃｅｖｅ￣ｈｉｃｌｅｄｅｓｉｇｎｃｒｉｔｅｒｉａ[Ｒ].ＮＡＳＡＳＰ￣８００１(ＲＥＶ)ꎬ１９６４. [３]ＥｒｉｃｓｓｏｎＬＥꎬＲｅｄｉｎｇＪＰ.Ａｎａｌｙｓｉｓｏｆｆｌｏｗｓｅｐａｒａｔｉｏｎｅｆｆｅｃｔｓｏｎｔｈｅｄｙｎａｍｉｃｓｏｆａｌａｒｇｅｓｐａｃｅｂｏｏｓｔｅｒ[Ｊ].Ｊｏｕｒ￣ｎａｌｏｆＳｐａｃｅｃｒａｆｔａｎｄＲｏｃｋｅｔｓꎬ１９６５ꎬ２(４):４８１￣４９０. [４]ＲｅｄｉｎｇＪＰꎬＥｒｉｃｓｓｏｎＬＥ.Ｅｆｆｅｃｔｏｆａｅｒｏｅｌａｓｔｉｃｃｏｎｓｉｄｅｒａ￣ｔｉｏｎｓｏｎｓｅａｓａｔ￣Ａｐａｙｌｏａｄｓｈｒｏｕｄｄｅｓｉｇｎ[Ｊ].ＪｏｕｒｎａｌｏｆＳｐａｃｅｃｒａｆｔａｎｄＲｏｃｋｅｔｓꎬ１９８１ꎬ１８(３):２４１￣２４７. [５]程镇煌.宇航飞行器跨音速气动弹性问题探讨[Ｊ].上海航天ꎬ１９９７(６):１６￣２１.ＣｈｅｎｇＺＨ.Ｄｉｓｃｕｓｓｉｏｎｏｎｔｒａｎｓｏｎｉｃｐｎｅｕｍａｔｉｃｅｌａｓｔｉｃｉｔｙｏｆｓｐａｃｅｃｒａｆｔ[Ｊ].ＡｅｒｏｓｐａｃｅＳｈａｎｇｈａｉꎬ１９９７(６):１６￣２１(ｉｎＣｈｉｎｅｓｅ).[６]倪嘉敏.我国运载火箭气动设计回顾[Ｃ].近代空气动力学研讨会论文集ꎬ２００５.ＮｉＪＭ.ＲｅｖｉｅｗｏｆａｅｒｏｄｙｎａｍｉｃｄｅｓｉｇｎｏｆＣｈｉｎａᶄｓｌａｕｎｃｈｖｅｈｉｃｌｅ[Ｃ].Ｐｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅｍｏｄｅｒｎａｅｒｏｄｙｎａｍｉｃｓｓｙｍｐｏｓｉｕｍꎬ２００５(ｉｎＣｈｉｎｅｓｅ).[７]臧涛成ꎬ胡焕性.大长细比弹箭弹性效应研究综述[Ｊ].弹道学报ꎬ１９９９ꎬ１１(３):８９￣９３ꎬ９６.ＺａｎｇＴＣꎬＨｕＨＸ.Ａｒｅｖｉｅｗｏｆｇｒｅａｔｓｌｅｎｄｅｒｎｅｓｓｒａｔｉｏｐｒｏｊｅｃｔｉｌｅｅｌａｓｔｉｃｅｆｆｅｃｔｒｅｓｅａｒｃｈ[Ｊ].ＪｏｕｒｎａｌｏｆＢａｌｌｉｓｔｉｃｓꎬ１９９９ꎬ１１(３):８９￣９３ꎬ９６(ｉｎＣｈｉｎｅｓｅ). [８]吴志刚ꎬ杨超.细长体弹箭的气动弹性问题与研究方法[Ｃ].第九届(２００５年)全国空气弹性学术交流会论文ꎬ２００５.ＷｕＺＧꎬＹａｎｇＣ.Ａｅｒｏｅｌａｓｔｉｃｉｔｙｐｒｏｂｌｅｍｓａｎｄｒｅｓｅａｒｃｈｍｅｔｈｏｄｓｏｆｓｌｅｎｄｅｒｂｏｄｙｍｉｓｓｉｌｅｓ[Ｃ].９ｔｈＮａｔｉｏｎａｌＳｙｍ￣ｐｏｓｉｕｍｏｎＡｅｒｏｅｌａｓｔｉｃｉｔｙꎬ２００５(ｉｎＣｈｉｎｅｓｅ). [９]张贺ꎬ黄晓鹏.弹性细长旋转弹箭运动稳定性问题的研究进展[Ｃ].第九届(２００５年)全国空气弹性学术交流会论文ꎬ２００５.ＺｈａｎｇＨꎬＨｕａｎｇＸＰ.Ｒｅｓｅａｒｃｈｐｒｏｇｒｅｓｓｏｎｔｈｅｓｔａｂｉｌｉｔｙｏｆｅｌａｓｔｉｃｓｌｅｎｄｅｒｒｏｔａｔｉｎｇｐｒｏｊｅｃｔｉｌｅｓ[Ｃ].９ｔｈＮａｔｉｏｎａｌＳｙｍｐｏｓｉｕｍｏｎＡｅｒｏｅｌａｓｔｉｃｉｔｙꎬ２００５(ｉｎＣｈｉｎｅｓｅ). [１０]ＣｏｌｅＨＡＪｒ.Ｄｙｎａｍｉｃｒｅｓｐｏｎｓｅｏｆｈａｍｍｅｒｈｅａｄｌａｕｎｃｈｖｅｈｉｃｌｅｓｔｏｔｒａｎｓｏｎｉｃｂｕｆｆｅｔｉｎｇ[Ｒ].ＮＡＳＡＴＮＤ￣１９８２ꎬ１９６３.[１１]ＢａｒｔｅｌｓＲＥꎬＷｉｅｓｅｍａｎＣＤꎬＭｉｎｅｃｋＲＥ.ＣｏｍｐｕｔａｔｉｏｎａｌａｅｒｏｅｌａｓｔｉｃａｎａｌｙｓｉｓｏｆｔｈｅＡｒｅｓｌａｕｎｃｈｖｅｈｉｃｌｅｄｕｒｉｎｇａｓ￣ｃｅｎｔ[Ｒ].ＡＩＡＡ２０１０￣４３７４ꎬ２０１０.[１２]ＡｚｅｖｅｄｏＪＬＦ.Ａｅｒｏｅｌａｓｔｉｃａｎａｌｙｓｉｓｏｆｌａｕｎｃｈｖｅｈｉｃｌｅｓｉｎｔｒａｎｓｏｎｉｃｆｌｉｇｈｔ[Ｊ].ＪｏｕｒｎａｌｏｆＳｐａｃｅｃｒａｆｔａｎｄＲｏｃｋｅｔｓꎬ１９８９ꎬ２６(１):１４￣２３.[１３]ＳｉｎｃｌａｉｒＡꎬＦｌｏｗｅｒｓＧ.Ｌｏｗ￣ｏｒｄｅｒａｅｒｏｅｌａｓｔｉｃｍｏｄｅｌｏｆｌａｕｎｃｈ￣ｖｅｈｉｃｌｅｄｙｎａｍｉｃｓ[Ｒ].ＡＩＡＡ２０１０￣７７２５ꎬ２０１０. [１４]ＤｏｔｓｏｎＫＷ.Ｔｒａｎｓｉｅｎｔｃｏｕｐｌｉｎｇｏｆｌａｕｎｃｈｖｅｈｉｃｌｅｂｅｎｄｉｎｇｒｅｓｐｏｎｓｅｓｗｉｔｈａｅｒｏｄｙｎａｍｉｃｆｌｏｗｓｔａｔｅｖａｒｉａｔｉｏｎｓ[Ｊ].ＪｏｕｒｎａｌｏｆＳｐａｃｅｃｒａｆｔａｎｄＲｏｃｋｅｔｓꎬ２００１ꎬ３８(１):９７￣１０４.[１５]ＣｏｌｅＳＲꎬＨｅｎｎｉｎｇＴＬ.Ｂｕｆｆｅｔｒｅｓｐｏｎｓｅｏｆａｈａｍｍｅｒｈｅａｄｌａｕｎｃｈｖｅｈｉｃｌｅｗｉｎｄ￣ｔｕｎｎｅｌｍｏｄｅｌ[Ｊ].ＪｏｕｒｎａｌｏｆＳｐａｃｅｃｒａｆｔａｎｄＲｏｃｋｅｔｓꎬ１９９２ꎬ２９(３):３７９￣３８５.[１６]崔尔杰.流固耦合力学研究与应用进展[Ｃ].钱学森科学贡献暨学术思想研讨会论文集ꎬ２００１.ＣｕｉＥＪ.Ｒｅｓｅａｒｃｈａｎｄａｐｐｌｉｃａｔｉｏｎｐｒｏｇｒｅｓｓｏｆｆｌｕｉｄ￣ｓｔｒｕｃ￣ｔｕｒｅｉｎｔｅｒａｃｔｉｏｎｍｅｃｈａｎｉｃｓ[Ｃ].ＰｒｏｃｅｅｄｉｎｇｓｏｆＳｅｍｉｎａｒｏｆＱｉａｎＸｕｅｓｅｎＳｃｉｅｎｔｉｆｉｃＣｏｎｔｒｉｂｕｔｉｏｎｓａｎｄＡｃａｄｅｍｉｃＴｈｏｕｇｈｔｓꎬ２００１(ｉｎＣｈｉｎｅｓｅ).[１７]冯明溪ꎬ王志安.火箭跨音速动导数和抖振实验[Ｊ].宇航学报ꎬ１９８７(１):５５￣６３.ＦｅｎｇＭＸꎬＷａｎｇＺＡ.Ｅｘｐｅｒｉｍｅｎｔｓｏｆｔｒａｎｓｏｎｉｃｄｅｒｉｖａ￣ｔｉｖｅｓａｎｄｂｕｆｆｅｔｉｎｇｏｆｒｏｃｋｅｔ[Ｊ].ＪｏｕｒｎａｌｏｆＡｓｔｒｏｎａｕｔｉｃｓꎬ１９８７(１):５５￣６３(ｉｎＣｈｉｎｅｓｅ).[１８]白葵ꎬ冯明溪.弹性模型实验技术[Ｊ].流体力学实验与测量ꎬ１９９９ꎬ１３(１):３８￣４２.ＢａｉＫꎬＦｅｎｇＭＸ.Ａｅｒｏｅｌａｓｔｉｃｍｏｄｅｌａｎｄｔｈｅｂｕｆｆｅｔｅｘ￣ｐｅｒｉｍｅｎｔａｌｔｅｃｈｎｉｑｕｅ[Ｊ].ＥｘｐｅｒｉｍｅｎｔｓａｎｄＭｅａｓｕｒｅｍｅｎｔｓｉｎＦｌｕｉｄＭｅｃｈａｎｉｃｓꎬ１９９９ꎬ１３(１):３８￣４２(ｉｎＣｈｉｎｅｓｅ). [１９]ＪｉＣꎬＲａｎＪＨꎬＬｉＦꎬｅｔａｌ.Ｔｈｅａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇ９３气体物理２０２３年㊀第８卷ｔｅｓｔｏｆｅｌａｓｔｉｃｌａｕｎｃｈｖｅｈｉｃｌｅｍｏｄｅｌｉｎｔｒａｎｓｏｎｉｃｆｌｏｗ[Ｃ].Ｐｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅ６４ｔｈＩｎｔｅｒｎａｔｉｏｎａｌＡｓｔｒｏｎａｕｔｉｃａｌＣｏｎ￣ｇｒｅｓｓꎬ２０１３.[２０]季辰ꎬ吴彦森ꎬ何岗ꎬ等.运载火箭气动阻尼风洞试验研究[Ｃ].第十二届全国空气弹性学术交流会论文集ꎬ２０１１.ＪｉＣꎬＷｕＹＳꎬＨｅＧꎬｅｔａｌ.Ｅｘｐｅｒｉｍｅｎｔａｌｓｔｕｄｙｏｎａｅｒｏ￣ｄｙｎａｍｉｃｄａｍｐｉｎｇｗｉｎｄｔｕｎｎｅｌｏｆｌａｕｎｃｈｖｅｈｉｃｌｅ[Ｃ].Ｐｒｏ￣ｃｅｅｄｉｎｇｓｏｆｔｈｅ１２ｔｈＮａｔｉｏｎａｌＡｅｒｏｅｌａｓｔｉｃｉｔｙＣｏｎｆｅｒｅｎｃｅꎬ２０１１(ｉｎＣｈｉｎｅｓｅ).[２１]刘子强ꎬ白葵ꎬ毛国良ꎬ等.锤头体弹性振动跨音速气动阻尼系数的确定[Ｊ].宇航学报ꎬ２００２ꎬ２３(６):１￣７.ＬｉｕＺＱꎬＢａｉＫꎬＭａｏＧＬꎬｅｔａｌ.Ｔｈｅｄｅｔｅｒｍｉｎａｔｉｏｎｏｆａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｎｈａｍｍｅｒｈｅａｄｌａｕｎｃｈｖｅｈｉｃｌｅｓａｔｔｒａｎｓｏｎｉｃｓｐｅｅｄｓ[Ｊ].ＪｏｕｒｎａｌｏｆＡｓｔｒｏｎａｕｔｉｃｓꎬ２００２ꎬ２３(６):１￣７(ｉｎＣｈｉｎｅｓｅ).[２２]冉景洪ꎬ刘子强ꎬ胡静ꎬ等.减阻杆气动阻尼研究[Ｊ].力学学报ꎬ２０１４ꎬ４６(４):６３６￣６４１.ＲａｎＪＨꎬＬｉｕＺＱꎬＨｕＪꎬｅｔａｌ.Ｒｅｓｅａｒｃｈｏｆａｅｒｏ￣ｄａｍｐｉｎｇｆｏｒｂｌｕｎｔｗｉｔｈｓｐｉｋｅ[Ｊ].ＣｈｉｎｅｓｅＪｏｕｒｎａｌｏｆＴｈｅ￣ｏｒｅｔｉｃａｌａｎｄＡｐｐｌｉｅｄＭｅｃｈａｎｉｃｓꎬ２０１４ꎬ４６(４):６３６￣６４１(ｉｎＣｈｉｎｅｓｅ).[２３]朱剑ꎬ冉景洪ꎬ吴彦森ꎬ等.捆绑式运载火箭的气动阻尼数值计算方法[Ｃ].第十三届全国空气弹性学术交流会论文集.哈尔滨:中国力学学会ꎬ中国空气动力学会ꎬ２０１３.ＺｈｕＪꎬＲａｎＪＨꎬＷｕＹＳꎬｅｔａｌ.Ｎｕｍｅｒｉｃａｌｃａｌｃｕｌａｔｉｏｎｍｅｔｈｏｄｆｏｒａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｏｆｂｕｎｄｌｅｌａｕｎｃｈｖｅｈｉｃｌｅｓ[Ｃ].Ｐｒｏｃｅｅｄｉｎｇｓｏｆｔｈｅ１３ｔｈＮａｔｉｏｎａｌＡｅｒｏｅｌａｓｔｉｃｉｔｙＣｏｎ￣ｆｅｒｅｎｃｅ.Ｈａｒｂｉｎ:ＣｈｉｎｅｓｅＳｏｃｉｅｔｙｏｆＴｈｅｏｒｅｔｉｃａｌａｎｄＡｐ￣ｐｌｉｅｄＭｅｃｈａｎｉｃｓꎬＣｈｉｎａＡｅｒｏｍｅｃｈａｎｉｃｓＳｏｃｉｅｔｙꎬ２０１３(ｉｎＣｈｉｎｅｓｅ).[２４]ＳｐａｌａｒｔＰＲ.Ｄｅｔａｃｈｅｄ￣ｅｄｄｙｓｉｍｕｌａｔｉｏｎ[Ｊ].ＡｎｎｕａｌＲｅｖｉｅｗｏｆＦｌｕｉｄＭｅｃｈａｎｉｃｓꎬ２００９ꎬ４１:１８１￣２０２. [２５]ＧｒｉｔｓｋｅｖｉｃｈＭＳꎬＧａｒｂａｒｕｋＡＶꎬＳｃｈüｔｚｅＪꎬｅｔａｌ.Ｄｅｖｅｌ￣ｏｐｍｅｎｔｏｆＤＤＥＳａｎｄＩＤＤＥＳｆｏｒｍｕｌａｔｉｏｎｓｆｏｒｔｈｅｋ￣ωｓｈｅａｒｓｔｒｅｓｓｔｒａｎｓｐｏｒｔｍｏｄｅｌ[Ｊ].ＦｌｏｗꎬＴｕｒｂｕｌｅｎｃｅａｎｄＣｏｍｂｕｓｔｉｏｎꎬ２０１２ꎬ８８(３):４３１￣４４９.[２６]季辰ꎬ吴彦森ꎬ侯英昱ꎬ等.捆绑式运载火箭跨声速气动阻尼特性试验研究[Ｊ].实验流体力学ꎬ２０２０ꎬ３４(６):２４￣３１.ＪｉＣꎬＷｕＹＳꎬＨｏｕＹＹꎬｅｔａｌ.Ｅｘｐｅｒｉｍｅｎｔａｌｓｔｕｄｙｏｆａｅｒｏｄｙｎａｍｉｃｄａｍｐｉｎｇｃｈａｒａｃｔｅｒｉｓｔｉｃｓｏｆａｌａｕｎｃｈｖｅｈｉｃｌｅｗｉｔｈｂｏｏｓｔｅｒｓｉｎｔｒａｎｓｏｎｉｃｆｌｏｗ[Ｊ].ＪｏｕｒｎａｌｏｆＥｘｐｅｒｉｍｅｎｔｓｉｎＦｌｕｉｄＭｅｃｈａｎｉｃｓꎬ２０２０ꎬ３４(６):２４￣３１(ｉｎＣｈｉｎｅｓｅ). [２７]ＡｌｌｅｎＣꎬＲｅｎｄａｌｌＴＣＳ.ＵｎｉｆｉｅｄａｐｐｒｏａｃｈｔｏＣＦＤ￣ＣＳＤｉｎｔｅｒｐｏｌａｔｉｏｎａｎｄｍｅｓｈｍｏｔｉｏｎｕｓｉｎｇｒａｄｉａｌｂａｓｉｓｆｕｎｃｔｉｏｎｓ[Ｒ].ＡＩＡＡ２００７￣３８０４ꎬ２００７.０４。

列车风洞试验综述1列车风洞模型试验系统1.1风洞的基本类型及基本原理当对列车的空气动力学特性进行试验研究时，直接而真实的方法是在线实车试验，但进行一次试验需要耗费大量的人力、物力、财力，组织一次试验很不容易，得到的数据有限，加之自然条件千变万化，如环境的风速和风向不可控制等，重复性难以保证，而且，实车试验需在列车制造出来后才能进行，用于研制新车代价太高，因此实车试验一般以验证、评估、考核试验为主，兼顾研究性试验。

于是，人们就想用模型试验来代替实车试验。

风洞是能人工产生和控制气流，以模拟飞行器或物体周围气体的流动，并可量度气流对物体的作用以及观察物理现象的一种管道状实验设备，它是进行空气动力实验最常用、最有效的工具。

风洞模型试验是研究列车气动特性中应用最广泛的手段之一。

它具有试验理论和试验手段成熟、测量精密，气流参数如速度、压力等易于控制，并且基本不受天气变化的影响等优点。

为了满足不同类型空气动力试验的要求，现代风洞的种类繁多。

风洞通常按照试验段气流的马赫数来分类，有低速风洞（Ma<0.3）、亚音速风洞（0.3<Ma<0.8）、跨音速风洞（0.8<Ma<1.5）、超音速风洞（1.5<Ma<4.5）、高超音速风洞（4.5<Ma<10）、极高速风洞（Ma>10）等。

列车模型风洞试验一般在低速风洞中进行。

低速风洞按通过试验段气流循环形式来分，有直流式和回流式两种基本类型。

按试验段结构不同，低速风洞又有“开口”和“闭口”之别。

直流式风洞的特点是把通过试验段的气流排在风洞外部，如图1。

回流式风洞的特点是通过试验段的气流经循环系统再返回试验段，如图2。

图1 直流式风洞图2回流式风洞对列车在空气中的等速直线运动，按照运动的相对性原理，在空气动力特性研究中，可以认为列车静止不动，与列车速度大小相同方向相反的空气流过列车，列车上承受的空气动力与类车运动在静止的空气中承受的空气动力完全相同。

大展弦比机翼跨声速静气动弹性风洞试验郭洪涛;陈德华;吕彬彬;余立;祖孝勇【摘要】基于静气动弹性风洞试验研究了某翼身组合体的跨声速静气动弹性效应.试验结果表明:在设计巡航点,静气动弹性对大展弦比超临界机翼的气动特性影响明显,可使机翼的升力系数降低21%、升阻比增加8%、焦点前移约1%bA;在超过巡航马赫数后,静气动弹性效应使得机翼气动特性有恶化的趋势.跨声速时,马赫数和速压对机翼的静气动弹性效应具有较大影响,且影响规律呈复杂非线性,难以依据现有理论分析准确预计.%Static aeroelasticity has great influence on the aerodynamic performance and safety of aircrafts. Based on a wing-body model for aerodynamic/structural similarity, a static aeroelastic wind tunnel test was carried out to research the transonic static aeroelastic effects. The results show that the static aeroelasticity has significant influence on the aerodynamic performance for the hypercritical high-aspect-ratio sweepback wing at its designed cruise point. The static aeroelastic effects lead to 21% decrement in the lift coefficient,8% increment in the lift-drag ratio,and the forward movement of the aerodynamic centre about 1%b A .With the Mach number higher than the designed cruise value,the static aeroelastic effects worsen the wing′s aerodynamic performance. In addition, in transonic range, the Mach number and dynamic pressure have significant influence on the static aeroelastic effects,and the influence has nonlinear characteristics.It is difficult to precisely predict the influence law and magnitude by current theory.【期刊名称】《空气动力学学报》【年(卷),期】2017(035)006【总页数】5页(P841-845)【关键词】静气动弹性;大展弦比机翼;跨声速流动;风洞试验;气动特性【作者】郭洪涛;陈德华;吕彬彬;余立;祖孝勇【作者单位】中国空气动力研究与发展中心空气动力学国家重点实验室,四川绵阳621000;中国空气动力研究与发展中心高速空气动力研究所,四川绵阳 621000;中国空气动力研究与发展中心空气动力学国家重点实验室,四川绵阳 621000;中国空气动力研究与发展中心高速空气动力研究所,四川绵阳 621000;中国空气动力研究与发展中心高速空气动力研究所,四川绵阳 621000;中国空气动力研究与发展中心高速空气动力研究所,四川绵阳 621000【正文语种】中文【中图分类】V224;V211.753静气动弹性表征了气动力和结构弹性力之间的耦合作用，对飞行器的气动性能与飞行安全有着重要影响。

第９卷㊀第１期２０２４年１月气体物理ＰＨＹＳＩＣＳＯＦＧＡＳＥＳＶｏｌ.９㊀Ｎｏ.１Ｊａｎ.２０２４㊀㊀ＤＯＩ:１０.１９５２７/ｊ.ｃｎｋｉ.２０９６￣１６４２.１０７１双喉道Ｌｕｄｗｉｅｇ管风洞启动过程及其有效运行时间延长李创创１ꎬ㊀李志远１ꎬ㊀张振辉２ꎬ㊀吴㊀杰１(１.华中科技大学航空航天学院ꎬ湖北武汉４３００７４ꎻ２.中国空气动力研究与发展中心ꎬ四川绵阳６２１０００)ＳｔａｒｔｉｎｇＰｒｏｃｅｓｓｏｆａＤｏｕｂｌｅ￣ＴｈｒｏａｔＬｕｄｗｉｅｇＴｕｂｅＴｕｎｎｅｌａｎｄｔｈｅＥｘｔｅｎｓｉｏｎｏｆＩｔｓＥｆｆｅｃｔｉｖｅＲｕｎｎｉｎｇＴｉｍｅＬＩＣｈｕａｎｇｃｈｕａｎｇ１ꎬ㊀ＬＩＺｈｉｙｕａｎ１ꎬ㊀ＺＨＡＮＧＺｈｅｎｈｕｉ２ꎬ㊀ＷｕＪｉｅ１(１.ＳｃｈｏｏｌｏｆＡｅｒｏｓｐａｃｅＥｎｇｉｎｅｅｒｉｎｇꎬＨｕａｚｈｏｎｇＵｎｉｖｅｒｓｉｔｙｏｆＳｃｉｅｎｃｅａｎｄＴｅｃｈｎｏｌｏｇｙꎬＷｕｈａｎ４３００７４ꎬＣｈｉｎａꎻ２.ＣｈｉｎａＡｅｒｏｄｙｎａｍｉｃｓＲｅｓｅａｒｃｈａｎｄＤｅｖｅｌｏｐｍｅｎｔＣｅｎｔｅｒꎬＭｉａｎｙａｎｇ６２１０００ꎬＣｈｉｎａ)摘㊀要:Ｌｕｄｗｉｅｇ管风洞是开展高超声速空气动力学实验基础研究的重要平台ꎮ但是ꎬ快开阀启动式高超声速Ｌｕｄｗｉｅｇ管风洞长期受快开阀影响ꎬ产生不同类型的来流扰动模态ꎮ双喉道气动布局可有效消除快开阀启动式高超声速Ｌｕｄｗｉｅｇ管风洞上游部件的扰动来源ꎬ但是会导致风洞有效运行时间大幅缩短ꎮ针对该问题ꎬ通过非定常数值模拟对双喉道气动布局高超声速Ｌｕｄｗｉｅｇ管风洞的启动特性进行研究ꎬ然后对第１喷管扩张段与稳定段进行了融合设计ꎬ研究了不同扩张角与稳定段组合对风洞启动时间以及流场品质的影响ꎮ结果表明ꎬ采用减小扩张角组合设计能够使双喉道气动布局高超声速Ｌｕｄｗｉｅｇ管风洞的有效运行时间提升近２０％ꎬ并且对下游实验段内的静态流场品质几乎无影响ꎬ有效提高了风洞的实验能力ꎮ同时ꎬ相较于较大的扩张角组合ꎬ较小的扩张角设计能够减少约１０％的总压损失ꎮ关键词:高超声速风洞ꎻＬｕｄｗｉｅｇ管风洞ꎻ双喉道布局ꎻ风洞启动特性ꎻ有效运行时间㊀㊀㊀收稿日期:２０２３￣０７￣０３ꎻ修回日期:２０２３￣０８￣１７第一作者简介:李创创(２０００ )㊀男ꎬ硕士ꎬ主要研究方向为高超声速空气动力学ꎮＥ￣ｍａｉｌ:２９３０３６９４００＠ｑｑ.ｃｏｍ通信作者简介:吴杰(１９８６ )㊀男ꎬ教授ꎬ主要研究方向为高速飞行器空气动力学㊁智能流动控制㊁高速风洞设计ꎮＥ￣ｍａｉｌ:中图分类号:Ｖ２１１.７５１㊀㊀文献标志码:ＡＡｂｓｔｒａｃｔ:ＴｈｅＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌｓｅｒｖｅｓａｓａｃｒｕｃｉａｌｔｅｓｔｉｎｇｇｒｏｕｎｄｆｏｒｆｕｎｄａｍｅｎｔａｌｓｔｕｄｉｅｓｏｆｈｙｐｅｒｓｏｎｉｃａｅｒｏｄｙ￣ｎａｍｉｃｓ.Ｈｏｗｅｖｅｒꎬｔｈｅｑｕｉｃｋ￣ｏｐｅｎｉｎｇｖａｌｖｅ￣ｓｔａｒｔｅｄｈｙｐｅｒｓｏｎｉｃＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌｈａｓｂｅｅｎａｆｆｅｃｔｅｄｂｙｔｈｅｑｕｉｃｋ￣ｏｐｅｎｉｎｇｖａｌｖｅｆｏｒａｌｏｎｇｔｉｍｅꎬｒｅｓｕｌｔｉｎｇｉｎｄｉｆｆｅｒｅｎｔｔｙｐｅｓｏｆｆｌｏｗｄｉｓｔｕｒｂａｎｃｅｍｏｄｅｓ.Ｔｈｅｄｏｕｂｌｅ￣ｔｈｒｏａｔａｅｒｏｄｙｎａｍｉｃｃｏｎｆｉｇｕ￣ｒａｔｉｏｎｃａｎｅｆｆｅｃｔｉｖｅｌｙｅｌｉｍｉｎａｔｅｔｈｅｄｉｓｔｕｒｂａｎｃｅｓｏｕｒｃｅｏｆｔｈｅｕｐｓｔｒｅａｍｃｏｍｐｏｎｅｎｔｓｏｆｔｈｅｑｕｉｃｋ￣ｏｐｅｎｉｎｇｖａｌｖｅ￣ｓｔａｒｔｅｄｈｙｐｅｒ￣ｓｏｎｉｃＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌꎬｂｕｔｉｔｗｉｌｌｌｅａｄｔｏａｓｉｇｎｉｆｉｃａｎｔｒｅｄｕｃｔｉｏｎｉｎｔｈｅｅｆｆｅｃｔｉｖｅｒｕｎｎｉｎｇｔｉｍｅｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌ.Ｉｎｏｒｄｅｒｔｏｓｏｌｖｅｔｈｉｓｐｒｏｂｌｅｍꎬｔｈｅｕｎｓｔｅａｄｙｎｕｍｅｒｉｃａｌｓｉｍｕｌａｔｉｏｎｗａｓｕｓｅｄｔｏｓｔｕｄｙｔｈｅｓｔａｒｔ￣ｕｐｃｈａｒａｃｔｅｒｉｓｔｉｃｓｏｆｔｈｅｈｙｐｅｒ￣ｓｏｎｉｃＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌｗｉｔｈｄｏｕｂｌｅ￣ｔｈｒｏａｔａｅｒｏｄｙｎａｍｉｃｃｏｎｆｉｇｕｒａｔｉｏｎ.Ｔｈｅｎꎬａｆｕｓｉｏｎｄｅｓｉｇｎｗａｓｃａｒｒｉｅｄｏｕｔｆｏｒｔｈｅｆｉｒｓｔｎｏｚｚｌｅｅｘｐａｎｓｉｏｎｓｅｃｔｉｏｎａｎｄｔｈｅｓｔａｂｌｅｓｅｃｔｉｏｎ.Ｔｈｅｅｆｆｅｃｔｓｏｆｄｉｆｆｅｒｅｎｔｃｏｍｂｉｎａｔｉｏｎｓｏｆｅｘｐａｎｓｉｏｎａｎｇｌｅｓａｎｄｓｔａｂｌｅｓｅｃ￣ｔｉｏｎｓｏｎｔｈｅｓｔａｒｔ￣ｕｐｔｉｍｅａｎｄｆｌｏｗｆｉｅｌｄｑｕａｌｉｔｙｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌｗｅｒｅｓｔｕｄｉｅｄ.ＴｈｅｒｅｓｕｌｔｓｓｈｏｗｔｈａｔｔｈｅｅｆｆｅｃｔｉｖｅｒｕｎｎｉｎｇｔｉｍｅｏｆｔｈｅｈｙｐｅｒｓｏｎｉｃＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌｗｉｔｈｄｏｕｂｌｅ￣ｔｈｒｏａｔａｅｒｏｄｙｎａｍｉｃｃｏｎｆｉｇｕｒａｔｉｏｎｃａｎｂｅｉｎｃｒｅａｓｅｄｂｙｎｅａｒｌｙ２０％ｂｙｕｓｉｎｇｔｈｅｃｏｍｂｉｎｅｄｄｅｓｉｇｎｏｆｒｅｄｕｃｉｎｇｔｈｅｅｘｐａｎｓｉｏｎａｎｇｌｅꎬａｎｄｔｈｅｓｔａｔｉｃｆｌｏｗｆｉｅｌｄｑｕａｌｉｔｙｉｎｔｈｅｄｏｗｎｓｔｒｅａｍｅｘ￣ｐｅｒｉｍｅｎｔａｌｓｅｃｔｉｏｎｉｓａｌｍｏｓｔｕｎａｆｆｅｃｔｅｄꎬｗｈｉｃｈｅｆｆｅｃｔｉｖｅｌｙｉｍｐｒｏｖｅｓｔｈｅｅｘｐｅｒｉｍｅｎｔａｌａｂｉｌｉｔｙｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌ.Ａｔｔｈｅｓａｍｅｔｉｍｅꎬｃｏｍｐａｒｅｄｗｉｔｈｔｈｅｌａｒｇｅｒｅｘｐａｎｓｉｏｎａｎｇｌｅｃｏｍｂｉｎａｔｉｏｎꎬｔｈｅｓｍａｌｌｅｒｅｘｐａｎｓｉｏｎａｎｇｌｅｄｅｓｉｇｎｃａｎｒｅｄｕｃｅｔｈｅｔｏｔａｌｐｒｅｓ￣ｓｕｒｅｌｏｓｓｂｙａｂｏｕｔ１０％.第１期李创创ꎬ等:双喉道Ｌｕｄｗｉｅｇ管风洞启动过程及其有效运行时间延长Ｋｅｙｗｏｒｄｓ:ｈｙｐｅｒｓｏｎｉｃｗｉｎｄｔｕｎｎｅｌꎻＬｕｄｗｉｅｇｔｕｂｅｔｕｎｎｅｌꎻｄｏｕｂｌｅ￣ｔｈｒｏａｔｃｏｎｆｉｇｕｒａｔｉｏｎꎻｓｔａｒｔ￣ｕｐｃｈａｒａｃｔｅｒｉｓｔｉｃｓｏｆｗｉｎｄｔｕｎ￣ｎｅｌꎻｅｆｆｅｃｔｉｖｅｒｕｎｎｉｎｇｔｉｍｅ引㊀言高超声速飞行器由于其高速度㊁强机动㊁超远程等特点ꎬ是当今世界航空航天强国的研究热点ꎬ具有重要的军事战略意义及广阔的民用应用前景[１]ꎮ高超声速飞行器趋向于临近空间发展ꎬ在大气层内飞行时间长ꎬ所涉及的飞行环境㊁速度以及气体流动特性十分复杂ꎬ给气动研究带来了巨大挑战ꎮ通过近几十年的研究ꎬ研究人员在高超声速流动领域积累了大量的理论基础与经验ꎬ但是高超声速空气动力学仍存在诸多难点问题ꎮ虽然计算空气动力学已成为飞行器精确㊁高效设计的预测手段ꎬ但在高超声速流动领域由于流动机理十分复杂ꎬ如高超声速湍流建模㊁高温下多物理场耦合与非化学平衡状态等ꎬ数值方法不能完全刻画其流场演化[２￣５]ꎮ对于一些特殊的流动现象ꎬ如层流㊁湍流边界层转捩以及激波￣边界层相互作用机理[６ꎬ７]ꎬ目前还无普适的数值模型能够进行预测ꎬ因此需要通过地面气动试验与真实飞行试验相互配合来研究ꎬ以降低飞行器设计中的不确定性ꎮ虽然飞行试验能够获得丰富且更真实的试验数据ꎬ但飞行试验成本昂贵且测试手段也受到较大约束ꎮ同时常规超声速及高超声速风洞设备建设成本昂贵ꎬ风洞的运行和操作也较为复杂ꎬ提高了科研机构进行基础研究的门槛[８ꎬ９]ꎮ随着我国国防事业的发展ꎬ飞行器设计的战术技术指标要求不断提高ꎬ意味着对风洞试验技术也有了更高的要求ꎬ风洞试验趋向于高精度㊁低成本和精细化等方向发展[１０]ꎮ在高超声速风洞中ꎬＬｕｄｗｉｅｇ管风洞由于其基于中心膨胀波的运行模式ꎬ具有来流Ｒｅｙｎｏｌｄｓ数高㊁流场重复性好㊁流场稳定等优点[１１]ꎮ最初ꎬＬｕｄｗｉｅｇ[１２]将该设施设计为在高Ｒｅｙｎｏｌｄｓ数下进行亚声速/跨声速试验的低成本替代方案ꎮ由于Ｌｕｄ￣ｗｉｅｇ管风洞具有能够产生低湍流均匀自由来流的优点ꎬ后来被用于高超声速试验[１３ꎬ１４]ꎮ典型Ｌｕｄｗｉｅｇ管风洞由长直储气段㊁Ｌａｖａｌ喷管㊁试验段和真空罐组成ꎬ通过这种设施可以产生相对较长的运行时间(０.１~０.２ｓ)㊁较大的试验截面和较高的Ｒｅｙｎｏｌｄｓ数(５ˑ１０６~５０ˑ１０６ｍ－１)的高超声速流场ꎮ为了进一步提高Ｌｕｄｗｉｅｇ管风洞的运行效率ꎬ德国宇航院Ｋｏｐｐｅｎｗａｌｌｎｅｒ等[１５]研制了适用于Ｌｕｄｗｉｅｇ管风洞的快开阀控制系统ꎬ大幅度提高了该类风洞的试验效率ꎬ并降低了风洞的运行成本ꎮ虽然Ｌｕｄ￣ｗｉｅｇ管风洞能够以较低的建设和运行成本实现高质量的高超声速流场ꎬ但其存在以下不足:首先ꎬＬｕｄｗｉｅｇ管风洞的设计Ｍａｃｈ数单一ꎬ向低宽Ｍａｃｈ数拓展存在困难[１６]ꎮ其次ꎬ由于气体重力作用ꎬ储气段内气体加热不均匀ꎬ从而产生较大的熵波扰动[１７ꎬ１８]ꎮ再者ꎬ快开阀对Ｌａｖａｌ喷管喉部上游流场的干扰会增加试验段自由流扰动幅值[１９]ꎻ同时ꎬ阀门形成的涡脱落也会在风洞流场中引入一定的涡波扰动[２０]ꎮ上述扰动源的引入会显著降低风洞的流场品质ꎮ为了将高Ｍａｃｈ数Ｌｕｄｗｉｅｇ管风洞向中低Ｍａｃｈ数扩展ꎬＷｕ等[２１]基于常规高超声速Ｌｕｄｗｉｅｇ管采用双喉道气动布局ꎬ即额外增加１个Ｌａｖａｌ喷管和用于连接两个喷管的稳定段ꎻ第１喷管和稳定段联合起到节流作用ꎬ可以低成本地实现将原有的Ｌｕｄｗｉｅｇ管风洞向多个运行Ｍａｃｈ数扩展ꎮ此外ꎬ采用双喉道布局Ｌｕｄｗｉｅｇ管风洞自由流扰动的特征表明ꎬ稳定段可以相当大程度地消除来自快开阀及其他上游组件产生的干扰ꎮＳｃｈｒｉｊｅｒ等[２２]也采用了相同设计ꎬ通过串联喷管的配置实现了风洞由高Ｍａｃｈ数流动到低Ｍａｃｈ数的转变ꎮ国内高亮杰等[２３ꎬ２４]详细分析了该气动布局风洞喷管的工作模态ꎬ并采用非定常数值模拟技术研究了风洞的启动特性ꎬ实现不同Ｍａｃｈ数下的尺寸匹配ꎮ然而ꎬ按照Ｗｕ等[１７]的设计思路ꎬＬｕｄｗｉｅｇ管风洞只能向低Ｍａｃｈ数工况拓展ꎬ传统的高超声速Ｌｕｄｗｉｅｇ管风洞要消除来自Ｌａｖａｌ喷管上游部件的干扰仍存在较大难度ꎮ同时ꎬＷｕ等[１７]的设计中第１喷管的喉道面积比第２喷管要小得多ꎬ导致两个喉道之间产生较大的总压损失ꎬ会影响风洞的气动性能ꎮ针对该问题ꎬＬｉ等[１６]设计了一种双喉道气动布局的Ｌｕｄｗｉｅｇ管风洞ꎬ使风洞的第２喉道截面积小于第１喉道ꎬ数值模拟与实验结果均验证了该方案可行ꎬ大幅提升了高超声速来流的流场品质ꎬ且总压损失较小ꎮ但是ꎬ由于第１Ｌａｖａｌ喷管和稳定段的引入ꎬ使该风洞的启动时间较长ꎬ因此须对其启动过程进行进一步优化分析ꎬ以期缩短采用该类气动布局风洞的启动时间ꎮ针对该问题ꎬ对该双喉道布局的Ｌｕｄｗｉｅｇ管风洞启动特性进行深入研究ꎬ并重点对缩短风洞启动过程进行优化ꎮ文章首９５气体物理２０２４年㊀第９卷先介绍采用该双喉道气动布局的Ｌｕｄｗｉｅｇ管风洞的工作原理ꎻ之后ꎬ对风洞的启动过程进行非定常数值模拟与分析ꎻ而后ꎬ将对第１段Ｌａｖａｌ喷管与稳定段进行融合设计ꎬ并分析其对风洞的启动过程和流场品质的影响ꎻ最后ꎬ根据数值模拟结果ꎬ确定进一步缩短采用该双喉道气动布局的Ｌｕｄｗｉｅｇ管风洞运行时间的有效设计ꎮ１㊀双喉道Ｌｕｄｗｉｅｇ管风洞设计１.１㊀双喉道Ｌｕｄｗｉｅｇ管风洞运行原理该双喉道Ｌｕｄｗｉｅｇ管风洞结构与Ｗｕ等[１７]以及Ｓｃｈｒｉｊｅｒ[２５]描述的相似ꎬ但由于第１喷管与第２喷管的喉道面积比大于１ꎬ导致其工作原理完全不同ꎮ采用第１喉道较大的布局方式可以消除风洞稳定运行时第１喉道处的激波结构和由双喉道面积比引入的压力损失ꎬ其气动布局如图１所示ꎮ快开阀开启瞬间ꎬ储气段内的高温高压气体首先通过第１喷管ꎬ并在第１喷管喉道处发生壅塞ꎮ随后气体通过稳定段ꎬ在第２喷管喉道处再次发生壅塞ꎬ而后流动沿第２喷管扩张至设计的超声速流动ꎮ此时第１和第２喷管喉道处的流动均为声速流动ꎬ膨胀波(图１中的ＯＡ和ＯＢ)向上游传播到储气段中ꎮ在该过程中ꎬ第１喷管和第２喷管中的流动符合一维等熵流动控制方程ꎮ稳定段内气体的总压和总温用Ｐ１和Ｔ表示ꎬ储气段中气体总压和总温用Ｐ０和Ｔ０表示ꎮ流场建立初期ꎬ流经第１和第２喉道的气体质量流量均可由下面的最大质量流量公式表示ꎬ其中Ａ∗对应喷管喉道截面积ꎮｍ＝Ｐ０ꎬ１Ｔ０ꎬ１Ａ∗γＲ２γ＋１æèçöø÷γ＋１２γ－２在风洞的整个流场中ꎬ忽略热传导ꎬ可认为是绝热过程(Ｔ０＝Ｔ１)ꎮ由于第１喉道面积大于第２喉道面积且Ｐ０大于Ｐ１ꎬ所以风洞启动过程中流经第１喉道的质量流量大于第２喉道ꎮ稳定段内的气体质量增加使压力Ｐ１上升ꎬ进而使通过第２喉道的质量流量增大ꎬ同时稳定段内的流动逐渐转变为亚声速流动ꎬ随后流经第１喉道的质量流量开始减小ꎮ一段时间后通过两喉道的质量流量会达到平衡ꎬ此时风洞处于稳定运行状态ꎬ声速点位于第２喉道ꎬ第２喉道为整个流道的几何喉道ꎮ当膨胀波系经储气段末端反射回来到达第１喉道时ꎬ风洞１个运行周期结束ꎬ快开阀随后关闭ꎮ该风洞布局中ꎬ稳定段内的主要扰动来源是阀门的涡脱落和第１喷管与稳定段交接处附近的流动分离ꎬ都属于涡波扰动ꎬ扰动形式较为单一ꎮ该气动布局总压损失较小ꎬ并且在稳定段中容易设计整流措施ꎬ因此通过该风洞实验平台可以获得更高质量的自由来流ꎮ但这种气动布局会使风洞的启动耗时较长ꎬ须对其启动过程进行研究ꎮ图１㊀新型双喉道Ｌｕｄｗｉｅｇ管风洞运行原理图Ｆｉｇ.１㊀Ｏｐｅｒａｔｉｏｎｐｒｉｎｃｉｐｌｅｏｆｎｅｗｄｏｕｂｌｅ￣ｔｈｒｏａｔＬｕｄｗｉｅｇｔｕｂｅｗｉｎｄｔｕｎｎｅｌ[１６]１.２㊀第１Ｌａｖａｌ喷管与稳定段设计风洞运行稳定时ꎬ第１喷管与稳定段内的流动均为亚声速流动ꎬ所以无须使用特征线法进行设计ꎮ同时在第１喷管喉道出口处预期会存在微弱的射流效应ꎬ在靠近壁面区域可能存在流动分离ꎬ因此第１喷管型线无须进行精细化设计ꎮ其中ꎬ第１喷管收缩段采用Ｗｉｔｏｓｚｙｎｓｋｉ方法ꎬ其型线表达式如下ｙ＝ｙ∗１－１－ｙ０ｙｉæèçöø÷２éëêêùûúú１－ｘ２Ｌ２ｉæèçöø÷２１＋ｘ２３Ｌ２ｉæèçöø÷３０６第１期李创创ꎬ等:双喉道Ｌｕｄｗｉｅｇ管风洞启动过程及其有效运行时间延长式中ꎬｙ为横坐标为ｘ处的截面半径ꎬｙ０为喉部半径ꎬｙｉ为收缩段入口半径ꎬＬｉ为收缩段长度ꎮ为了方便实现第１喷管与稳定段的融合设计ꎬ喷管扩张段出口采用５次曲线并引出与之相切的直线与稳定段相连ꎬ保证了第１喉道出口处的光滑过渡ꎬ避免在该位置处产生激波干扰ꎮ第１喷管与稳定段的整体示意图如图２所示ꎬ扩张段５次曲线公式如下ｙ－ｙ０ｙｏｕｔ－ｙ０＝１＋１０ｘＬｅｘｐæèçöø÷３＋１５ｘＬｅｘｐæèçöø÷４＋６ｘＬｅｘｐæèçöø÷５式中ꎬｙｏｕｔ为扩张段出口半径ꎬＬｅｘｐ为扩张段长度ꎮ图２㊀第１Ｌａｖａｌ喷管示意图Ｆｉｇ.２㊀ＳｋｅｔｃｈｏｆｔｈｅｆｉｒｓｔＬａｖａｌｎｏｚｚｌｅ第１喷管的收缩段和扩张段的出口半径与稳定段保持一致ꎬ稳定段和扩张段总长保持不变ꎬ因此喷管扩张角θ的大小决定了稳定段和扩张段的长度ꎮ稳定段的设计依赖于其内的流动状态ꎬ为了使经过第１喷管加速后的流场便于通过阻尼材料进行整流来获得高质量流场ꎬ稳定段内的流速通常被限制在Ｍａ＝０.１ꎮ根据Ｍａｃｈ￣面积关系式ꎬ可由第１喉道面积获得稳定段半径ꎬ而第１喉道面积受限于第１喉道与第２喉道面积比ꎮ本次计算的风洞气动外形尺寸参考华中科技大学Φ０.２５ｍ口径Ｍａｃｈ数６Ｌｕｄｗｉｅｇ管风洞[２６]ꎬ在本次计算中第１喉道与第２喉道截面面积比为２.２６ꎬ稳定段出口到第１喉道的距离为７２０ｍｍꎬ第１喷管扩张角为１７ʎꎬ第２喷管长度为１４００ｍｍꎮ１.３㊀数值设置使用ＳＵ２代码对该风洞的启动过程进行非定常Ｒｅｙｎｏｌｄｓ平均Ｎａｖｉｅｒ￣Ｓｔｏｋｅｓ(ＵＲＡＮＳ)模拟ꎬＳＵ２代码能够预测高超声速流经复杂几何形状周围的黏性流动[２７ꎬ２８]ꎮ在本计算中ꎬ采用的湍流模式为标准ＭｅｎｔｅｒＳＳＴ两方程湍流模型[２９ꎬ３０]ꎬ对流项格式采用２阶迎风格式[３１ꎬ３２]ꎮ出于稳定性考虑ꎬ采用双时间格式ꎬ并将物理时间步长限制为Δｔ＝１０－６ｓꎮ每个物理时间步长的子迭代次数为４０ꎬ确保残差能够减少至少４个数量级ꎮ为了得到精确的结果ꎬ流场采用结构化网格ꎬｙ＋均取１左右ꎬ部分网格与边界条件如图３所示ꎮ储气段总长为２２ｍꎬ图中未完全画出ꎮ第２Ｌａｖａｌ喷管沿流向和垂向的网格数目为４８０ˑ６０ꎬ其余部分的网格生成采用了类似的网格密度ꎬ并对喷管壁面处的网格进行加密处理ꎮ流场的边界条件定义如下:风洞出口定义为压力出口ꎬ所有的壁面指定为绝热无滑移壁面边界条件ꎬ风洞轴线采用轴对称边界条件ꎮ图３㊀网格划分与边界条件设置Ｆｉｇ.３㊀Ｍｅｓｈｔｏｐｏｌｏｇｙａｎｄｂｏｕｎｄａｒｙｃｏｎｄｉｔｉｏｎｓｏｆｆｌｏｗｆｉｅｌｄ风洞启动时ꎬ快开阀可在５ｍｓ内完成启闭ꎬ快开阀开启过程对风洞启动过程的影响几乎可以忽略不计ꎮ因此在快开阀与Ｌａｖａｌ喷管喉部的接触面上布置了１个数值膜片ꎬ以初始化流场并进行时间瞬态模拟ꎮ储气段的初始条件为总压１ＭＰａ和总温４３４Ｋꎬ出口设置为压力出口条件ꎬ出口压力为１００Ｐａꎬ单位Ｒｅｙｎｏｌｄｓ数为Ｒｅ/ｌ＝１.０６ˑ１０７ｍ－１ꎮ赵家权等[３３]通过将ＳＵ２计算的结果与其他求解器进行比较ꎬ检验了其求解高超声速流场的能力ꎬ结果表明ＳＵ２和其他求解器吻合良好ꎬ可用于高超声速管道流动模拟ꎮ２㊀数值结果与分析２.１㊀双喉道Ｌｕｄｗｉｅｇ管风洞的启动过程分析目前对该气动布局风洞启动过程研究较少且不够深入ꎬ因此着重对风洞的启动过程进行分析ꎮ为了对风洞的启动过程有进一步的理解ꎬ并确定该气动布局的可用性ꎬ以及更好地对第１喷管与稳定段进行融合设计ꎬ对该风洞的启动过程进行了非定常数值模拟ꎮ图４展示了阀门开启后１~３０ｍｓ时刻从快开阀到试验段的Ｍａｃｈ数分布情况ꎮ阀门开启后ꎬ流动在初始高压比作用下首先在第１喷管喉道处达１６气体物理２０２４年㊀第９卷到声速ꎬ然后超声速流动迅速占据整个稳定段并在第２喉道入口处形成壅塞流动ꎮ此时ꎬ第１和第２喉道处的流动均为声速流动ꎬ流动经过第２喷管膨胀加速至高超声速ꎬ并在试验段内形成不稳定的高超声速流动ꎮ在ｔ＝２ｍｓ时ꎬ为了匹配两喉道截面处的流量ꎬ稳定段内的压力开始增加ꎬ第２喉道入口处的流动开始由超声速流动转变成亚声速流动并在收缩段内形成一个亚声速流动区ꎻ该亚声速流动区逐渐沿风洞壁面向上游发展ꎬ在ｔ＝５ｍｓ时亚声速流动占据整个稳定段ꎮ之后第１喷管处开始出现微弱的射流现象ꎬ并随时间逐渐增强然后减弱ꎬ最后在ｔ＝３０ｍｓ时完全消失ꎮ此时风洞的流场已经建立ꎬ第１喷管与稳定段内的流动完全为亚声速流动ꎮ图４㊀风洞启动过程中不同时刻的Ｍａｃｈ数云图Ｆｉｇ.４㊀Ｍａｃｈｎｕｍｂｅｒｃｌｏｕｄｃｈａｒｔａｔｄｉｆｆｅｒｅｎｔｔｉｍｅｓｄｕｒｉｎｇｔｈｅｓｔａｒｔ￣ｕｐｐｒｏｃｅｓｓｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌ密度梯度云图反映了风洞启动过程中流场内部的激波结构(如图５所示)ꎮ风洞启动瞬间ꎬ超声速气流通过第１喉道并在第１喷管与稳定段的交接处产生斜激波ꎮ由于风洞流场的对称性ꎬ斜激波会在风洞轴线上汇聚ꎬ之后在风洞壁面与轴线间来回反射ꎮ由于能量的耗散ꎬ激波强度沿流向逐渐衰减ꎮ随着稳定段内压力升高ꎬ第２喉道入口处的激波结构开始消失并逐渐向上游过渡ꎮｔ＝５ｍｓ时稳定段内的流动比较稳定ꎬ第１喷管处于欠膨胀状态ꎬ但却呈现典型过膨胀激波特征ꎮ激波后的逆压梯度导致边界层与喷管壁面分离ꎬ在靠近壁面两侧形成λ形激波ꎮ随着稳定段与储气段出口压力比的进一步升高ꎬ激波逐渐向第１喷管喉道处靠近ꎬ并伴随出现一系列的激波串ꎬ激波串的强度随时间逐渐减弱并消失ꎮ在ｔ＝３０ｍｓ时稳定段内流动完全为亚声速流动ꎬ因此不存在激波结构ꎮ图５㊀风洞启动过程中不同时刻的密度梯度云图Ｆｉｇ.５㊀Ｄｅｎｓｉｔｙｇｒａｄｉｅｎｔｃｌｏｕｄｃｈａｒｔａｔｄｉｆｆｅｒｅｎｔｔｉｍｅｓｄｕｒｉｎｇｔｈｅｓｔａｒｔ￣ｕｐｐｒｏｃｅｓｓｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌ为了更清楚地了解稳定段内的流动情况ꎬ图６给出了不同时刻稳定段内沿轴向的动量及流线云图ꎮｔ＝１ｍｓ时ꎬ在第２喉道入口处ꎬ由于超声速气流流动发生壅塞而形成流动分离区ꎮ该流动分离区内的涡卷成股旋转着沿风洞壁面向上游流去ꎬ并最终汇聚在第１喷管与稳定段的交接位置处ꎬ在后续过程中导致了激波串的出现ꎮ即使流动达到稳定后ꎬ稳定段内也存在较大的流动分离区ꎬ第１喉道的核心流动区域从壁面分离ꎬ该处的流动以恒定的面积向下游流去ꎮ由于流动分离区集中在第１喷管与稳定段交接位置处ꎬ该区域内的流动可能会对稳定段内的流动产生新的干扰ꎬ该扰动可能通过第２喉道进而对下游试验段内的流场产生影响ꎮ流动分离区的存在也会增加试验段内的总压损失ꎮ图６㊀风洞启动过程中不同时刻沿轴线方向动量云图Ｆｉｇ.６㊀Ｍｏｍｅｎｔｕｍｃｌｏｕｄｃｈａｒｔａｌｏｎｇｔｈｅａｘｉｓｄｉｒｅｃｔｉｏｎａｔｄｉｆｆｅｒｅｎｔｔｉｍｅｓｄｕｒｉｎｇｔｈｅｓｔａｒｔ￣ｕｐｐｒｏｃｅｓｓｏｆｔｈｅｗｉｎｄｔｕｎｎｅｌ为了分析该风洞启动时的流动演化过程ꎬ在图４所示的４个特殊流向位置提取Ｍａｃｈ数和静压以及总压随时间变化的数据进行分析(如图７所示)ꎮ其中位置１~４分别位于储气段出口㊁稳定段入口㊁稳定段出口和试验段内位置处ꎬ数据提取点距离风洞轴线法向１０ｍｍ处ꎮ２６第１期李创创ꎬ等:双喉道Ｌｕｄｗｉｅｇ管风洞启动过程及其有效运行时间延长(ａ)ｓｌｏｔ１㊀㊀㊀(ｂ)ｓｌｏｔ２(ｃ)ｓｌｏｔ３㊀㊀５．００(ｄ)ｓｌｏｔ４图７㊀沿双喉道Ｌｕｄｗｉｅｇ管风洞不同位置的压力和Ｍａｃｈ数随时间的变化Ｆｉｇ.７㊀ＶａｒｉａｔｉｏｎｏｆｐｒｅｓｓｕｒｅａｎｄＭａｃｈｎｕｍｂｅｒｗｉｔｈｔｉｍｅａｔｄｉｆｆｅｒｅｎｔｌｏｃａｔｉｏｎｓａｌｏｎｇｔｈｅｄｏｕｂｌｅ￣ｔｈｒｏａｔＬｕｄｗｉｅｇｔｕｂｅ储气段出口处(如图７(ａ)所示)的Ｍａｃｈ数在２０ｍｓ前后先增大后减小ꎬ在４０~１００ｍｓ之间基本保持不变ꎬ对应于双喉道Ｌｕｄｗｉｅｇ管风洞的有效运行时间ꎮ储气段出口处压力在启动过程中先下降ꎬ直到第１喉道不存在声速流动ꎻ此时第１喉道处激波串完全消失ꎬ储气段压力才开始逐渐上升ꎬ通过第１喉道质量流量开始减小ꎮ然后储气段压力维持稳定并存在９％的总压损失ꎬ１００ｍｓ后膨胀波再次到达快开阀ꎬ风洞一个运行周期结束ꎮ在稳定段入口处(如图７(ｂ)所示)ꎬＭａｃｈ数一开始由于流动膨胀而增加ꎬ随后稳定段内Ｍａｃｈ数下降到亚声速流动ꎮ风洞运行时ꎬ由于稳定段入口Ｍａｃｈ数仍然较高ꎬ稳定段入口静压值与总压值相比较小ꎮ值得注意的是ꎬ在１０ｍｓ前ꎬ稳定段入口Ｍａｃｈ数和压力存在两次急剧的波动ꎮ从前面的分析可以看出ꎬ第１次波动是由流动在第２喉道处发生壅塞㊁亚声速区向上游过渡导致的ꎬ而第２次波动是受此处激波串的出现和消失的影响ꎮ在稳定段出口(如图７(ｃ)所示)观察到同样现象ꎬ但两次Ｍａｃｈ数波动的幅值明显减小ꎮ随后Ｍａｃｈ数又出现明显的上升和下降趋势ꎮ该趋势是由第１喷管处的射流效应导致的ꎬ从侧面反映了当激波串逐渐消失时ꎬ随后射流效应的影响区域逐渐扩大然后又逐渐减小ꎮ从位置１ꎬ２ꎬ３中的压力变化可以看出ꎬ稳定段内的质量流量匹配过程以第１喉道不存在声速流动为节点ꎬ可以分为两个阶段:在０~１８ｍｓ时间段内通过第１喉道的质量流量基本不变ꎬ而通过第２喉道的质量流量逐渐增大ꎬ该阶段稳定段内的流场结构复杂ꎬ流动较为混乱ꎬ是流场结构趋于稳定的过程ꎻ在１８~３４ｍｓ时间段内ꎬ通过第１喉道的质量流量逐渐减小而通过第２喉道的质量流量逐渐增大ꎬ该阶段稳定段内流动为亚声速流动ꎬ无激波结构存在ꎬ是储气段与稳定段内压力的匹配过程ꎮ图７(ｄ)描绘了试验段内Ｍａｃｈ数和压力的演化过程ꎬ其中前３０ｍｓ是流场演化所需的时间ꎻ３０~１００ｍｓꎬ风洞处于稳定运行状态ꎬＭａｃｈ数和压力基本保持不变ꎮ风洞稳定运行时ꎬ从储气段到试验段中的总压基本保持不变ꎬ并维持在０.９ＭＰａ左右ꎮ由此说明ꎬ该气动布局带来的额外总压损失基本可以忽略不计ꎬ总压损失主要来源于风洞启动时产生３６气体物理２０２４年㊀第９卷的总压损失ꎮ从图７(ｃ)㊁(ｄ)可以看出ꎬ压力在稳定段内达到稳定大约需要３４ｍｓꎬ这意味着该双喉道Ｌｕｄｗｉｅｇ管风洞的有效运行时间为７０ｍｓꎮ值得注意的是ꎬ静压达到稳定所需的时间比Ｍａｃｈ数要长ꎮ因此ꎬ应以试验段静压幅值波动小于１％为判断标准来评价后续风洞试验段的有效运行时间ꎮ为了说明试验段内距离风洞轴线法向１０ｍｍ处的提取点数据能够用来评估风洞的有效运行时间ꎬ选取距离风洞轴线法向１０ꎬ２０ꎬ３０ꎬ４０ｍｍ处数据进行分析(如图８所示)ꎮ可以看到不同位置点所体现的风洞有效运行时间基本相同ꎬ采用距离风洞轴线１０ｍｍ位置处的提取点可以用来评估风洞的有效运行时间ꎮ图８㊀试验段内不同法向位置静压分布Ｆｉｇ.８㊀Ｓｔａｔｉｃｐｒｅｓｓｕｒｅｄｉｓｔｒｉｂｕｔｉｏｎａｔｄｉｆｆｅｒｅｎｔｎｏｒｍａｌｌｏｃａｔｉｏｎｓｉｎｔｈｅｔｅｓｔｓｅｃｔｉｏｎ为了进一步了解风洞稳定运行时稳定段内的流动情况ꎬ图９给出了风洞启动后６０ｍｓꎬ距离稳定段入口处０ꎬ０.２ꎬ０.４ꎬ０.６ｍ位置处Ｍａｃｈ数沿轴线法向方向的分布情况ꎮ可以看到ꎬ稳定段内依然存在微弱的射流效应ꎬ流动沿着稳定段发展ꎬ流速降低的同时流动会更加均匀ꎮ流动在稳定段出口处的Ｍａ<０.１ꎬ能够满足稳定段的设计要求ꎬ一般来说稳定段越长其匀流效果越好ꎬ但会带来更长的风洞启动时间ꎮ２.２　扩张段与稳定段优化设计对图６的分析发现ꎬ流动分离区主要集中在扩张段与稳定段的交接位置处ꎬ因此可以通过减小扩张段的扩张角并保持稳定段与扩张段总长不变ꎬ使喷管壁面占据部分分离区来减少稳定段内流场的演化时间ꎬ进而减少风洞的启动时间ꎮ图９㊀稳定段不同位置处Ｍａｃｈ数沿轴线法向分布Ｆｉｇ.９㊀Ｍａｃｈｎｕｍｂｅｒｄｉｓｔｒｉｂｕｔｉｏｎａｌｏｎｇｔｈｅｎｏｒｍａｌｄｉｒｅｃｔｉｏｎｏｆｔｈｅａｘｉｓａｔｄｉｆｆｅｒｅｎｔｌｏｃａｔｉｏｎｓｉｎｔｈｅｓｔａｂｌｅｓｅｃｔｉｏｎ为了研究不同扩张角与稳定段组合对风洞启动时间以及稳定段内流场品质的影响ꎬ对扩张角为１７ʎꎬ６.７ʎꎬ４.１ʎ和３.２ʎ的风洞进行非定常数值模拟ꎮ随着扩张角减小ꎬ扩张段变长而稳定段相对变短ꎬ图１０给出了风洞启动后１ｍｓ和１５ｍｓ时不同扩张角和稳定段组合下沿轴向密度梯度云图ꎮ在ｔ＝１ｍｓ时ꎬ随着扩张角减小ꎬ由于第１喷管与稳定段交接处的拐角增大ꎬ形成的斜激波有明显减弱ꎬ同时稳定段变短会使稳定段内的激波结构更少ꎮ在ｔ＝１５ｍｓ时ꎬ随着扩张角减小ꎬ在第１喉道下游附近的激波串数量逐渐减少ꎬ说明较小扩张角稳定段内的流场演化速度更快ꎬ即压力匹配的第１阶段完成更快ꎮ对不同扩张角与稳定段组合下沿风洞轴线动量以及流线分布云图进行分析(见图１１)ꎮｔ＝１ｍｓ时ꎬ第２喉道入口处的流动分离区并没有随扩张角的减小产生明显变化ꎮ在ｔ＝１５ｍｓ时ꎬ对于不同扩张角与稳定段组合ꎬ流动分离区主要分布在稳定段内动量较低的区域ꎬ且随扩张角减小而减小ꎻ而较高动量区域为流动的核心区域ꎬ其分布基本一致ꎮ这说明改变扩张角会影响稳定段内分离区大小ꎬ而基本不会对稳定段内核心区域内流动产生较大影响ꎮ(ａ)ｔ＝１ｍｓ４６第１期李创创ꎬ等:双喉道Ｌｕｄｗｉｅｇ管风洞启动过程及其有效运行时间延长(ｂ)ｔ＝１５ｍｓ图１０㊀不同扩张角下沿轴线方向上密度梯度分布云图Ｆｉｇ.１０㊀Ｄｅｎｓｉｔｙｇｒａｄｉｅｎｔｃｌｏｕｄｃｈａｒｔａｌｏｎｇｔｈｅａｘｉｓｄｉｒｅｃｔｉｏｎａｔｄｉｆｆｅｒｅｎｔｅｘｐａｎｓｉｏｎａｎｇｌｅｓ(ａ)ｔ＝１ｍｓ(ｂ)ｔ＝１５ｍｓ图１１㊀不同扩张角下沿轴线方向上动量分布云图Ｆｉｇ.１１㊀Ｍｏｍｅｎｔｕｍｃｌｏｕｄｃｈａｒｔａｌｏｎｇｔｈｅａｘｉｓｄｉｒｅｃｔｉｏｎａｔｄｉｆｆｅｒｅｎｔｅｘｐａｎｓｉｏｎａｎｇｌｅｓ为了分析不同扩张角与稳定段组合对风洞启动时间的影响ꎬ取不同扩张角下位置１ꎬ３ꎬ４处的Ｍａｃｈ数和压力分布随时间变化曲线如图１２~１４所示ꎮ从图１２(ａ)㊁(ｂ)可以看出ꎬ储气段入口处Ｍａｃｈ数在１０~２０ｍｓ之间的峰值随扩张角减小而降低ꎬ达到稳定的时刻逐渐提前ꎮ储气段入口处的压力变化反映了稳定段内质量流量匹配的第１阶段与第２阶段所需时间随扩张角减小均缩短ꎬ说明较小扩张角与稳定段组合加速了风洞的启动过程ꎮ稳定段出口处的Ｍａｃｈ数和压力(如图１３(ａ)㊁(ｂ)所示)达到稳定的时间随扩张角减小逐渐提前ꎮ风洞运行时ꎬＭａｃｈ数会随着扩张角减小而有略微增大但压力基本保持不变ꎬ可能是由于稳定段越短匀流效果越差ꎬ同时受第１喉道射流效应的影响越大ꎮ从图１４(ａ)㊁(ｂ)两图可以得出结论ꎬ不同扩张角与稳定段组合几乎不会对风洞稳定运行时试验段内Ｍａｃｈ数和压力产生影响ꎬ较小扩张角组合可以明显提高风洞的有效运行时间ꎮ由于风洞的有效运行时间应以试验段内静压幅值波动小于１％作为判断标准ꎬ依照此判断标准ꎬ扩张角为１７ʎꎬ６.７ʎꎬ４.１ʎ和３.２ʎ对应的启动时刻为３４ꎬ２９ꎬ２２ꎬ１８ｍｓꎬ对应的结束时刻均为１０４ｍｓꎬ因此相应的风洞有效运行时间分别为７０ꎬ７５ꎬ８２ꎬ８６ｍｓꎮ通过采用减小扩张角并保持扩张段与稳定段总长不变的方式ꎬ风洞有效运行时间增加将近１６ｍｓꎬ相较之前提升了约２０％ꎬ大大提高了风洞的试验能力ꎬ这对于高超声速风洞设施来说十分关键ꎮ(ａ)Ｍａｃｈｎｕｍｂｅｒ(ｂ)Ｓｔａｔｉｃｐｒｅｓｓｕｒｅ图１２㊀位置１处不同扩张角下Ｍａｃｈ数和压力随时间的变化Ｆｉｇ.１２㊀ＶａｒｉａｔｉｏｎｏｆＭａｃｈｎｕｍｂｅｒａｎｄｐｒｅｓｓｕｒｅｗｉｔｈｔｉｍｅａｔｄｉｆｆｅｒｅｎｔｅｘｐａｎｓｉｏｎａｎｇｌｅｓａｔｓｌｏｔ１ｐｏｓｉｔｉｏｎ５６。

RAE2822翼型跨声速绕流的CFD计算专业：工程力学学号：******姓名：**指导老师：***********2015/4/12问题描述：模型：RAE2822翼型来流条件：Ma∞＝0.729，α=2.31°计算状态：（1）无粘绕流（欧拉方程、无附面层网格）。

（2）至少3种不同湍流模型计算粘性绕流（同一带附面层网格，y＋≈30）。

（3）采用S-A湍流模型，建立4种不同y＋的网格计算（y＋<1, y＋≈10, y＋≈30, y ＋≈50）。

（4）对y＋≈50的网格，额外采用流场求解网格自适应功能（基于压力梯度）进行计算。

第一章物理模型1.1 RAE2822翼型将RAE2822翼型点导入CATIA，扫描得到翼型曲线如下。

图1.11.2 计算参数来流条件：Ma∞＝0.729，α=2.31来流压力和温度：P=101325Pa，T=300K第二章网格划分及求解器设置2.1网格划分2.1.1 无附面层网格划分本文使用Gridgen软件作为网格生成器，并采用结构网格进行划分。

远场边界上下、左侧距离约为弦长10倍，向后远场范围设置为翼型弦长15倍。

在翼型上下弧线分别布节点100，减小两端节点间距，在翼面附近进行适当加密。

分块生成domain。

图2.1 无附面层网格2.1.2 带附面层网格划分对于考虑粘性的计算模型。

需要建立附面层，附面层第一层厚度通过Y+公式计算。

对于给定y+，可以得到对应第一层厚度如下：在翼型上下弧线及尾迹上分别布置100节点，并在翼型两端进行适当加密。

使用extrude生成附面层网格，共15层，初始厚度设置为0.1259mm，增长率为1.1。

按无附面层网格划分方式设置外围流场边界。

布置节点，生成结构网格。

生成体。

图 2.2 带附面层网格2.2 求解器及边界条件设置将流动分析求解程序设置为2D FLUENT。

边界条件设置如下：外围设置为压力远场，上下翼面设置为wall。

输出能被FLUENT读入的文件。

翼型气动特性数值模拟研究翼型是飞行器的重要组成部分，其气动特性（如升力、阻力、升力系数、升阻比等）对飞行器的性能有很大影响。

由于实验设备和费用的限制，气动试验成本高昂，因此数值模拟成为了研究翼型气动特性的主要方法之一。

数值模拟方法数值模拟方法主要包括计算流体力学（CFD）方法和边界元方法。

其中，CFD方法是一种利用数值计算方法处理流体动力学问题的方法，可根据所建立的数学模型，通过计算机模拟流体的运动状态，获得流体介质的相应物理量。

而边界元方法则是一种计算机辅助工具，针对问题内部的微观变化关系较弥散的情况下，仅需检查问题外缘的变化，即可通过边界元法反映问题内部变化。

两种方法的原理和适用范围存在区别。

本文主要讨论CFD方法，根据不同模型和假设，CFD方法分为欧拉方程模型、纳维—斯托克斯方程模型等。

其中，普遍认为海拔高度2000米，马赫数0.3的常温常压环境下，采用欧拉方程模型就能较为精确地预测翼型的气动特性。

欧拉方程模型及其应用欧拉方程模型的基本假设是流体为理想气体，连续性方程为无穷小量，流体的运动状态由欧拉方程控制。

其中，欧拉方程考虑了三个物理量：密度（rho）、速度（v）、热力学气压（p），并描述了它们之间的关系。

欧拉方程模型的适用范围很广，可以处理多种气流复杂情况，可以在空气、液体（如水）及其它流体的流动中预测相关的力学变量，有效地用于翼型气动特性数值模拟。

实例分析以NACA 0012翼型为例，它是由美国航空航天局设计的一支标准组合翼，被广泛应用于飞行器领域。

研究采用Ansys Fluent 15.0数值模拟软件，通过对NACA 0012翼型的气动特性的分析，验证了欧拉方程模型在预测翼型的气动特性方面的有效性。

翼型模型的几何尺寸定义采用了标准的NACA 4位数型号，其的绘制遵守了标准的绘制规则。

通常，翼型的比尺寸Re数（不能大于100万）是气动特性数值模拟的一个关键因素，它决定了模拟结果的准确度。

低速风洞内部流场数值模拟代燚;陈作钢;马宁;任泽斌【摘要】以上海交通大学在建的多功能低速风洞为研究对象,针对其设计上的独特性,采用CFD方法首次实现对风洞内部流场的整体数值模拟,风扇段流场利用MRF 模型处理,数值模拟中计及了阻尼网的整流效果.同经验设计法的气动设计结果进行对比发现,数值模拟得到的流量与设计流量差距很小,压力损失在一定程度上较为吻合.数值模拟揭示了风洞内部存在的流动分离现象,为提高风洞的气动性能提供了改进方向.对风洞同轴度下降引起的流场品质变化进行了分析,结果表明微小同轴度偏差对流场品质的影响不大.【期刊名称】《空气动力学学报》【年(卷),期】2014(032)002【总页数】7页(P203-208,213)【关键词】低速风洞的气动设计;CFD;阻尼网;流动分离;数值风洞【作者】代燚;陈作钢;马宁;任泽斌【作者单位】上海交通大学船舶海洋与建筑工程学院,上海 200240;上海交通大学海洋工程国家重点实验室,上海200030;上海交通大学海洋工程国家重点实验室,上海200030;中国空气动力研究与发展中心,绵阳621000【正文语种】中文【中图分类】V2110 引言风洞是一种能在管道内产生人工可控气流，利用相似原理进行各种空气动力学实验的设备。

低速风洞是指试验段气流速度小于130m／s（马赫数小于0.4）的风洞［1］。

从20世纪四十年代开始，国外开始建造低速风洞，迄今低速风洞已经得到了很大的发展。

计算机信息技术的迅速发展推动了计算机与风洞的一体化，风洞的试验能力得到了极大提升，获取的信息量增多，试验数据质量明显提高，试验费用下降［2］。

技术进步导致试验市场对于风洞性能的要求逐步提高，那些性能差、效率低的风洞必将被逐步淘汰。

风洞的气动设计主要依靠经验公式。

随着计算流体力学（Computational Fluid Dynamics，CFD）的快速发展，CFD方法已开始用于风洞试验结果对比及流动机理解释等方面，数值风洞的研究与应用也取得了很大的进展［3－4］，在结构风工程研究中具有战略意义，其长期目标是用计算机部分甚至全部代替风洞。

2017年发达国家风洞试验研究进展-力学论文-物理论文——文章均为WORD文档，下载后可直接编辑使用亦可打印——摘要：从三个方面综述了2017年国外发达国家风洞试验发展动态, 包括国家风洞试验设备改造情况、风洞试验技术发展情况以及大型风洞开展的试验研究项目。

在此基础上, 给出了风洞试验未来发展趋势的分析研究结果。

关键词：风洞试验; 风洞测试技术; 风洞设备;0、引言大型风洞试验设备被视为国家的战略资源。

风洞试验是开展先进飞行器预研、型号设计/评估和CFD工具验证的重要手段。

透过2017年度国外航空航天发达国家风洞试验设备、试验技术和风洞试验情况的发展动态, 我们可以从一个侧面了解和认识国外发达国家风洞试验和飞行器研究的现状, 分析发展趋势, 从而为我国风洞设备建设和型号研制提供参考。

1、加强核心风洞设备改造, 提高试验模拟能力根据美国国会颁布的81-415公共法案, 美国国家风洞试验设备主要集中于建设NASA (美国国家航空航天局) 和军方的AEDC (阿诺德工程发展综合体) 。

进入21世纪以来, 美国国家风洞试验设备已完成去产能工作, 国家资源向国家基本核心风洞集中, 提高资金利用率和风洞更新改造的科学化管理水平。

根据NASA航空评估和试验能力项目(AETC) , 2017年, 美国NASA完成了对兰利、格林和艾姆斯三个研究中心12座核心风洞设备2016财年的评估工作, 掌握了主要风洞设备现状、可靠性以及满足未来五年试验的能力, 更新了设备管理数据库, 为科学管理风洞维修改造奠定了基础。

NASA格林中心IRT结冰风洞采用组合使用标准喷嘴和Mod 1喷嘴的方式, 调试完成IRT水滴分布, 满足FAA 25部附录O冻雨(FZDZ) 模拟MVD40的试验模拟要求。

世界最大的全尺寸风洞美国NASA的NFAC (国家空气动力设施) 在经历了2003年关停、2008年交由军方AEDC管理运营后, 美国国防部投资、由Jacobs工程/宇航试验联盟(ATA) 负责对其进行了全面恢复和升级改造工作。

基于风洞试验模型的跨声速颤振研究侯良学;张戈;刘南;王冬;钱卫;杨希明【摘要】飞行器跨声速工况下颤振边界快速下降,是结构设计和强度校核重点关注的状态之一.目前工程中采用基于偶极子格网法的线性分析手段无法准确预测跨声速颤振边界,风洞试验仍然是研究飞行器跨声速颤振特性的重要手段;以两套颤振试验标准模型为研究对象,在FL-3风洞中开展了风洞跨声速试验研究,采用PEAK-HOLD亚临界响应分析方法预测了模型的跨声速颤振边界,并利用ZAERO和CFD/CSD耦合两种数值计算方法预测了试验模型的颤振边界.结果表明:PEAK-HOLD亚临界响应分析方法预测颤振边界具有较好的趋势性,颤振边界合理可靠;ZAERO线性方法对跨声速颤振边界的预测精度较低,而基于CFD/CSD耦合的非线性方法得到的跨声速颤振边界与试验吻合较好,相互验证了风洞试验和数值计算的可靠性.【期刊名称】《振动与冲击》【年(卷),期】2019(038)013【总页数】7页(P236-241,260)【关键词】跨声速颤振;颤振标模;亚临界响应;风洞试验;数值计算【作者】侯良学;张戈;刘南;王冬;钱卫;杨希明【作者单位】中国航空工业空气动力研究院,沈阳 110034;高速高雷诺数气动力航空科技重点试验室,沈阳110034;中国航空工业空气动力研究院,沈阳 110034;高速高雷诺数气动力航空科技重点试验室,沈阳110034;中国航空工业空气动力研究院,沈阳 110034;高速高雷诺数气动力航空科技重点试验室,沈阳110034;中国航空工业空气动力研究院,沈阳 110034;高速高雷诺数气动力航空科技重点试验室,沈阳110034;大连理工大学航空航天学院,辽宁大连116024;中国航空工业空气动力研究院,沈阳 110034;高速高雷诺数气动力航空科技重点试验室,沈阳110034【正文语种】中文【中图分类】V215.3在跨声速区，由于压缩性和气动非线性的影响，飞行器的颤振速压边界一般降低较大，颤振边界随马赫数的变化曲线会呈现出所谓的跨声速“凹坑”现象，最小临界颤振速压值通常出现在跨声速范围。

2020年8月Aug. ,2020第36卷第4期Vol. 36, No. 4滨州学院学报Journal of Binzhou University【航空科学与工程研究】跨声速层流机翼变保真度分层优化设计张同任1,陈永彬2,魏德宸彳(1.滨州学院飞行学院；2.滨州学院航空工程技术研究院，山东滨州256603)摘要：层流减阻作为一项革新性技术，为显著提高飞机气动性能提供了可能。

针对跨声速自然层流机翼设计方法开展研究，建立了基于变保真度流动分析的分层优化设计模型来提高机 翼的优化效率。

将求解欧拉方程耦合线性稳定性理论的方法作为流动转按判断的依据，采用建 立的变保真度分层优化方法耦合多岛并行遗传算法开展跨声速层流机翼优化设计研究。

优化 后，机翼表面层流区域由初始外形的23.05%增大到优化后外形的45. 50%；表面摩擦阻力系数 明显降低，由初始外形的0.004 74降低到优化后外形的0. 003 94。

另外，采用变保真度分层优 化方法只需较少的时间便得到高精度模型上的收敛解，对比采用非并行遗传算法，计算效率提高 了 9倍以上，验证了变保真度多层优化设计模型的高效性和在求解实际工程问题中的潜力。

关键词：跨声速层流机翼；自然层流技术；多岛并行遗传算法；变保真度分层优化方法中图分类号：V 211. 3 文献标识码：A EOI ：10.13486/j. cnki. 1673 - 2618. 2020. 04. 0010引言近年来，民用航空对环境的影响越来越受到公众的关注，进一步减少二氧化碳排放和燃油的消耗是航 空领域关注的热点和研究方向，对提高民航客机的性能指标具有重要指导意义⑴幻o 飞机飞行效率的提升 有望使大中型客机的油耗减小20%以上，其中，7%的贡献来自发动机工作效率的提升,4. 5%来自新型材 料的使用,8.5%来自飞机外形布局优化、空气动力效率的提升可见，飞机气动布局设计仍将是未来 民航客机商业竞争的一个焦点。

小攻角时风力机翼型边界层特性的数值模拟李德顺;李秋燕;李仁年;李银然;王成泽【摘要】通过求解三维不可压雷诺时均N-S方程,当雷诺数为3×106,攻角为0°时,研究了NACA0012翼型绕流的边界层厚度分布和边界层内流体切向速度的变化规律.研究发现,沿着翼型弦线方向,从前缘到后缘,边界层的名义厚度、位移厚度、动量损失厚度以及能量损失厚度均呈增大趋势,且数值均很小,四种厚度的最大值分别占翼型弦长的1.25％、0.36％、0.17％和0.29％;在边界层内流体黏性的影响明显,流体切向速度与其势流解的比值沿着翼型吸力面外法线方向先迅速增大,之后增长率逐渐减小,当法向高度大于边界层名义厚度后基本保持不变,呈现出典型的边界层速度剪切特性.%By means of solving three-dimensional Reynolds time-averaged equation of incompressible flow with Reynolds number of 3 × 106,the distribution of thickness of the boundary layer and the variation pattern of the tangential velocity in the boundary layer on NACA0012 airfoil with attack angle of 0 deg was studied.It was found by the investigation that the nominal thickness,displacement thickness,the momentum thickness,and energy thickness of the boundary layer would all exhibit a increasing trend from leading edge to trailing edge of airfoil and the four thicknesses would be very small.Their maximum value will be1.25％,0.36％,0.17％,and 0.29％ of airfoil chord,respectively.Besides,the effect of the fluid viscosity will be obvious in the boundary layer,and the ratio of the tangential velocity to that of potential flow will increase quickly at first along the external normal of airfoil suction side,then its increase rate will gradually reduce,and at last will remain basically unchanged whenthe normal distance was greater than the nominal thickness of the boundary layer,showing the typical velocity shear characteristics.【期刊名称】《兰州理工大学学报》【年(卷),期】2018(044)001【总页数】5页(P47-51)【关键词】风力机;翼型;边界层;数值模拟【作者】李德顺;李秋燕;李仁年;李银然;王成泽【作者单位】兰州理工大学能源与动力工程学院,甘肃兰州730050;兰州理工大学甘肃省流体机械及系统重点实验室,甘肃兰州 730050;兰州理工大学甘肃省风力机工程技术研究中心,甘肃兰州 730050;兰州理工大学能源与动力工程学院,甘肃兰州730050;兰州理工大学能源与动力工程学院,甘肃兰州730050;兰州理工大学甘肃省流体机械及系统重点实验室,甘肃兰州 730050;兰州理工大学甘肃省风力机工程技术研究中心,甘肃兰州 730050;兰州理工大学能源与动力工程学院,甘肃兰州730050;兰州理工大学甘肃省流体机械及系统重点实验室,甘肃兰州 730050;兰州理工大学甘肃省风力机工程技术研究中心,甘肃兰州 730050;兰州理工大学能源与动力工程学院,甘肃兰州730050【正文语种】中文【中图分类】TK83翼型绕流是常见的绕流问题之一,像换热器中的扰流子、船用舵、螺旋桨的桨叶、机翼等,均为翼型绕流问题.由于实际流体都具有黏性,当发生翼型绕流时,在翼型表面会形成一层很薄的剪切层,即边界层.对于风力机,流体绕翼型边界层的流动特性直接影响风轮的性能和载荷,因此,研究翼型绕流的边界层问题非常必要.国内,边界层问题的研究主要集中于飞机机翼绕流方面[1-2],而且以跨声速绕流问题居多,主要通过Euler方程与边界层积分方程耦合方法进行求解[3-4]；水力机械的边界层问题方面,刘彦[5]将两相流动与边界层理论相结合,分析了含沙水流绕流过流部件时边界层内的流动规律,得出一种定常固液两相流动边界层的数学模型.朱玉才等[6]通过研究得到离心泵叶片压力面的固液两相流的边界层特性；关于风力机翼型边界层问题,杨从新等[7]用直接模拟方法对风力机专用翼型开展了三维数值模拟研究,发现该方法对翼型的气动性能在线性区有很高的预测精度.李银然等[8]发现计算模型的维数对于翼型的气动性能有一定的影响,用三维模型计算得到的翼型的气动性能要比二维计算的结果更加准确.李隆键等[9]通过数值模拟发现分离流动的三维特性是影响翼型气动性能的重要因素,二维模型并不适用于翼型气动性能的计算.曲佳佳[10]的相关研究表明三维旋转效应对于风力机叶片气动性能的影响不能忽略,按动量叶素理论计算的叶片表面的载荷分布与二维结果存在一定差异.戴韧等[11]研究发现边界层理论和方法在风力机的气动性能分析领域有很好的应用前景.国外意大利帕多瓦大学Castelli等[12]通过数值方法研究了NACA0012翼型边界层的转捩特性,美国亚利桑那大学Balzer等[13]通过直接数值模拟的方法,研究了低雷诺数条件下翼型吸力面的分离特性以及分离控制方法.本文通过求解三维不可压雷诺时均N-S方程,研究了NACA0012翼型的边界层特性,得到边界层的名义厚度、位移厚度、动量损失厚度和能量损失厚度沿翼型弦向的分布规律,以及边界层内流体切向速度的分布情况.1 研究对象选取NACA0012对称翼型作为研究对象,该翼型的相对厚度为12%,弦长c为1 m,展向长度为0.25 m其几何外形如图1所示.图1 NACA0012 翼型的几何外形Fig.1 Geometry of airfoil NACA00122 数值计算方法2.1 数学模型在风力机翼型绕流问题中,对于不可压缩黏性流体的三维定常流动,直角坐标系下的微分方程如下：连续性方程：(1)动量方程：2ui(2)式中：下标i,j为直角坐标系中x、y、z三个方向的坐标分量,i,j=1,2,3；p为流体的平均压力；ρ为流体的密度；为拉普拉斯算子.2.2 计算域及网格划分采用C型计算域,如图2所示.入口边界IN距翼型尾缘3c,出口边界OUT距翼型尾缘4c,上下两个入口边界IN1和IN2均距翼型弦线3c.计算网格采用结构化网格,沿翼型弦向布置111个节点,展向布置99个节点,翼型表面采用C型贴体网格,并进行局部加密处理,首层网格高度d0为3×10-4m(见图3).整个计算区域的网格总数约6 000 000个.图2 计算域Fig.2 Computational domain图3 翼型周围网格Fig.3 Grids around airfoil2.3 湍流模型和边界条件设定本文参考1945年NACA LTT Wind Tunnel针对NACA0012翼型的一个风洞实验工况,雷诺数R e=3×106,攻角α=0°,湍流强度为0.03%,翼型表面光滑.采用Fluent软件进行模拟计算,求解器采用基于压力的三维隐式求解器；分别采用无黏的Inviscid 模型和k-ω SST模型求解理想流体的势流流动和实际流体的黏性流动；离散格式均采用二阶迎风格式；速度和压力耦合方式采用SIMPLEC；来流为水平均匀来流.图2中,计算域的IN、IN1和IN2边界均为速度进口,给定来流沿水平方向,速度大小为53.67 m/s；OUT边界为压力出口；翼型表面固定无滑移；三维计算域的前后两界面均采用对称边界.另外,边界层外边界的边界条件为(3)式中:δ为边界层名义厚度;ue(x)为边界层外边界上势流的切向速度.3 边界层的确定及数值方法验证3.1 边界层的确定方法根据普朗特边界层理论,在势流的近似计算中可略去边界层的厚度,得到沿物体表面的流速和压力分布,并将其认为是边界层外边界上的速度和压力分布.本文采用Inviscid模型求解无黏势流流动(来流为理想气体,速度分布均匀、无旋),得到翼型吸力面上沿弦向的流体的切向速度ue(x)；基于RANS方程求解流体绕翼型的黏性流动,将吸力面上各点外法线上流体切向速度为0.99ue(x)的点连接而得到的曲线认为是翼型边界层的外边界,外边界与翼型吸力面之间的距离就是边界层的名义厚度.即根据式(4～6)的积分关系得到NACA0012翼型的位移厚度δ*、动量损失厚度θ以及能量损失厚度δ**沿翼型弦向的分布情况.3.2 黏性求解数值模拟方法的验证为确保NACA0012翼型绕流黏性解的准确性,进而得到更为准确的边界层的相关特性.首先将通过数值模拟所得的翼型的升、阻力系数CL和CD与风洞实验数据,以及通过翼型设计软件Profili计算所得到的数据进行了对比.数值模拟计算采用Fluent软件,将来流设为实际定常流体,湍流模型选择k-ω SST模型,通过模拟分别得到0°～19°内不同攻角下NACA0012翼型的升、阻力系数CL和CD,对比结果如图4所示.图4 NACA0012翼型升、阻力系数的对比 Fig.4 Comparison of lift and drag coefficient of NACA0012 airfoil由图4可以看出,NACA0012翼型的升、阻力系数CL和CD的数值模拟结果、实验数据以及Profili软件的计算结果在α=6°之前吻合得很好,当α>6°,三种形式下的CL和CD值出现了一定差距,并且随着攻角α的增大这种差距也逐渐增大.这是由于随着攻角的增大,在攻角为6°左右,在翼型的吸力面逐步发生了流动分离.由于本文针对小攻角(0°攻角)下NACA0012翼型的边界层特性开展研究工作,因此所采用的数值计算方法满足精度要求.4 计算结果及其分析4.1 无黏流动求解结果对0°攻角下NACA0012翼型绕流流动采用无黏的Inviscid模型进行模拟求解,得到翼型吸力面上流体的切向速度ue(x)沿弦向的分布,如图5所示.图5 翼型吸力面势流切向速度沿弦向的分布Fig.5 Distribution of tangential velocity of potential flow on airfoil suction side由图5可以看出,在NACA0012翼型的吸力面,当0<x/c<0.14时,流体的切向速度u迅速增大,约在x/c=0.14附近达到最大值,umax约为60 m/s；当0.14<x/c<0.99时,流体切向速度u逐渐减小；当0.99<x/c<1时,流体切向速度u 迅速减小.表1为不同弦长位置处翼型吸力面上流体的切向速度.表1 Re=3×106、α=0°时不同弦长位置NACA0012翼型吸力面势流的切向速度Tab.1 Tangential velocity of potential flow at different chord position on suction side of NACA0012 airfoil under condition of Re=3×106,α=0°x/cue(x)/(m·s-1)0.99ue(x)/(m·s-1)0.010437.798837.42080.101959.824159.22590.207459.142058.55060.305 358.797458.20950.403856.716956.14970.502555.336454.78300.601653.882 253.34340.700851.970551.45080.800550.463649.95890.902647.707147.230 00.991343.178342.74654.2 边界层厚度分布根据3.1节边界层的确定方法,研究了NACA0012翼型绕流的边界层特性,得到该翼型吸力面边界层外边界的分布、边界层的名义厚度δ、位移厚度δ*、动量损失厚度θ和能量损失厚度δ**沿弦向的分布情况.由图6a可以看出,在翼型的吸力面,沿着翼型弦向,当0<x/c<0.1时,边界层外边界几乎紧贴着壁面,边界层名义厚度接近于0;当x/c>0.1后,边界层外边界逐渐脱离壁面,并且与壁面间的距离逐渐增大,说明边界层的名义厚度沿着翼型弦向逐渐增大.由图6b可知,随着弦向位置的向后推移,边界层名义厚度δ逐渐增大,且基本呈现出线性增长的规律.另外,由图6可知,边界层的名义厚度是一个很小的值,在翼型前缘处几乎为0,在翼型尾缘附近(x/c=0.991 3)边界层厚度达到最大,约为0.001 3 m,占弦长1.25% .图6c中边界层的位移厚度δ*、动量损失厚度θ图6 边界层的外边界及四种厚度的分布规律Fig.6 Outer boundary of boundary layer and distribution of four kinds of thickness以及能量损失厚度δ**的值也是沿着翼型弦向不断增长的,增长趋势也基本一致,且三种厚度的值都非常小,其最大值分别占翼型弦长的0.36%、0.17%和0.29%.另外,对比NACA0012翼型边界层四种厚度的分布规律,可以发现,δ>δ*>δ**>θ.4.3 边界层内流体切向速度分布通过后处理软件CFD-POST得到翼型弦线方向上x/c=0.010 4、0.101 9、0.502 5、0.991 3四个位置处吸力面外法线上流体的切向速度u,并绘制出边界层内这四个位置流体切向速度的分布曲线,如图7所示.图7 翼型吸力面边界层内切向速度分布Fig.7 Distribution of tangential velocity of boundary layer on airfoil suction side图7中横坐标为相对速度,纵坐标为高度,图7a和7b是翼型前缘(x/c=0.010 4和0.101 9)处相对速度u/ue的值与法向高度y的关系,图7c和图7d分别是翼型中部(x/c=0.502 5)和尾缘(x/c=0.991 3)处u/ue的值与y的关系.可明显看出,在边界层内沿着翼型弦向,从前缘到尾缘,相对速度u/ue与高度y的关系曲线相似,u/ue的值均从壁面处开始先迅速增大,之后增长率有所减小,当高度y大于边界层名义厚度δ之后保持最大值不变.这表明在边界层内流体具有较大的速度梯度,这是由于流体黏滞力牵制的结果.总体来看,数值模拟方法得到的边界层内切向速度分布曲线表现出典型的边界层速度剪切的特点.5 结论1) 边界层名义厚度沿翼型弦向呈增长趋势,翼型尾缘处厚度最大.2) 边界层名义厚度相比翼型的弦长要小很多,翼型尾缘处名义厚度的最大值为0.001 3 m,占弦长的1.25%.3) 边界层的位移厚度、动量损失厚度以及能量损失厚度的值均非常小,最大值同样出现在尾缘处,分别占翼型弦长的0.36%、0.17%和0.29%.4) 在边界层内流体黏性剪切作用影响明显,沿着翼型吸力面外法线方向,流体切向速度与其势流解的比值先迅速增大,之后增长率有所减小,当法向高度大于边界层名义厚度之后基本保持不变,呈现出典型的边界层速度剪切特性.参考文献：[1] 卞荫贵.边界层研究进展 [J].力学进展,1980,2(增刊):7-26.[2] 包涵龄.机翼边界层计算及其在飞机设计中的应用 [J].航空学报,1987,8(11):553-560.[3] 蔡罕龙,李锋.Euler方程与边界层积分方程耦合求解跨音速翼型绕流 [J].航空学报,1991,12(5):221-227.[4] 李凤蔚,鄂秦.跨声速翼型绕流的Euler/边界层方程干扰数值解 [J].空气动力学学报,1993,11(1):49-56.[5] 刘彦.固液两相流动边界层的数学模型及流场计算 [D].兰州:甘肃工业大学,2001.[6] 朱玉才,吴玉林,曲衍国,等.离心泵叶片压力面固液两相流的边界层分析 [J].工程热物理学报,2005,26(3):429-431.[7] 杨从新,金开,王秀勇.一种风力机专用翼型气动性能的三维数值模拟 [J].兰州理工大学学报,2012,38(6):39-43.[8] 李银然,李仁年,王秀勇,等.计算模型维数对风力机二维翼型气动性能预测的影响[J].农业机械学报,2011,42(2):115-119.[9] 李隆键,张义华,唐胜利.风力机翼型边界层分离流动三维特性的数值模拟 [J].重庆理工大学学报,2010,24(9):41-46.[10] 曲佳佳.风力机叶片气动载荷的计算方法研究 [D].北京:中国科学院研究生院,2014.[11] 戴韧,金科逸,张国强,等.风力机翼型绕流分析边界层方法与应用:二维流动 [J].工程热物理学报,2010,31(11):1855-1858.[12] CASTELLI M R, GARBO F, BENINI E. Numerical investigation of laminar to turbulent boundary layer transition on a NACA0012 airfoil for vertical-axis wind turbine applications[J].Wind Engineering,2012,35(6):661-686. [13] BALZER W,FASEL H F.Direct numerical simulation of laminar boundary-layer separation and separation control on the suction side of an airfoil at low Reynolds number conditions[C]//Fluid Dynamics Conference and Exhibit.Arizona:University of Arizona Press,2010:1-20.。

用ANSYS实现二维翼型风洞试验数值模拟.
张鹏;马玉清
【期刊名称】《科学技术与工程》
【年(卷),期】2008(008)002
【摘要】基于有限元的ANsYs软件,可以分析二维翼型流场,解决作用于气动翼(叶)型上的升力和阻力,并可以得到流场中翼型表面的压力与速度分布.与此同时,通过对ANSYs的二次开发,编制相关的三分力系数求解模块,得到各系数曲线.现以标模翼型为例,利用此模块得到的三分力系数与NF-3风洞试验值作对比,验证了用ANSYS 实现二维翼型风洞数值模拟的可行性.
【总页数】5页(P570-574)
【作者】张鹏;马玉清
【作者单位】西北工业大学国家翼型叶栅重点实验室,西安710072;西北工业大学国家翼型叶栅重点实验室,西安710072
【正文语种】中文
【中图分类】V211.4
【相关文献】
1.基于Euler方程的二维翼型大迎角数值模拟 [J], 陈哲
2.水平轴风机叶片截面与二维翼型的数值模拟 [J], 孔秀;李龙
3.翼型低速风洞试验的数值模拟 [J], 傅建明
4.二维多段翼型网格生成与数值模拟研究 [J], 胡偶;赵宁
5.二维翼型的积冰数值模拟精度提升方法 [J], 胡义明
因版权原因，仅展示原文概要，查看原文内容请购买。