水面艦船流噪聲數(shù)值預(yù)報(bào)研究與應(yīng)用

張維俊, 楊瓊方

(1.海軍裝備部駐武漢軍事代表局 ,湖北 武漢 430064) (2.海軍工程大學(xué) 動(dòng)力工程學(xué)院，湖北 武漢 430033)

提高艦船全航速范圍內(nèi)的噪聲性能和艦船具有盡可能高的空化初始航速的隱蔽性,需要應(yīng)用低噪聲船型.在該應(yīng)用背景下,將認(rèn)識(shí)艦船在3類噪聲源貢獻(xiàn)程度最小的流噪聲從定性層面,轉(zhuǎn)為對(duì)給定航速下流噪聲定量總聲壓級(jí)以及流噪聲等效聲中心確定,是值得開展的研究方向.

文中以某型船為對(duì)象,采用計(jì)算流體力學(xué)(computational fluid dynamics,CFD)和計(jì)算聲學(xué)(computational acoustics,CA)相結(jié)合的數(shù)值預(yù)報(bào)方法,對(duì)中高航速下的流噪聲進(jìn)行了數(shù)值預(yù)報(bào),求取了不同航速下的等效聲中心.研究結(jié)果可為實(shí)船海試噪聲測(cè)量過程中測(cè)點(diǎn)的布置和距離換算提供參考.

1 流噪聲數(shù)值預(yù)報(bào)方法

目前在用的噪聲定性和定量數(shù)值計(jì)算方法主要有4種:同時(shí)求解流場(chǎng)與聲場(chǎng)控制方程的直接計(jì)算法(computational aere-acoustics,CAA),CFD與聲波動(dòng)方程耦合求解的混合法(hybrid method),基于聲類比方程的積分模型法(integral method)和基于局部湍流尺度的聲源強(qiáng)度統(tǒng)計(jì)法[1].聲源強(qiáng)度統(tǒng)計(jì)法能快速地表明物體表面流動(dòng)對(duì)噪聲起主要作用的部位,并不計(jì)算接收點(diǎn)的聲壓信息,不能預(yù)報(bào)遠(yuǎn)場(chǎng)噪聲,計(jì)算精度最低,只能用于聲場(chǎng)的定性分析.CAA與CFD計(jì)算中的直接數(shù)值模擬(direct numerical simulation,DNS)一樣,因計(jì)算量巨大,現(xiàn)階段尚無法在工程中得到應(yīng)用.目前在氣動(dòng)噪聲預(yù)報(bào)中應(yīng)用較廣泛的是混合法和積分法,可根據(jù)分析對(duì)象的不同和獲得脈動(dòng)壓力場(chǎng)途徑的不同而具體選擇.

文獻(xiàn)[1]對(duì)脈動(dòng)流體輻射噪聲的數(shù)值預(yù)報(bào)方法及相應(yīng)的理論給予了較為全面地論述.當(dāng)流體域中存在剛體和靜止壁面時(shí),如果使用自由空間的格林函數(shù),就可以得到萊特希爾方程的Curle解[2-3],從而考慮固體壁面的影響,將流體流經(jīng)固體壁面產(chǎn)生的噪聲等同于壁面上的等價(jià)力偶極源和壁面周圍的湍流四極源輻射發(fā)聲.當(dāng)流體域中存在旋轉(zhuǎn)壁面時(shí),文獻(xiàn)[4]中FW-H方程進(jìn)一步將流體與旋轉(zhuǎn)壁面的相互作用產(chǎn)生的輻射噪聲等同于3類噪聲源輻射發(fā)聲:單極源產(chǎn)生的厚度噪聲、旋轉(zhuǎn)偶極源產(chǎn)生的負(fù)載噪聲和湍流四極源產(chǎn)生的渦流噪聲.FW-H方程是以控制面上的速度和壓力分布對(duì)Lighthill波動(dòng)方程進(jìn)行的數(shù)學(xué)描述,允許區(qū)分3個(gè)源項(xiàng):兩個(gè)面積分項(xiàng)和一個(gè)體積分項(xiàng).但計(jì)算費(fèi)時(shí)且較難實(shí)現(xiàn)的體積分仍然包含在內(nèi).在體積源項(xiàng)可以忽略不計(jì)的情況下,通過FW-H方程求解聲場(chǎng)時(shí)僅需要包含面積分.這一點(diǎn)正是FW-H方法相對(duì)于Lighthill聲類比公式的強(qiáng)大之處.FW-H方程也是目前主要采用的噪聲預(yù)報(bào)方法的基礎(chǔ).

