LES和无限元耦合方法预报螺旋桨均匀流噪声

王超，张立新，郑小龙，魏胜任

(1.哈尔滨工程大学船舶工程学院，黑龙江 哈尔滨150001;2.武汉第二船舶设计研究所，湖北 武汉430074)

螺旋桨噪声是船舶三大噪声之一，并在其中占有极其重要的位置，主要由旋转噪声、涡流噪声和空泡噪声等组成。螺旋桨的旋转噪声是一种周期性的噪声，它与螺旋桨桨叶负荷和厚度有关。涡流噪声不只是水流对螺旋桨的冲击，还包括由于桨叶叶梢和螺旋桨毂涡流破裂而产生的噪声。旋转噪声和涡流噪声都与声学理论中的单极子、偶极子和四极子声源机理相联系。对于均匀来流下的螺旋桨无空化噪声，主要由单极子噪声和偶极子噪声所组成［1］。关于螺旋桨噪声的理论预报大致分为两种:一种是利用势流理论程序预报，一种就是利用粘流软件进行数值模拟。朱锡清等［2-4］结合非定常升力面理论和声类比方法对船舶螺旋桨的低频线谱噪声和宽带谱噪声进行了研究;Seol［5-6］等采用面元法求解流场和FWH方法时域内求解声传播的耦合方法，计算了非均匀进流条件下DTRC4119无侧斜标准桨无空化和空化状态下的总声级中厚度噪声和负载噪声各自所占的比重以及声指向性，对非均匀流螺旋桨噪声做了全面的分析;杨琼方［7］采用流场大涡模拟(large eddy simulation，LES)和声场边界元数值声学的弱耦合方法在频域内进行了螺旋桨噪声分析，并对不同频率桨叶表面激励对声辐射的影响进行分析。这些计算结果大多没有表现出叶频处声压峰值等信息。对这一问题，杨琼方在文献［8］中做了很好的完善。本文采用FLUENT软件进行了大涡模拟(LES)非定常运算，对比了敞水性能和桨叶表面压力分布，并对螺旋桨周围流场的脉动压力幅值进行了分析，且利用CFD和声学无限元方法相结合的耦合方法预报了螺旋桨噪声特性。

1 流-声耦合方法理论

1.1 大涡模拟SGS模型

SGS模型在LES方法中占有十分重要的地位，本文应用Smagorinsky-Lilly模型来模拟亚格子应力:［9-10］

式中:μt是亚格子尺度的湍动粘度，在文献［11］中推荐用下式计算:

1.2 Lighthill声学类比理论

Lighthill考虑的模型为:在无限大的均匀、静态声介质中包含一个有限的湍流运动区域V，因此与流动有关的声源都集中在该区域内。在区域V外，远离湍流区域的流体中密度的波动和声波相似，因此整理连续方程和动量方程，并简化得到远离湍流区域流体中的匀质声学波动方程:

Lighthill声类比理论是从流体力学基本方程纳维-斯托克斯(Navier-Stokes)方程导出的:

式中:c0是等熵条件下的声速值;ρ＇=ρ-ρ0为有声扰动时的密度分量，ρ与ρ0分别是扰动与未扰动时的流体密度;Tij是 Lighthill应力张量［12］:

式中，σij表示雷诺应力张量的粘性部分。

对于低速等熵流动，粘性应力张量σij对Lighthill应力张量Tij的贡献远小于雷诺应力项ρuiuj，可以忽略不计，同时可以得到

因此，可得到Lighthill应力张量的近似式Tij=ρuiuj。

1.3 ACTRAN 声学原理

ACTRAN/LA基于 Lighthill方法，并结合 Curle’s理论:

1)Curle’s方程的体积分作为有限元区域的体源;

2)Curle’s方程的面积分作为边界条件;

3)自由场的格林函数作为其他的边界条件。

对式(6)在Ω边界上积分，乘以测试函数δρ，并应用分布积分产生弱变分形式，在面积分上应用应力张量，可得到下式:

1.4 螺旋桨噪声模拟基本流程

ACTRAN软件处理流致噪声问题时，CFD计算与声学计算是耦合的，即首先进行CFD仿真，提取出湍流信息，然后再利用Lighthill声类比方法分析声场。对于声学分析中，只要满足每波长6网格的规则即可，使用积分法将流场信息加载到声学网格上，因此不需要对声源区的网格做特别的优化。图1为流-声耦合方法的计算流程图。

图1 流－声耦合计算流程图Fig.1 Flowchart of fluid-acoustics interaction calculation

2 计算前处理

2.1 几何模型建立和计算域的划分

文章研究螺旋桨为DTRC4119桨，其具有比较详细的流场数据，被ITTC选为考证数值方法预报精度的标准桨，数据取自文献［13］，几何参数如表1所示。

表1 螺旋桨几何参数Table 1 Geometry coefficients of propeller

