飞机尾涡演变及快速预测的大涡模拟研究1)

林孟达 崔桂香 张兆顺 许春晓 黄伟希

(清华大学航天航空学院，北京100841)

创刊60周年专栏

飞机尾涡演变及快速预测的大涡模拟研究1)

林孟达 崔桂香2)张兆顺 许春晓 黄伟希

(清华大学航天航空学院，北京100841)

随着我国人民生活水平的提高，航空运输的重要性与日俱增，航班延误问题也日益严重.尾流间隔(保障后机不受前机尾流影响的最小安全间隔)是制约机场效率的关键因素.针对这一工程应用问题，采用大涡模拟方法研究飞机尾涡在大气中的演变特性.研究工作首先发展了飞机尾涡演变的大涡模拟方法，将自适应网格技术应用于飞机尾涡演变的大涡模拟，大幅减少所需的网格量，提高计算效率.提出了升力面尾涡生成方法，在不增加计算量的情况下实现了尾涡卷起过程和远场衰减的组合模拟.在系列算例分析研究基础上，创建了基于大涡模拟计算结果的尾流间隔快速预测系统.该系统可以根据实时大气风场和进出港的前后飞机机型，快速预测并输出所需的尾流间隔.经过与场地测试数据比较表明，在北京市2014年的平均风速条件下，本系统预测的尾流间隔可在现有标准基础上缩减7%∼50%，能够有效提高机场容量.

大涡模拟,飞机尾涡,动态尾流间隔,快速预测

引言

自20世纪中叶以来，随着商业客机制造技术的日益成熟，航空运输的重要性与日俱增.机场的起降能力限制了航班数量的进一步增加，因此，机场容量成为制约航空运输发展的瓶颈.在飞机起降过程中，前机的尾流(尾涡)会对后机安全产生威胁，因此，在同一跑道起降的飞机必须保证一个最小安全间隔，在民用航空管理领域称为尾流间隔.尾流间隔限制了机场单位时间起降的航班数量，制约了机场容量的提升.现有的尾流间隔标准有较大的缩减空间.尾流的消散速度与气象条件有关，当前标准是根据不利于尾涡消散的气象条件制定的，因此在多数情况下比较保守.相反，在极个别特殊情况下，气象条件极其不利于尾涡消散，尾涡可能有超常的持续时间，现行标准并不足以保证安全.因此，研究气象条件对尾流消散速率的影响，具有非常重要的意义.随着中国人民生活水平的提高，我国的航空运输业迅速发展，航班数量增加的压力日益增大，航班延误也日趋严重.然而，城市用地的高昂成本制约了机场的增加和扩建.研究尾涡演变特征与规律，建立尾流间隔预测系统，对提高现有机场的容量，缓解航空压力具有非常重要的现实意义.相比美国[1-2]，欧盟[3-4]和日本[5-6]都已经开始研发动态尾流间隔系统，我国这方面的工作起步较晚，正加速相关研究.

图1 飞机尾涡的形成过程示意图Fig.1 The generation process of aircraft wake vortex

飞机尾涡为一对反向旋转的柱涡，其形成过程如图1所示.在飞行过程中，从机翼后缘脱落的涡面，经过近场演变阶段，或卷起阶段[7]，形成尾涡.飞机尾涡常用三个基本参数来描述：初始总环量Γ0，初始涡核半径rc0，以及初始涡间距b0.Γ0与b0可以通过椭圆翼假设近似求得[8-9]

其中M，B为飞机质量和翼展，g为重力加速度.ρa为当地空气密度，VA为飞机飞行空速.此外，根据Γ0与b0，可以定义尾涡演变的特征速度w0和特征时间t0

式中，w0为根据毕奥萨伐尔定律计算所得尾涡在互诱导下的初始下降速度，t0为尾涡以速度w0下降距离b0所需的时间.单个尾涡的环量和切向速度分布可以用Burnham-Hallock(BH)涡模型来近似描述[10]

其中r为空间点到涡心的距离或半径.

尾涡的强度衰减通常分为两个阶段，扩散阶段与快速衰减阶段[11].在扩散阶段，尾涡环量衰减主要依靠涡量的径向扩散，速度较慢；进入快速衰减阶段后，环量衰减速度突然加快.Crow[12]在1970年采用运动学方法分析了一对柱涡在互诱导作用下的最优失稳模态，指出尾涡的长短波不稳定性是尾涡快速衰减的机制之一.尾涡的演变与大气环境密切相关.Proctor等[13]使用大涡模拟方法研究了尾涡在稳定和中性大气中的演变，并给出了以大气湍流耗散率ε和Brunt-Vaisälä(BV)频率N(与温度垂直梯度有关，为大气稳定性参数)为参数的尾涡演变经验模型.Proctor的研究表明，较强的大气湍流耗散率和 BV频率将使尾涡更快进入快速衰减.Holzäpfel等[14]采用数值模拟方法研究了尾涡在不稳定分层下的演变特性，发现尾涡在此情况下的衰减速度远大于中性和稳定状况，这是由于不稳定分层往往对应着较强的大气湍流以及大尺度的上升下降气流.以现有研究为基础，出现了一系列的尾涡演变预测模型.Proctor等在经验模型[13]基础上，发展了TDAWP(TASS driven algorithms for wake prediction)模型[15]，Holzäpfel等[11]在2003年提出了P2P模型(probabilistic two-phase wake vortex decay and transport model)，随后又经过了几次改进[16-17]，现已在部分尾涡快速预测系统中实现了试运行[3-6].但我国尚缺少相关研究和适用于本国机场的快速预测模型.

