顧信忠，李舜酩

?

圓柱繞流氣動噪聲數(shù)值分析

顧信忠1,2，李舜酩1

(1. 南京航空航天大學(xué)能源與動力學(xué)院，江蘇南京210016; 2. 南京航空航天大學(xué)金城學(xué)院車輛系，江蘇南京211156)

為了快速預(yù)測剛性圓柱繞流的氣動噪聲，研究了一種將離散渦方法(DVM)和渦聲理論結(jié)合起來計算低馬赫數(shù)、高雷諾數(shù)流場氣動噪聲的方法。首先用Oseen粘性渦模型改進(jìn)了離散渦方法并模擬了圓柱繞流，分析結(jié)果與實(shí)際情況相符。根據(jù)流場計算的結(jié)果，應(yīng)用渦聲理論進(jìn)一步計算了遠(yuǎn)場的聲壓。測點(diǎn)的總聲壓級與實(shí)驗(yàn)值及其他數(shù)值計算結(jié)果都比較吻合。最后繪制了聲場的指向性特性曲線，表明圓柱繞流聲場明顯的偶極子特性。

氣動噪聲；圓柱繞流；離散渦方法；渦聲理論；Oseen渦

0 引言

圓柱體繞流產(chǎn)生聲音是一種非常普遍的現(xiàn)象。產(chǎn)生的機(jī)理包括：流體與圓柱壁面的相互作用使圓柱產(chǎn)生振動，振動的圓柱排開流體就形成了單極子聲源；當(dāng)旋渦從圓柱表面交替脫落時，會在圓柱表面產(chǎn)生正負(fù)環(huán)量，柱面附近區(qū)域因此產(chǎn)生周期性壓力波動，從而形成了偶極子聲源；圓柱表面邊界層中的渦旋和脫落的渦旋會因湍流應(yīng)力而不斷衰減，此湍流應(yīng)力就是四極子聲源[1]。

圓柱繞流產(chǎn)聲問題在航空航天、航海、風(fēng)工程等實(shí)際工程中非常具有代表性，因此眾多的學(xué)者對此進(jìn)行過深入研究，所采用的方法大致可以分為聲類比方法和直接數(shù)值模擬方法兩大類。

文獻(xiàn)[2]采用有限體積法分別求解二維、三維可壓縮平均雷諾數(shù)N-S方程，對馬赫數(shù)為0.2、雷諾數(shù)為100～5000000的圓柱繞流進(jìn)行了數(shù)值分析，并用Lighthill聲類比方法獲得了遠(yuǎn)場聲信息。文獻(xiàn)[3]用大渦模擬方法模擬了雷諾數(shù)為50000～90000的圓柱繞流湍流流場，提取聲源的強(qiáng)度和分布信息后，用FW-H方程的積分外推法計算了遠(yuǎn)場氣動噪聲，研究了圓柱繞流氣動噪聲的聲場特性。文獻(xiàn)[4]則采用高精度的空間和時間離散算法直接求解N-S方程（DNS方法），計算結(jié)果可以清晰顯示聲波的產(chǎn)生和傳播過程而不受馬赫數(shù)和雷諾數(shù)的限制。文獻(xiàn)[5]同樣采用DNS方法分析了馬赫數(shù)為0.2、雷諾數(shù)為150時圓柱繞流的氣動噪聲，著重研究了渦與聲壓之間的內(nèi)在關(guān)系。

上述兩類方法都需要采用非常細(xì)密的網(wǎng)格和高精度算法對流場進(jìn)行準(zhǔn)確的計算，因此需要大量的存儲空間和很長的計算時間。本文則是基于離散渦方法和渦聲理論，采用Oseen粘性渦模型，研究了可快速預(yù)測繞流產(chǎn)生氣動噪聲的方法。

1 流場模擬

1.1 離散渦方法

離散渦方法通過計算渦量，求解Navier-Stokes方程，以拉格朗日方式跟蹤脫落渦的運(yùn)動軌跡，恰能滿足流場渦聲計算的需要，因此本文采用離散渦方法模擬流場。

離散渦方法的基本思路是將流場分成有旋和無旋兩部分，有旋部分主要集中在柱面的邊界層和分離后的剪切層，并認(rèn)為兩者都是薄薄的連續(xù)渦層，可用離散的渦元來模擬，并將其置于無旋流場中，通過計算離散渦元之間相互作用和演化實(shí)現(xiàn)對流場的快速模擬。當(dāng)柱面的軸向尺寸遠(yuǎn)大于其斷面尺寸時，則可用二維離散渦方法模擬流場[6-7]，如圖1所示。

