串聯(lián)圓柱體繞流氣動(dòng)噪聲三維數(shù)值仿真*

□ 寧方立 □ 王善景 □ 馬 堯 □ 郭琪磊

西北工業(yè)大學(xué) 機(jī)電學(xué)院 西安 710072

飛機(jī)噪聲問(wèn)題已成為航空界研究的前沿領(lǐng)域和技術(shù)難點(diǎn)之一［1］。飛機(jī)噪聲主要包括發(fā)動(dòng)機(jī)噪聲和機(jī)體噪聲，機(jī)體噪聲包括起落架氣動(dòng)噪聲和增升裝置氣動(dòng)噪聲，尤其是在飛機(jī)的起飛和著陸階段，飛機(jī)的機(jī)體噪聲與發(fā)動(dòng)機(jī)噪聲已經(jīng)處在同一水平。Chow L.C.等［2］對(duì)空客A340進(jìn)行試驗(yàn)發(fā)現(xiàn)，起落架噪聲比襟翼噪聲高6dB，所以研究飛機(jī)的起落架噪聲對(duì)降低飛機(jī)總體噪聲具有重要意義。由于起落架的結(jié)構(gòu)復(fù)雜，對(duì)起落架的氣動(dòng)噪聲直接仿真計(jì)算十分困難。起落架的很多部件都具有圓柱體形狀，如：機(jī)輪、支撐柱、減震器、軟管等，可以將這些結(jié)構(gòu)體簡(jiǎn)化為圓柱體。串聯(lián)圓柱相對(duì)單個(gè)圓柱來(lái)說(shuō)，流動(dòng)更加復(fù)雜，呈現(xiàn)更加明顯的三維特性，因此通過(guò)對(duì)串聯(lián)圓柱體的研究，將為低噪聲起落架的噪聲預(yù)測(cè)和制造提供理論性的指導(dǎo)。同時(shí)，物體繞流也是一個(gè)廣泛存在于航空航天、船舶、機(jī)械等工程中的實(shí)際問(wèn)題。

國(guó)內(nèi)外對(duì)圓柱體繞流研究主要包括數(shù)值仿真和試驗(yàn)兩類方法。龍雙麗等［3］對(duì)Re=90 000的二維圓柱繞流氣動(dòng)噪聲進(jìn)行了數(shù)值仿真。趙良舉等［4］對(duì)二維串聯(lián)圓柱體繞流氣動(dòng)噪聲進(jìn)行了數(shù)值模擬。劉敏等［5］使用大渦模擬（Large Eddy Simulation：LES）和 Farassat-1A方程對(duì)串聯(lián)圓柱體進(jìn)行了流場(chǎng)和聲場(chǎng)模擬，研究在不同間距比下流場(chǎng)和遠(yuǎn)場(chǎng)聲場(chǎng)的聲壓頻譜圖的變化。Cox J.S.等［6］采用基于 Lighthill的聲類比方法與雷諾時(shí)均Navier-Stokes（RANS）方程相結(jié)合，對(duì)單個(gè)圓柱體的流場(chǎng)及遠(yuǎn)場(chǎng)輻射噪聲進(jìn)行仿真。Lockard D.P.等［7］應(yīng)用基于有限體積法求解三維RANS的CFL3D軟件，模擬了圓柱體間的流動(dòng)對(duì)上、下游圓柱體的不同影響。Brès G.A. 等［8］使 用 離 散 波 爾 茲 曼 方 法 （Lattice Boltzmann Method：LBM）與FW-H方程相結(jié)合的方法，對(duì)串聯(lián)圓柱體的遠(yuǎn)場(chǎng)氣動(dòng)噪聲進(jìn)行了預(yù)測(cè)。 UZUN A.等［9］使用延遲分離渦方法 （Delayed Detached Eddy Simulation：DDES）計(jì)算了串聯(lián)圓柱體的流場(chǎng)結(jié)果。

本文基于LES和FW-H方程相結(jié)合的方法，對(duì)串聯(lián)圓柱體繞流氣動(dòng)噪聲進(jìn)行三維數(shù)值仿真。計(jì)算分為兩步：首先基于LES得到圓柱繞流的非定常湍流流場(chǎng)分布；然后求解FW-H方程模擬遠(yuǎn)場(chǎng)氣動(dòng)噪聲。