当尾涡与地面距离小于1.5b0时，尾涡的演变将受到地面的影响，此阶段可称为近地演变阶段[18].此阶段往往为飞机进近的最后阶段，一旦发生尾涡相遇事故，留给飞行员的反应时间很短.因此，尾涡在近地层中的演变预测尤为重要.Robert等[18]使用镜像涡和二次涡模型描述了尾涡的近地轨迹，并给出了一种轨迹预测方法.关于尾涡强度在近地阶段的衰减，Luton等[19]在1996年采用三维直接数值模拟方法研究了尾涡在壁面附近的运动，发现尾涡与它的镜像涡之间也会产生Crow不稳定性的发展，并发现了涡和壁面互相作用产生的二次涡结构.1999年，Proctor等[20]使用三维大涡模拟方法研究了地面对尾涡环量的影响，发现在尾涡距地面0.6b0高度以下时，尾涡环量会在地面影响下急剧下降.Proctor[21]在2000年的数值研究中指出，尾涡进入快速衰减的时间为尾涡到达轨迹最低点的时间之后0.25t0，而且在无量纲大气湍流耗散率较弱时(ε*<0.12)，尾涡的环量衰减速度对大气湍流强度不敏感.Proctor据此给出了一个描述尾涡环量近地衰减的经验模型.虽然尾涡在近地演变中更早进入快速衰减，尾涡的寿命却有可能比在高空更长，这是由于两个尾涡的间距在地面效应下增大，成为两个孤立涡.尾涡之间相互作用的减小阻止了Crow不稳定性的发展，导致孤立涡往往具有较长的寿命.这种长时间存在的孤立涡大大增加了后机遭遇尾流的风险.Holzäpfel等[22]在2007年通过对观测数据的统计，给出了一个尾涡近地衰减模型，提出了尾涡在快速衰减阶段的衰减速度与侧风大小有关.此外，Stenhan与Holzäpfel等在2013至2016年[23-27]的一系列数值研究表明，地面的障碍物会加快尾涡的衰减，并提出了在地面铺设挡板以加速尾涡消散的设想.尾涡近地演变研究的最终目标是构建尾涡近地快速预测模型.近地层中的多种因素，包括剪切、风速、大气湍流和浮力等等都被证实对尾涡演变有影响，如何将这些因素体现到预测模型中需要进一步的研究.研究它们的影响机制和评估这些因素的影响程度，对构建近地快速预测模型有重要的指导意义.

随着计算技术的进步，利用计算流体力学方法进行数值模拟，已成为尾涡研究的重要手段.与场地观测方法相比，数值模拟可以控制变量研究大气参数的影响，能够得到全流场的时空信息，可用于研究尾涡演变的机制.30年来，尾涡数值模拟的技术已经不断发展.在飞机参考系中，尾涡流场是沿空间发展的，尾涡在机翼后缘延伸距离可长达数千米，要求计算域有很长的轴向长度，在现有计算条件下较难实现.因此通常考察一个相对地面静止并垂直于飞行路径的截面内尾涡的时间演变.早期的研究者多使用二维直接数值模拟方法[9,28-30]，这些研究揭示了尾涡移动轨迹在地面效应、剪切、浮力作用下的特点.由于尾涡的长短波不稳定性发生在三维空间中，从1999年起，三维大涡模拟逐渐成为主要的研究方法.

表1列出了部分尾涡数值研究的计算域与网格分辨率信息.可见随着计算技术的发展和计算精度要求的提高，尾涡数值模拟的计算量急剧增大，由于尾涡核心区域剪切强烈，需要较高网格分辨率，而尾涡的移动导致常规网格需要全场加密.此外，为了分辨长、短波不稳定性，尾涡的轴向也需要长计算域和高分辨率，这使得尾涡数值模拟的网格总量非常巨大.因此需要研究更为高效准确的数值模拟方案.

表1 部分尾涡数值模拟研究的计算域与分辨率设置Table 1 The simulationg domain and resolution of part of previous wake vortex numerical simulation

作为尾涡研究成果的实际应用，美国航空航天局在1997年提出了动态尾流间隔的概念[1].动态尾流间隔系统可以根据当前的气象条件和机型实时计算需要的尾流间隔，可以有效地提高机场容量.现在已经进入试运行的动态尾流间隔系统有美国航空航天局设计的 AVOSS[37]，德国宇航中心开发的WSVBS[3-4,38]和日本宇宙航空研究开发机构开发的WVAS[5-6].日本的仿真结果[6]显示，在机场应用WVAS系统，可以在49%的情况下缩减尾流间隔，机场的航班延误可以相应减少20%∼70%.这表明动态尾流间隔系统在提升机场容量方面有良好的应用前景.现有动态尾流间隔系统的基本工作流程大致相同，可以用图2表示：将实时测量的气象参数和前后机的机型参数输入快速预测系统中，计算出所需的尾流间隔，为即将进入下滑道的飞机设置安全间隔.有些动态尾流间隔系统还包括尾涡监测系统，通过实测的尾涡位置来检验预测结果，以进一步保证安全.动态尾流间隔系统的核心是尾流间隔快速预测模型，此模型要求在输入的气象和机型参数基础上快速预测出尾涡的演变，并确定需要的尾流间隔.虽然使用计算流体力学方法可以预测尾涡演变过程，但是现有计算技术和计算资源远远达不到实时预测的要求，因此现有系统多采用经验模型预测.目前用于高空演变快速预测的模型较为成熟；关于尾涡的近地快速预测，除了上文提到的经验模型如Proctor模型[21]和Holzäpfel模型[22]外，还有Visscher等[39]2013年提出的DVM(deterministic wake vortex model)模型.由于近地大气运动复杂，构建尾涡的近地快速预测模型较为困难.近地的温度分层、侧风、逆风以及它们的剪切率都对尾涡演变有影响，经验模型却只能考虑有限的参数.例如P2P模型[22]中考虑了10米高度侧风大小，而TDAWP模型[13]考虑了侧风剖面的一阶导数和二阶导数.考虑参数过少将降低预测精度，因此必须通过增大安全余量的方法来应对，这会限制动态尾流间隔系统的效能.而引入更多的参数，例如逆风及其剪切，将使经验模型变得过于复杂而难以构建.因此，尾涡近地预测模型的建立方法还需要进一步研究.

图2 动态尾流间隔系统工作原理示意图Fig.2 The overview of dynamic wake separation system

在我国，吴永华和胡以华等对激光雷达尾涡测量方法进行了相关研究[40-43].关于飞机遭遇尾流的安全问题，安崇君等[44]通过对一次尾涡事故中黑匣子数据的分析，分析了尾涡对后机的滚转力矩和滚转率的影响.何德富[45]研究了飞机横穿尾涡和从正后方进入尾涡两种情况下所受的安全威胁，论证了遵守一定尾流安全间隔的必要性.韩红蓉等[46]从飞行控制角度分析了飞机进入尾涡后的响应，综合分析了阻尼特性、反应时间和操纵品质等，建立了飞机的安全分析模型并预测尾流间隔.魏志强等[47-48]以美国航空航天局的TDAWP模型[13]为基础研究了动态尾流间隔计算方法.相对于欧美等发达国家，我国在飞机尾涡演变领域的研究起步较晚，缺少原创性工作.因此，开展相关研究，建立我国自己的动态尾流间隔预测系统，对提高我国机场容量，减少航班延误，提高我国民用航空的安全和管理水平，具有重大的意义.

