三维孔腔流激噪声的大涡模拟与声学类比预报与验证研究

张 楠，沈泓萃，朱锡清，姚惠之，谢 华

（中国船舶科学研究中心，江苏 无锡 214082）

1 前 言

流声耦合是流体力学与声学的交叉学科，是国际瞩目的重点与难点。流激噪声或称流动发声的预报是流声耦合领域的重要组成部分，具有相当的技术难度。这是因为水是可压缩性极小的流体介质，要研究的问题一般是低马赫数的流动，而水中运载体的雷诺数都很大，声转换效率相当低，因此水中噪声预报是十分复杂和困难的。研究水中流激噪声的预报方法，探讨涡旋与声辐射之间的相互关系，具有重要意义。

近年来，大涡模拟（LES）方法的逐渐成熟以及计算机的快速发展，使得计算一些试验状态下的流动发声成为可能。动态Smagorinsky模型（Germano 1991，Lilly 1992）[1-2]对于大涡模拟的可靠性与准确性做出了重大的改善。 Seror等人 （2000，2001）[3-4]，He等人 （2002，2004）[5-6]的工作都研究了动态Smagorinsky模型在预报流动辐射噪声方面的可靠性。Wang，Freund与Lele（2006）[7]指出对于实际应用的雷诺数范围，在可预见的将来，LES是一种最有希望取得重大突破的用于流激噪声预报的CFD模拟方法。Wang与Moin（2000）[8]利用LES结合Lighthill理论模拟了翼型支杆的随边噪声，弦长雷诺数为Re=2.1×106，其计算结果与 Blake（1975）的试验结果进行了对比，达到了合理的吻合程度。Wang（2004）[9]通过计算弦长雷诺数为Re=1.5×105的风扇叶片模型非定常压力的时空特性，进一步验证了LES方法。这些工作都为预报复杂流动的辐射噪声奠定了基础。

孔腔流动从属于与自持振荡密切相关的一类基本流动范畴。孔腔流动能够引起结构振动与疲劳，激发很强的噪声，其中包含了剪切层不稳定性，湍流与结构及噪声之间的相互作用等流固耦合、流声耦合复杂现象，给理论分析、数值模拟和试验测量带来了很大的挑战。在过去的五十多年里，这一类型的流动被广泛的研究过，最初的工作是Roshko在1955年所做的试验测量与分析[10]。Rossiter（1964）[11]给出了辨识流动模态振荡共振频率的公式。近年来，人们又做了大量的试验和计算研究。

Ask等人（2003）[12]采用不可压缩CFD代码CALC-BFC计算了长深比L/D=4的二维长方形孔腔的绕流，采用Curle方程计算了远场9个接收点处的辐射噪声。来流马赫数M∞=0.15。计算充分说明了壁面脉动压力是重要的偶极子声源。计算得到的声压谱趋势与参考结果完全一致，幅值最大误差为6dB。

Ching（2004）[13]采用MILES方法，利用非结构化网格研究了低马赫数下的二维“L”型孔腔发声问题（M∞=0.0788～0.147）。计算域中含85300非结构化网格。在M∞=0.147时，对于流体动力振荡引起的声谐振峰，频率的计算误差为1%，幅值的计算误差为1dB，对于Hemholtz共振引起的声谐振，频率的计算误差为2.2%，幅值的计算误差为17dB。对于其他谐振峰，计算得到的声压谱幅值都明显大于试验结果。

在空气动力学领域，M 219孔腔是一种被广泛研究的经典模型。M219孔腔长（L）、宽（W）、深（D）之比为，L:W:D=5:1:1，来流马赫数M∞=0.85，是接近于跨音速状态下的一种高雷诺数流动问题。研究者做了大量的脉动压力与噪声试验，并采用这些详细的试验结果来校验各种数值计算方法的可信度。