图2 螺旋桨几何模型和计算域Fig.2 Geometry model of propeller and Computational domain

在建模过程中使用的是直角坐标系O-XYZ，X轴方向代表来流方向，它沿着螺旋桨的旋转轴指向下游，Y轴沿螺旋桨某一叶片的母线，Z轴服从右手定则。为了方便对螺旋桨进行结构化网格划分，只需要建立一个桨叶的几何体和附近小域即可。创建完成的螺旋桨三维模型如图2(a)所示。

计算域采用与螺旋桨同轴的圆柱流域，并被划分为两个区域，内域包含螺旋桨。上游速度入口位置为桨前2倍直径，下游自由出口设定在桨后5倍直径，外边界半径大小为3倍直径，如图2(b)所示。

2.2 网格划分及参数设定

网格划分是CFD模拟过程中比较耗时的环节，也是直接影响模拟精度和效率的因素之一。本文在划分网格时采用全结构化网格。根据螺旋桨周期性特点，首先划分单个桨叶流域内的网格，然后进行旋转，得到全流域内的网格。划分单桨叶流域网格时，采用H型网格，在桨叶位置处单独布置C型网格与桨叶形状匹配，此外在桨叶位置处向叶面外一侧开O网(如图3所示)，此种方式可以达到在桨叶表面进行局部加密，同时还可控制网格数量，本文边界层第一层网格尺寸取 0.000 3，Y+控制在 10～100。

对于流场的关键区域(如桨叶随边、导边、叶根与桨毂连接处、叶梢等)进行加密，以便捕捉到重要的流场信息;而对于距离螺旋桨较远的区域网格将其密度适当降低，便于控制总网格数，这样，在网格模型总节点数一定的情况下可以提高计算精度，还可以避免流场变化平缓区域的计算资源浪费。图4(a)是随边叶根处O型网格局部放大图。桨叶网格如图4(b)所示，整体计算域网格约为200万。

图3 桨叶外O网Fig.3 O-block out of propeller

螺旋桨壁面定义为不可滑移壁面条件，外壁面不考虑粘性作用，采用滑移壁面。计算时采用单一旋转参考坐标系模型，将旋转参考坐标系固定于螺旋桨中心，并以600 r/min的角速度进行旋转，外域采用绝对静止坐标系，两个域之间利用INTERFACE边界进行连接，流场信息通过插值进行传递。

图4 螺旋桨网格拓扑结构Fig.4 Mesh topology of propeller

3 大涡模拟计算结果分析

计算螺旋桨流场，首先采用k-ε湍流模型进行定常运算，待获得稳定流场后，改用LES湍流模型进行非定常运算。计算时采用有限体积法进行离散，扩散项采用中心差分格式，压力速度耦合采用SIMPLEC算法，连续性曲线小于0.000 01时认为计算收敛，时间步长设置为 0.002 5。

3.1 敞水性能计算结果验证

为便于计算结果表述，定义如下参数:

式中:J表示进速系数，KT表示推力系数，KQ表示扭矩系数。

对进速系数在0.5～0.9进行大涡模拟，对于非定常运算，流场并不是恒定不变，存在微小波动，待计算稳定后，取整周期读取其平均值作为计算结果，对比其敞水性能。

图5为计算结果与试验值［14］的比较，从中可以看出螺旋桨水动力性能计算值与试验结果吻合较好，KT曲线最大误差为1.5%，KQ曲线最大误差1.417%，这说明本方法计算螺旋桨水动力性能具有很好的准确性。

图5 螺旋桨敞水性能曲线Fig.5 Open water performance of propeller

3.2 桨叶表面压力系数对比分析

在研究螺旋桨的诱导脉动压力、空泡和噪声等问题时，求解出流场中的桨叶压力分布是十分重要的。对比螺旋桨0.3R、0.7R叶剖面处的压力系数CP，以进一步验证计算方案的可靠性，如图6所示。CP=，其中 (P-P0)为相对压力，相对进流速度

图6DTRC4119螺旋桨表面压力系数分布(J=0.833)Fig.6 Pressure coefficient distribution of propeller DTRC4119(J=0.833)

由图6可以看出，通过大涡模拟计算出的结果与试验值之间误差较小，0.7R处压力分布与试验值吻合较好，0.3R处相对吻合稍差，这可能是因桨毂简化为圆柱体而省略了毂帽部分，而0.3R较贴近桨毂，使得该处表面压力分布存在些许偏差。导边(x/c=0.0)与随边(x/c=1.0)处偏差比叶面中部大。这是由于导边与随边处流动梯度变化大，压力变化明显，造成模拟计算结果偏差较大，文献［16］也揭示了相同的情况。

3.3 脉动压力计算分析

为进一步获取流场信息，对桨叶周围流场中特征点的脉动压力进行监控分析，特征点的分布情况如图7所示。