基于飞机尾涡演变的研究现状，本文主要开展了以下研究工作：首先发展了飞机尾涡的数值模型，以提高数值模拟的效率和准确性.将自适应网格技术引入飞机尾涡的数值模拟以大幅减少计算网格，并采用升力面法，模拟尾涡卷起过程对远场衰减的影响，并通过与观测实验的对比，验证了所建数值模型的准确性.在此基础上，应用大量大涡模拟研究成果，创建了近地尾流间隔快速预测系统.此系统可根据大气参数快速预测前机的尾流影响区域，以确定后机所需的安全间隔.全文内容安排如下：第二章介绍尾涡演变的数学模型和数值方法，并通过与现有文献和实验观测数据对比，验证模型的准确性；第三章简要介绍近地尾涡快速预测系统算例结果；主要结论和工作展望将在第四章表述.

1 飞机尾涡数值研究方法及验证

1.1 大涡模拟控制方程

飞机尾涡流场的最大速度通常不超过30m/s，属于低速流动，可应用不可压缩Navie-Stokes(N-S)方程组求解尾涡流场.对于温差产生的浮力效应，本文采用Boussinesq假设，通过动量方程中的浮力项体现.综上所述，本文求解N-S方程过滤后的大涡模拟方程组为

其中xi(i=1,2,3)或x,y,z表示3个方向的坐标，ui或u,v,w代表对应的速度分量;p,θ表示压强和温度，θ0为参考温度;ρ,ν,κ分别表示空气的密度、运动黏度和热扩散率，fi表示流体受到的体积力分量;顶标横线表示过滤后的的物理量;为亚格子应力，而为亚格子热通量.本文使用Meneveau等[49]提出的沿流体迹线平均的Lagrange动力模式来封闭方程组.

1.2 自适应网格与离散方法

本文在动网格上离散和求解控制方程组，首先介绍本文使用的自适应网格算法.本文将自适应网格应用在尾涡数值模拟中，在尾涡中心区域使用较密的网格以保障计算精度，在远离尾涡中心的区域采用较稀疏的网格以减少计算量(图3).随着尾涡的移动，网格的加密区域通过自适应算法做相应的调整.

图3 尾涡中心附近的局部加密Fig.3 The fi ner grid space near vortex cores

本文的计算程序ATTLES在TTLES程序的基础上发展而来.TTLES程序最初用于城市大气环境大涡模拟[50-52].本文在TTLES基础上增加了自适应网格模块，并应用直角正交网格.为了使网格的加密区域始终跟随尾涡中心，各网格平面需要沿其法线方向运动，而网格运动速度需要通过自适应算法来确定.本文所用的网格自适应算法在Gno ff o[53]1983年提出的弹簧比拟法基础上进行改进.Gno ff o的自适应网格算法适用于定常流动，基本思想是假定网格节点间由弹簧相连，当地的流场梯度越大，则弹簧刚度越大.当系统达到受力平衡时，即得到了与流场相适应的网格，本文将Gno ff o的方法推广到了非定常流动.图4显示了本文采用的自适应网格算法原理，以x方向为例，Ki+1/2表示连接编号为i与i+1的网格面的弹簧刚度，(ug)i表示网格面i的运动速度.首先，根据流场梯度计算弹簧刚度值

F(f)为一定义域和值域均为[0,1]的增函数，A为大于1的常数，代表了x方向最大网格与最小网格的尺寸比，在本文中取A=7.fi为网格面i上标准化的流场梯度，取值在0∼1之间

ai为网格面i内的最大流场梯度，可由下式计算得到

amin与amax分别为ai(i=1,2,···,Nx)的最小值与最大值，Nx为x方向的网格数量.关于式(7)中F(f)的形式，本文根据文献[54]取为

其中B为一可调节正数.计算出Ki的值后，可以列出网格面的静力平衡方程组

在本文的计算中，最小网格间距出现在尾涡中心区域，为了保证尾涡中心的分辨率，∆xmin应当预先给定.据此，联立式(7)、式(10)和式(12)可迭代解得B的值.最后，求解三对角方程组(11)得到目标网格位置本文令网格在Nu个时间步内匀速运动到目标网格.网格速度由下式确定

应用此方法时，网格并不能准确适应当前的流场，而是存在时长为Nu∆t的滞后，为了保证计算准确性，Nu∆t应尽可能小.然而Nu∆t过小会导致网格速度过大，影响计算的稳定性.在本文的计算中，Nu∆t取为0.1s.至此，x方向的网格面移动速度已经确定，可以利用同样的方法确定y，z方向的网格移动速度(vg)i,(wg)i.

在动网格上，本文采用有限体积方法对控制方程(4)∼(6)进行离散，使用动量插值的SIMPLE方法求解控制方程.在时间方向，本文应用四阶龙格库塔积分.为实现大规模计算，本文计算程序采用分区并行方法.

1.3 尾涡生成模型

在国内外现有尾涡数值模拟的研究论文中绝大多数直接利用一对尾涡的模型，如前文式(1)∼式(3)的BH模型，不考虑尾涡从涡片卷起过程.本文采用升力面模型来初始化尾涡.图5(a)显示了升力面模型.x，y，z分别为飞机的飞行方向(即尾涡轴向)，翼展方向和垂直方向.在一个厚为h，宽为飞机翼展B，x方向无限长的矩形区域内，施加一个椭圆分布的，竖直向下的，作用时长为τ体积力，其形式为

图5 升力面模型和尾涡模拟流程Fig.5 The lift-drag model and the fl ow chart of wake vortex simulation

此体积力的物理意义为飞机机翼对空气的反作用力.与第一章中尾涡的初始环量和间距的确定过程类似，式(14)同样应用了椭圆翼假设.h和τ应尽可能小，为了兼顾数值稳定性，本文中取h=0.014b0，τ=0.003b0.图5(b)显示用升力面法的尾涡数值模拟过程：在力fz作用下，计算域中将形成一个涡层；此涡层在自诱导下卷起，最终形成一对尾涡；此后尾涡进入远场衰减阶段.尾涡的卷起过程在图6中显示，其中云图为使用特征时间t0无量纲化的轴向涡量分量.为了合理模拟飞行阻力的影响，沿飞行方向的力fx也被施加在这个矩形区域中，其大小为