瑞典国防研究署（SDRA）的Peng（2005）[14]对于M 219孔腔的非定常流动进行了数值模拟，计算工作被SDRA当作一种基准式的计算（baseline-type computation），为随后所有该类型的计算提供可资参考的技术细节。研究者首先在跨音速风洞中测量了脉动压力，给出了脉动压力的频谱，试验清晰反映了不同频率下的谐振峰。作者使用结构化的NS求解器EURANUS，采用非定常RANS方法结合SA模型，数值模拟了二维孔腔的脉动压力。计算中采用了4套网格，网格数分别为：19884（C）、20930（NM）、33036（M）、51410（F），分别进行了计算。研究表明，只有正确求解了来流边界层和孔口剪切层以及孔中的流谱，才能获得对压力谱和声谐振的准确预报。M网格和F网格计算得到的幅值与试验值接近程度超过前两套网格，特别是对于F网格，其预报的第一阶模态的幅值与频率都与试验十分吻合。

Larcheveque（2003）[15]采用大涡模拟方法，研究了M219孔腔的流动与声辐射，网格数为3.2×106。与试验结果相比，其预报的声谐振峰幅值的误差分别为-3dB、-14dB、2dB、-4dB，只有第四阶模态频率的计算误差在5%之内。数值模拟计算表明最重要的谐振峰为第三阶模态，而试验测量表明第二阶模态是最主要的。

Mendonca（2003，2004）[16-17]使用离散涡模拟（DES）方法研究了 M 219 孔腔流动与发声，计算中采用了两套网格1.1×106，2.8×106。从脉动压力均方根的计算结果来看，精细的网格并没有明显提高计算精度。Ashworth（2005）[18]比较了非定常RANS（URANS）方法与DES方法在计算M219孔腔绕流时的差异，网格数为1.68×106。计算结果表明，URANS丢失了第二阶模态，但是DES方法捕捉到了这一模态。

Lai等人（2006）[19]针对M 219孔腔，采用大涡模拟和FW-H声类比技术研究了抑制流激噪声的渗流壁面方法。来流马赫数M＝0.85，雷诺数Re＝1.36×106。网格数为441万。Lai等人在孔腔底部和后壁上采用了渗流介质（porous media），流入孔腔的流体可以通过渗流墙流入或流出，从而孔腔之内高度非定常的流动得到调节，涡强显著减小。壁面渗流减小了流动与壁面的撞击以及波反射，改变了流动振荡的相位，有效降低了近场模态振荡的幅值，因而可以明显降低孔腔内的各种脉动，从而使辐射噪声减小。

Chen（2007）[20]采用CFD数值模拟结合气动声学预报方法研究了三维M219孔腔流动的近场与远场流/声特性。研究者旨在通过商用代码SBLI和SotonCAA来发展一种有效的数值方法来完成从近场流动到远场噪声的完整的数值模拟。对于M219孔腔的流动和噪声进行了大涡模拟（LES）和离散涡模拟（DES）并与试验结果进行了对比。二维LES算例的网格数为7.2×104，三维LES算例的网格数为5.74×106。 二维DES算例的网格数为8.09×104与 8.17×104，三维DES算例的网格数为 1.05×106与4.08×106。与试验结果相比，可以看出三个特点：（1）LES预报精度高于DES；（2）三维算例预报精度高于二维算例；（3）数值计算低估了频率而高估了噪声。

本文采用大涡模拟方法求解非定常流动，进而计算出流体动力声源。在此基础上，采用基于Kirchhoff积分（源于电磁场求解的Kirchhoff-Helmholtz定理）的FW-H声学类比方法计算出远场辐射噪声。基于Kirchhoff积分的FW-H方法可以克服传统声学类比方法计算体积分的困难，具有很高的灵活性和实用性。目前，在国际气动声学界此方法得到一定程度的运用，而在国内，特别是水动力学界，此方法的应用还非常稀少。本文即采用上述方法求解了水中孔腔的流动辐射噪声。

2 流激噪声计算方法

2.1 大涡模拟（LES）方法

通过将非定常的NS方程进行滤波，得到大涡模拟的控制方程。滤波过程有效地过滤掉了那些尺度小于滤波宽度（或网格尺度）的小涡。

滤波函数G（x,x′）取为：

则滤波后的连续性方程和NS方程可以表示为：

其中，σij为分子粘性引起的应力张量，τij为亚格子应力

本文采用动态Smagorinsky模型，其推导过程简述如下：

Germano（1991）在推导动力模型的过程中定义了两种滤波算子。

