李曉東，姜蔚婉，2，屈 崑，*，蔡晉生

(1. 西北工業(yè)大學(xué) 航空學(xué)院，西安 710072；2. 中國航空工業(yè)集團(tuán) 成都飛機(jī)設(shè)計研究所，成都 610091)

0 引 言

內(nèi)埋式武器在投放分離過程中，彈艙將演變?yōu)椋◣撻T）大尺度空腔。大尺度空腔是航空航天武器裝備必不可少且普遍采用的一種典型結(jié)構(gòu)形式，主要用于作戰(zhàn)武器裝載、動力燃油存儲和起落架存放等，具有尺度大、結(jié)構(gòu)復(fù)雜等特點(diǎn)。在動載荷激勵下，空腔薄壁結(jié)構(gòu)容易產(chǎn)生嚴(yán)重的振動和噪聲，不僅影響武器裝備的作戰(zhàn)效能，甚至危及武器裝備安全。因此，空腔繞流是空氣動力學(xué)領(lǐng)域的一個重要且典型問題。

針對空腔中流動和噪聲的復(fù)雜現(xiàn)象與物理機(jī)理，國內(nèi)外科研工作者通過風(fēng)洞試驗(yàn)和計算手段對干凈彈艙（空腔）模型[1-4]進(jìn)行了廣泛深入的研究。具體包括分析流致振蕩現(xiàn)象的成因[5-6]，對空腔類型[7]、來流馬赫數(shù)、邊界層厚度等[3,8]參數(shù)進(jìn)行詳細(xì)研究，并探討了多種不同的控制措施[2,9]，對空腔設(shè)計具有重要的指導(dǎo)意義。Rossiter[10]提出了預(yù)測空腔噪聲模態(tài)頻率的半經(jīng)驗(yàn)?zāi)Ｐ?；Heller 等[11]改進(jìn)了此模型，將其推廣到高馬赫數(shù)（Ma= 2、3）的可壓縮流動中。

空腔是彈艙的一種理想模型，實(shí)際中艙門的存在會使得空腔內(nèi)非定常流動更加復(fù)雜，因而研究學(xué)者們對帶艙門的空腔流動開展了研究和分析。首先是對帶有靜態(tài)艙門的空腔流動的研究：Blair 和Stallings[12]采用風(fēng)洞試驗(yàn)研究了不同艙門開度和艙門高度對導(dǎo)彈類物體在超聲速空腔中投放的影響，結(jié)果表明艙門開度對彈體的氣動載荷有顯著影響，并且降低艙門高度會減少彈體的氣動載荷；Lawson 和Barakos[13]采用脫體渦模擬（DES）方法對艙門開啟90°的空腔標(biāo)模M219 和1303 無人機(jī)模型進(jìn)行了模擬，并分析了艙門打開對彈體拋投的影響。Casper 等[14]采用試驗(yàn)的方式模擬了帶有艙門的復(fù)雜空腔中彈體拋投的流固耦合問題，發(fā)現(xiàn)增加空腔復(fù)雜度會引發(fā)彈體展向共振。其次在動態(tài)艙門的研究上，Sheta 等[15]采用RANS/LES混合方法，模擬了超聲速（Ma= 1.44）空腔艙門從5°打開到35°過程中的流動，發(fā)現(xiàn)在艙門運(yùn)動到30°時，剪切層和艙內(nèi)聲波產(chǎn)生較強(qiáng)的相互作用，腔內(nèi)壁面出現(xiàn)強(qiáng)的壓力脈動；Loupy 等[16]對比分析了固定開度艙門和動態(tài)艙門對馬赫數(shù)為0.85 的空腔流動、載荷和氣動噪聲的影響，發(fā)現(xiàn)艙門打開過程中的過渡階段，艙門氣動載荷增大，流動的脈動和總聲壓級增強(qiáng)，對空腔結(jié)構(gòu)產(chǎn)生顯著影響。

國內(nèi)多位學(xué)者采用RANS 方法模擬分析靜態(tài)艙門[17]和艙門開啟過程中[18]空腔動態(tài)氣動特性和載荷的變化。吳繼飛等[19-20]對Ma= 0.6 空腔艙門的開閉問題開展了研究，設(shè)計了不同的運(yùn)動方式，發(fā)現(xiàn)艙門運(yùn)動時氣動載荷發(fā)生明顯變化，艙門打開30°時氣動載荷最大，并分析了腔內(nèi)和艙門的總聲壓級和聲壓譜特性。由于帶艙門的空腔流動非定常效應(yīng)顯著，且流動的湍流結(jié)構(gòu)復(fù)雜，所以傳統(tǒng)的RANS方法捕捉復(fù)雜流動結(jié)構(gòu)比較困難。閆盼盼等[21]采用基于k?ω的分離渦模擬方法（IDDES）計算了艙門不同開啟姿態(tài)對內(nèi)埋武器投放的影響，發(fā)現(xiàn)艙門會改變超聲速流場中的波系結(jié)構(gòu)，對彈體偏航影響較大。

考慮運(yùn)動艙門時，空腔中氣動與噪聲特性會變得更加復(fù)雜，預(yù)測難度增加，而國內(nèi)對帶運(yùn)動艙門的復(fù)雜空腔流動與噪聲的研究較少。為了研究艙門運(yùn)動情況下的空腔復(fù)雜流動現(xiàn)象，深入分析靜、動態(tài)艙門對空腔復(fù)雜流動的影響，本文采用嵌套重疊網(wǎng)格方法對艙門和空腔分別生成網(wǎng)格， 采用改進(jìn)的脫體渦模擬方法對亞聲速（Ma= 0.6）空腔流場進(jìn)行數(shù)值模擬，從流場特性、湍流結(jié)構(gòu)、脈動特性上分析艙門對湍流流場和腔內(nèi)噪聲的影響。