其中LD是飞机的升阻比，在本文中，LD根据Keye[55]的研究结果取为12，模拟客机进场的情形.本节介绍的升力面尾涡初始化模型可以模拟尾涡的卷起过程，进而可以在模拟中反映卷起过程对远场衰减的影响.此外，此模型并不在传统涡模型法基础上增加计算量.它将飞机对空气的作用简化为平面分布的升力和阻力，因此并不需要模拟整机绕流[36]，尾涡的卷起阶段和远场衰减可在同一算例中完成.

图6 尾涡生成过程的轴向涡量分布Fig.6 The vortex generation process,shown by the contour of axial vorticity

1.4 模型验证

自适应网格方法虽然能有效减少计算网格总量，但是，不均匀网格的存在和网格的移动均有可能影响计算的数值精度，本节通过5个不同气象条件下的算例来检验自适应网格法的准确性.在本节算例中，飞机参数取为Γ0=446m2/s，b0=47.4m以模拟一架空中客车A330-300的尾涡.根据式(2)计算特征速度和时间分别为w0=1.50m/s和t0=31.7s.空气的运动黏度系数取ν=1.8×10−5m2/s，使用环量定义的雷诺数为ReΓ=Γ0/ν≈2.5×107.

表2列出了验证算例的气象参数.在尾涡研究中，为扩大研究结果的应用范围，可将背景大气参数(如大气湍流耗散率ε，BV 频率N=[g/θ0(dθ/dz)]1/2等)用机型的特征参数无量纲化

而环量、长度、速度和时间变量分别用特征量Γ0,b0,w0和t0无量纲化.现有研究表明，采用这种无量纲化的计算和分析方法，不同机型，不同初始强度尾涡的演变具有相似性[11,21]，可以有效扩展研究结果的适用范围.如无特殊说明，后文将沿用此无量纲方案.本节选取无量纲耗散率为0.01,0.05和0.23三个算例来验证自适应网格在不同背景湍流下的表现.此外，本文还模拟了两个稳定分层条件下的算例S06和S10.本文选取了无量纲BV频率为0.6和1.0的两种情况.算例计算域取为Lx×Ly×Lz=8.0b0×6.3b0×6.3b0，尾涡轴向计算域取为8b0以分辨长波不稳定性.计算网格数为Nx×Ny×Nz=380×144×120.在尾涡轴向采用了均匀网格间距∆X=0.021b0(1m),在垂直于尾涡轴线的平面内用自适应网格算法，涡心附近分辨率为∆y=∆z=0.010b0(0.47m).远离涡心区域分辨率为0.07b0(3.3m).网格分辨率在轴向与文献[33]相同，而在垂直于涡轴平面内，涡心附近分辨率则超过了文献[33](0.021b0).此算例网格总数为660万左右，约为文献[33]的11%，可见采用自适应网格法可以大幅减少计算量.

表2 自适应网格验证算例中的背景气象参数Table 2 The ambient atmospheric parameters in the validation of sela adaptive grid

对于中性算例(N01∼N23)，计算域的3个方向均为周期边界条件.对于稳定算例，由于垂直温度梯度的存在，在计算域的顶部和底部边界采用零梯度边界条件，而其他边界仍应用周期边界条件.计算中取固定时间步长∆t=t0/5000(1/150s).

本节算例的初始流场由 3个部分叠加而成.背景湍流场，背景温度场及尾涡场.本文采用了Rogallo[56]提出的方法生成了满足目标能谱的均匀湍流场，目标能谱为文献[57]提出的改进的冯·卡门谱

其中kp为能谱峰值波数，本文取为2π/90m−1，与对照算例相同[33].kkol=(ε/ν3)1/4是 Kolmogorov 波数.K0可以由目标湍动能耗散率根据下式确定

算例的初始温度场取为

其中θ0=290K为参考温度，温度的垂直梯度可根据BV频率求得

图7通过λ2=−0.1的等值面[58]显示了算例N05和S10网格适应流场的情况，可以看到，随着尾涡演变，网格的加密部分都可以很好地跟随尾涡中心移动，并且加密区域面积逐渐扩大以适应尾涡形状变化.

图7 网格对流场演变的自适应过程Fig.7 The adaption of the grid to the fl ow fi eld

图8(a)显示了全部5个算例尾涡环量随时间衰减的情况，尾涡环量采用其初始值(t=0)无量纲化.高度使用b0无量纲化而时间使用t0无量纲化.图中可以清楚地看到尾涡环量的两阶段衰减特性.在扩散阶段，所有算例的环量的衰减速度几乎相同.用表示尾涡进入第二阶段衰减(快速衰减)的无量纲时间，与气象条件有关.对比算例N01(ε*=0.01,N*=0.0)，N05(ε*=0.05,N*=0.0)和 N23(ε*=0.23,N*=0.0)可以看出，随着湍流耗散率的增大而减小，表明耗散率的增大将使尾涡更快进入快速衰减.对比算例 N05(ε*=0.05,N*=0.0)，S06(ε*=0.05,N*=0.6)和 S10(ε*=0.05,N*=1.0)可以看出，t2*随着BV频率的增大而减小.可见温度分层对尾涡的快速衰减也有促进作用.

图8 尾涡环量，涡心高度随时间的演化Fig.8 The temporal evolution of circulation and vortex altitude

图8(b)显示了算例 N01(ε*=0.01,N*=0.0),N05(ε*=0.05,N*=0.0)和 S10(ε*=0.05,N*=1.0)的结果.对比文献[33]在2011年相同气象条件下的固定网格大涡模拟结果可以看出，在初始的扩散阶段，本文结果的衰减速度略小于固定网格结果，这是由于文献[33]使用了耗散更强的标准Smagorinsky亚格子模式.此外，本文的涡核分辨率(0.01b0)高于文献[33]的分辨率(0.021b0)，已有研究表明较低的网格分辨率会导致涡核附近的数值伪扩散，从而加快扩散阶段的衰减速度，但是对尾涡进入快速衰减的时间影响不大[34].关于进入快速衰减阶段的时间t2*，本文自适应网格方法与文献[33]得到的结果几乎相同.图9显示了本文得到的t2*与Sapkaya模型[59]及P2P模型[11]的对比.Sapkaya于2000年提出了中性大气中t2*的经验模型，由于Sapkaya模型是基于长波失稳下尾涡连接的理论提出的，而在算例 N23中尾涡并未发生连接 (后文将详细讨论)，因此算例N23的结果与 Sapkaya模型有一定偏差.Holzäpfel于2003年通过对大涡模拟结果的拟合得出P2P模型[11]，给出了BV频率对的影响，可看出本文的结果与P2P模型基本符合.综上所述，虽然本文使用的网格量只有文献[33]的11%，对环量衰减的模拟精度却达到了几乎相同的效果.