混合法和積分法預(yù)報(bào)噪聲的核心都在于非定常壓力場(chǎng)的求解.針對(duì)分析對(duì)象、計(jì)算資源和對(duì)求解精度要求的差異,通常選擇非定常雷諾時(shí)均模擬(unsteady Reynolds averaged Navier-Stokes simulation,URANS)、分離渦模擬(detached eddy simulation,DES)、尺度適應(yīng)模擬(scale-adaptive simulation, SAS)、大渦模擬(large eddy simulation, LES)以及雷諾時(shí)均模擬和LES相結(jié)合的混合模擬方法中的一種來進(jìn)行.

文獻(xiàn)[5]全面地論述了從近場(chǎng)CFD計(jì)算到遠(yuǎn)場(chǎng)聲預(yù)報(bào)過程中采用的面積分方法的數(shù)值模型,包括Kirchhoff方法和FW-H方程方法,給出了Farassat Kirchhoff公式、擴(kuò)展的Kirchhoff公式以及剛體壁面的FW-H方程和虛擬滲透面的Porous FW-H方程,并且比較了兩種積分方法的差異.Kirchhoff公式因不需要體積分項(xiàng),可以避免傳統(tǒng)的聲類比方法中體積分計(jì)算的相關(guān)問題.但該方法的一個(gè)缺點(diǎn)是Kirchhoff面必須選擇在線性流動(dòng)區(qū)域,以使其輸入面壓力、法向和時(shí)間導(dǎo)數(shù)滿足均相波動(dòng)方程.然而,線性區(qū)域的位置在定義時(shí)卻并不容易實(shí)現(xiàn),并與所分析的問題直接相關(guān).通常的CFD計(jì)算并不能滿足Kirchhoff面對(duì)遠(yuǎn)離噪聲源的線性區(qū)域要求,所以,Kirchhoff面的放置通常需要折中考慮.Porous FW-H方程可以較好的解決上述問題,其既能利用Kirchhoff方法相對(duì)于傳統(tǒng)聲類比方法的優(yōu)點(diǎn),同時(shí)又能避免公式中積分面的放置位置所帶來的靈敏性問題.當(dāng)控制面(可以是實(shí)際的剛體壁面,也可以是虛擬的滲透面)合適地置于遠(yuǎn)離物體且包含所有聲源時(shí),只需計(jì)算面積分項(xiàng)即可,類似于Kirchhoff方法中的處理思路.在該情況下,積分面的位置只需要出于計(jì)算便利考慮,而不要求必須置于線性流動(dòng)區(qū).因該公式同時(shí)融合了FW-H和Kirchhoff方法的優(yōu)點(diǎn),通常也被稱之為Kirchhoff-FWH公式(KFWH).除了上述Kirchhoff面方法和FW-H方法以外,以有限元(finite element method, FEM)數(shù)值聲學(xué)方法為主的體積分方法(volume integral method)在湍流噪聲預(yù)報(bào)中也得到了一定的應(yīng)用.文獻(xiàn)[6]采用該方法對(duì)簡(jiǎn)單的二維長方體和三維立方體湍流噪聲進(jìn)行了預(yù)報(bào).其CFD計(jì)算時(shí)同時(shí)采用了SAS和LES模擬,聲場(chǎng)模擬得到100Hz內(nèi)單根主要譜線頻率與實(shí)測(cè)值相比誤差小于3Hz,總聲壓級(jí)誤差小于3dB.當(dāng)前,由比利時(shí)聲學(xué)專家開發(fā)的噪聲預(yù)報(bào)商用軟件LMS Virtual-Lab中的Acoustics模塊中的FEM聲學(xué)子模塊即采用了該方法,其在汽車行業(yè)噪聲與振動(dòng)方面的應(yīng)用正日益擴(kuò)大.