1 數(shù)值方法

針對空腔流動強(qiáng)脈動、強(qiáng)湍流的特點(diǎn)，采用RANS/LES 混合方法進(jìn)行計算，本文采用改進(jìn)的脫體渦湍流模擬方法（Improved Delayed DES，IDDES）[22]，解析空腔內(nèi)分離區(qū)域的湍流結(jié)構(gòu)和脈動特性。基于以上的數(shù)值方法，本文發(fā)展了一套基于有限體積方法的多塊結(jié)構(gòu)網(wǎng)格求解器，并通過多種復(fù)雜流動問題驗(yàn)證其準(zhǔn)確性和可靠性[23-26]。

1.1 高精度低耗散格式

式中，

1.2 空腔標(biāo)模（M219）驗(yàn)證

為驗(yàn)證高精度數(shù)值格式在空腔問題中的有效性，采用本文提出的格式對馬赫數(shù)0.85 空腔標(biāo)模M219進(jìn)行模擬分析，并與參考文獻(xiàn)[32]的LES 計算結(jié)果進(jìn)行對比。網(wǎng)格單元總數(shù)量約800 萬，空腔內(nèi)部與附近的網(wǎng)格數(shù)量約380 萬。對于馬赫數(shù)為0.85 的流動，時間步長為 Δt=2×10?5s， 約為Tref=L/U∞的百分之一，足以捕捉第一階Rossiter 頻率。圖1 給出了來流馬赫數(shù)為0.85 時，中截面上時均流向速度、縱向速度和湍動能的分布，結(jié)果與參考文獻(xiàn)符合得較好。由此可見，本文采用的計算方法可以比較準(zhǔn)確地捕捉空腔復(fù)雜湍流的時均統(tǒng)計特性。

圖1 空腔標(biāo)模M219 在Ma = 0.85 時均速度與湍動能分布Fig. 1 Wall-normal profiles of the mean velocity and turbulence kinetic energy in canonical cavity M219 at Ma = 0.85

2 動態(tài)嵌套重疊網(wǎng)格方法

針對運(yùn)動物體的網(wǎng)格生成問題，發(fā)展了快速的動態(tài)嵌套重疊網(wǎng)格方法[33]。首先通過采用分層嵌套重疊網(wǎng)格策略（如圖2）對網(wǎng)格進(jìn)行劃分，然后采用自適應(yīng)樹結(jié)構(gòu)完成網(wǎng)格切割和貢獻(xiàn)單元搜索。引入壁面距離概念，設(shè)計了統(tǒng)一的分層嵌套重疊策略，同時實(shí)現(xiàn)了重疊網(wǎng)格切割和嵌套網(wǎng)格切割兩個功能，完成重疊最小化的功能，使得重疊網(wǎng)格組裝質(zhì)量得到了改善。

圖2 分層嵌套重疊網(wǎng)格策略Fig. 2 Hierarchically embedded strategy of overset grids

根據(jù)運(yùn)動規(guī)律的不同，可將所有運(yùn)動物體分為若干動態(tài)組，而所有靜態(tài)物體歸為唯一的靜態(tài)組。各組分別構(gòu)造獨(dú)立的二叉樹和八叉樹輔助網(wǎng)格，并且每個組內(nèi)的流場網(wǎng)格和二/八叉樹網(wǎng)格的坐標(biāo)均保持不變。即，各個動態(tài)分組仍然用各自原有的局部坐標(biāo)系，但在每個時間步上記錄下各動態(tài)組運(yùn)動的位移。而在更新切割信息時，將所有被切割的坐標(biāo)點(diǎn)變換到切割物體所在的局部坐標(biāo)內(nèi)，從而完成運(yùn)動問題的網(wǎng)格切割。對于同一組內(nèi)的切割，由于相對位置不變，只需進(jìn)行一次切割，后續(xù)不需再改變。此方法使得二叉樹和八叉樹在運(yùn)動過程中不必重新構(gòu)造，只需進(jìn)行坐標(biāo)變換操作，大大節(jié)省了計算時間。表1 給出了艙門運(yùn)動時網(wǎng)格組裝時間的對比。對于總計2 100 多萬的網(wǎng)格單元來說，兩種動態(tài)嵌套重疊網(wǎng)格組裝方法均可滿足該類非定常流動計算的效率要求。

表1 考慮運(yùn)動艙門的空腔動態(tài)重疊網(wǎng)格組裝的耗時統(tǒng)計Table 1 Assembly time of overset grids for a cavity with moving doors

3 空腔標(biāo)模C201 與艙門模型

3.1 空腔標(biāo)模C201 與艙門模型

空腔標(biāo)模C201 由中國空氣動力研究與發(fā)展中心設(shè)計，空腔底部、前后壁分別布置了靜壓測點(diǎn)和脈動壓力測點(diǎn)（見圖3）。本文采用的計算模型如圖4（a）所示，與試驗(yàn)的模型外形一致?？涨婚L度（L）為200 mm，深度（D）為33.3 mm，寬度（W）為66.7 mm，空腔長深比L/D為6，寬深比W/D為2；與空腔相連壁面的長度和寬度分別是Wa和La。

圖3 空腔標(biāo)模C201 中截面的壓力測點(diǎn)Fig. 3 Distribution of pressure probes in the central section of canonical cavity C201