另一个是测试滤波算子（test filter operator），假设测试滤波宽度比网格滤波大，也就是说测试滤波相当于较粗网格的滤波。Germano等人经研究槽道湍流发现较优的测试滤波宽度与网格滤波宽度的比为 2，亦即

考虑可解的湍流应力Lij。Lij表征介于网格滤波宽度与测试滤宽度之间的雷诺应力的贡献，亦即小的可解尺度的贡献。

（14）式就相当于在最小二乘意义下的最优值。这个办法提供了一个自洽（self－consistent）的动态模型，直接通过大涡模拟本身的计算结果，可以在每个空间网格点、在每一时间步上计算模型参数C。

2.2 Ffowcs Williams-Hawkings（FW-H）方程

Ffowcs Williams与Hawkings（1969）对于Lighthill理论进行了推广，考虑了流体中任意移动壁面的发声问题。他们使用了广义函数的数学理论，依照Lighthill的推导步骤，将连续性方程和动量方程整理为含有两个面源和一个体源的有源波动方程，即FW-H方程：

其中，ui：xi方向的流体速度分量；un：垂直于物体表面f＝（0）的流体速度分量；vi：xi方向的物体表面速度分量；vn：物体表面法向速度分量；δ（f）：Dirac delta函数；H（f）：Heaviside阶跃函数。

p′为远场声压。f=0表示物体表面，f>0表示外部无界的自由空间。ni为壁面外法线，指向流体内部。c0为远场声速，ρ0为远场密度，Tij为Lighthill应力张量。

利用自由空间格林函数δ（g）/4πr，再结合Kirchhoff可渗透表面积分法可以得到远场解。远场辐射声压中包含厚度噪声pT′（x,t），负荷噪声pL′（x,t）以及四极子噪声pQ′（x,t）。

详细的公式推导以及公式中符号含义和公式表达的物理意义见文献[21]。

2.3 数值计算方法

时间项采用二阶隐式格式离散，动量方程采用有界中心差分格式离散，压力速度耦合采用PISO算法。计算中时间步长Δt=1×10-5s。壁面y+≈1～5。采用代数多重网格法加速收敛。

3 计算结果分析

3.1 计算模型

本文计算了流经孔腔的非定常流动以及孔腔的流动辐射噪声。孔腔开设在载体上，载体长4m，宽0.8m，厚0.3m，首尾及两侧光顺过渡。孔腔编号、几何特征以及在算例中所用的网格数见表1。五个孔腔在载体上的位置以及相互之间的大小关系见图1。计算中水速取V=5m/s，与试验保持一致，基于载体长度的雷诺数为Re=2×107。试验中每次只开设一个孔腔。

表1 孔腔编号与尺寸Tab.1 Serial number and size of cavity

3.2 非定常流动计算结果分析

Colonius（2001）[22]指出三维孔腔绕流与发声问题的复杂程度远远超过二维孔腔。本文所研究的正是三维孔腔，下面的计算结果足以说明这一问题的复杂程度。流动充分发展之后，某时刻下孔腔内的流谱见图2～6，图中左侧为孔腔中纵剖面内的流线分布，右侧为孔腔中的三维流线。载体开设4#孔腔后表面涡量分布见图7。

从图2～6可知，孔腔内的流动是十分复杂的非定常三维涡旋流动，孔内充满了多种形式的旋涡，它们的翻卷和扰动影响十分显著。在孔口处，有着明显的内外流动交换。这些旋涡产生、输运并相互干扰，构成了主要的发声根源。从图7可知，孔腔极大地影响了后方及其周围的涡量分布，其影响范围和增大与减小的趋势十分复杂。我们下面分析一下引起远场辐射噪声的主要成因。

Rossiter（1964）将孔腔中的自持振荡（self-sustaining oscillation、self-controlled oscillation）归结为孔腔流动发声的根源。孔穴自持振荡可以分为三种类型：流体—动力（fluid-dynamic）振荡、流体—共振（fluid-resonant）振荡；流体—弹性（fluid-elastic）振荡。下面我们就上述三种类型逐一分析，找出本研究中孔腔辐射噪声的来源。