與積分法對(duì)應(yīng)的混合法噪聲預(yù)報(bào)中,仍然是以求解Farassat 1A公式為主,且通常是在頻域內(nèi)計(jì)算非均相Helmholtz波動(dòng)方程.此時(shí)必須要指定合適的邊界條件,如在考慮散射面影響時(shí)的Dirichelet，Neumann或者是Robin邊界條件,或者是計(jì)算輻射噪聲時(shí)的Sommerfeld邊界條件,LMS Virtual-Lab中的Acoustics模塊中的BEM聲學(xué)子模塊即采用了該方法.該方法中噪聲源項(xiàng)的獲得與積分法一樣,仍然是由CFD計(jì)算或者是實(shí)驗(yàn)測(cè)量得到,只需將得到的脈動(dòng)壓力的時(shí)域信息轉(zhuǎn)換為頻域信息即可.除該方法外,湍流的波數(shù)-頻率譜頻域內(nèi)預(yù)報(bào)方法也在氣動(dòng)噪聲預(yù)報(bào)中得到了一定的應(yīng)用.

綜上所述,當(dāng)前在流噪聲工程數(shù)值預(yù)報(bào)中采用的方法主要有基于FW-H方程的積分法和混合法兩種.具體是在時(shí)域還是頻域內(nèi)數(shù)值求解流體噪聲,取決于聲源量向聲網(wǎng)格節(jié)點(diǎn)傳遞時(shí)是時(shí)域量還是頻域量[7-8].文中船體定常流場(chǎng)采用RANS方法進(jìn)行計(jì)算,湍流模型取為剪切應(yīng)力輸運(yùn)模型,非定常流場(chǎng)計(jì)算采用LES方法,亞格子尺度模型取為DSM模型,并將流場(chǎng)定常RANS計(jì)算作為假定初值,以加快計(jì)算收斂.聲場(chǎng)計(jì)算采用頻域內(nèi)分步耦合的方法進(jìn)行分析.

2 流噪聲混合法預(yù)報(bào)數(shù)值模型

LES是一種直接模擬大尺度脈動(dòng)渦而將小尺度渦對(duì)大尺度渦運(yùn)動(dòng)的影響通過亞格子尺度模型(sub-grid scale, SGS)來模擬的湍流數(shù)值模擬方法,不再如RANS模擬中那樣需要采用湍流模型的假定.在計(jì)算源流場(chǎng)的脈動(dòng)時(shí),大渦模擬方法在精度上是最能給予保證的,但其對(duì)空間網(wǎng)格尺度和時(shí)間尺度的要求也最高.LES的計(jì)算精度主要取決于區(qū)分大尺度渦和小尺度渦的濾波函數(shù)和模擬小尺度渦運(yùn)動(dòng)的SGS模型[9].參照文獻(xiàn)[10]中對(duì)螺旋槳緊急倒車流動(dòng)的LES模擬設(shè)置,文中采用動(dòng)力學(xué)Smagorinsky-Lilly模型來模擬亞格子應(yīng)力:

(1)

式中：C為Smagorinsky系數(shù)，計(jì)算時(shí)取C=0.1;Δ為局部網(wǎng)格尺度,取為單元體積的1/3.