1.2 渦模型選擇

當(dāng)兩個點(diǎn)渦接近時產(chǎn)生誘導(dǎo)速度非常大，與實(shí)際情況不符，因此很多文獻(xiàn)[8-9]采用Rankin渦模型進(jìn)行光滑處理。但是渦聲與渦量的變化率有關(guān)，而Rankin渦的渦量是恒定的，因此本文采用的是Oseen粘性渦模型，其周向速度為[10]

其中：為渦的環(huán)量；為離渦核中心的距離；為流體的運(yùn)動粘性系數(shù)；為計算時間。

1.3 離散渦方法計算流程

1.3.1 邊界層渦量的確定

將圓柱面用若干個等長首尾相連的直線段來離散，直線段的端點(diǎn)布置適當(dāng)強(qiáng)度的渦元，稱為渦點(diǎn)，直線段的中點(diǎn)稱為控制點(diǎn)。根椐Dirichet邊界條件——物面的切向速度為零，要求任意控制點(diǎn)處的流速與渦元產(chǎn)生的誘導(dǎo)速度相應(yīng)[11]，即

其中：(s,s)為相關(guān)系數(shù)，表示s處單位強(qiáng)度點(diǎn)渦在s控制點(diǎn)處產(chǎn)生的誘導(dǎo)速度；(s)為s渦點(diǎn)處渦的強(qiáng)度；u、v分別為來流的速度分量；β為控制點(diǎn)s的位矢與水平方向的夾角。寫成矩陣形式為

(3)

求解該線性方程組即可得到各渦點(diǎn)處渦元的強(qiáng)度。

1.3.2 渦的脫落

邊界層中的渦會不斷從柱面分離進(jìn)入尾流，本文采用固定分離點(diǎn)的方式產(chǎn)生脫落渦，即預(yù)先在物面上確定幾個渦點(diǎn)產(chǎn)生脫落渦。根據(jù)庫塔(Kutta)條件可確定初生脫落渦的環(huán)量為[12]

其中：u為相應(yīng)渦點(diǎn)處邊界層外的流速；為計算的時間步長。

初生脫落渦的平均對流速度為[13]

1.3.3 渦的對流

已知二維不可壓縮無粘流場中所有渦元（包括邊界層中的渦元和脫落渦）的強(qiáng)度和位置，則可用Biot-Savart公式計算出流場中任意點(diǎn)處的誘導(dǎo)速度為

式中：為渦元在復(fù)平面中的位矢；x、y分別為其水平和豎直分量。

將誘導(dǎo)速度與來流速度疊加就可以得到流場中任意位置的實(shí)際速度為

假設(shè)在時刻渦元的位置為()、()，則可得到渦元對流的計算公式為

(8)

其中：和為渦元因粘性擴(kuò)散產(chǎn)生的運(yùn)動，可采用隨機(jī)走位估算[14]：

式中：和分別是(0，1)和(0，2π)區(qū)間內(nèi)兩個相互獨(dú)立的隨機(jī)數(shù)。97

1.3.4 渦的合并與鏡像

實(shí)際情況下，當(dāng)兩渦元之間距離足夠小時，就會因渦核被破壞而發(fā)生合并。于是Spalart提出了渦的融合機(jī)制[15]，即認(rèn)為當(dāng)兩個渦元滿足一定條件時，就被合并成一個渦元。

式中：0和0都是控制參數(shù)；0控制圓柱附近的渦元數(shù)目；0控制流場中渦元總數(shù)。合并后新渦元的位置和環(huán)量分別為：

， (11)

當(dāng)渦元運(yùn)動到圓柱內(nèi)部時，文中采用鏡像的方法將其強(qiáng)制移到圓柱外部。

1.3.5 作用在圓柱上的力

流體流經(jīng)圓柱表面時會產(chǎn)生相互作用，圓柱上的作用力可用布拉休斯公式計算，將該作用力在豎直方向和水平方向進(jìn)行分解，并除以0.5∞2進(jìn)行無量綱化，即可得到圓柱上的升力和阻力系數(shù)[12]：