1 計(jì)算方法

1.1 LES理論

LES的主要思想是：直接模擬大尺度湍流運(yùn)動(dòng)，而利用亞格子模型模擬小尺度的湍流流動(dòng)對(duì)大尺度湍流流動(dòng)的影響。LES在計(jì)算時(shí)間和成本方面優(yōu)于DNS（Direct Numerical Simulation：DNS）， 而在計(jì)算精度方面優(yōu)于RANS。LES計(jì)算方法能夠獲得比RANS更多的湍流信息以及比DNS更有效的快速計(jì)算。LES方法是在現(xiàn)有計(jì)算條件下較為精確的方法。

湍流運(yùn)動(dòng)黏性控制方程：

式中：ρ為流體密度；u為流動(dòng)速度；P為壓力項(xiàng)；S為拉伸率張量，Sij=（?ui/?xj+?uj/?xi）/2 ；γ 為黏性系數(shù)。 下標(biāo)i、j取值范圍是（1，2，3）。

將式（1）中 ui分為大尺度渦 （以上標(biāo)“-”表示）和的小尺度渦（以上標(biāo)“’”表示），即 ui=ui+ui′。

過(guò)濾函數(shù)選擇帽型函數(shù)（Top-hat）：

式中：Δ為網(wǎng)格平均尺度。

將過(guò)濾函數(shù)作用于N-S方程的各項(xiàng)，得到過(guò)濾后的不可壓縮N-S方程：

式中：t為時(shí)間變量。

本文選擇的亞網(wǎng)格尺度模型為Smagorinsky-Lilly模型。

1.2 遠(yuǎn)場(chǎng)聲場(chǎng)計(jì)算

通過(guò)連續(xù)性方程和N-S方程可以得到FW-H方程［10］，F(xiàn)W-H 方程右端含有三項(xiàng)，如式（5），右端第一項(xiàng)代表湍流應(yīng)力，第二項(xiàng)為施加在某些面上非穩(wěn)定力的散度，第三項(xiàng)包括進(jìn)入到流體中的非穩(wěn)定質(zhì)量流，這三項(xiàng)噪聲源項(xiàng)分別帶有四極子、偶極子、單極子的特性。

式中：f為積分表面；ui為xi方向的速度分量；un為在f=0面上的法向速度；vn為物面速度的法向分量；δ（f）為 Dirac 函數(shù)；p′為遠(yuǎn)場(chǎng)聲壓，p′=p-p0；a0為遠(yuǎn)場(chǎng)的聲速；H（f）為 Heaviside 廣義函數(shù)，H=（f）T為L(zhǎng)ighthill應(yīng)力張ij量，為表面的載荷，

對(duì)于串聯(lián)圓柱體氣動(dòng)噪聲計(jì)算，定義積分表面為圓柱體表面，是不可滲透面。在低馬赫數(shù)產(chǎn)生的聲場(chǎng)中，將四極子噪聲源忽略。

2 試驗(yàn)與計(jì)算模型

2.1 基準(zhǔn)實(shí)驗(yàn)

NASA Langley Research Center在BART和QFF中對(duì)串聯(lián)圓柱體進(jìn)行了試驗(yàn)研究［7］，這是串聯(lián)圓柱體氣動(dòng)噪聲的基準(zhǔn)試驗(yàn)。BART中，Ma=0.128，Re=166 000，展向長(zhǎng)度S=12.4D，圓柱體間距L=3.7D，這是串聯(lián)圓柱體氣動(dòng)噪聲的臨界間距。QFF裝置在消聲實(shí)驗(yàn)室中進(jìn)行，該裝置中圓柱體直徑D=0.057 15 m，S=16D，L=3.7D，Re=166 000。

2.2 模型與網(wǎng)格

計(jì)算串聯(lián)圓柱體模型與QFF試驗(yàn)一致，截面如圖1所示。圓柱方位角θ以上游駐點(diǎn)為0°，逆時(shí)針?lè)较驗(yàn)檎较?。XY平面網(wǎng)格如圖2所示，壁面網(wǎng)格的y+≤1。計(jì)算區(qū)域選擇長(zhǎng)為20D、寬為10D、高為16D的長(zhǎng)方體模型。上游圓柱體圓心距入口5D，距出口15D，入口為速度入口，出口為壓力出口，圓柱體表面為壁面，圓柱體末端所在平面為壁面，其余邊界為壓力出口。