图9 尾涡快速衰减起始时间与气象条件关系Fig.9 The relationship between the onset time of rapid decayand the atmospheric condition

图8(c)中显示了尾涡高度下降的过程.在演变的初始阶段，尾涡在自诱导下以速度w0下降.在算例 N01(ε*=0.01,N*=0.0)和 N05(ε*=0.05,N*=0.0)中，尾涡的下降速度随时间变化不大，而其他算例中尾涡的下降速度随时间减小.在算例N23(ε*=0.23,N*=0.0)中，尾涡强度的衰减导致了互诱导效应的减弱，进而导致了下降速度的减慢.在算例S06(ε*=0.05,N*=0.6)和 S10(ε*=0.05,N*=1.0)中，除了尾涡强度衰减的因素外，浮力效应[28]也是引起尾涡下降速度减慢的重要原因.

图10通过λ2=−0.1的等值面显示出了尾涡形态的演变过程.算例N05是一个典型的尾涡连接失稳的例子(图10(a)).可以看到，尾涡经历了正弦长波失稳过程，在t*=5.7时刻，两个尾涡发生连接，导致尾涡的快速衰减，在连接点附近产生了螺旋结构[34,42].图8(a)对照可以看出，尾涡快速衰减开始的时间与连接时间吻合.涡环形成后，继续在自诱导下发生形变 (图 10(a)，t*=8.0).而算例 S10(图10(b))则表现出短波失稳的特性.在快速衰减的开始时刻t*=1.5，可以看到大量肋状二次涡结构.这些结构是短波不稳定发展的产物[28]，进而促使尾涡更早地进入快速衰减[35].图11显示了算例N23中尾涡的衰减过程.在较强大气湍流作用下，尾涡的长短波不稳定性都快速增长.由于尾涡强度的快速衰减，尾涡在连接成涡环之前已经消散，因此未能观察到涡环形成.

图10 λ2=−0.1等值面.算例N05和算例S10Fig.10 Iso-surface of λ2= −0.1.Case N05 and Case S10

图11 λ2=−0.1等值面，算例N23Fig.11 Iso-surface of λ2= −0.1.Case N23

Crow指出，尾涡长波正弦失稳的振幅将按指数规律增长[12]

其增长率α≈0.81.本文通过下列方法计算了算例N01和N05的长波失稳振幅，以考察其增长率

其中Ry和Rz分别为水平方向和垂直方向尾涡中心线偏离原平衡位置的最大距离，如图12所示.长波失稳的波长λ取为轴向计算域长度8b0.RA随时间增长如图13所示，可以看到，N01和N05算例中均能观察到很长的指数增长阶段，增长率分别为0.74与0.71，略小于Crow给出的理论值.由于Crow的结论是在理想流体假设下得出的[12]，而真实情况下尾涡环量将逐渐衰减并导致尾涡的互诱导作用减弱，因此增长率略低于理论值是合理的.

图12 尾涡振幅计算示意图Fig.12 The calculation of wake vortex amplitude

图13 长波失稳振幅随时间的增长Fig.13 The growth of long wave instability

以上算例结果证明了自适应网格方法应用于尾涡演变的数值模拟是可行的.与固定均匀网格大涡模拟相比，应用自适应网格大幅减少了计算量，且得到同样可靠的结果：自适应网格方法可以准确模拟出尾涡的两阶段衰减以及尾涡演变中的长短波失稳过程，可以分辨出尾涡演变中产生的二次结构，长波失稳的增长率与理论值吻合.

1.5 激光雷达观测验证

为了进一步验证数值方法的可靠性，本节模拟了一个实际观测算例.Idaho Falls B-757 Run 9(IDF9)是美国联邦航空管理局(FAA)主持的一系列尾涡观测实验之一[60]，该观测实验使用激光雷达测量尾涡在近地阶段的演变，主要目标是评估波音757和767飞机的尾涡特性.本文选取这个实例进行验证，是由于此例提供了详细的信息，包括实时气象条件以及尾涡环量与位置演变数据.表3给出了IDF9算例的飞机参数，算例的环量雷诺数为ReΓ=Γ0/ν≈2.0×107.大气的实时侧风剖面和计算采用的剖面在图14中给出.可见本文生成的风廓线与实测结果基本一致.

图14 IDF9的实测风廓线与本文计算用的背景风廓线Fig.14 The real time wind pro fi le and the ambient wind fi eld used in this paper,Case IDF9

表3 IDF9算例的机型参数Table 3 The aircraft parameters in case IDF9

本文所提模型的计算域为Lx×Ly×Lz=240m×300m×105m，计算网格为366×326×120.在yz平面内，涡心附近网格分辨率为0.0074b0(0.22m).翼展方向和飞行方向采用周期边界条件，地面应用无滑移条件，而计算域顶部边界应用零梯度条件.

图15(a)显示了尾涡强度的时间演变.由图15可看出，数值模拟的结果与实际观测数据符合良好.其中上游涡在70s左右进入快速衰减，这与尾涡高度到达最低点的时间一致(图15(b))，可见此快速衰减是由地面效应导致的.下游涡的衰减速度明显快于上游涡，这种不对称演变与平均侧风剪切有关，将在下一节中详细论证.图15(b)和图15(c)显示了尾涡的位置随时间的变化.其中尾涡高度的模拟结果与实测数据有一定偏差，这是由于尾涡高度对风剖面的形状非常敏感[61-62].尾涡的水平位置和实测符合较好.可见，本文数值模型可以较为准确地预测尾涡的强度和水平位置，而对尾涡高度的预测精度略低.总体看来，本文提出的尾涡数值方法，可以用于尾涡近地快速预测模型的建立.

