開式方腔可壓縮流動(dòng)及噪聲LES計(jì)算分析

馮 峰，王 強(qiáng)

(中國(guó)航天空氣動(dòng)力技術(shù)研究院，空氣動(dòng)力科學(xué)中心， 北京 100074)

1 引言

戰(zhàn)斗機(jī)和轟炸機(jī)內(nèi)埋武器艙可壓縮空腔流動(dòng)噪聲可超過(guò)160 dB，能夠致使機(jī)載電子設(shè)備失靈、艙內(nèi)結(jié)構(gòu)疲勞損傷，甚至誘發(fā)武器在分離過(guò)程中非常規(guī)運(yùn)動(dòng)，導(dǎo)致投放精度降低[1]。此外，空腔流強(qiáng)噪聲輻射還會(huì)影響商用航空運(yùn)輸?shù)脑肼曔m航性及軍用戰(zhàn)斗機(jī)、轟炸機(jī)飛行隱蔽性[2]。20世紀(jì)50年代，Krishnamurty[3]最早開始研究空腔流動(dòng)現(xiàn)象，并指出空腔流中存在聲-流耦合共振的聲反饋環(huán)機(jī)制。當(dāng)高速氣流掠過(guò)空腔上方表面后，形成高度有序的剪切層渦列結(jié)構(gòu)，該渦列與空腔后壁撞擊產(chǎn)生反射聲波擾動(dòng)空腔前緣初始剪切層，當(dāng)擾動(dòng)波與反射壓力波頻率相匹配時(shí)，空腔內(nèi)出現(xiàn)聲場(chǎng)-流場(chǎng)相互耦合的聲反饋環(huán)流動(dòng)振蕩現(xiàn)象。之后，Rossiter[4]發(fā)展了半經(jīng)驗(yàn)公式對(duì)聲-流共振現(xiàn)象提供了更完整的描述，并可預(yù)測(cè)空腔流振蕩頻率，在空腔流研究中得到廣泛應(yīng)用[5]。

針對(duì)超聲速飛行條件下戰(zhàn)斗機(jī)內(nèi)埋武器分離投放可靠性問(wèn)題，高保真數(shù)值技術(shù)在探索超聲速空腔流機(jī)制方面展現(xiàn)出良好的前景。其中LES(large eddy simulation，大渦模擬)方法可以精確求解空腔內(nèi)各類渦結(jié)構(gòu)及其產(chǎn)生的噪聲特性，具有很強(qiáng)的預(yù)測(cè)分析空腔流噪聲機(jī)制的能力。Sinha等[6]和Smith等[7]較早采用三維LES方法分別計(jì)算了Mach 2.0和Mach 1.5超聲速開式空腔流動(dòng)噪聲。其中，Smith等[7]的LES結(jié)果顯示壓力頻譜峰值與實(shí)驗(yàn)結(jié)果一致。Rizzetta和Visbal[8]采用了基于動(dòng)態(tài)Smagorinsky亞格子模型的LES技術(shù)求解了長(zhǎng)高比L/D=5，Mach數(shù)1.19，Reynolds數(shù) 200 000的超聲速方腔流動(dòng)，并應(yīng)用了高精度緊致格式和大量計(jì)算網(wǎng)格，獲得的方腔底部多個(gè)位置上壓力譜與實(shí)驗(yàn)結(jié)果在幅值和特征頻率上趨勢(shì)一致。但上述研究缺乏細(xì)致的流場(chǎng)動(dòng)力學(xué)分析，近期在深化空腔流機(jī)理認(rèn)知以改進(jìn)噪聲控制的需求牽引下，人們?cè)俅伍_始關(guān)注展向受限條件下空腔流動(dòng)基本特性。Brès和Colonius[9]使用DNS(direct numerical simulation，直接數(shù)值模擬)方法和整體穩(wěn)定性理論分析了低Reynolds數(shù)1500、亞聲速0.325空腔流動(dòng)展向穩(wěn)定性，并發(fā)現(xiàn)三維空腔模態(tài)具有一個(gè)與方腔寬度相當(dāng)?shù)恼瓜虿ㄩL(zhǎng)，其頻率比二維Rossiter模態(tài)低一個(gè)量級(jí)。Sun等[10]完成了一個(gè)有限寬度超聲速開式方腔流動(dòng)LES計(jì)算分析，研究了側(cè)壁效應(yīng)、自由來(lái)流Mach數(shù)、Reynolds數(shù)等對(duì)空腔流動(dòng)的影響。整體而言，目前對(duì)可壓縮空腔流動(dòng)及噪聲特性及機(jī)理細(xì)節(jié)認(rèn)識(shí)仍不足。本文采用基于計(jì)算氣動(dòng)聲學(xué)技術(shù)的LES方法計(jì)算分析三維Mach 1.4和Mach 0.6超/亞聲速來(lái)流條件下可壓縮開式方腔流動(dòng)。針對(duì)空腔流動(dòng)中剪切層演化、多尺度渦結(jié)構(gòu)、聲反饋、可壓縮效應(yīng)等基本動(dòng)力學(xué)過(guò)程開展分析，以深化可壓縮空腔流動(dòng)及噪聲現(xiàn)象認(rèn)識(shí)和機(jī)制理解。