首先，计算中采用的是刚性壁面，孔腔没有变形，也不存在结构振动。试验模型也进行了加固和刚性处理，通过加速度传感器测量分析，可知也不存在结构振动。所以，那种由于结构变形和结构振动而引起的流体—弹性振荡在本研究中没有出现，因而流体—弹性振荡的影响可以不予考虑。

其次，我们知道，当流体动力振荡与孔腔的某阶声驻波共振的固有频率相吻合时，将发出极强的单频音（tone），产生声驻波耦合共振，这种流/声共振将产生强烈的辐射噪声。Rossiter（1964）给出了辨识流体动力模态振荡频率的半经验公式：

其中，U0为来流速度，L为孔腔沿流向的长度，M为马赫数，n为模态阶数，fn为模态频率。

第一阶孔腔声驻波共振的固有频率为：

其中，C为环境声速，L为孔腔沿流向的长度，H为孔腔深度。

要注意的是水中声速约为1500m/s，水中航行体的运动都是低马赫数或极低马赫数下的运动（一般M≈0.001～0.01），再考虑到水中航行体上常见孔腔的尺寸，可知上述流体动力振荡的模态频率和声驻波的固有频率之间相差很远，一般而言，不太可能出现上述的耦合共振，以2#孔腔为例，则可以算得第一阶fA=14.2Hz，fB=1421.8Hz，二者在数量级上相差甚远，其他各阶模态也有很大差异，不太可能出现耦合共振。无论是从计算得到的孔腔辐射噪声频谱上，还是从试验得到的孔腔辐射噪声频谱上，都没有出现声谐振尖峰，说明耦合共振现象不存在，因而流体-共振振荡的影响可以不予考虑。

综上所述，可以肯定地讲，本研究中孔腔发出的辐射噪声完全是由于流体动力振荡引起的。孔腔开口处存在自由湍流剪切层，孔腔周围为不连续的壁面湍流边界层，孔腔内为非定常的涡旋运动，它们三者相互作用，共同组成辐射噪声的声源。归根到底，还在于三者所包含的脉动剪应力、脉动压力与脉动速度，它们构成了最根本的偶极子声源与四极子声源。

3.3 流动辐射噪声计算结果分析

本节分析了计算得到的孔腔流激辐射噪声。在试验中，单水听器位于开孔部位的正下方。文中所使用的试验结果已经过自由场修正和距离修正，测量到的流动辐射噪声声级统一进行归一化处理，即测量距离为1m，分析带宽1Hz，参考声压Pref=1μPa的声源级。详细的处理过程和结果见文献[23]，本文不再赘述。计算与试验状态保持一致。计算与试验结果的比较见图8～12。

从图8～12可知，计算反映出的孔腔辐射噪声频谱趋势与量级与试验吻合很好，各孔腔发声的大小关系与试验略有差异。在2000～10kHz孔腔主要发声频段，1#孔声压谱级计算误差为-4.4～1.7dB；2#孔声压谱级计算误差为-0.6～3.3dB；3#孔声压谱级计算误差为-0.9～4.4dB；4#孔声压谱级计算误差为-0.8～4.6dB；5#孔声压谱级计算误差为-3.5～0.8dB。可见在此频率范围内计算与试验结果十分吻合，表明建立的计算方法是可靠的，可以用于复杂流动结构诱发辐射噪声的预报研究。在较低频段，计算误差增大，但还在可接受的范围内。

4 结 语

流声耦合是流体力学与声学的交叉学科，该领域的研究需要深厚的数学基础、繁琐的数值计算和精细的试验验证，其难度是不言而喻的。孔腔流动与流激噪声的计算与试验研究是流声耦合领域中的重要研究课题，其深刻的流/声内涵以及计算与试验上的难度已经得到国际上的普遍公认，具有相当的技术难度。本文在前期工作的基础上，首先利用大涡模拟方法预报了三维孔腔在水中绕流的流动特征，然后利用FW-H声学类比方法预报由此产生的远场辐射噪声，分析了噪声产生的流动根源，并与相应的试验结果进行了对比，证明本文方法的计算流程切实可行，计算量级可靠，可用于水中流激噪声的预报分析问题，并应在将来的工作中做进一步的分析和验证。