參照文獻(xiàn)[1]的數(shù)據(jù)選取艙門，艙門厚度1.5 mm，艙門四周采用弧形物面進(jìn)行過渡。此構(gòu)型與空腔標(biāo)模M219 的艙門外形有一定相似性。在艙門與空腔、艙門與艙門之間保留1 mm（約為艙門深度的3%）縫隙，這樣可以保證縫隙處保留足夠的網(wǎng)格用于數(shù)值模擬。艙門三維模型如圖4（b）所示。

圖4 空腔標(biāo)模C201 與艙門計算模型Fig. 4 Geometry of canonical cavity C201 and doors

3.2 計算網(wǎng)格

整個模型網(wǎng)格采用一個網(wǎng)格層，包含三個網(wǎng)格簇：空腔網(wǎng)格簇（其中包含背景網(wǎng)格）和兩個艙門網(wǎng)格簇。圖5 給出了空腔的網(wǎng)格，由于網(wǎng)格在展向（z方向）是對稱的，所以圖中顯示了一半網(wǎng)格。邊界條件設(shè)置如下：在圖5 給出的模型中，紅色為無滑移物面邊界條件；綠色表示的是空腔前后的自由端，采用對稱邊界條件；展向的兩個邊界均采用對稱邊界條件；最上面的外場邊界設(shè)置為入口/出口邊界；前后外場分別為入口/出口邊界?？涨痪W(wǎng)格采用結(jié)構(gòu)網(wǎng)格，共計約1 600 萬網(wǎng)格單元。

圖5 空腔標(biāo)模C201 計算模型與網(wǎng)格（每4 個點(diǎn)顯示一個網(wǎng)格點(diǎn)）Fig. 5 Computational domain and grids for canonical cavity C201

艙門采用沿著物面法向?qū)⑽锩婢W(wǎng)格向外推進(jìn)的方法生成貼體網(wǎng)格，然后在外圍采用矩形外場，使得外場網(wǎng)格的周向尺寸得到控制。單個艙門的網(wǎng)格單元共計約454 萬，圖6（a）給出了艙門的物面網(wǎng)格，艙y+=1門的前向剖視圖如圖6（b）所示。艙門壁面的第一層網(wǎng)格厚度由 給出。艙門開度的定義如圖6（c）所示。

圖6 艙門計算網(wǎng)格與艙門打開角度（開度）定義Fig. 6 Computational grids for cavity doors and the definition of the door-opening angle

剪切層流動是空腔復(fù)雜流場的重要因素之一，所以需要保留高質(zhì)量的空腔剪切層網(wǎng)格。通過改進(jìn)空腔壁面距離的計算方法，使得在四種不同艙門打開角度下均能保留空腔剪切層網(wǎng)格（見圖7），從而保證剪切層流動計算的精細(xì)和準(zhǔn)確。

圖7 四種不同艙門開度在 x/L=0.5處yOz 視圖上網(wǎng)格的組裝結(jié)果Fig. 7 The grid slice of assembled overset grids in yOz planes for cavities with different door-opening angles

4 空腔標(biāo)模C201 數(shù)值模擬

對無艙門的空腔標(biāo)模C201 亞聲速流動（Ma=0.6）進(jìn)行模擬，來流條件為T∞=286.66 K，Re=12.14×106。 時 間步 長 采用 Δt=1.0×10?6s，約為 參 考 時 間Tref=L/U∞的千分之一，足以分辨流動的細(xì)節(jié)。采用雙時間步推進(jìn)求解非定常流場，每步偽時間推進(jìn)采用多重網(wǎng)格加速收斂，并采用MPI 并行技術(shù)提高計算效率。非定常計算采用IDDES 模擬方法，在廣州天河二號超級計算機(jī)上進(jìn)行。采用SA 湍流模型計算的穩(wěn)態(tài)流場作為初場，在非定常流場充分發(fā)展之后（約15Tref～20Tref），開始非定常流場的數(shù)據(jù)采集和時均場的累計。時均場的累計大約為 50Tref，約為0.05 s，足以分辨第一階Rossiter 頻率。

4.1 流場時均特性分析

圖13（a）給出了中截面時均流場的壓強(qiáng)系數(shù)分布和流線形狀，在空腔后部存在一個較小的第二環(huán)流。在較低速來流條件下，剪切層厚度增長較大，大范圍沖擊到后壁、角落和底部。沖擊后壁與角落的氣流較為靠近剪切層核心區(qū)，流速較高，形成局部高壓環(huán)流；而沖擊底部的氣流較為遠(yuǎn)離剪切層核心區(qū)，流速較低，被角落的高壓環(huán)流所阻礙，提前回流形成了主環(huán)流。

圖8 對比了三種不同網(wǎng)格尺度下腔內(nèi)底部壓力系數(shù)的分布與試驗(yàn)值的比較。隨著網(wǎng)格的加密，計算結(jié)果趨于穩(wěn)定。盡管其值高于試驗(yàn)結(jié)果，但是計算結(jié)果可以較好地預(yù)測壓力系數(shù)的變化。從圖中可以看出，底部壓強(qiáng)分布表現(xiàn)為前段平緩下降、中段增長和后段快速上升，最后氣流沖擊后壁形成高壓區(qū)域。

圖8 三種網(wǎng)格對腔內(nèi)底部壓力系數(shù)時均分布的影響Fig. 8 Grid convergence of the mean pressure coefficients at cavity floor compared to experimental distribution

4.2 總聲壓級

總聲壓級（overall sound pressure level，OASPL）可以通過湍流模擬得到的壓力脈動來獲得，其定義如下：

圖9 空腔標(biāo)模C201 腔內(nèi)底部總聲壓級Fig. 9 OASPL distribution at the floor of canonical cavity C201 without doors