2 數(shù)值方法

2.1 控制方程及SGS模型

為實(shí)現(xiàn)LES方法，采用無(wú)量綱化三維可壓縮Favre濾波Navier-Stokes方程

(1)

(2)

(3)

(4)

(5)

(6)

2.2 計(jì)算格式

控制方程采用基于計(jì)算氣動(dòng)聲學(xué)的高精度數(shù)值格式進(jìn)行離散，空間格式使用Tam等[11-12]提出的低色散低耗散誤差7點(diǎn)4階DRP格式，時(shí)間推進(jìn)采用Berland等[13]提出的優(yōu)化低存儲(chǔ)6步4階Runge-Kutta格式。由于DRP數(shù)值格式不能精確離散高波數(shù)短波，為避免數(shù)值高頻波導(dǎo)致Gibbs振蕩，引入Tam[11-12]構(gòu)造的人工選擇性阻尼項(xiàng)對(duì)非物理數(shù)值短波進(jìn)行濾波。選擇性阻尼項(xiàng)、黏性項(xiàng)及SGS項(xiàng)在計(jì)算中均在Runge-Kutta最后一個(gè)子步引入。對(duì)于含強(qiáng)湍流脈動(dòng)和激波間斷的超聲速空腔流場(chǎng)，采用Bogey等[14]提出的自適應(yīng)空間濾波法對(duì)守恒變量濾波以捕捉激波并維持計(jì)算穩(wěn)健性。

2.3 模型設(shè)置

計(jì)算模型為一個(gè)有限寬度的開式方腔，如圖1所示，方腔長(zhǎng)深比L/D=6，寬深比W/D=2，D、L、W分別為方腔深度、長(zhǎng)度和寬度。以方腔深度為特征長(zhǎng)度的Reynolds數(shù)Re=10 000，來(lái)流邊界層厚度為δ=0.167。

圖1 方腔模型及計(jì)算網(wǎng)格示意圖

亞/超聲速方腔流動(dòng)使用相同Cartesian網(wǎng)格，900萬(wàn)個(gè)網(wǎng)格點(diǎn)求解空腔流動(dòng)及噪聲?？涨粌?nèi)設(shè)置199×103×91個(gè)網(wǎng)格點(diǎn)，計(jì)算域范圍為x1=0～6，x2=-1～0，x3=-1～1。空腔外計(jì)算域339×141×183個(gè)網(wǎng)格點(diǎn)，計(jì)算域范圍為x1=-3～13，x2=-1～9，x3=-2.5～2.5。壁面附近對(duì)網(wǎng)格進(jìn)行加密，各方向最小網(wǎng)格尺度為Δx1=0.002，Δx2=0.001，Δx3=0.002。超聲速方腔算例進(jìn)行500 000次迭代，時(shí)間步長(zhǎng)Δt=0.000 5，流場(chǎng)統(tǒng)計(jì)在第100 000步后啟動(dòng)。亞聲速方腔算例進(jìn)行250 000次迭代，時(shí)間步長(zhǎng)Δt=0.000 5。為加速計(jì)算，LES程序采用MPI方法，使用84個(gè)進(jìn)程，每網(wǎng)格點(diǎn)每次迭代約花費(fèi)0.5 μs CPU計(jì)算時(shí)間，一次計(jì)算花費(fèi)約480 h。