在得到聲源節(jié)點(diǎn)時(shí)域脈動(dòng)量后,需要將其轉(zhuǎn)移到聲網(wǎng)格節(jié)點(diǎn)上,繼而求解聲傳播.網(wǎng)格節(jié)點(diǎn)變量傳遞的插值方式會(huì)對(duì)計(jì)算結(jié)果產(chǎn)生一定的影響.采用邊界元積分法在頻域內(nèi)計(jì)算輻射聲場(chǎng)時(shí),聲源節(jié)點(diǎn)和聲節(jié)點(diǎn)可以是一一對(duì)應(yīng)的變量守恒傳遞(一對(duì)一),也可以是在適當(dāng)稀化聲網(wǎng)格節(jié)點(diǎn)的基礎(chǔ)上,采用節(jié)點(diǎn)加權(quán)的方式進(jìn)行變量傳遞(多對(duì)一),這樣可以適當(dāng)減小計(jì)算量[11].為盡量提高計(jì)算精度,文中采用“一對(duì)一”的節(jié)點(diǎn)變量傳遞方式.

對(duì)于平動(dòng)剛體輻射發(fā)聲,非均相波動(dòng)方程可表述為:

(2)

Tij≈ρvivj

(3)

式中：G為聲格林函數(shù)；Ω為積分面.當(dāng)忽略四極源項(xiàng)的體積分時(shí),等式右邊僅剩下面積分項(xiàng).

采用邊界元數(shù)值聲學(xué)方法頻域內(nèi)計(jì)算時(shí),Pij即對(duì)應(yīng)為船體表面的脈動(dòng)壓力,可采用直接邊界元積分公式,其一般形式為[7,11]:

(4)

3 某型船流噪聲預(yù)報(bào)與分析

3.1 船體非定常流動(dòng)的大渦模擬

分析對(duì)象取為某排水型船,按照?qǐng)D紙,由三維CAD軟件建立船體幾何模型.船體流場(chǎng)計(jì)算域范圍設(shè)置為:進(jìn)流方向2倍船長,出流方向3倍船長,側(cè)邊界1.5倍船長.整個(gè)計(jì)算域采用六面體全結(jié)構(gòu)化網(wǎng)格進(jìn)行空間離散,對(duì)船體首、尾和自由液面區(qū)域均進(jìn)行局部加密.數(shù)值計(jì)算區(qū)域及船體壁面網(wǎng)格如圖1.

(5)

式中：Awp為水線面積；Iy為水線面積對(duì)y軸的面積慣性矩.船體縱傾及重心升沉具體修正過程如下:

1) 在給定航速下,定常計(jì)算船體處于正浮狀態(tài)下的流場(chǎng);

2) 根據(jù)定常計(jì)算得到的縱向合力Fz修正船體吃水;根據(jù)縱傾力矩My修正縱傾角Δα,并按縱傾角調(diào)整船體姿態(tài)及船體附近的離散網(wǎng)格;

3) 按照調(diào)整后的參數(shù)和網(wǎng)格重新計(jì)算;

4) 迭代2)和3)兩個(gè)步驟,直至力和力矩達(dá)到平衡.

圖1 船體壁面和計(jì)算域離散網(wǎng)格Fig.1 Hexed mesh on hull surface and wall in computation domain

因船體結(jié)構(gòu)左右對(duì)稱,數(shù)值計(jì)算時(shí)取半船體流場(chǎng)區(qū)域,通過對(duì)稱邊界條件來代替另一半船.首先采用RANS方法求解船體定常流場(chǎng),作為非定常流場(chǎng)求解的迭代初始條件,以加快計(jì)算迭代.非定常計(jì)算時(shí)取船尾底板不同監(jiān)控點(diǎn)以分析其流場(chǎng)脈動(dòng)特征,測(cè)點(diǎn)布置如圖2.

圖2 船尾脈動(dòng)流場(chǎng)監(jiān)控點(diǎn)布置Fig.2 Moitoring points of fluctuating pressure in ship wake flow