[1]Germano M,Piomelli U,Moin P,Cabot W H.A dynamic subgrid-scale eddy viscosity model[J].Phys.Fluids,1991,A3(7):1760-1765.

[2]Lilly D K.A proposed modification of the Germano Subgrid-scale closure method[J].Phys.Fluids,1992,A4(3):633-635.

[3]Seror C,Sagaut P,Bailly C,Juv’e D.Subgrid-scale contribution to noise production in decaying isotropic turbulence[J].AIAA J,2000，38(10):1795-1803.

[4]Seror C,Sagaut P,Bailly C,Juv’e D.On the radiated noise computed by large eddy simulation[J].Phys.Fluids,2001,13(2):476-487.

[5]He G W,Rubinstein R,Wang L P.Effects of subgrid scale modeling on time correlations in large eddy simulation[J].Phys.Fluids,2002,14(7):2186-2193.

[6]He G W,Wang M,Lele S K.On the computation of space-time correlations by large-eddy simulation[J].Phys.Fluids,2004,16(11):3859-3867.

[7]Wang Meng,Freund J B,Lele S K.Computational prediction of flow-generated sound[J].Annu.Rev.Fluid Mech.,2006,38:483-512.

[8]Wang M,Moin P.Computation of trailing-edge flow and noise using large eddy simulation[J].AIAA J,2000,38(12):2201-2209.

[9]Wang M,Moreau S,Iaccarino G,Roger M.LES prediction of pressure fluctuations on a low-speed airfoil.Annu[Z].Res.Briefs,Cent.Turbul.Res.,Stanford Univ./NASA Ames,Stanford,CA,2004:183-193.

[10]Roshko A.Some measurements of flows in a rectangular cutout[K].NACA TN 3488,1955.

[11]Rossiter J E.Wind tunnel experiments on the flow over rectangular cavities at subsonic and transonic speeds[R].Rep.Mem.3438,Aeronautical Research Council,1964.

[12]Jonas Ask,Lars Davidson,Hans Enwald,Johan Larsson.An acoustic analogy applied to incompressible flow fields[C]//Computational Aeroacoustics:From Acoustic Sources Modeling to Far-Field Radiated Noise Prediction Colloquium EUROMECH 449,December 9-12,2003.Chamonix,France,2003.

[13]Loh Ching Y.Computation of low speed cavity noise[C]//42nd Aerospace Sciences Meeting and Exhibit sponsored by the American Institute of Aeronautics and Astronautics Reno.Nevada,January 5-8,2004.(AIAA-2004-0680)

[14]Peng Shiahui.Unsteady RANS simulation of turbulent cavity flow:Summary of 2D baseline computations[R].Swedish Defence Research Agency.FOI-R-1915-SE,Technical report,December 2005.

[15]Larcheveque L,Sagaut P,Le T H.Large-eddy simulations of flows in weapon bays[G].AIAA Paper 2003-0778,2003.

[16]Mendonca F,Allen R,Charentenay J,Kirkham D.CFD prediction of narrowband cavity acoustics at M=0.85[G].AIAA Paper 2003-3303,2003.

[17]Allen R,Mendonca F.DES validations of cavity acoustics over the subsonic to supersonic range[G].AIAA Paper 2004-2862,2004.

[18]Ashworth R M.Prediction of acoustic resonance phenomena for weapon bays using detached eddy simulation[J].Aeronautical Journal,2005:631-638.

[19]Lai H,Luo K H.Large-eddy simulation and control of cavity aeroacoustics[C]//Conference on Turbulence and Interactions TI2006,May 29-June 2,2006.Porquerolles,France,2006.

[20]Chen Xiaoxian,Sandham Neil D,Zhang Xin.Cavity Flow Noise Predictions[R].University of Southampton,Report No.AFM-07/05,2007.

[21]Brentner K S.An efficient and robust method for prediction helicopter high-speed impulsive noise[J].JSV,1997,203(1):87-100.

[22]Tim Colonius.An overview of simulation,modeling,and active control of flow/Acoustic resonance in open cavities[G].AIAA 2001-0076,2001.

[23]熊紫英,孙红星.孔腔模型流动噪声测量试验报告[R].无锡:中国船舶科学研究中心科技报告，2009.