2.3 流場(chǎng)計(jì)算

基于FLUENT軟件計(jì)算非定常湍流流動(dòng)，選擇基于壓力的求解器進(jìn)行計(jì)算。湍流數(shù)值模擬方法選擇LES方法，壓力和速度的耦合采用PISO算法，壓力差值算法為PRESTO！。時(shí)間和空間的參數(shù)采用二階精度計(jì)算，時(shí)間步長(zhǎng)為 5×10-6s。

3 計(jì)算結(jié)果與分析

計(jì)算0.1 s后，上下游圓柱體的升力系數(shù)和阻力系數(shù)呈現(xiàn)出穩(wěn)定周期性特征，以此開始計(jì)算流場(chǎng)和聲場(chǎng)數(shù)據(jù)。

3.1 圓柱表面壓力系數(shù)

壓力系數(shù)為：

式中：p0為遠(yuǎn)場(chǎng)靜壓；U0為初速度。

▲圖1 模型和坐標(biāo)系統(tǒng)示意圖

▲圖2 計(jì)算區(qū)域網(wǎng)格

▲圖3 上游和下游圓柱體的壓力系數(shù)在t=0.24 s時(shí)的瞬態(tài)分布

▲圖4 中間截面處CP的平均分布

▲圖5 特征平面上的時(shí)均速度分布

從圖3可以看出，Cp=1的區(qū)域 （流動(dòng)趨于停滯狀態(tài)）主要位于兩個(gè)圓柱體的迎風(fēng)面上，下游圓柱體的停滯區(qū)域較少且分布規(guī)律性較差，因此形成更大的脈動(dòng)力，產(chǎn)生較大的聲壓波動(dòng)，所以下游圓柱體能夠產(chǎn)生較大的遠(yuǎn)場(chǎng)氣動(dòng)噪聲。

圖4是在展向中間截面處的時(shí)均壓力系數(shù)分布。通過(guò)與BART和QFF得到的試驗(yàn)數(shù)據(jù)對(duì)比，仿真結(jié)果與試驗(yàn)值相當(dāng)吻合。本文的仿真結(jié)果能夠準(zhǔn)確計(jì)算出圓柱體表面的壓力分布情況。

3.2 時(shí)均速度的分布

圖5為平面上的時(shí)均速度分布，平面A位于串聯(lián)圓柱體展向長(zhǎng)度的中間位置，平面B位于y=0處。圓柱體前段速度較小是因?yàn)榱鲃?dòng)處于停滯狀態(tài)，圓柱體后端速度較小原因是由于圓柱體后端形成大量的渦，圓柱體兩側(cè)會(huì)形成加速區(qū)域。在上游圓柱體方位角θ約為100°和260°位置時(shí)速度變化最大，邊界層在圓柱體的這個(gè)位置開始分離，形成漩渦進(jìn)而撞擊到下游圓柱體上，使得下游圓柱體周圍的流動(dòng)更加復(fù)雜。

3.3 遠(yuǎn)場(chǎng)噪聲分布

功率譜密度（Power Spectral Density，PSD）是描述脈動(dòng)壓力能量隨頻率的分布。聲壓級(jí) （Sound Pressure Level，SPL）為將待測(cè)聲壓有效值p′與參考聲壓pref的比值取常用對(duì)數(shù)，再乘以20。OASPL（Overall Sound Pressure Level）為測(cè)量點(diǎn)的總聲壓級(jí)。

圖6是QFF試驗(yàn)中使用麥克風(fēng)對(duì)串聯(lián)圓柱體遠(yuǎn)場(chǎng)氣動(dòng)噪聲測(cè)量的點(diǎn)，3個(gè)點(diǎn)的坐標(biāo)分別為A（-8.33D，27.815D， 8D）、B （9.11D，32.49D，8D）、C （26.55D，27.815D，8D），D為圓柱體直徑。

