雷安鹏
（南昌航空大学 飞行器工程学院，江西 南昌 330063）
桨涡干扰（Blade-Vortex Interaction，BVI）噪声在旋翼噪声中占主导地位，且在低速下降飞行情况下具有较强的指向性，严重影响了直升机在人口密集区域的广泛使用[1]。BVI 噪声的数值直接计算可以为研究BVI 的产生机理和传播过程提供有效途径，从而对BVI 噪声进行控制。
伴随着现代计算机性能的提升和计算手段的创新，使用数值计算方法来预测气动噪声以及分析气动噪声产生机理渐渐成为了一种趋势。计算气动声学主要围绕着以下两种方法[2]进行发展，其中一种为混合法，该方法先是通过计算流体力学（Computational Fluid Dynamics，CFD）方法来获得较为准确的非定常流场的信息，然后使用声比拟方法对这些流场的信息进行求解从而得到噪声数据；另一种则为直接法，使用该方法直接对非定常、可压缩的Navier-Stokes（N-S）方程进行求解，便可以同时获取流场信息以及声场数据。混合法一般用于远场噪声的预测，Cao Y[3]和朱正[4]等人采用混合法对BVI 远场噪声进行了计算，他们的结果显示混合法能够较为准确地预测BVI 远场噪声。虽然混合法不需要占用大量的计算资源就能够得到较为准确的远场噪声数据，但当所研究的问题涉及到噪声生成机理的研究时，使用混合法就无法将噪声的生成与传播分开进行研究，相比之下直接法能够更好对噪声复杂的生成机理进行研究[5]。
虽然直接法有着混合法所不具备的优势，但因与主流能量相比，流场声波以及远场声场的能量在梯度上有很大的不同，故其需采用低耗散与低弥散的差分格式，才能对计算域内声波的传播过程进行较为准确的数值模拟[2，6-9]，即需要直接计算方法使用更高阶的计算格式对N-S 方程进行求解[10]。
与其他传统的噪声计算方法不同，格子玻尔兹曼方法（Lattice Boltzmann Method，LBM）较为容易实现并行运算，导致越来越多的研究人员习惯使用LBM 对相关气动声学的问题进行研究。针对LBM 的大量研究结果表明，因LBM 在空间上的计算精度较低，但因其时间离散方程使用的是具有低耗散、高精度的显式方程[11]，故其在用于求解非定常、低速不可压的湍流问题时却展现出非常好的效果[6]，并由此发展出了基于LBM 计算气动噪声的直接法[12-15]。
直接法计算气动噪声主要可以通过两种方式得以实现，一种为直接数值模拟法（Direct Numerical Simulation，DNS），另一种为大涡模拟法（Large Eddy Simulation，LES）。其中DNS 方法需要通过大量的网格来对各个的流动与噪声进行模拟，因此其需要占用大量的计算资源来求解远场噪声。辐射噪声主要是由较大的含能涡而引起，LES 方法通常对较大的含能涡进行直接求解，而选用亚格子模型来模拟较小的含能涡，不需要占用太多的计算资源就能较为精准地捕捉到近场的湍流特征与非定常的壁面压力[16]，因此使用LES 方法能够对气动噪声进行较为高效的直接数值计算。LES方法在求解高频噪声时有一定的局限性，因高频部分噪声的能量很小，故在空间与时间上的要求相当高；但BVI 噪声属于窄带中低频噪声[17]，其声场比较容易从动力学仿真结果中得到，其预测声场所需的时间和空间精度与流场所需的精度没有实质上的差别[18]，这就为使用LES 方法对BVI 近场以及远场噪声进行直接计算提供了可能性。
目前，国内外的研究大都使用混合法对BVI 远场气动噪声进行预测，很少有使用直接法的。但因直接法在研究BVI 气动噪声传播过程与产生机理等方面具有独特的优势，故本文在低马赫数条件下，使用LBM-LES 方法对平行BVI 气动噪声进行了直接数值模拟，对LES 方法的网格需求与适用性进行了分析，为研究BVI 噪声的传播特征以及生成机理提供较为适合的方法。
1 计算方法
1.1 LBM 方法
LBM 离散方程如下：
式中：f 表示粒子分布函数；ri代表节点i 的位置矢量；eα{α=0，1，…，Q－1}为离散速度；M 表示转换矩阵；表示碰撞矩阵；m 代表粒子分布函数f 通过转换矩阵M 映射到矩空间的量。
平衡态粒子分布函数的表达式如下：
式中，cs代表格子声速；ωα表示权重系数，ρ 代表流体密度；u 表示流体速度。
粒子分布函数统计可得到宏观变量密度和速度：
1.2 LES 方法
LES 方法的基本思想是对大涡进行直接数值求解，对小的涡则使用亚格子模型进行近似求解。BVI 的流动问题属于多尺度流动问题，但使用Smagorinsky 亚格子模型对多尺度流动问题进行数值模拟在技术上较难实施，因此本文LES 方法考虑使用动态Smagorinsky 亚格子模型[9]来求解小的涡。
动态Smagorinsky 亚格子模型中涡粘度的表达式如下：
1.3 远场声压计算方法
本文远场声压p′可以由下式直接得到
式中ρ′为密度脉动值，可由LES 方法得出。
2 网格分布及边界条件
2.1 网格分布
为了深入地研究不同网格密度与计算结果的相关性，如图1 所示，根据网格尺度的不同将整个计算区域大致分成四块区域：网格尺度很小的近壁面区域（Region 0）、网格尺度较小的翼型流动近场区域（Region 1）、网格尺度较大的声波辐射场区域（Region 2）还有网格尺度很大的声波吸收区域（Region 3）。其中，声波吸收区域之所以使用尺度很大的网格，是为了增大亚格子模型中的涡粘度系数，减少边界处虚假声波对BVI 噪声的影响，从而达到提高计算效率与噪声计算精确度目的[20]。
图1 计算区域的划分及网格分布示意图
2.2 边界条件
数值模拟的区域大小选取为80C×80C，C 代表翼型的弦长，将上、下两个边界y=±40C 均设定为压力出口边界，左边界x=-40C 与右边界x=±40C 则分别被设定速度进口边界、压力出口边界。为了避免压力波在边界处的反射，对声场数据造成过大的影响，所以采用求解局部一维无粘方程对每个边界都进行求解，从而实现无反射边界[21-22]。
3 计算结果与分析
本文算例来流速度U∞=12.2m/s，翼型弦长C=0.46m；旋涡模型为Lamb-Oseen，环量Γ=0.15U∞C，涡核半径rC=0.05C，涡的初始位置（xv，yv）=（-4C，0），其中xv为涡相对于翼型前缘的位置，yv为涡与翼型之间的垂直距离。
表1 为采用不同加密方式的5 组网格。当BVI 强度最大的时候，Mesh 5 近壁面网格y+<6，这说明Mesh 5 近壁面网格的尺度能够达到直接求解翼型壁面的标准，所以不需要使用壁面函数，但Mesh 1-4 近壁面网格的尺度都过大，因此都需要使用壁面函数对近壁区域的流动进行求解[23]。
表1 5 组不同区域网格尺度划分情况
3.1 网格尺度分析
图2 为观测点A（0，20C）在5 组网格中的声压时间历程，观察Mesh 4 与Mesh 5，发现两者除了近壁面区域网格的尺度不一样之外，其他三个区域网格的尺度都一样。观察Mesh 4 与Mesh 5 的声压时间历程曲线可以发现，两者声压随时间变化的整体趋势大致一样，且各个时间点的声压值相差不大。这就表明在BVI 噪声的数值模拟中，可以使用壁面函数近似的方法取代直接求解法，这样能够节省一定的计算资源。
图2 观测点A 的声压时间历程
根据mesh1，2，4 的计算数据可以发现，远场声压与壁面网格的密度相关性较小，而近场网格密度与声压相关性相对较大。这说明在计算BVI 噪声时，计算结果的精确性主要取决于远离壁面的近场流动，而不是近壁面的流动。
将Mesh 2 与Mesh 3 进行对比可以看出，两者除了近场区域的网格尺度不同之外，其他区域的网格尺度都相同。观察图2 可知，Mesh 2 和Mesh 3 近场采用不同尺度网格所得到的声压时间历程曲线基本重合，由此可以说明LES 方法在近场网格尺为1.56250×10-3m 的条件下可以较为准确地对近场声源信息进行捕捉。
将Mesh 3 与Mesh 4 进行对比可知，两者除了远场区域的网格尺度不同之外，其他区域的网格尺度都相同，但使用两种网格计算得到的远场声压的结果非常接近，这说明远场声压与远场区域尺度并无太大关联。近壁面区域与近场区域网格的尺度都要远小于远场区域网格的尺度，一方面表明使用LES 方法能在涡粘度较大的情况下，计算含能较高的BVI 噪声的传播过程；另一方面也表明近场区域的网格尺度决定了BVI 噪声直接数值计算结果的准确度。
3.2 声场计算
直接法的优势之一就是能够直接得到声场数据，图3 是不同时刻的膨胀度云图（Θ＝▽·V），从图中可以看出，本文基于LES 的直接法可以得到BVI 噪声的传播过程，可以判断噪声源的辐射特性、噪声源的位置，以及噪声的指向性；且可以看出BVI 噪声呈现出明显的脉冲噪声特性，向远场传播过程中，没有发现从边界来的虚假声波与BVI 噪声干扰。
图3 BVI 噪声向远场传播过程中的膨胀度云图
BVI 噪声向下远场传播过程中的声压曲线如图4 所示。图中BVI 噪声脉冲声压呈现出正弦分布，声压在向远场传播的过程中渐渐变小，但仍然保持正弦分布，这表明基于LES 的直接法能够用于BVI 噪声传播特征的研究和分析。
图4 BVI 噪声向下远场传播过程中声压曲线
4 结论
本文在低马赫数条件下，使用LBM-LES 方法对平行BVI 气动噪声进行了直接数值模拟，得到以下结论：
（1）当平行BVI 气动噪声可通过使用LES 方法的直接计算得到声场的信息和声波的传播过程；本文算例马赫数较低，马赫数增大，BVI 气动脉冲噪声会更强，更容易捕捉，对LES 方法尺度和耗散特性要求更容易满足，本文基于LBM-LES 的直接计算法能够对BVI 噪声的生成机理、传播规律及噪声特性进行研究。
（2）在BVI 气动噪声的直接计算中采用分区域不同尺度的网格相较于不分区域尺度相同的网格，其计算效率明显有所提高；因近壁面区域网格的尺度对噪声数值模拟结果的影响可以忽略不计，故LES 方法使用壁面函数对近壁面区域的流动直接进行近似求解，也能够获得比较准确的数值解；处于远场区域的网格其尺度可以很大，因此BVI 气动噪声的直接数值模拟所需要的计算资源相较于其他方法有所减少。
虽然目前本文的研究工作局限于二维，但随着计算机性能提高以及计算方法的创新，相信使用LBM-LES方法对直升机旋翼BVI 气动噪声进行三维直接数值计算将会变得可能。