圖10 空腔標(biāo)模C201 底部聲壓級成分合成示意圖Fig. 10 Sketch of the development of OASPL in canonical cavity C201

4.3 脈動壓力功率譜

通過脈動壓力功率譜分析空腔底部的頻譜特性，對脈動壓力p′采用快速傅里葉分析得到功率譜密度（PSD）， 然后將PSD 用SPL[13]呈現(xiàn)，進(jìn)行功率譜的分析。圖11 給出了空腔標(biāo)模C201 腔內(nèi)底部第3 和第11 監(jiān)測點(diǎn)脈動壓強(qiáng)功率譜，圖中灰色條帶是根據(jù)Heller 公式計算的峰值頻率帶。計算結(jié)果與試驗(yàn)結(jié)果吻合，主導(dǎo)峰值頻率都是第二峰值，因此在工程上需要特別注意第二峰值頻率的影響。

圖11 空腔標(biāo)模C201 腔內(nèi)底部p3 和p11 測點(diǎn)脈動壓強(qiáng)功率譜Fig. 11 Comparison of SPL of fluctuating pressure at probes p3 and p11 between numerical simulation and experiment

5 靜態(tài)艙門空腔模擬

在滿足準(zhǔn)定常假設(shè)時，可以對靜態(tài)艙門直接模擬分析。采用發(fā)展的嵌套重疊網(wǎng)格方法和IDDES 湍流模擬方法，對艙門開度分別為30°、60°、90°和120°的四種工況進(jìn)行準(zhǔn)定常流動模擬，分析其流動時均特性和流場動態(tài)特性。從時均流場（中截面時均流場和底部壓力分布）、流場的動態(tài)特性（總聲壓級、脈動壓強(qiáng)功率譜）和瞬時場湍流結(jié)構(gòu)三個方面分析艙門對空腔流場產(chǎn)生的影響。本節(jié)采用的數(shù)值方法、計算方式與干凈空腔模擬方法一致。

5.1 時均場特性分析

將空腔底部的靜壓測點(diǎn)時均壓力與干凈空腔的結(jié)果進(jìn)行對比分析，圖12 給出了不同狀態(tài)下底部壓力系數(shù)的分布。從腔內(nèi)底部壓力分布來看，艙門打開30°時，底部壓力系數(shù)的分布與干凈空腔的差別較大，前部壓強(qiáng)有所升高，后部壓強(qiáng)明顯降低。在艙門打開60°、90°、120°時，腔內(nèi)前半部分底部壓力系數(shù)基本與干凈空腔的計算結(jié)果一致，后半部分壓力系數(shù)隨著艙門打開角度的增加逐漸趨于干凈空腔。開啟角度為90°和120°時，時均壓強(qiáng)系數(shù)分布與干凈空腔基本一致。

圖12 四種艙門開度（30°、60°、90°、120°）對空腔底部壓力系數(shù)的影響Fig. 12 Impacts of the door-opening angle on mean pressure coefficients at the cavity floor

圖13 給出了干凈空腔和不同艙門開度的中截面時均流場壓力系數(shù)云圖和流線分布情況。從圖中可以看出在艙門打開角度為30°時，空腔前緣附近的平板上，壓力系數(shù)會升高，這是由于艙門的阻礙使得流動減速造成的。這種壓強(qiáng)升高作用也使得空腔前部的負(fù)壓區(qū)域的壓強(qiáng)系數(shù)相對干凈空腔要偏高。同時，對于30°開啟角度，外流與空腔附近的流動受到艙門的阻擋，只能通過中截面附近的開口相互作用。特別是阻礙外流對剪切層的能量供給，使得后壁附近的壓強(qiáng)相對干凈空腔明顯降低，導(dǎo)致30°時，空腔底面后部壓強(qiáng)系數(shù)水平整體降低。此外，后壁附近的高壓也受到艙門阻礙，在縱向方向只能在中截面開口處泄壓，后部流線有更加明顯地向上彎曲。

圖13 干凈空腔和四種艙門開度下中截面平均流場壓力系數(shù)與流線分布Fig. 13 Impacts of door-opening angles on the streamlines and mean pressure contours in the middle section of cavity C201

艙門開度為60°、90°和120°時，空腔前緣平板同一位置的壓力系數(shù)逐漸降低，且后壁高壓區(qū)和后壁拐角處低壓區(qū)的范圍較艙門開度為30°時逐漸變大，并與干凈空腔的結(jié)果趨于一致。這與“艙門打開角度越大，對流動的限制明顯越小”的物理原因是一致的。

另外，相比于干凈空腔，隨著艙門打開角的不斷增大，主渦位置沿流向不斷地后移并靠近干凈空腔主渦位置。艙門打開120°時，湍流結(jié)構(gòu)基本上與干凈空腔的一致，但是在前壁附近的流場結(jié)構(gòu)與干凈空腔有所區(qū)別。由此可知，隨著空腔打開角度的增加，腔內(nèi)時均流場的形態(tài)、渦的結(jié)構(gòu)與干凈空腔的結(jié)果逐漸趨于一致。

5.2 總聲壓級

圖14 給出了四種艙門開度下空腔底部總聲壓級的分布。從圖中可以看出，艙門開度30°時，底部各個位置的總聲壓級較干凈空腔均明顯下降，但變化趨勢與干凈空腔的結(jié)果基本一致。艙門打開90°時，總聲壓級與120°的變化趨勢一致，前者的值略高1～2 dB。艙門開度為60°時，總聲壓級腔內(nèi)前半部分的分布特點(diǎn)與90°、120°的結(jié)果基本一致；而在腔內(nèi)后半部分，其總聲壓級明顯有所改變。