图15 IDF9算例，尾涡的时间演变Fig.15 The temporal evolution of wake vortex,Case IDF9

2 尾涡近地演变快速预测

图16 基于尾涡演变数值模拟结果的尾流间隔快速预测系统Fig.16 The fast time wake separation prediction system based on the LES database

基于上节研究成果，创建了基于大涡模拟的近地尾流间隔快速预测系统.图16示意尾流间隔快速预测系统的基本工作流程.系统的核心是尾涡演变的大量数值模拟结果，包括预先计算的一系列不同气象条件下的尾涡演变算例结果.系统运行步骤如下：首先，输入气象参数和机型参数，通过气象数据前处理，计算前机(产生尾涡的飞机)的特征参数(特征速度w0和尾涡初始间距b0)，并利用前机特征参数无量纲化气象数据；然后计算尾流间隔；最后，将预测结果可视化输出.

根据预测得到的尾涡演变特征，可以作出飞机后方的尾流影响区域.根据得到的尾流危险区最终确定所需的最小进场间隔.图17显示了最小安全间隔的确定方法.图17(a)显示了危险区在侧风作用下离开后机下滑道的情况，图17(b)则显示了后机由于尾涡强度衰减而脱离危险区的情况.

图17 安全走廊与安全间隔的确定示意图Fig.17 The safety corridor and the determination of safety separation

图18给出了一个尾涡间隔快速预测系统的运行示例，系统预测了一架空中客车A380客机作为前机时，后方的尾流危险区以及建议的最小尾流间隔.这里没有选取特定的后机，而是同时输出了后机为波音 B747-400(Heavy，重型)、波音 B737-800(Medium，中型)和塞斯纳C650(Light，轻型)时的尾流间隔结果.输入的机型和大气参数如表4所示.输入大气参数波动范围取为平均值的正负34%(可进一步考核与优化),最终获得的安全间隔分别为：3.90km(后机为重型)，4.92km(后机为中性)和 5.65km(后机为轻型).

图18 示例：系统可视化输出的空中客车A380后方的尾流危险区及最小尾流间隔Fig.18 Example:The visualization of the prediction result.The dangerous area and the minimum wake separation

表4 示例：计算空中客车A380后方飞机的最小尾流间隔.机型参数及气象参数Table 4 The calculation of minimum wake separation when approaching after an A380.The aircraft parameters and the atmospheric parameters

表5是另外5个预测示例.根据《北京统计年鉴2015》[63]，北京市 2014年的平均风速为2.1m/s(约4kn).表中示例1∼4为10m高全风速为4kn的情况，由于我国现行规定雷达间隔为6km，当预测的尾流间隔小于6km时，需按6km执行.从表中的预测结果可见，与现行规定相比，在北京市全年平均风速条件下，使用本系统在绝大多数情况下可以缩减尾流间隔.尤其是在侧风情况下(示例1,3)，中型机和轻型机跟随重型机进场的间隔可以比现有标准缩减40%∼50%.若没有雷达间隔的限制，缩减比例还可以进一步提高.可见，如果此系统得到应用，对机场容量的提升作用将非常明显.示例2比较特殊，当前机为A380，后机为轻型机时，预测的尾流间隔比现行标准大4%(用负值表示)，可见随着超重新机型(A380)的出现，当前标准在特殊气象条件下不足以保证安全，需要做一定的补充或修正.示例5是一个较为极端的天气，10m高空全风速仅为1kn(0.514444m/s)，几乎为静风条件.后机为中型机时，系统预测间隔大于现行标准.说明本系统还可以在极端不利于尾涡消散的天气条件下提供预警，以保障安全.

表5 本文间隔预测系统运行示例Table 5 Some of the predictions made by the current system

3 总结

本文采用大涡模拟方法研究了飞机尾涡的演变特性，并在此基础上创建了近地尾流间隔快速预测系统.其主要贡献如下：

(1)发展了飞机尾涡的数值预测模型.首先将自适应网格技术引入尾涡大涡模拟计算，有效减少了网格总量，提高了计算效率.同时提出了升力面尾涡生成模型，在不增加计算量的前提下，引入尾涡卷起过程的计算，进而体现了卷起过程对远场衰减的影响.

(2)创建了基于大涡模拟的尾流间隔快速预测系统.此系统可通过输入大气和机型参数，计算并输出所需的最小尾流间隔.经过场地观测数据评估，在北京市全年平均风速条件下，使用本系统在绝大多数情况下可以缩减尾流间隔.在部分算例中，间隔缩减比例可达40%∼50%.

进一步工作需要结合大量的场地实际观测数据，精细考核、优化和完善快速尾涡预测平台，建立我国可实际应用的动态尾流间隔系统.

1 Perry RB,Hinton DA,Stuever RA.NASA wake vortex research for aircraft spacing.AIAA Paper,1997,57

2 O’Connor CJ,Rutishauser DK.Enhanced airport capacity through safe,dynamic reductions in aircraft separation:NASA’s Aircraft Vortex Spacing System(AVOSS).NASA,2001-211052

3 Holzäpfel F,Gerz T,Frech M,et al.The wake vortex prediction and monitoring system WSVBS Part I:Design.Air Traffic Control Quarterly,2009,17(4):301

4 Gerz T,Holzäpfel F,Gerling W,et al.The wake vortex prediction and monitoring system WSVBS Part II:Performance and ATC integration at Frankfurt airport.Air Traffic Control Quarterly,2009,17(4):323-346

5 Matayoshi N.Dynamic wake vortex separation according to weather Condition//2013 Aviation Technology,Integration,and Operations Conference.2013:4425

6 Matayoshi N,Yoshikawa E.Dynamic wake vortex separation combining with AMAN/DMAN concept//15th AIAA Aviation Technology,Integration,and Operations Conference.2015:3397

7 Gerz T,Holzäpfel F,Darracq D.Commercial aircraft wake vortices.Progress in Aerospace Sciences,2002,38(3):181-208

8 KaramchetiK.PrinciplesofIdeal- fl uidAerodynamics.KriegerPublishing Company,1980

9 Proctor FH.Numerical simulation of wake vortices measured during the Idaho Falls and Memphis field programs.AIAA Paper,1996,96-2496

10 Burnham DC,Hallock JN.Chicago monostatic acoustic vortex sensing system.Volume IV.Wake vortex decay.Department of Transportation,Rep.No.DOT/FAA/RD-79-103 IV,July 1982,206 pp.[Available from NTIS]