計(jì)算得到在不同航速下船體波形如圖3.可以看出,隨著航速增加,尾波峰逐漸后移,且波峰和波谷幅值變大,對(duì)輻射噪聲的有效貢獻(xiàn)區(qū)域增加.興波凱爾文角減小,首波峰和首二次波峰之間的過渡區(qū)變窄.非定常計(jì)算時(shí)間步取為Δt=2.5×10-5s,對(duì)應(yīng)流噪聲有效分析頻率為2 kHz.將時(shí)域上船體壁面壓力分布進(jìn)行FFT變化后,得到14 kn航速時(shí)頻域內(nèi)特征頻率下船體壁面壓力分布如圖4，可以看出球首部位以及球首后方對(duì)整個(gè)頻段的聲貢獻(xiàn)都最強(qiáng).

a)14 kn

b)18 kn

c)22 kn圖3 不同航速下船體波形Fig.3 Wave patterns in different ship speeds

a) 31.5 Hz

b) 80 HZ

c) 4 000 Hz圖4 船體表面特征頻率下脈動(dòng)源強(qiáng)分布Fig.4 Fluctuating source intesity at characteristic frequencies on hull

3.2 流噪聲計(jì)算

光體聲邊界元網(wǎng)格如圖5.將脈動(dòng)源流場(chǎng)傳遞到聲網(wǎng)格后,考慮流場(chǎng)介質(zhì)可壓縮性影響,在一個(gè)聲波長范圍內(nèi)至少布置8個(gè)聲網(wǎng)格單元,以計(jì)算聲傳播.與流場(chǎng)計(jì)算相對(duì)應(yīng),光體縱中剖面采用聲對(duì)稱面處理,自由液面采用聲反對(duì)稱面處理,以考慮聲反射影響.計(jì)算得到在航速14 kn時(shí)特征頻率下船體壁面自噪聲場(chǎng)如圖6.圖中可以看出,聲貢獻(xiàn)量大區(qū)域集中于艏波首波峰、二次波峰和艉波峰以及球首尾流區(qū),且球首尾流區(qū)對(duì)于全頻段都有明顯的貢獻(xiàn),與非定常計(jì)算脈動(dòng)聲源對(duì)應(yīng).

圖5 船體聲邊界元網(wǎng)格Fig.5 Acoustic boundary element nodes on hull

a) 25 Hz

b) 1 000 Hz

c) 2 000 Hz圖6 特征頻率下船體壁面噪聲源強(qiáng)分布Fig.6 Acoustic source intensity at characteristic frequencies on hull

非定常計(jì)算時(shí),仍然取空間特征點(diǎn)脈動(dòng)壓力譜進(jìn)一步計(jì)算脈動(dòng)流場(chǎng)的輻射聲場(chǎng),布置空間測(cè)點(diǎn)如圖7,同時(shí)對(duì)應(yīng)有近場(chǎng)和遠(yuǎn)場(chǎng)測(cè)點(diǎn),以校驗(yàn)大幾何尺度下聲輻射是否滿足球面衰減規(guī)律.計(jì)算得各測(cè)點(diǎn)在2 kHz頻段內(nèi)1/3oct中心頻率輻射聲壓譜級(jí)曲線如圖8.可知,由船首到船尾的測(cè)點(diǎn),因受聲指向性影響,隨著離船體縱中剖面距離的增大,頻段22 Hz～2 kHz內(nèi)總聲壓級(jí)衰減規(guī)律并不相同,如圖9.僅有F點(diǎn)和O點(diǎn)處衰減值近似于球面衰減規(guī)律,表明只有在聲源集中處或者是聲主要貢獻(xiàn)量區(qū)域,才能近似應(yīng)用球面衰減規(guī)律得出離物體1 m處的聲源級(jí)值.依據(jù)聲壓級(jí)衰減量,可將該航速下船體流噪聲的等效聲源中心取于船體尾波處,則船尾處在20 Hz～2 kHz頻段內(nèi)總聲源級(jí)可達(dá)133 dB,船首二次波峰處該頻段內(nèi)總聲源級(jí)可達(dá)131 dB;采用同樣的計(jì)算方法,將非定常計(jì)算時(shí)間步進(jìn)一步減小至有效頻率為20 kHz時(shí),船尾總聲源級(jí)達(dá)143.3 dB.