2.4 Rossiter空腔頻率預(yù)測(cè)模型

Rossiter[4]于1964年提出了預(yù)測(cè)空腔流噪聲特征頻率的半經(jīng)驗(yàn)公式，或稱Rossiter模態(tài)預(yù)測(cè)，其形式如下

(7)

式中:fm為第m個(gè)模態(tài)的頻率，u、Ma為自由來(lái)流速度和Mach數(shù)，γ為比熱比，α為相位移。κν為大尺度渦結(jié)構(gòu)對(duì)流速度與自由來(lái)流速度的比值，是依賴于空腔幾何形狀和測(cè)試條件常數(shù)。本文LES計(jì)算中在方腔底部(3，-1，0)和中心(3，-0.5，0)設(shè)置兩個(gè)測(cè)壓點(diǎn)，并將計(jì)算壓力脈動(dòng)值與Rossiter預(yù)測(cè)結(jié)果對(duì)比，相互驗(yàn)證并確認(rèn)空腔噪聲特征頻率。

3 驗(yàn)證與校核

將計(jì)算統(tǒng)計(jì)得到的空腔流時(shí)均流向速度和湍流脈動(dòng)速度與Dudley和Ukeiley[15]的實(shí)驗(yàn)結(jié)果和DES數(shù)值結(jié)果進(jìn)行對(duì)比，以驗(yàn)證LES結(jié)果的準(zhǔn)確性。圖2展示了相同來(lái)流條件下本文LES計(jì)算和Dudley和Ukeiley的實(shí)驗(yàn)測(cè)試的方腔x3=0截面上平均流向速度，兩圖等值云圖范圍均在u/U∞=-0.25～1之間，增量均為0.25。圖中LES和實(shí)驗(yàn)獲得的超聲速空腔流流向速度沿空間分布完全一致，均展現(xiàn)了明顯的剪切流動(dòng)特征，即方腔上方剪切層沿下游逐漸增厚，受方腔后壁遮擋而浸入腔內(nèi)的動(dòng)力學(xué)現(xiàn)象，空腔前部及底部平均流向速度相對(duì)較低。

圖2 時(shí)均流向速度場(chǎng)云圖

圖3展示了方腔x3=0截面上流向脈動(dòng)速度均方根，LES結(jié)果與文獻(xiàn)中DES結(jié)果[15]采用等值云圖范圍一致，均為〈u1rms〉/U∞=0～0.25，增量0.05。計(jì)算與實(shí)驗(yàn)結(jié)果均顯示了湍流脈動(dòng)速度沿剪切層向空腔后部逐漸增強(qiáng)，在方腔后部剪切層區(qū)域湍流脈動(dòng)達(dá)到最大值的。受方腔后緣誘導(dǎo)，后緣表面也形成一個(gè)小范圍的強(qiáng)湍流區(qū)。值得注意的是，在靠近空腔后底部小的湍流區(qū)域，文獻(xiàn)DES獲得的湍流脈動(dòng)速度強(qiáng)度高于本文LES計(jì)算統(tǒng)計(jì)結(jié)果。這是由于DES計(jì)算的底壁面邊界層較薄，易于在腔內(nèi)剪切層于后壁回流誘導(dǎo)的逆壓梯度下，形成較大范圍湍流分流區(qū)，致使后底壁區(qū)呈現(xiàn)較LES結(jié)果高的湍流脈動(dòng)強(qiáng)度。