本文認(rèn)為導(dǎo)致噪聲升高的原因在于艙門阻礙了腔內(nèi)噪聲向外傳播。從圖14 可以看出，艙門分別開啟60°、90°、120°時，三種條件下的聲壓級分布幾乎一致，因此可知噪聲升高與艙門開啟角度無關(guān)，只與艙門的存在有關(guān)。艙門的存在使得空腔內(nèi)部的脈動能量無法以聲波的形式沿展向傳播，只能沿著兩艙門開口方向傳播，第二階模態(tài)的能量無法得到有效釋放，這樣就導(dǎo)致腔內(nèi)噪聲強(qiáng)度升高。本文中艙門開度為90°時，總聲壓級的計算結(jié)果與空腔標(biāo)模M219 在馬赫數(shù)0.85 下的試驗(yàn)結(jié)果[1]分布特點(diǎn)相似，呈現(xiàn)出“W”形態(tài)的分布。

圖14 干凈彈艙與四種艙門開度下的彈艙底部總聲壓級分布Fig. 14 Impacts of door-opening angles on OASPL at the cavity floor

5.3 脈動壓強(qiáng)功率譜

空腔中的噪聲主要由寬頻白噪聲和窄頻噪聲（Rossiter 模態(tài)）兩部分組成[13]。白噪聲是由流動結(jié)構(gòu)產(chǎn)生，其來源為自由來流的脈動、剪切層和湍流等。第二部分也稱為Rossiter 模態(tài)，是不同流場結(jié)構(gòu)直接相互作用產(chǎn)生的，其特點(diǎn)為幅值大且集中在某幾個頻率。圖15 給出了四種不同艙門開啟角度下，腔內(nèi)9 個脈動測點(diǎn)的脈動壓強(qiáng)功率譜曲線，柱狀條帶是Heller 公式計算的峰值頻率。從圖中可以看出：

圖15 不同艙門開度下腔內(nèi)9 個脈動測點(diǎn)脈動壓強(qiáng)功率譜Fig. 15 SPL of fluctuating pressure at nine probes on cavity floors for cases with different door-opening angles

1）艙門開度為30°時，二階Rossiter 所占優(yōu)勢顯著減弱，并且該模態(tài)頻率偏高，其他高頻（第三、第四）模態(tài)已經(jīng)不再明顯。這是由于艙門打開較小時，空腔流動中剪切層與回傳壓縮波的相互作用受到限制，所以流場結(jié)構(gòu)相互作用的強(qiáng)度減弱，致使聲壓幅值降低。

2）60°、90°和120°艙門開度下，各脈動壓力測點(diǎn)的Rossiter 模態(tài)明顯，且模態(tài)頻率與經(jīng)驗(yàn)公式基本一致。這說明此時流場中的結(jié)構(gòu)與干凈空腔一致。二階Rossiter 模態(tài)占主導(dǎo)，幅值明顯比其他模態(tài)高10 dB以上。

3）圖中虛線、實(shí)線分別表示腔內(nèi)前、后半部分采樣點(diǎn)的結(jié)果。從圖中可以看出，四種不同艙門開啟角度下，后半部分的脈動強(qiáng)度明顯高于前半部分。

為研究艙門對腔內(nèi)脈動特性的影響，進(jìn)一步選取后壁（p11）的脈動壓強(qiáng)測點(diǎn)上的脈動壓強(qiáng)功率譜進(jìn)行比較分析，如圖16 所示。各圖中黑色曲線為干凈空腔試驗(yàn)曲線，褐色曲線為干凈空腔計算曲線，其他彩色曲線為各角度情況下的計算曲線。

圖16 不同艙門開度腔內(nèi)后壁測點(diǎn)p11 脈動壓強(qiáng)功率譜Fig. 16 SPL of fluctuating pressure at probe p11 on cavity floors for cases with different door-opening angles

對于后壁測點(diǎn)p11 來說，當(dāng)開啟角為30°時，功率水平明顯低于干凈空腔，其主峰頻率比干凈空腔的主峰頻率偏高。這可能是由于小角度開啟情況下，空腔構(gòu)型更接近封閉腔體，其展向或者縱向的共振成分加強(qiáng)。而在其他開啟角度情況下，變化基本一致，峰值頻率大致和干凈空腔一致，但是第二模態(tài)峰值都顯著升高，三、四模態(tài)峰值沒有升高甚至幅值降低。這也印證了第二階模態(tài)被強(qiáng)化的觀點(diǎn)。同時這也說明，當(dāng)打開角度較大時，在各處，艙門對主模態(tài)的影響是一致的，是通過改變流場的主結(jié)構(gòu)來影響噪聲的產(chǎn)生。在前壁和腔內(nèi)底部測點(diǎn)的脈動壓強(qiáng)功率譜中，出現(xiàn)了相似的變化和影響，不在此贅述。

整體來看，出現(xiàn)不同的變化的原因是，在艙門打開小角度時，空腔內(nèi)的法向（y向）和展向（z向）流動均受到限制，流場脈動會受到抑制，噪聲也會減小。艙門打開角度較大（60°以后）時，艙門的限制主要在展向，對法向的限制已經(jīng)不明顯；但是此時艙門的存在會增加空腔流動的復(fù)雜性，阻斷空腔近壁面噪聲沿展向向外傳播，從而會增強(qiáng)腔內(nèi)的噪聲。

5.4 湍流結(jié)構(gòu)