式中：Г、u、v、β分別為第渦點(diǎn)處的渦的環(huán)量、柱面的切向速度、法向速度和方位角；∞為平均來流速度；為圓柱直徑。

1.4 流場仿真結(jié)果

本文以實(shí)驗(yàn)驗(yàn)證的模型[16]為研究對象，圓柱直徑為19 mm，空氣流速為69.2 m/s，此時雷諾數(shù)為Re≈90000，馬赫數(shù)約為0.2。圓柱面上均勻布置了40個渦點(diǎn)，其中5個為預(yù)定的脫落點(diǎn)，如圖1所示。設(shè)置計算時間步長為10-5 s，圖2是=0.02 s時的渦元分布圖，此時流場中有3600多個脫落渦元，可見從圓柱脫落的渦并沒有形成穩(wěn)定的卡門渦街，而是擴(kuò)散到尾流區(qū)很寬的范圍內(nèi)形成了紊流，與實(shí)際情況相符。

升力系數(shù)和阻力系數(shù)如圖3和圖4所示，因?yàn)樵谟嬎氵^程中存在大量的渦元合并與鏡像，使得升力系數(shù)的變化范圍非常大，流場穩(wěn)定后升力系數(shù)在0值的上下等幅振蕩。阻力系數(shù)的時間平均值約為1.52，比實(shí)驗(yàn)值1.32[16]大15.2%，比文獻(xiàn)[17]中的大渦模擬結(jié)果1.47大3.4%，如表1所示。說明文中采用的離散渦方法能夠較為準(zhǔn)確地模擬圓柱繞流。

表1 圓柱繞流平均阻力系數(shù)(Re≈90000，馬赫數(shù)≈0.2)

2 聲場計算

2.1 渦聲理論

根據(jù)氣動聲學(xué)理論[18]，在大雷諾數(shù)、低馬赫數(shù)無粘三維流場中，運(yùn)動渦元產(chǎn)生的遠(yuǎn)場聲壓為

式中：0為流體密度；為渦量；為渦的運(yùn)動速度；為聲源發(fā)射聲波的時刻；為緊致格林函數(shù)；1為觀測點(diǎn)的位矢；2為聲源即渦的位矢。

若柱體的軸向尺寸遠(yuǎn)大于其斷面尺寸，可在二維中計算渦聲[18]，此時緊致格林函數(shù)為

其中：2為基爾霍夫矢量；0為聲波傳播速度；H(.)為Heaviside廣義函數(shù)。將式(15)代入式(14)即可得到二維流場的聲壓為

(16)

分別將=1、2代入上式，就可以得到柱體上脈動升力和阻力產(chǎn)生的遠(yuǎn)場聲壓。若聲場中的柱體為圓柱，則基爾霍夫矢量為

2.2 聲場計算結(jié)果

2.2.1 聲壓

為了與實(shí)驗(yàn)值比較，與文獻(xiàn)[16]相同，在遠(yuǎn)場取三個測量點(diǎn)，坐標(biāo)分別為：測量點(diǎn)1(0.665，0)、測量點(diǎn)2(0，-0.665)、測量點(diǎn)3(0，2.432)，如圖5所示，各測點(diǎn)的聲壓如圖6所示。

對測點(diǎn)的聲壓信號進(jìn)行快速傅里葉變換，得到頻域下的聲壓規(guī)律，利用式(18)計算各測點(diǎn)的聲壓級為

(a)

(b)

圖6 各測點(diǎn)聲壓

Fig.6 The time-histories of sound pressures at different testing points

測量點(diǎn)聲壓級的頻譜特性曲線如圖7所示，噪聲頻率分布范圍很寬，且測點(diǎn)的聲壓級在750 Hz左右有明顯的峰值，與實(shí)際流場中渦從柱面脫落的周期相應(yīng)。

還可進(jìn)一步計算各測點(diǎn)的總聲壓級：

得到測量點(diǎn)1總聲壓級為118.52 dB，而實(shí)驗(yàn)測量值為117 dB[16]，誤差僅為1.3%，比大渦模擬的結(jié)果115.72[3]高2.4%，比115.16 dB[17]高2.9%；測量點(diǎn)3 總聲壓級的仿真值為113.74 dB，實(shí)驗(yàn)測量值為100 dB，誤差為13.7%，比大渦模擬的結(jié)果104.52 dB[3]高8.8%，比103.89 dB[17]高9.5%，如表2所示。