圖4展示了方腔內(nèi)7個(gè)流向站位上LES與實(shí)驗(yàn)測(cè)量獲得的時(shí)均流向速度型和流向脈動(dòng)速度均方根。如圖所示，二者很好地吻合，表明LES計(jì)算準(zhǔn)確地捕捉到了剪切層和方腔流動(dòng)對(duì)流場(chǎng)平均速度及湍流脈動(dòng)的影響。在空腔下游第六和第七個(gè)站位處，數(shù)值和實(shí)驗(yàn)的時(shí)均流向速度和流向脈動(dòng)速度均方根在靠近空腔底壁面時(shí)出現(xiàn)最大約10%的相對(duì)誤差。對(duì)比本文LES及文獻(xiàn)DES結(jié)果[15]，發(fā)現(xiàn)在空腔后底壁面區(qū)域，LES計(jì)算湍流脈動(dòng)速度強(qiáng)度略低于實(shí)驗(yàn)測(cè)試結(jié)果，而DES結(jié)果則略高于實(shí)驗(yàn)結(jié)果。可見未來(lái)仍需不斷深入探索以完善更精確的空腔流動(dòng)數(shù)值分析方法。

圖3 時(shí)均流向脈動(dòng)速度均方根云圖

圖4 空腔內(nèi)不同流向站位時(shí)均流向速度和 流向脈動(dòng)速度均方根示意圖

4 LES結(jié)果與分析

4.1 瞬態(tài)流場(chǎng)

LES方法可精確計(jì)算空腔流動(dòng)渦、聲、激波等全部流場(chǎng)細(xì)節(jié)，便于對(duì)復(fù)雜振蕩空腔流動(dòng)的演化特征和機(jī)制開展分析探索。如圖5，使用Q準(zhǔn)則等值面(Q=3)和壓力云圖展示了Mach 1.4超聲速開式空腔瞬態(tài)流場(chǎng)。來(lái)流進(jìn)入計(jì)算域后，在壁面誘導(dǎo)的氣流壓縮效應(yīng)下計(jì)算域前緣形成一道激波。氣流經(jīng)過(guò)空腔前緣后，與腔內(nèi)氣體混合使得邊界層向剪切層轉(zhuǎn)化，且由于流動(dòng)不穩(wěn)定性剪切層卷起形成展向渦結(jié)構(gòu)。渦結(jié)構(gòu)在向下游運(yùn)動(dòng)過(guò)程中沿展向失穩(wěn)，產(chǎn)生大量Λ渦及其他更小尺度的結(jié)構(gòu)，在高速來(lái)流裹挾下繼續(xù)向下游對(duì)流。在空腔后緣，剪切層被撞擊產(chǎn)生大量小尺度湍流渦結(jié)構(gòu)，一部分侵入空腔后部與腔內(nèi)流體急劇混合，另一部分溢出空腔沿壁面繼續(xù)向下游流出計(jì)算域。剪切層撞擊空腔后緣還形成向上游傳播的聲波，聲波輻射出空腔后被超聲速來(lái)流的對(duì)流作用下向下游遠(yuǎn)場(chǎng)傳播，在腔內(nèi)聲波穿過(guò)亞聲速流動(dòng)輻射至空腔前壁面。

圖5 瞬態(tài)Q準(zhǔn)則等值面(Q=3)及壓力云圖

圖6展示了空腔中心x1x2截面(x3=0)上4個(gè)連續(xù)時(shí)刻密度梯度幅值云圖，由圖可見復(fù)雜超聲速空腔流場(chǎng)動(dòng)力學(xué)演化。首先，圖中清晰地展示了剪切層大尺度結(jié)構(gòu)的運(yùn)動(dòng)發(fā)展過(guò)程。在t=139時(shí)刻大尺度結(jié)構(gòu)形成并向下游對(duì)流，t=139.5時(shí)刻該大尺度結(jié)構(gòu)向方腔后緣靠近，同時(shí)新的渦結(jié)構(gòu)卷起，t=140時(shí)刻原大尺度結(jié)構(gòu)繼續(xù)逼近方腔后緣且頭部結(jié)構(gòu)開始變形，新的渦結(jié)構(gòu)顯著增強(qiáng)，t=140.5原大尺度結(jié)構(gòu)與方腔后緣撞擊，新大尺度渦結(jié)構(gòu)形成并向下游對(duì)流，進(jìn)而完成了一次大尺度渦結(jié)構(gòu)演化周期。同時(shí)，如圖6(c)所示，可觀察到剪切層附近的亞聲速流動(dòng)區(qū)域內(nèi)形成向上游輻射的聲波，其波長(zhǎng)與大尺度渦結(jié)構(gòu)相當(dāng)。這證實(shí)了空腔振蕩流動(dòng)中剪切層與空腔后壁面撞擊產(chǎn)生向空腔前方傳播的聲波，并擾動(dòng)空腔前緣剪切層，形成自持振蕩的聲反饋環(huán)現(xiàn)象。此外，超聲速來(lái)流在開式空腔誘導(dǎo)下在空腔前緣也形成一道激波，如圖6(a)所示，由于與剪切層相互作用，其始終處于振蕩狀態(tài)，并反向作用影響剪切層的發(fā)展演化。整體而言，三維超聲速開式空腔流動(dòng)的聲反饋環(huán)特性較為顯著。