11 Holzäpfel F.Probabilistic two-phase wake vortex decay and transport model.Journal of Aircraft,2003,40(2):323-331

12 Crow SC.Stability theory for a pair of trailing vortices.AIAA Journal,1970,8(12):2172-2179

13 Proctor FH,Switzer GF.Numerical simulation of aircraft trailing Vortices//The Ninth Conference on Aviation,Range and Aerospace Meteorology,11–15 September 2000,Paper 7.12

14 Holzäpfel F,Gerz T,et al.Wake vortices in convective boundary layer and their in fl uence on following aircraft.Journal of Aircraft,2000,37(6):1001-1007

15 Proctor FH,Hamilton DW,Switzer GF.TASS driven algorithms for wake prediction.44th AIAA Aerospace Sciences Meeting and Exhibit,AIAA.2006,1073

16 Holzäpfel F.Probabilistic two-phase aircraft wake-vortex model:further development and assessment.Journal of Aircraft,2006,43(3):700-708

17 Frech M,Holzäpfel F.Skill of an aircraft wake-vortex model using weather prediction and observation.Journal of Aircraft,2008,45(2):461-470

18 Robins RE,Delisi DP,Greene GC.Algorithm for prediction of trailing vortex evolution.Journal of Aircraft,2001,38(5):911-917

19 Luton JA,Ragab SA.The three-dimensional interaction of a vortex pair with a wall.Physics of Fluids,1997,9(10):2967-2980

20 Proctor FH,Han J.Numerical study of wake vortex interaction with the ground using the terminal area simulation system//37th Aerospace Sciences Meeting&Exhibit,1999,AIAA-99-0754

21 Proctor FH,Hamilton DW,Han J.Wake vortex transport and decay in ground e ff ect:Vortex linking with the ground.AIAA Paper,2000,2000-0757

22 Holzäpfel F,Steen M.Aircraft wake-vortex evolution in ground proximity:analysis and parameterization.AIAA Journal,2007,45(1):218-227

23 Stephan A,Holzäpfel F,Misaka T.Aircraft wake-vortex decay in ground proximity—Physical mechanisms and arti fi cial enhancement.Journal of Aircraft,2013,50(4):1250-1260

24 Stephan A,Holzäpfel F,Misaka T.Hybrid simulation of wakevortex evolution during landing on fl at terrain and with plate line.International Journal of Heat and Fluid Flow,2014,49:18-27

25 Holzäpfel F,Stephan A,Körner S,et al.Wake vortex evolution during approach and landing with and without plate lines.AIAA Paper,2014,925:2014

26 Holzäpfel F,Stephan A,Tchipev N,et al.Impact of wind and obstaclesonwakevortexevolutioningroundproximity//6thAIAAAtmospheric and Space Environments Conference,Atlanta(June 2014).2014

27 Holzäpfel F,Stephan A,Heel T,et al.Enhanced wake vortex decay in ground proximity triggered by plate lines.Aircraft Engineering and Aerospace Technology:An International Journal,2016,88(2):206-214

28 Proctor FH,Hinton DA,Han J,et al.Two dimensional wake vortex simulations in the atmosphere:Preliminary sensitivity studies.AIAA Paper,1997,97-0056

29 Robins RE,Delisi DP.Numerical study of vertical shear and strati fication e ff ects on the evolution of a vortex pair.AIAA Journal,1990,28(4):661-669

30 Corjon A,Poinsot T.Behavior of wake vortices near ground.AIAA Journal,1997,35(5):849-855

31 Shen SH,Ding F,Han J,et al.Numerical modeling studies of wake vortices:real case simulations//37th Aerospace Sciences Meeting&Exhibit,11-14 January,1999,Reno,NV,AIAA Paper No.99-0755

32 Holzäpfel F,Gerz T,Baumann R.The turbulent decay of trailing vortex pairs in stably strati fi ed environments.Aerospace Science and Technology,2001,5(2):95-108.

33 Hennemann I,Holzäpfel F.Large-eddy simulation of aircraft wake vortex deformation and topology//Proceedings of the Institution of Mechanical Engineers,Part G:Journal of Aerospace Engineering,2011:0954410011402257

34 Misaka T,Holzäpfel F,Hennemann I,et al.Vortex bursting and tracer transport of a counter-rotating vortex pair.Physics of Fluids,2012,24(2):025104

35 De Visscher I,Bricteux L,Winckelmans G.Aircraft vortices in stably strati fi ed and weakly turbulent atmospheres:Simulation and modeling.AIAA Journal,2013,51(3):551-566

36 Misaka T,Holzäpfel F,Gerz T.Large-eddy simulation of aircraft wake evolution from roll-up until vortex decay.AIAA Journal,2015,53(9):2646-2670

37 AhmadNN,VanValkenburgRL,PruisM.NASAAVOSSFast-Time Wake prediction models:User’s guide.National Aeronautics and Space Administration,NASA/TM-2014-218152,2014

38 Gerz T,Holzäpfel F,Bryant W,et al.Research towards a wakevortex advisory system for optimal aircraft spacing.Comptes Rendus Physique,2005,6(4):501-523

39 De Visscher I,Lon fi ls T,Winckelmans G.Fast-time modeling of ground e ff ects on wake vortex transport and decay.Journal of Aircraft,2013,50(5):1514-1525

40 吴永华,胡以华,戴定川等.基于 1.5µm多普勒激光雷达的飞机尾涡探测技术研究.光子学报,2011,40(6):811-817(Wu Yonghua,Hu Yihua,Dai Dingchuan,et al.Research on the technique of aircraft wake vortex detection based on 1.5µm doppler lidar.Acta Photonica Sinica,2011,40(6):811-817(in Chinese))

41 胡以华,吴永华.飞机尾涡特性分析与激光探测技术研究.红外与激光工程,2011,40(6):1063-1069(Hu Yihua,Wu Yonghua.Research on aircraft evolution and the Lidar measurement technology.Infrared and Laser Engineerin,2011,40(6):1063-1069(in Chinese))

42 Moet H,Laporte F,Chevalier G,et al.Wave propagation in vortices and vortex bursting.Physics of Fluids,2005,17(5):054109