為了進(jìn)一步分析不同艙門開度下空腔流動的特點(diǎn)，通過Q等值面比較瞬時場中湍流結(jié)構(gòu)的異同。正的Q值區(qū)域表示渦張量比應(yīng)變率張量大的流動區(qū)域，是與旋渦有關(guān)的部分，是常被用作判斷旋渦區(qū)域的標(biāo)準(zhǔn)[34]。圖17 給出了干凈空腔和四種不同艙門開度下，瞬時流場中的Q= 3 000 等值面，采用來流方向的無量綱速度u進(jìn)行染色。

圖17 瞬時流場的Q = 3 000 等值面（用u 進(jìn) 行染色）Fig. 17 Instantaneous iso-surfaces (extracted by Q = 3 000) colored by u for clean cavity and the cavity with stationary doors

從瞬時圖中可以看出，由于艙門的存在，空腔內(nèi)的流動會受到不同程度的影響，具體表現(xiàn)為：

1）30°的湍流結(jié)構(gòu)明顯受到艙門的限制，從而使得流場的湍流結(jié)構(gòu)較干凈空腔要小，艙內(nèi)渦結(jié)構(gòu)很稀疏，尤其是空腔前部區(qū)域。

2）在艙門打開60°時，渦結(jié)構(gòu)變得密集，可以明顯看到艙內(nèi)的流動結(jié)構(gòu)沿法向（y向）向外運(yùn)動，但是空腔前部仍然比較稀疏。

3）到90°和120°時，渦結(jié)構(gòu)的豐富程度變得相似。在艙門開度為120°時，艙門仍然會限制空腔展向的流動。這從側(cè)面反映了艙門的存在會阻礙空腔周圍噪聲的傳播，從而使得腔內(nèi)噪聲一定程度地增大。

4）相比干凈空腔的瞬時Q準(zhǔn)則等值面圖，可以看出，有艙門時空腔上方的剪切層渦結(jié)構(gòu)運(yùn)動高度更高，這仍然是剪切層在展向受艙門限制的結(jié)果。而隨著艙門開度的不斷增大，展向的空腔流動會受到不同程度的限制。

5）流動經(jīng)過前緣平板后形成的剪切層，在流場結(jié)構(gòu)中會出現(xiàn)大尺度的剪切層渦結(jié)構(gòu)。艙門對腔內(nèi)前半部分湍流結(jié)構(gòu)的相互作用較弱；對后腔內(nèi)半部分不同流動結(jié)構(gòu)之間相互作用的影響較大，從而影響空腔流動的脈動和噪聲特性。

6 運(yùn)動艙門空腔模擬

為進(jìn)一步模擬真實(shí)艙門打開過程對空腔流動和聲場的影響，本節(jié)采用前文的空腔和艙門，在考慮艙門運(yùn)動的情況下，對空腔的湍流流動進(jìn)行模擬。對考慮運(yùn)動的非穩(wěn)態(tài)非定常流動問題，采用經(jīng)驗(yàn)?zāi)B(tài)分解（empirical mode decomposition，EMD）方法提取非穩(wěn)定特征，定量分析其噪聲特性的變化。

6.1 艙門運(yùn)動設(shè)計

艙門運(yùn)動規(guī)律的選取由條件相似準(zhǔn)則確定，而斯特勞哈爾數(shù)（St）是非定常流動中須考慮的相似準(zhǔn)則，物理上代表非定常運(yùn)動與定常運(yùn)動的慣性力之比，其定義為：

式中，LD為艙門長度， ΔT為艙門每次開啟的總時間，V∞表示來流速度。戰(zhàn)斗機(jī)武器艙艙門開啟過程中，需1 ～2 s 打開100°[35]。為使得艙門運(yùn)動規(guī)律對流動的影響與真實(shí)飛行情況一致，需要保證St數(shù)相同。因此，在來流條件相同的情況下，當(dāng)模型縮比為1∶10 時，在縮比模型模擬過程中，艙門需要0.1～0.2s左右完成一次動作。本文選取0.1 s，艙門運(yùn)動的角速度為900 rad/s；從艙門完全關(guān)閉到打開120°，模擬的時間步長取 Δt=2×10?5。

6.2 非穩(wěn)態(tài)非定常分析方法

首先對艙門關(guān)閉狀態(tài)進(jìn)行非定常計算，保證流動充分發(fā)展；然后艙門開始運(yùn)動，進(jìn)行完整的非穩(wěn)態(tài)非定常流動模擬。圖18 給出了艙門從0°打開到100°過程中，空腔內(nèi)動態(tài)測點(diǎn)的壓力系數(shù)變化。從圖中可以看出，整個流動大致分為三個階段：初期的小幅穩(wěn)定階段（≤30°）、過渡階段（30°～60°）、充分發(fā)展階段（≥60°）。同時空腔前半?yún)^(qū)域各動態(tài)采集點(diǎn)（p3～p7）的壓力系數(shù)呈對稱擺動趨勢，而空腔后半?yún)^(qū)域動態(tài)采集點(diǎn)的壓力系數(shù)呈增長趨勢。以p11 測點(diǎn)為例，圖19展示了其壓力值的變化情況。

圖18 艙門打開過程中各動態(tài)采集點(diǎn)的壓力系數(shù)變化Fig. 18 Variations of pressure coefficients at different probes during the opening of cavity doors

圖19 艙門開啟過程中彈艙后壁測點(diǎn)p11 處壓力變化與EMD 模態(tài)曲線Fig. 19 Variation of pressure and its EMD modes at probe p11 during the opening of cavity doors