圖6 密度梯度幅值(紋影圖)演化云圖

4.2 空腔壓力脈動(dòng)

在空腔底部壁面中心P1點(diǎn)(3，-1，0)和空腔空間中心P2點(diǎn)(3，-0.5，0)兩個(gè)典型位置設(shè)置測(cè)壓點(diǎn)以觀測(cè)腔內(nèi)壓力變化。如圖7(a)所示，P2點(diǎn)壓力脈動(dòng)幅值略高于P1點(diǎn)，這是由于湍流在空腔空間中的脈動(dòng)強(qiáng)度強(qiáng)于空腔底壁所致。使用了Fourier變換計(jì)算兩點(diǎn)上壓力頻譜，由于壓力信號(hào)計(jì)算及存儲(chǔ)時(shí)間有限，對(duì)原始?jí)毫γ}動(dòng)信號(hào)做Hamming窗濾波處理以方便識(shí)別主頻信息。

圖7 空腔內(nèi)P1、P2兩測(cè)點(diǎn)的壓力脈動(dòng)和頻譜曲線

從圖7(b)可見，壓力脈動(dòng)頻譜具有多個(gè)峰值，其中前3個(gè)最優(yōu)峰值頻率被標(biāo)記為f1=0.88、f2=1.88和f3=2.84(圓頻率)。兩觀測(cè)點(diǎn)上前2個(gè)主導(dǎo)頻率一致，但第3個(gè)主導(dǎo)頻率在空腔中心位置較為顯著，展示了腔內(nèi)復(fù)雜的非線性渦系運(yùn)動(dòng)特性。結(jié)合第4.1節(jié)流場(chǎng)演化特征，分析可知主頻f2與聲反饋導(dǎo)致的剪切層渦運(yùn)動(dòng)相關(guān)，而頻率f1與空腔后部運(yùn)動(dòng)較緩慢的大尺度渦運(yùn)動(dòng)相關(guān)，據(jù)信f3應(yīng)是2種運(yùn)動(dòng)交互作用衍生頻率。值得說(shuō)明的是，相似的脈動(dòng)主頻在Dudley和Ukeiley[15]和Sheta等[16]的工作中也被發(fā)現(xiàn)過(guò)，表明超聲速空腔流存在固有動(dòng)力學(xué)特性。

此外，如圖7(b)所示，虛線為Rossiter公式[4]預(yù)測(cè)的空腔Rossiter模態(tài)，其中設(shè)定常參數(shù)κν=0.66，相位移由經(jīng)驗(yàn)公式α=0.062×L/D計(jì)算得到，本文方腔模型α=0.397。LES計(jì)算得到的3個(gè)壓力脈動(dòng)主導(dǎo)頻率與Rossiter預(yù)測(cè)頻率非常吻合，證明Rossiter公式對(duì)超聲速空腔流動(dòng)共振模態(tài)也具有預(yù)測(cè)能力。

4.3 壓縮性效應(yīng)