由圖7可以看出，本文得到的遠(yuǎn)場(chǎng)噪聲與基準(zhǔn)試驗(yàn)結(jié)果能較好吻合，3個(gè)測(cè)量點(diǎn)的最大聲壓級(jí)均出現(xiàn)在基頻（fs=185 Hz）處，在倍頻處也出現(xiàn)了局部峰值。而基準(zhǔn)試驗(yàn)［10］得到的渦脫落頻率是 178～180 Hz，模擬值和基準(zhǔn)試驗(yàn)值峰值位置較吻合，誤差在3.3%。PSD最大值出現(xiàn)在渦脫落頻率位置處，說(shuō)明該噪聲是由于有規(guī)律的渦脫落引起的脈動(dòng)力所引起的。噪聲在很寬的頻率上都有分布，串聯(lián)圓柱產(chǎn)生的氣動(dòng)噪聲屬于寬頻噪聲。

3.4 遠(yuǎn)場(chǎng)輻射噪聲指向性

▲圖6 麥克風(fēng)測(cè)量位置示意

▲圖7 遠(yuǎn)場(chǎng)輻射噪聲頻譜圖

▲圖8 遠(yuǎn)場(chǎng)噪聲輻射指向性曲線

▲圖9 渦量的等值面

如圖8所示，在離相對(duì)圓點(diǎn) O'（1.85D，0，8D）為1 m的位置，共取36個(gè)測(cè)量點(diǎn)，指向性角度即為圖1中角度θ。串聯(lián)圓柱體的輻射噪聲指向性具有典型的偶極子聲源特性，在圖中，上下部分對(duì)稱，而右邊數(shù)值大于左邊，因?yàn)樵谄叫杏趤?lái)流的方向上，圓柱體后方的總聲壓級(jí)大于圓柱體前方相應(yīng)位置的總聲壓級(jí)，圓柱體后方的湍流應(yīng)力大于圓柱體前方的湍流應(yīng)力。

4 不同直徑的串聯(lián)圓柱繞流比較

下游圓柱體是整個(gè)串聯(lián)圓柱體的主要噪聲源，改變下游圓柱體的直徑尺寸，可以更好地研究串聯(lián)圓柱體氣動(dòng)噪聲的相關(guān)規(guī)律，對(duì)下游圓柱體直徑分別為0.5D、1.0D、1.5D進(jìn)行計(jì)算，渦量瞬態(tài)等值面、遠(yuǎn)場(chǎng)頻譜特性和噪聲總體聲壓級(jí)分別如圖9、圖10和表1所示。從圖9可以看出，下游圓柱體直徑為0.5D時(shí)，由于圓柱體間距離較遠(yuǎn)，一部分渦沒(méi)有附著在下游圓柱體上，減小了脈動(dòng)力；下游圓柱體直徑為1.5D時(shí)，圓柱體間距太小，上游形成的一部分渦沒(méi)有充分發(fā)展就打在下游圓柱體表面。在下游圓柱體的直徑為1.0D時(shí)，圓柱體間形成相對(duì)最為復(fù)雜的漩渦結(jié)構(gòu)。

在圖10中，當(dāng)下游圓柱體直徑為0.5D、1.0D、1.5D時(shí)，PSD 的最大值出現(xiàn)在 193 Hz、185 Hz、81 Hz處。 在低頻處，1.5D的圓柱體產(chǎn)生的噪聲最大；在高頻處，0.5D的圓柱體產(chǎn)生的噪聲較大。當(dāng)下游圓柱體直徑為1.5D時(shí)，PSD的幅值明顯降低，最高點(diǎn)位置對(duì)應(yīng)的頻率是81 Hz，下游圓柱體直徑較大，渦脫落頻率較小。

表1 各部件產(chǎn)生總聲壓級(jí)對(duì)比/dB

由表1得出，下游圓柱體是主要噪聲元，對(duì)總體氣動(dòng)噪聲起了決定性的貢獻(xiàn)。在下游圓柱體直徑為1.5D時(shí)，遠(yuǎn)場(chǎng)噪聲明顯減小。下游圓柱體直徑0.5D時(shí)，由于圓柱體間距離較遠(yuǎn)，一部分渦沒(méi)有附著在下游圓柱體上，減小了脈動(dòng)力，也使總體氣動(dòng)噪聲略微降低。

▲圖10 下游圓柱體在不同直徑下頻譜特性

5 結(jié)論