由于艙門運(yùn)動使得空腔內(nèi)流動呈現(xiàn)非穩(wěn)態(tài)的變化，即壓力脈動（p′）需要使用局部平均值而不是簡單的代數(shù)平均。為了準(zhǔn)確描述艙門運(yùn)動過程的非穩(wěn)態(tài)變化，本文提出采用希爾伯特-黃轉(zhuǎn)換（Hilbert-Huang Transform） 中 的 經(jīng) 驗(yàn) 模 態(tài) 分 解（ Empirical Mode Decomposition，EMD）進(jìn)行分析。該方法的核心在于尋找當(dāng)前數(shù)據(jù)中的所有局部極大值以及局部極小值，然后分別采用三次樣條構(gòu)建極大值包絡(luò)線和極小值包絡(luò)線，兩者的平均值即為EMD 模態(tài)；再對扣除EMD 模態(tài)的剩余數(shù)據(jù)進(jìn)行分析，直到剩余數(shù)據(jù)為單調(diào)函數(shù)為止。

為證明該方法的有效性，采用艙門打開90°的準(zhǔn)定常流動進(jìn)行分析，并用代數(shù)平均的計算結(jié)果進(jìn)行對比驗(yàn)證。以測壓點(diǎn)p11 為例，圖20 給出了不同EMD 模態(tài)的分布曲線。分析結(jié)果表明，低階EMD 模態(tài)分析包含較多高頻信息，與原始數(shù)據(jù)分布具有一致性，無法直接用于捕捉非穩(wěn)態(tài)的信息；高階的EMD 模態(tài)多數(shù)與低頻流場信息有關(guān)，可用于描述艙門運(yùn)動引起的非穩(wěn)態(tài)流場的變化。從圖20 中，EMD5 和EMD6模態(tài)的分布差異可以觀察到此特征。而更高階的EMD 模態(tài)（如EMD9）基本與代數(shù)平均值一致。從而說明EMD 方法具有與代數(shù)平均計算方式相同的功能。

圖20 艙門打開90°時p11 點(diǎn)壓力變化與不同EMD 模態(tài)曲線Fig. 20 Time histories of pressure and its EMD modes at probe p11 in the cavity with vertical doors

圖21 中給出了采用不同EMD 模態(tài)作為時均值時，p11 點(diǎn)總聲壓級的變化。結(jié)果表明，采用低階EMD 模態(tài)所得的OASPL 遠(yuǎn)低于采用代數(shù)平均值方式計算的值。而提高采用的EMD 模態(tài)，OASPL 會逐漸收斂于代數(shù)平均計算的值（相差在0.5%以內(nèi)）。這說明采用EMD 模態(tài)進(jìn)行OASPL 計算是有效可行的。盡管EMD 模態(tài)的選取無統(tǒng)一標(biāo)準(zhǔn)，但可以采用OASPL 的收斂性判斷選取的EMD 模態(tài)。

圖21 艙門打開90°時考慮不同EMD 模態(tài)的p11 點(diǎn)OASPL 計算結(jié)果Fig. 21 OASPL computed using different EMD modes at probe p11 in the cavity with vertical doors

圖22 給出了EMD 模態(tài)對p11 測點(diǎn)功率譜的影響，主要特點(diǎn)為低頻區(qū)域有差別，高頻區(qū)域無明顯差異。由于EMD 本身就是對高頻濾波的結(jié)果，不同的EMD 曲線代表不同的低頻模態(tài)，所以會對低階頻率產(chǎn)生影響。計算結(jié)果中的典型Rossiter 頻率與采用平均值計算的方式一致，可見采用EMD 模態(tài)對壓力功率譜的計算方式是可行的。

6.3 運(yùn)動艙門的噪聲分析

采用經(jīng)驗(yàn)?zāi)B(tài)分解方法，對亞聲速空腔艙門開啟過程中的非穩(wěn)態(tài)動態(tài)特性進(jìn)行分析。圖23 給出了空腔后壁p11 測點(diǎn)的OASPL 隨著選取不同EMD 模態(tài)所產(chǎn)生的變化。當(dāng)取到EMD5 模態(tài)時，計算得到的OASPL 開始趨于穩(wěn)定，收斂值為艙門運(yùn)動過程中的OASPL。值得注意的是，該值低于靜態(tài)艙門的計算值。這是因?yàn)椋谂撻T的開啟過程中，當(dāng)艙門打開角度小時，流動的脈動較小，艙門打開后期的脈動變大（見圖19），空腔內(nèi)流動與聲場逐漸產(chǎn)生相互作用并增強(qiáng)，進(jìn)而產(chǎn)生了流動的振蕩現(xiàn)象。圖24 給出了采用艙門開啟到90°附近的數(shù)據(jù)進(jìn)行分段動態(tài)分析的OASPL 計算結(jié)果。此時的結(jié)果與靜態(tài)艙門計算值相當(dāng)，說明對靜態(tài)艙門打開90°的準(zhǔn)定常計算跟真實(shí)運(yùn)動過程中的狀態(tài)吻合，這對工程問題的分析具有指導(dǎo)意義。

圖23 艙門開啟過程中不同EMD 模態(tài)的p11 點(diǎn)OASPL 計算結(jié)果Fig. 23 OASPL computed using different EMD modes at probe p11 during the opening of cavity doors

圖24 艙門開啟到90°時不同EMD 模態(tài)的p11 點(diǎn)的OASPL 分析結(jié)果Fig. 24 OASPL computed using different EMD modes at probe p11 when cavity doors are opened to 90°

圖25 給出了采用不同EMD 模態(tài)時，對應(yīng)的脈動壓力功率譜分析結(jié)果。艙門運(yùn)動過程中，二階Rossiter 模態(tài)仍為主導(dǎo)模態(tài)，幅值明顯高于其他模態(tài)10 dB 以上，但低于靜態(tài)艙門打開90°的分析結(jié)果，此處與總聲壓級降低的原因一致。一階Rossiter 模態(tài)強(qiáng)度降低，而三、四階Rossiter 模態(tài)略有升高。