图7 特征点位置Fig.7 Positions of the feature points

图8是经过FFT变换后的各特征点的脉动压力幅值，由于螺旋桨旋转作用，能够明显看出桨叶的叶频(blade passing frequency，BPF)及其各倍叶频。螺旋桨转速为10 r/s，叶数为3，故BPF=30，同时可以看出一阶叶频脉动压力幅值远高于其他倍叶频处，这也是研究脉动压力时一般只对一阶脉动压力进行研究的原因。通过对比桨叶径向P1、P2、P3三点的脉动特性，得知P1和P2之间压力幅值衰减速度明显大于P2和P3之间的衰减速度，这是由于螺旋桨径向作用范围有限，距离桨叶较近范围，螺旋桨影响程度随距离增加迅速减弱，向外延伸，流场逐渐进入稳定状态。P4和P5两点BPF幅值达到了300 Pa以上，明显高于径向上的P1点，说明该处受螺旋桨旋转作用影响很大。P6和P7两点较P4和P5两点衰减相当显著，尤其是P7的BPF压力幅值已降低到了90 Pa左右，说明沿桨叶轴向向后移动，桨叶尾流影响作用逐渐减小。对比P4～P7，可以看出，靠近桨毂的P5、P7两点压力幅值要小于各自同轴向位置的P4、P6两点，说明P5、P7两点受螺旋桨旋转作用影响相对也略小。

图8 特征点的脉动压力幅值Fig.8 Pulsating pressure amplitude of feature points

4 螺旋桨噪声数值计算

通过大涡模拟获得稳定的非定常流场，创建声学网格，包括声源区和声传播区，其中声源区取自CFD计算区域，略小于CFD计算区域。声场的外边界是一层无限元边界面，无限元可以传播声而不反射声，因此较好地解决了声音反射为计算带来误差的问题。同时为了有利于传播声波，声学网格尺寸满足每波长至少6个网格节点。CFD节点和声学网格节点之间通过保守整合法进行信息传递以保证计算精度。

取3个距离桨中心5D距离的特征点分别位于轴向桨前和桨后、径向，3点的频谱曲线如图9所示，由于时间步长设置为0.002 5，对应有效频率上限为200 Hz，参考声压为 1 μPa。

图9 特征点的声压频谱曲线Fig.9 Sound pressure spectrum at feature points

由图9可以看出，螺旋桨各特征点的频谱曲线能够反映出叶频BPF信息，且在各阶叶频处都存在声压峰值。径向上桨前桨后两点的声压频谱曲线几乎重合。径向与轴向相比，各阶叶频处声压峰值相差不大，但在一阶叶频内的较低频段，径向位置声压明显小于轴向位置。

对应频谱曲线中特殊频率，给出各频率时通过螺旋桨轴线的截面声压云图，如图10所示。图10(a)以5 Hz声压云图为例，展示了在一阶叶频内频段的声压云图，该频段偶极子噪声尤为显著，声压分布呈现出横8字形特征，螺旋桨轴向声压明显高于径向声压。图10(b)～(d)为螺旋桨各阶叶频 BPF=30、60、90 Hz处的声压分布云图，观察螺旋桨周围，可以发现声压分布同样呈8字形，这也体现出螺旋桨周围流场的脉动压力幅值大小，轴向位置受桨叶旋转影响很小，压力幅值明显小于桨叶周围其他位置，此规律在一阶叶频BPF处表现尤其明显。同时可以看出在二阶叶频和三阶叶频处8字形明显减小，桨叶后方声压高于桨叶前方，说明后方受桨叶旋转作用影响要强于前方。总体来说，可以看出螺旋桨噪声主要由单极子噪声组成，低频段内比重较大的偶极子噪声在各阶叶频处相对于单极子噪声要小很多，在云图中已经无法体现，这与图9频谱曲线表述是一致的。

图10 螺旋桨噪声的声压云图Fig.10 Sound pressure contours of propeller noise

在通过螺旋桨轴线平面上距离桨中心5D距离取32个测点平均分布于螺旋桨周围，根据总声压级计算:

式中:SLi是第i个1/3oct中心频率点处声压级。

计算得到各点处总声压级，如图11所示。可见，螺旋桨周围5D处各点声压级大小相近，符合均匀来流单极子声源的特性，各点最大声压级为112.3 dB，最小声压级为111.3 dB，误差可能是流体网格和声学网格插值时产生的，其差值仅为1 dB。

图11 螺旋桨声指向性(距离中心5D)Fig.11 Noise directivity of propeller at distance 5D

5 结论

本文针对NREL5MW风力机设计了Spar式风力机平台，并采用频域水动力分析方法对Spar式风力机平台的水动力性能进行数值模拟，分析平台有无纵荡板、平台的重本文采用流-声耦合方法对螺旋桨均匀流噪声进行了分析，进而得出以下结论:

1)通过分析敞水性能、桨叶表面压力分布，验证LES湍流模型结合全结构化网格技术的计算方案具有较好的准确性;

2)通过对桨叶周围流场的脉动压力计算，可以显示出各阶叶频处的峰值，均匀来流时的一阶叶频(BPF)脉动压力幅值明显高于其他叶频;

3)通过计算得到的桨叶周围特征点的声压频谱曲线，可以看出，在各阶叶频处存在声压峰值，且在一阶叶频处达到最大值，总声压级大小和该值大小相近;

4)通过计算得到的声压云图分析，验证了均匀流螺旋桨噪声主要由单极子噪声和偶极子噪声组成，其中单极子噪声占主要部分。

通过本文大涡模拟计算螺旋桨水动力性能，结合Actran声学软件对螺旋桨均匀流噪声进行数值模拟，为进一步非均匀流噪声的数值计算打下了基础。

［1］ROSS D.Mechanics of underwater noise［M］.New York:Pergamon Press，1976:113-119.

［2］朱锡清，李亚，孙红星.船舶螺旋桨叶片与艉部湍流场互作用噪声的预报研究［J］.声学技术，2006，25(4):361-364.ZHU Xiqing，LI Ya，SUN Hongxing.Ship propeller blades and stern turbulence field interaction noise prediction research［J］.Journal of Acoustic Technology，2006，25(4):361-364.

［3］朱锡清，唐登海，孙红星，等.船舶螺旋桨低频噪声研究［J］.水动力学研究与进展(A辑)，2000(1):74-81.ZHU Xiqing，TANG Denghai，SUN Hongxing，et al.Ship low-frequency noise of propeller study［J］.Water Dynamics Research and Progress(A)，2000(1):74-81.

［4］熊紫英，朱锡清，刘小龙，等.船尾伴流场-导管-螺旋桨互作用噪声预报研究［J］.声学学报，2009(2):117-123.XIONG Ziying，ZHU Xiqing，LIU Xiaolong，et al.Stern wake field-catheter-interaction noise prediction research［J］.Acta Acustica，2009(11):117-123.

［5］SEOL H，JUNG B，SUH J C，et al.Prediction of non-cavitation underwater propeller noise［J］.Journal of Sound and Vibration，2002，257(1):131-156.

［6］SEOL H，SUH J C，LEE S.Development of hybrid method for the prediction of underwater propeller noise［J］.Journal of Sound and Vibration，2005，288:345-360.

［7］杨琼方，王永生，曾文德，等.大侧斜螺旋桨负载噪声的边界元数值声学方法频域内计算分析［J］.兵工学报，2011(9):1118-1125.YANG Qiongfang，WANG Yongsheng，ZENG Wende，et al.Big side of the noise of the propeller load boundary element numerical acoustic method in the frequency domain analysis［J］.Acta Armamentarii，2011(9):1118-1125.

［8］杨琼方，王永生，魏应三，等.基于URANS模拟的伴流场中大侧斜桨无空化噪声的时域和频域预报［C］//第十三届船舶水下噪声学术讨论会论文集.北京，2011:325-338.

［9］KOBAYASHI T.Large eddy simulation for engineering applications［J］.Fluid Dynamics Research，2006，38:84-107.

［10］HOLM D D，GEURTS B.Commentator errors in large-eddy simulation［J］.Journal of Physics A:Mathematical and General，2006，39:2213-2229.

［11］FELTEN F，FAUTRELLE Y，TERRAIL Y D，et al.Numerical modelling of electromagnetically-driven turbulent flows using LES methods［J］.Applied Mathematical Modelling，2004，28(1):15-27.

［12］ESCOBAR M.Finite element simulation of flow-induced noise using Lighthill＇s acoustic analogy［D］.Bayen:University Erlangen Nurnberg，2007:36-57.

［13］KOYAMA K.Comparative calculations of propeller by surface panel method-workshop organized by 20th ITTC propulsion committee［C］//Papers of Ship Research Institute，(s.l.)，1993:179-186.

［14］谭廷寿.非均匀流场中螺旋桨性能预报和理论设计研究［D］.武汉:武汉理工大学，2003:53-66.TAN Tingshou.Theory of propeller operating in non-uniform flow field performance prediction and design research［D］.Wuhan:Wuhan University of Science and Engineering，2003:53-66.

［15］STAINER M J.The application of a RANS code to model propeller DTRC 4119［R］.Grenoble:ITTC report，1998.

［16］王超.螺旋桨水动力性能、空泡及噪声性能的数值预报研究［D］.哈尔滨:哈尔滨工程大学，2010:78-79.WANG Chao.Propeller hydrodynamic performance and cavitation and noise performance of the numerical prediction research［D］.Harbin:Harbin Engineering University，2010:78-79.