為明確壓縮效應(yīng)對(duì)空腔流動(dòng)的影響，對(duì)Mach 0.6亞聲速空腔流進(jìn)行LES計(jì)算及對(duì)比分析。值得說(shuō)明的是，除來(lái)流速度外，該算例其余來(lái)流參數(shù)以及計(jì)算網(wǎng)格均與Mach 1.4超聲速空腔流一致。圖8展示了亞聲速空腔流動(dòng)密度梯度幅值云圖隨時(shí)間的演化。由圖8(a)可見，該流動(dòng)仍然受剪切層渦結(jié)構(gòu)運(yùn)動(dòng)主導(dǎo)，在空腔前緣位置卷起一列規(guī)則小渦。約在x1=2～5站位，細(xì)小渦積聚形成大尺度渦團(tuán)結(jié)構(gòu)。如圖8(b)，包含大量小渦的大尺度渦團(tuán)開始向空腔后緣撞擊，破碎成大量的小湍流渦，同時(shí)形成向上游傳播的聲波。聲波在空腔中沿腔內(nèi)向上游傳播，則與初始剪切層耦合形成聲-流共振狀態(tài)。聲波向腔外輻射，形成上游遠(yuǎn)場(chǎng)占優(yōu)的強(qiáng)聲場(chǎng)，這與超聲速空腔流中聲波無(wú)法向上游傳播不同。另外，由于沒有激波出現(xiàn)，缺少了激波與剪切層相互干擾現(xiàn)象。

圖8 亞聲速(Ma=0.6)空腔流動(dòng)密度梯度幅值云圖

圖9 亞/超聲速空腔流壓力脈動(dòng)及其頻譜曲線

5 結(jié)論

1) 超聲速(Mach 1.4)空腔流動(dòng)LES數(shù)值結(jié)果與文獻(xiàn)中實(shí)驗(yàn)及數(shù)值結(jié)果對(duì)比，表明基于低色散低耗散計(jì)算氣動(dòng)聲學(xué)技術(shù)的LES方法對(duì)空腔流計(jì)算的精確性、魯棒性，解決了常規(guī)CFD方法或商用軟件計(jì)算精度低不易捕捉噪聲及邊界非物理聲反射等問(wèn)題。

2) 超聲速來(lái)流條件下，剪切層大尺度渦結(jié)構(gòu)與空腔后壁發(fā)生撞擊并產(chǎn)生經(jīng)腔內(nèi)向前傳播的反射聲波，該聲波與空腔前緣初始剪切層相互作用，導(dǎo)致空腔內(nèi)形成聲-流耦合自持振蕩的聲反饋環(huán)現(xiàn)象?；谠撜J(rèn)識(shí)，可設(shè)置專門裝置在超聲速戰(zhàn)斗機(jī)武器艙等工程應(yīng)用中破壞來(lái)流剪切層內(nèi)大尺度渦結(jié)構(gòu)的一致性，以削弱聲-流耦合現(xiàn)象，降低腔內(nèi)噪聲水平。

3) 超聲速空腔流展現(xiàn)出復(fù)雜非線性渦系運(yùn)動(dòng)導(dǎo)致多個(gè)特征頻率。觀測(cè)發(fā)現(xiàn)前3個(gè)頻率較為顯著，與流場(chǎng)演化歷程結(jié)合分析可知第2個(gè)主頻f2最優(yōu)，其與聲反饋導(dǎo)致的剪切層渦運(yùn)動(dòng)相關(guān)，而第1個(gè)主頻f1與空腔后部較緩慢的大尺度渦運(yùn)動(dòng)相關(guān)，第3個(gè)主頻f3是前2種運(yùn)動(dòng)相互作用衍生。

4) 當(dāng)來(lái)流壓縮效應(yīng)減弱時(shí)，亞聲速(Mach 0.6)空腔流動(dòng)剪切層渦團(tuán)與后緣撞擊形成的聲-流共振現(xiàn)象與超聲速情形基本一致，但亞聲速空腔流動(dòng)的壓力脈動(dòng)強(qiáng)度明顯減弱，且流動(dòng)振蕩主導(dǎo)頻率減小，表明腔內(nèi)渦運(yùn)動(dòng)緩慢。