圖25 艙門開啟過程中不同EMD 模態(tài)的p11 點(diǎn)的功率譜比較Fig. 25 SPL computed using different EMD modes at probe p11 during the opening of cavity doors

采用EMD 方法對空腔內(nèi)全部脈動測點(diǎn)進(jìn)行OASPL 分析，如圖26 所示。計算結(jié)果表明，艙門開啟過程中，OASPL 分布狀態(tài)仍呈現(xiàn)“W”形狀，但其值比靜態(tài)大開度艙門整體降低5 dB 左右。這是由于艙門開啟過程中，存在由小幅脈動到強(qiáng)脈動的發(fā)展過程。小角度的艙門開啟過程中，流場的脈動較小，還沒有完全發(fā)展；強(qiáng)脈動的完全發(fā)展需要一定的時間。在艙門的整個運(yùn)動過程中流動呈現(xiàn)過渡發(fā)展的特點(diǎn)，因此出現(xiàn)OASPL 的降低。

圖26 艙門開啟過程中空腔底部總聲壓級分布Fig. 26 OASPL at cavity floor during the opening of cavity doors

6.4 流動分析

圖27 給出了艙門開啟過程中，艙門打開到不同角度時，物面壓力系數(shù)的分布情況。從圖中可以看出，壓力分布具有一定對稱性，兩側(cè)艙門承受的壓力相當(dāng)。在艙門打開小角度時，艙門限制了空腔內(nèi)部與外部流場的摻混與交換，壓力變化較小。隨著艙門的逐漸打開，空腔與外流場發(fā)生混合，空腔內(nèi)壓力變化增大。艙門打開角度較大（60°以后）時，空腔內(nèi)流動出現(xiàn)大幅度振蕩現(xiàn)象，流場與聲場出現(xiàn)強(qiáng)相互作用，從而造成流場脈動增強(qiáng)，噪聲特性變得顯著。與前文OASPL 分析和功率譜分析結(jié)果一致。

圖27 艙門開啟過程中物面壓力系數(shù)變化Fig. 27 Instantaneous wall pressure coefficients during the opening of cavity doors

7 結(jié) 論

針對考慮艙門的多尺度、強(qiáng)湍流、高雷諾數(shù)、強(qiáng)脈動的空腔流動問題，采用本文發(fā)展的嵌套重疊網(wǎng)格方法對帶運(yùn)動艙門的空腔生成網(wǎng)格，用IDDES 湍流模擬方法對亞聲速（Ma= 0.6）空腔流動進(jìn)行模擬，分析了靜、動態(tài)艙門對空腔復(fù)雜流動和噪聲特性的影響。具體結(jié)論如下：

1）發(fā)展的動態(tài)嵌套重疊網(wǎng)格方法，可自動高效生成大幅度運(yùn)動物體的計算網(wǎng)格，用于運(yùn)動艙門的空腔復(fù)雜流場模擬。

2）采用四階對稱重構(gòu)和三階MUSCL 混合的格式，在不明顯增加計算量的前提下，明顯降低了數(shù)值耗散，提升了計算的準(zhǔn)確度，可以捕捉空腔中的湍流結(jié)構(gòu)。對標(biāo)模M219 的計算結(jié)果與文獻(xiàn)給的結(jié)果吻合。對空腔標(biāo)模C201 腔內(nèi)噪聲的預(yù)測和對脈動壓力功率譜的分析結(jié)果與試驗(yàn)數(shù)據(jù)一致。

3）亞聲速空腔標(biāo)模C201 腔內(nèi)底部壓強(qiáng)分布表現(xiàn)出三段分布：前段平緩負(fù)壓區(qū)；中段壓力增長過渡區(qū)；后段高度增加區(qū)，形成后壁的高壓區(qū)域。腔內(nèi)主導(dǎo)峰值頻率為第二峰值。

4）與干凈空腔相比，當(dāng)艙門為靜態(tài)時，艙門小開度（30°）時，時均流場中截面上腔內(nèi)的流態(tài)變得更加復(fù)雜，并且腔內(nèi)底部壓力系數(shù)的值均會明顯減小。艙門開度較大（60°以后）時，時均流場中截面的腔內(nèi)渦結(jié)構(gòu)與干凈空腔一致，底部壓力分布趨于干凈空腔的計算結(jié)果；但是艙門的存在，會使得底部后半部分的壓力較干凈空腔出現(xiàn)一定程度上的降低。艙門開度為90°時，腔內(nèi)總聲壓級較干凈空腔有所升高，最大升高值近9 dB，呈現(xiàn)出“W”形態(tài)的分布；二階Rossiter模態(tài)占主導(dǎo)。艙門與空腔后半部分的不同流動結(jié)構(gòu)之間相互作用，從而對空腔流動的脈動和噪聲特性產(chǎn)生影響。

5）針對艙門運(yùn)動的非穩(wěn)態(tài)非定常特點(diǎn)，本文采用經(jīng)驗(yàn)?zāi)B(tài)分解方法對噪聲特性進(jìn)行分析，計算結(jié)果表明該方法有效可行。艙門開啟過程中，整體的OASPL有所下降，主導(dǎo)頻率依然為二階Rossiter 模態(tài)。對于小開度的艙門，靜態(tài)艙門模擬與運(yùn)動艙門結(jié)果相差較大，因而此類情況需要考慮運(yùn)動艙門的影響。在艙門開度較大（90°）時，兩者的分析結(jié)論一致，說明可采用準(zhǔn)定常來分析此狀態(tài)的流場與噪聲特性，指導(dǎo)內(nèi)埋彈艙的設(shè)計。