表2 圓柱繞流總聲壓級(Re≈90000，馬赫數(shù)≈0.2)

Tabie 2 Overall sound pressure levels at different testing points (Re≈90000, Mach number≈0.2)

研究方法測量點(diǎn)1 測量點(diǎn)3 總聲壓 誤差總聲壓 誤差 本文中方法118.52/ 113.74/ 二維大渦模擬[3]115.722.4% 104.528.8% 二維大渦模擬[17]115.162.9% 103.899.5% 實(shí)驗(yàn)結(jié)果[16]1171.3% 10013.7%

測點(diǎn)1、3的總聲壓級均大于實(shí)測值及大渦模擬的結(jié)果，可能產(chǎn)生誤差的原因包括以下幾個方面：第一，在仿真過程中用大量的離散渦元模擬邊界層和剪切層中的連續(xù)渦層，使得聲壓信號中的高頻成分偏高；第二，為了編程方便，文中將邊界層的離散渦元直接布置在物面上，而當(dāng)聲源靠近物面時，對外輻射的聲壓較大；第三，在計算渦元對流時采用的是精度較低的一階差分格式，可能對計算結(jié)果存在一定的影響?？傮w來看，仿真預(yù)測結(jié)果是較為可信的。

2.2.2 指向性特性

圖8是聲場在=0.019 s、=0.0195 s和=0.02 s時的輻射特性曲線，聲壓最大輻射值都在垂直于來流的軸線上，指向性與偶極子聲源的指向性相近，說明圓柱繞流的氣動噪聲主要是由偶極子聲源引起的。在平行于來流方向仍有一定的聲輻射，這是由作用在圓柱上的阻力周期性波動以及渦元粘性衰減的四極子聲源共同作用的結(jié)果。

=0.019

=0.0195

3 結(jié)論

(1) 采用離散渦方法對二維圓柱繞流流場的模擬結(jié)果顯示，流場圖與實(shí)際情況相符，升力系數(shù)和阻力系數(shù)的計算值都與實(shí)驗(yàn)結(jié)果較為吻合，模擬結(jié)果是可信的。

(2) 圓柱繞流產(chǎn)生的氣動噪聲的峰值頻率與渦的脫落頻率相應(yīng)，總聲壓級與實(shí)驗(yàn)值及大渦模擬仿真結(jié)果比較吻合，聲場的輻射特性與實(shí)際情況相符，聲場的計算結(jié)果是比較準(zhǔn)確的。

(3) 由于采用大量的離散渦元模擬邊界層和剪切層的連續(xù)渦層，獲得的聲信號中的高頻成分偏高。

(4) 采用離散渦方法和渦聲理論相結(jié)合的方法，對聲場進(jìn)行快速預(yù)測的方法是可行的。

[1] 張強(qiáng). 氣動聲學(xué)基礎(chǔ)[M]. 北京: 國防工業(yè)出版社, 2012.

ZHANG Qiang. Foundation of Acoustics[M]. Beijing: National Defense Industry Press, 2012.

[2] COX J S, Brentner K S, Rumsey C L. Computation of vortex shedding and radiated sound for a circular cylinder: subcritical to transonic Reynolds numbers[J]. Theoret. Comput. Fluid Dyn, 1988, 12(4): 233-253.

[3] 龍雙麗, 聶宏, 許鑫. 不同雷諾數(shù)下圓柱繞流氣動噪聲數(shù)值模擬[J]. 聲學(xué)技術(shù), 2011, 30(2): 111-116.

LONG Shuangli, NIE Hong, XU Xin. Numerical simulation of noise induced by flow around a cylinder at different Reynolds number[J]. Technical Acoustics, 2011, 30(2): 111-116.

[4] Mitchell B E, Lele S K, Moin P. Direct computation of the sound generated by vortex pairing in an axisymmetric jet[J]. J. Fluid Mech, 1999, 383(6): 113-142.

[5] Inoue O, Hatakeyama N. Sound generation by a two-dimensional circular cylinder in a uniform flow[J]. Journal of Fluid Mechanics, 2002, 471(01): 285-314.

[6] MUSTTO A A, Hirata M H, Bodstein G C R. Discrete vortex method simulation of the flow around a circular cylinder with and without rotation[J]. AIAA paper, 1998, 98-2409: 59-69.