43 徐世龙,胡以华,赵楠翔.基于激光雷达回波的飞机尾涡参量提取.光子学报,2013,42(1):54-58(Xu Shilong,Hu Yihua,Zhao Nanxiang.Extrication of wake vortex parameters based on lidar echo.Acta Photonica Sinica,2013,42(1):54-58(in Chinese))

44 安崇君,刘长卫.飞机尾涡对飞行安全的影响.飞行力学,2000,18(4):69-72(An Congjun,Liu Changwei.The in fl uence of aircraft wake vortex on air safety.Flight Dynamics,2000,18(4):69-72(in Chinese))

45 何德富.飞机翼尖尾涡对后面飞机飞行安全影响及安全措施.中国民航飞行学院学报,2005,16(1):12-19(He Defu.Safety e ff ect and precautions for fl owing aircraft by leading aircraft wake vortices.Journal of Civil Aviation Flight University of China,2005,16(1):12-19(in Chinese))

46 韩红蓉,李娜,魏志强.飞机遭遇尾涡的安全性分析.交通运输工程学报,2012,12(1):45-49(Han Hongrong,Li Na,Wei Zhiqiang.Safety analysis of aircraft encountering wake vortex.Journal of Traffic and Transportation Engineering,2012,12(1):45-49(in Chi-nese))

47 魏志强，徐肖豪.飞机尾涡流场的建模与仿真计算研究.交通运输系统工程与信息,2010,10(4):186-191(Wei Zhiqiang,Xu Xiaohao.Modeling and simulating of fl ow fi eld for aircraft wake vortex.Journal of Traffic and Transportation Engineering,2010,10(4):186-191(in Chinese))

48 魏志强.动态尾涡间隔计算工具的开发与实现.航空计算技术,2010,40(5):87-90(Wei Zhiqiang.Development and implementation of dynamic calculating tool for spacing of wake vortex.Aeronautical Computing Technique,2010,40(5):87-90(in Chinese))

49 Meneveau C,Lund TS,Cabot WH.A Lagrangian dynamic subgridscale model of turbulence.Journal of Fluid Mechanics,1996,319:353-385

50 Shi RF,Cui GX,Wang ZS,et al.Large eddy simulation of wind fi eld and plume dispersion in building array.Atmospheric Environment,2008,42(6):1083-1097

51 Liu YS,Cui GX,Wang ZS,et al.Large eddy simulation of wind field and pollutant dispersion in downtown Macao.Atmospheric Environment,2011,45(17):2849-2859

52 Wang BB,Wang ZS,Cui GX,et al.Study on the dynamic characteristics of fl ow over building cluster at high Reynolds number by large eddy simulation.Science China Physics,Mechanics&Astronomy,2014,57(6):1144-1159

53 Gno ff o PA.A fi nite-volume,adaptive grid algorithm applied to planetary entry fl ow fi elds.AIAA Journal,1983,21(9):1249-1254

54 Nakahashi K,Deiwert GS.Self-adaptive-grid method with application to airfoil fl ow.AIAA Journal,1987,25(4):513-520

55 Keye S.Fluid-structure coupled analysis of a transport aircraft and flight-test validation.Journal of Aircraft,2011,48(2):381-390

56 Rogallo RS.Numerical experiments in homogeneous turbulence.National Aeronautics and Space Administration,1981.TM-81315

57 Bechara W,Bailly C,Lafon P,et al.Stochastic approach to noise modeling for free turbulent fl ows.AIAA Journal,1994,32(3):455-463

58 Jeong J,Hussain F.On the identi fi cation of a vortex.Journal of Fluid Mechanics,1995,285:69-94

59 Sarpkaya T.New model for vortex decay in the atmosphere.Journal of Aircraft,2000,37(1):53-61

60 Garodz LJ,Clawson KL.Vortex wake characteristics of B757-200 and B767-200 aircraft using the tower fl y by technique.NOAA Technical Memorandum ERL-ARL-199,1993

61 Proctor F H.The NASA-Langley wake vortex modelling e ff ort in support of an operational aircraft spacing system.36th Aerospace Sciences Meeting&Exhibit,January 12-15,1998,Reno,NV.AIAA 98-0589

62 Holzäpfel F.E ff ects of environmental and aircraft parameters on wake vortex behavior.Journal of Aircraft,2014,51(5):1490-1500

63 北京市统计局.2015北京统计年鉴.北京：中国统计出版社，2015(Beijing Statistics Bureau.Beijing Statistical Yearbook 2015.Beijing：China Statistics Press,2015(in Chinese))

LARGE EDDY SIMULATION ON THE EVOLUTION AND THE FAST-TIME PREDICTION OF AIRCRAFT WAKE VORTICES1)

Lin Mengda Cui Guixiang2)Zhang Zhaoshun Xu Chunxiao Huang Weixi
(School of Aerospace，Tsinghua University，Beijing100841，China)

As the commercial air transport increasing in China,the problem of fl ight delays is increasingly serious.Wake separation(the minimon separation between the leading and following aircraft to avoid wake encounter)limits the capacity of airports.Aimed at this problem,this paper study the evolution of wake vortex with large eddy simulation(LES).The self-adaptive grid method is applied to the LES of wake vortex to improve the computation efficiency and a lift-drag model is applied to the wake vortex generation process to simulate the roll up phase.Based on the LES wake vortex evolution database,a fast-time wake separation prediction system is established.Given the real time ambient wind field and the aircraft parameters,the prediction system can output the suggested wake separation.The results show that under the average wind condition in Beijing Capital International Airport in 2014,the current wake separation can be reduced by 7%∼50%with the established system,thus the airport capacity can be considerably improved.

large eddy simulation,aircraft wake vortex,dynamic wake separation,fast-time prediction

V21

A doi：10.6052/0459-1879-17-198

2017–05–24 收稿，2017–10–18 录用,2017–10–19 网络版发表.

1)中国商飞与波音节能减排中心资助项目.

2)崔桂香，教授，主要研究方向：复杂湍流数值研究.E-mail:cgx@tsinghua.edu.cn

林孟达,崔桂香,张兆顺,许春晓,黄伟希.飞机尾涡演变及快速预测的大涡模拟研究.力学学报,2017,49(6):1185-1200

Lin Mengda,Cui Guixiang,Zhang Zhaoshun,Xu Chunxiao,Huang Weixi.Large eddy simulation on the evolution and the fast-time prediction of aircraft wake vortices.Chinese Journal of Theoretical and Applied Mechanics,2017,49(6):1185-1200