（1）仿真計(jì)算結(jié)果與基準(zhǔn)試驗(yàn)結(jié)果吻合，證明了采用LES和FW-H聲類比相結(jié)合的方法對(duì)串聯(lián)圓柱體繞流氣動(dòng)噪聲預(yù)測(cè)的準(zhǔn)確性，采用該方法可對(duì)簡(jiǎn)單起落架等繞流結(jié)構(gòu)進(jìn)行遠(yuǎn)場(chǎng)噪聲預(yù)測(cè)。

（2）串聯(lián)圓柱體的最大噪聲位置出現(xiàn)在渦脫落頻率處，頻譜能量主要位于渦脫落頻率的基頻處以及倍頻處，串聯(lián)圓柱體的氣動(dòng)噪聲屬于寬頻噪聲，輻射聲場(chǎng)的指向性具有明顯的偶極子輻射特性，下游圓柱體是最重要的噪聲貢獻(xiàn)源。

（3）在臨界間距比下，增大下游圓柱體直徑，產(chǎn)生較小的總聲壓級(jí)，噪聲集中在低頻區(qū)域。減小下游圓柱體直徑，總聲壓級(jí)略微降低，噪聲集中在高頻區(qū)域。為有效降低串聯(lián)圓柱體遠(yuǎn)場(chǎng)氣動(dòng)噪聲，在一定間距下，可增大圓柱體直徑尺寸減小氣動(dòng)噪聲。

根據(jù)以上研究結(jié)果，為進(jìn)一步研究串聯(lián)圓柱體在不同的形狀和尺寸下的氣動(dòng)噪聲特性提供了良好基礎(chǔ)。為起落架的氣動(dòng)噪聲的預(yù)測(cè)研究提供前期指導(dǎo)，以便在設(shè)計(jì)和制造階段，優(yōu)化起落架主要噪聲源的結(jié)構(gòu)、尺寸，制造出噪聲更低的飛機(jī)起落架。

［1］ 喬渭陽(yáng)，許開富，武兆偉，等.大型客機(jī)起飛著陸過(guò)程噪聲輻射特性對(duì)比分析［J］.航空學(xué)報(bào)， 2008， 29（5）：534-541.

［2］ Chow L C,Mau K,Remy H.Landing Gears and High Lift Devices Airframe Noise Research ［C］.8th AIAA/CEAS Aeroacoustics Conference & Exhibit, Breckenridge,Colorado,2002.

［3］ 龍雙麗，聶宏，許鑫.不同雷諾數(shù)下圓柱繞流氣動(dòng)噪聲數(shù)值模擬［J］.聲學(xué)技術(shù)， 2011， 30（2）： 111-116.

［4］ 趙良舉，楊南奇，吳朵，等.橫掠二維串列雙圓柱體繞流氣動(dòng)噪聲的數(shù)值模擬 ［J］.重慶大學(xué)學(xué)報(bào)， 2009， 32（8）：943-949.

［5］ 劉敏，劉飛，胡亞濤，等.三維串列雙圓柱繞流氣動(dòng)流場(chǎng)及聲場(chǎng)模擬 ［J］.工程熱物理學(xué)報(bào)，2008，29 （3）：403-406.

［6］ Cox Jared S,Brentner Kenneth S,Rumsey Christopher L.Computation of Vortex Shedding and Radiated Sound for a CircularCylinder:Subcriticalto TranscriticalReynolds Numbers ［J］.Theoretical and Computational Dynamics,1998,12:233-253.

［7］ Lockard D P,Khorrami M R,Choudhari M M.Tandem Cylinder Noise Predictions ［C］. 13th AIAA/CEAS Aeroacoustics Conference （28th AIAA Aeroacoustic Conference）,Rome,Italy，2007.

［8］ Bres G A,F D,Wessels M,et al.Flow and Noise Predictions forthe Tandem CylinderAeroacoustic Benchmark ［J］.Physics of Fluids,2012,24：1-25.

［9］ Uzun A,Hussaini M Y.An Application of Delayed Detached Eddy Simulation to Tandem Cylinder Flow Field Prediction［J］.Computers&Fluids,2012,60：71-85.

［10］ Ffowcs Williams J E,Hawkings D L.Sound Generation by Turbulence and Surfaces in Arbitrary Motion ［J］.Mathematical and Physical Sciences,1969,264:321-342.