[7] Etoh F, Kamemoto K, Matsumoto H, et al. Numerical simulation of flow around a rotary oscillating foil with constant amplitude angle by use of the vortex method[C]//Proc. 11th Symp. on CFD, Tokyo, 1997: 385-386.

[8] 潘巖松. 高層建筑二維流場的離散渦方法數(shù)值模擬[D]. 武漢: 華中科技大學(xué), 2004.

PAN Yansong. Numerical simulation of two-dimensional flow past a tall building using a discrete vortex method[D]. Wuhan: Huazhong University of Science & Technology, 2004.

[9] 徐曉亮, 黃海明, 章梓茂. 燒蝕環(huán)境下的圓柱繞流計算模型[J]. 空氣動力學(xué)學(xué)報, 2011, 29(2): 170-176.

XU Xiaoliang, Huang Haiming, Zhang Zimao. Numerical models for flow past a cylinder in the ablation environment[J]. Acta Aerodynam Ica Sinica, 2011, 29(2): 170-176.

[10] Schaefer J W, Eskinazi S. An analysis of the vortex street generated in a viscous fluid[J]. J. Fluid Mech.,1959, 6(2): 241-260.

[11] Lewis R I, Vortex element methods for fluid dynamic analysis of engineering systems[M]. Cambridge University Press, 2005.

[12] Sarpkaya T, Shoaff R L. An inviscid model of two dimensional vortex shedding for transient and asymptotically-steady separated flow over a cylinder [J]. AIAA Paper, 1979, 79-0281: 1-9.

[13] 童秉綱, 夏南, 李潛.物體繞流的離散渦方法[J]. 力學(xué)進(jìn)展, 1985, 15(3): 318-328.

TONG Bingang, XIA Nan, LI Qian. Discrete Vortex Method in external flow problems [J].Advances in Mechanics, 1985, 15(3): 318-328.

[14] 陳偉, 宗智. 二維圓柱繞流的離散渦數(shù)值模擬[J]. 艦船科學(xué)技術(shù), 2010, 32(5): 111-115.

CHEN Wei, ZONG Zhi. Numerical simulation of two- dimensional flow around circular cylinder using discrete vortex method[J]. Ship SciencE and Technology, 2010, 32(5): 111-115.

[15] SPALARTR S, LEONARD A, BAGANOFF D. Numerical simulation of separated flows[D]. Stanford University, California, 1983.

[16] Revell J D, Prydz R A, Hays A P. Experimental study of aerody-namic noise vs drag relationships for circular cylinders[J]. Anz Journal of Surgery, 1978, 16(6): 889-897.

[17] YANG Dangguo, LI Jianqiang, LIU Jun. Analysis on physical mechanism of sound generation inside cavities based on acoustic analogy method[J]. Journal of Fluid Dynamics, 2013, 3(1): 23-31.

[18] HOWE M S. Theory of vortex sound [M]. Cambridge Texts in Applied Mathematics, 2002.

Numerical analysis of noise induced by flow around a cylinder

GU Xin-zhong1,2, LI Shun-ming1

(1. College of Energy and Power Engineering,Nanjing University of Aeronautics &Astronautics, Nanjing 210016,Jiangsu, China；2. Department of Vehicle Engineering, Nanhang Jincheng College, Nanjing211156, Jiangsu,China)

The method of combination of Discrete Vortex Method (DVM) and the theory of vortex sound are used to predict two-dimensional far field noise induced by flow around a cylinder. Firstly, the field of flow around a cylinder is simulated using the Oseen vortex model. The simulation results are in agreement with experimental data. Based on the simulation results, the characteristic of the acoustic field is studied. The overall sound pressure levels of receivers also accord wellwith both experimental data and other numerical simulation results. The polar diagrams of the aeroacoutic pressure are plotted. And the dipolar nature of the generated sound is confirmed.

aerodynamic noise; flow around a cylinder; Discrete Vortex Method(DVM); theory of vortex sound; oseen vortex

TB533

A

1000-3630(2016)-02-0095-06

10.16300/j.cnki.1000-3630.2016.02.003

2015-04-10;

2015-06-10

顧信忠(1980－), 男, 江蘇儀征人, 博士研究生, 研究方向?yàn)闅鈩釉肼曨A(yù)測與控制技術(shù)。

顧信忠, E-mail: guxingzhong007@126.com