其它兩個(gè)航速下的流噪聲計(jì)算方法相同,不再重復(fù)敘述.通過計(jì)算得到:隨著航速增加,等效聲中心后移,總聲級(jí)增加.航速達(dá)到22 kn時(shí),等效聲中心位置移到P點(diǎn)處.

圖7 流噪聲測(cè)點(diǎn)布置Fig.7 Hudrophones for flow noise exanination

圖8 測(cè)點(diǎn)流噪聲譜曲線Fig.8 Noise spectrum of flow noise on different observers

圖9 測(cè)點(diǎn)流噪聲沿橫向距離增加20倍時(shí)聲壓衰減量Fig.9 Reduced total sound pressure level of flow noise when the perpendicular distance is enlarged 20 times

4 結(jié)論

采用CFD結(jié)合CA的分步耦合計(jì)算方法,對(duì)船體湍流脈動(dòng)壓力場(chǎng)及其輻射聲場(chǎng)進(jìn)行了計(jì)算和分析,得出了船體壁面聲貢獻(xiàn)量主要區(qū)域、流噪聲空間測(cè)點(diǎn)的聲壓頻譜及在船體正橫方向的聲壓衰減規(guī)律,定量求出了給定航速下的流噪聲等效聲中心.航速為14 kn時(shí),船尾波峰位置20 Hz～2 kHz頻段內(nèi)總聲源級(jí)達(dá)133 dB,當(dāng)計(jì)算有效頻段擴(kuò)展到20 kHz時(shí),船尾總聲源級(jí)達(dá)143.3 dB.光體流噪聲主要來源于興波引起的渦量,且主要集中在100 Hz～20 kHz頻段.球首尾流區(qū)和船體尾渦區(qū)對(duì)聲輻射量貢獻(xiàn)明顯,特別是球首尾流區(qū),對(duì)全頻段輻射聲都有明顯貢獻(xiàn),設(shè)計(jì)時(shí)需重點(diǎn)關(guān)注.文中方法也為水面艦艇流噪聲的研究探索了一種新途徑.

參考文獻(xiàn)(References)

[ 1 ] Wang M, Freund J B, Lele S K. Computation prediction of flow-generated sound[J].AnnualReviewofFluidMechanics,2006,38: 483-512.

[ 2 ] Lighthill M J.On sound generated aerodynamically, I:general theory[C]∥ProceedingsoftheRoyalSociety.1952, A211: 564-587.

[ 3 ] Curle N. The influence of solid boundaries upon aerodynamic Sound[C]∥ProcedingsoftheRoyalSociety.1955, A231: 505-514.

[ 4 ] Ffowcs Williams J E, Hawkings D L. Sound generation by turbulence and surfaces in arbitrary motions[J].PhilosophicalTransactionsoftheRoyalSociety,1969,A264:321-342.

[ 5 ] Lyrintzis A S. Surface integral methods in computational aeroacoutics-from the (CFD) near-field to the (acoustic) far-field[J].InternationalJournalofAeroacoustics,2003, 2(2): 95-128.

[ 6 ] Escobar M. Finite element simulation of flow-induced noise using lighthill′s acoustic analogy[D].Erlangen:UniversityErlangen-Nurnberg,2007.

[ 7 ] Yang Qiongfang, Wang Yongsheng,Zhang Mingmin. Scale effects on non-cavitation hydrodynamics and noise of highly-skewed propeller in wake flow[J].JournalofSoutheastUniversity：EnglishEdition, 2013, 29(2): 162-169.

[ 8 ] Kato C,Yamade Y,Wang H,et al.Numerical prediction of sound generated from flows with a low Mach number[J].Computers&Fluids,2007,36:53-68.

[ 9 ] Kobayashi T. Large eddy simulation for engineering applications[J].FluidDynamicsResearch,2006, 38: 84-107.

[10] Vysohlid M. Large eddy simulation of crashback in marine propellers[D].Twin cities:University of Minnesota,2007.

[11] LMS International. Numerical acoustics[R].Belgium:LMS Virtual Lab Theoretical Manual,2006.