郭 裕，姚世衛(wèi)，李邦明，茍金瀾

(武漢第二船舶設(shè)計(jì)研究所 熱能動力技術(shù)重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430205)

0 引 言

空腔是廣泛存在于船舶系統(tǒng)和管道運(yùn)輸?shù)膬σ合到y(tǒng)中的結(jié)構(gòu)形式，尤其在船舶系統(tǒng)中，一些有特殊用途的船體開口，都可以看作有流體流經(jīng)的空腔結(jié)構(gòu)[1]。以艦船通海口為例，空腔的非定常流動會對艦船造成嚴(yán)重危害，所以有必要深入了解空腔流的速度場和壓力脈動。因此，空腔流的研究具有重要的實(shí)際意義。20世紀(jì)50年代起，國外開始對空腔流的壓力振蕩和聲輻射進(jìn)行了大量的實(shí)驗(yàn)研究和數(shù)值計(jì)算研究[2–7]，國內(nèi)空腔流的研究起步較晚，且多集中在超聲速和跨聲速空腔流方面。例如Rossiter等[7]從實(shí)驗(yàn)研究出發(fā)，提出空腔振蕩的流聲共振反饋模型，并給出預(yù)估振蕩的半經(jīng)驗(yàn)公式。目前空腔的研究主要集中于可壓縮流動方面，對極低馬赫數(shù)或不可壓縮空腔流關(guān)注較少，極低馬赫數(shù)范圍或者不可壓縮空腔流條件下，空腔流可分為剪切層模式和尾流模式[8]。剪切層模式的空腔流與自由剪切層流動極為類似，腔內(nèi)渦的基本形態(tài)比較穩(wěn)定，腔內(nèi)流體與腔外主流之間的質(zhì)量交換不明顯。而空腔流的尾流模式，腔內(nèi)的流體不穩(wěn)定，在腔的前緣處會形成一個(gè)旋渦，此旋渦不斷變大直到脫離前緣，在主流的作用下向腔后緣運(yùn)動，最后從尾緣拋射向腔外，自由流周期性進(jìn)出空腔，腔外流體受腔內(nèi)流體的影響極為明顯。根據(jù)聲速附近可壓縮空腔流的機(jī)理，空腔跨度對振蕩頻率有明顯影響，而深度幾乎不影響振蕩頻率，可壓縮流空腔振蕩機(jī)理是否可以推廣至不可壓縮流中有待研究，所以本文針對不可壓縮空腔流，研究尾流模式振蕩頻率的影響因素。

大多數(shù)研究者通過求解非定常可壓縮N-S方程模擬空腔流，對于湍流項(xiàng)的處理，最初研究者采用BL代數(shù)湍流模型[9]，后來也有人采用一方程和二方程湍流模型，但是結(jié)果并不理想。之后大渦模擬（LES）、直接數(shù)值模擬（DNS）也被運(yùn)用到空腔流的數(shù)值模擬中，但這2種方法需要占用巨大的計(jì)算機(jī)內(nèi)存。而DES對壁面邊界層完全采用RANS湍流模型模擬，利用RANS可以有效處理壁面邊界層流動的優(yōu)點(diǎn)，不需要耗費(fèi)很大的計(jì)算機(jī)資源；在遠(yuǎn)離壁面的分離區(qū)域內(nèi)，對小尺度渦采用亞格子模型模擬，對大尺度渦進(jìn)行直接模擬，利用LES可以比較有效地模擬分離流動的優(yōu)點(diǎn)。通過結(jié)合RANS和LES各自的優(yōu)點(diǎn)，DES可以比較快速而有效的模擬空腔流。

綜上所述，以DES為代表的RANS/LES混合方法在空腔流中具有一定的優(yōu)勢，因此，采用基于k-ω-SST模型的IDDES方法模擬空腔的振蕩特性。首先根據(jù)肥皂膜水洞實(shí)驗(yàn)結(jié)果，數(shù)值模擬同一工況的空腔流，將仿真結(jié)果與實(shí)驗(yàn)結(jié)果進(jìn)行對比，驗(yàn)證數(shù)值計(jì)算方法的準(zhǔn)確性。在此基礎(chǔ)上對目標(biāo)模型進(jìn)行研究，得到空腔流一個(gè)周期內(nèi)渦的演化過程，著重分析了空腔的跨度和深度對振蕩頻率的影響。

1 數(shù)值方法及幾何模型

基于有限體積法求解RANS方程組，湍流模型采用基于兩方程k-ω-SST的IDDES[10], 以降低分離區(qū)域湍流粘性，解析流動小尺度結(jié)構(gòu)。以k-ω-SST模型構(gòu)造IDDES方法僅需要在湍動能輸運(yùn)方程引入混合長度尺度lIDDES，如下式所示：

其中，混合長度尺度lIDDES定義為，由RANS的長度尺度和LES的長度尺度混合構(gòu)成，總體上在邊界層外的流動分離區(qū)域，即IDDES方法表現(xiàn)為LES的特點(diǎn)。求解采用SIMPLE算法，其中壓力項(xiàng)用二階格式離散，動量方程和湍流方程用二階迎風(fēng)格式離散。

本文計(jì)算5類不同類型的空腔，空腔的具體參數(shù)如表1所示，其中模型4跨度為，空腔跨度為時(shí)模型就等價(jià)于后臺階模型，模型4為臺階高度20 mm的后臺階。空腔流在腔口附近流場的物理量變化極為劇烈，對腔口和壁面附近網(wǎng)格進(jìn)行加密處理，最小網(wǎng)格為0.01 mm。入口條件為速度入口，速度v=1.5，2，3，4 m/s，速度方向垂直于入口邊界；出口邊界條件為壓力邊界條件，流場參考壓強(qiáng)為0.1 MPa，腔體結(jié)構(gòu)為無滑移壁面條件，時(shí)間步長為0.001 s，總計(jì)算時(shí)間為5 s。在腔口距離導(dǎo)邊20 mm位置設(shè)置監(jiān)測點(diǎn)，監(jiān)測法向速度。

表 1 空腔模型參數(shù)Tab. 1 Parameters of the cavity model

2 數(shù)值方法的驗(yàn)證

為考察IDDES數(shù)值計(jì)算方法的準(zhǔn)確性，對跨度為80 mm，深度為15 mm的空腔進(jìn)行數(shù)值模擬，模型的外輪廓采用翼型曲線進(jìn)行改進(jìn)，來流速度為1.5 m/s，流體為肥皂液，粘性為2×10–5m2/s。

圖1（a）給出了基于肥皂膜水洞裝置的實(shí)驗(yàn)結(jié)果，圖1（b）為數(shù)值模擬的瞬時(shí)渦量圖。圖1（a）中明顯可以觀察到腔口位置的自由剪切層，起始階段自由剪切層擾動不明顯，自由剪切層與腔口位置基本重合，圖1（b）中自由剪切層起始階段與實(shí)驗(yàn)結(jié)果一致，自由剪切層位置與腔口基本重合；圖1（a）在腔口中間靠后位置自由剪切層明顯向腔外突出，腔內(nèi)靠近后壁位置出現(xiàn)明顯的渦結(jié)構(gòu)，圖1（b）數(shù)值仿真結(jié)果靠近腔后壁位置也出現(xiàn)了相同的渦結(jié)構(gòu)，導(dǎo)致腔口自由剪切層向腔外突出。數(shù)值計(jì)算與實(shí)驗(yàn)結(jié)果吻合良好，驗(yàn)證了IDDES模型在計(jì)算空腔流方面的準(zhǔn)確性，為下一步使用IDDES模型開展空腔流數(shù)值模擬研究打下了堅(jiān)實(shí)的基礎(chǔ)。

圖 1 空腔流場結(jié)構(gòu)Fig. 1 Flow field structure of cavity

3 結(jié)果分析

3.1 流場結(jié)構(gòu)

圖2給出了模型1（80:20）一個(gè)周期內(nèi)不同時(shí)刻的渦量云圖，流速v=1.5 m/s。來流到達(dá)腔口位置失穩(wěn)產(chǎn)生渦結(jié)構(gòu)（本文將此渦命名為腔口渦），受前一周期渦擠壓，腔口渦緊貼空腔前壁發(fā)展，如圖2（a）；腔口渦尺度逐漸增長，與空腔深度相當(dāng)時(shí)，受腔底限制，腔口渦轉(zhuǎn)而橫向發(fā)展并變得扁平，如圖2（b）；腔口渦與空腔前壁底角組成一個(gè)封閉區(qū)域，此封閉區(qū)域內(nèi)出現(xiàn)另一小渦（本文將此渦命名為腔角渦），如圖2（c）；腔角渦不斷增長并擠壓腔口渦，如圖2（d）；直至腔口渦脫離空腔導(dǎo)邊，之后腔角渦跟隨腔口渦向空腔后壁運(yùn)動，如圖2（e）；最終腔口渦撞擊空腔后壁，和腔角渦一起被甩至腔外，如圖2（f）。

根據(jù)一個(gè)周期內(nèi)渦發(fā)展規(guī)律，發(fā)現(xiàn)空腔流尾流模式振蕩主要是由腔口渦導(dǎo)致。當(dāng)來流到達(dá)空腔開口時(shí)，在導(dǎo)邊處發(fā)生流動分離，產(chǎn)生自由剪切層，剪切層下側(cè)腔內(nèi)流速較低，上側(cè)流速接近于來流速度。由于存在較大速度梯度，剪切層卷起形成腔口渦，腔口渦增長并脫離空腔前緣，最終撞擊至后緣角被甩出腔外，產(chǎn)生流激振蕩。

圖3給出了模型1（80:20）在5T/6和T時(shí)刻的速度矢量圖，來流流速v=1.5 m/s。圖3（a）腔口渦影響范圍超過空腔深度，腔內(nèi)流體不穩(wěn)定，腔外流體回流至腔內(nèi)；圖3（b）腔口渦從尾緣拋射向腔外，腔外流體受腔內(nèi)流體的影響極其明顯，腔內(nèi)流體與腔外流體有明顯的質(zhì)量交換，這是典型的尾流模式特征。

圖 2 空腔渦量云圖（80:20，v=1.5 m/s）Fig. 2 Vorticity contours of cavity（80:20，v=1.5 m/s）

圖 3 空腔速度矢量圖（80:20，v=1.5 m/s）Fig. 3 Velocity vector of cavity （80:20，v=1.5 m/s）

圖4給出了不同跨度空腔模型在T/2時(shí)刻的渦量云圖，來流流速v=1.5 m/s。圖4（a）為模型1（80:20）的渦量云圖，圖4（b）為模型2（100:20）的渦量云圖，圖4（c）為模型3（120:20）的渦量云圖，圖4（d）為模型4（step flow）的渦量云圖。渦量云圖顯示不同跨度空腔結(jié)構(gòu)在T/2時(shí)刻流動特征相似，都存在腔口渦和腔角渦，腔口渦的尺度大于腔角渦，腔角渦位于腔口渦與空腔前壁底角構(gòu)成的封閉空間內(nèi)。圖5給出了不同深度空腔模型在T/2時(shí)刻的渦量云圖，來流流速v=1.5 m/s。圖5（a）為模型1（80:20）的渦量云圖，圖5（b）為模型5（80:15）的渦量云圖，渦量云圖顯示不同深度空腔結(jié)構(gòu)在T/2時(shí)刻流動特征基本一致，也出現(xiàn)腔口渦和腔角渦，且出現(xiàn)位置和相對尺度基本一致，可見不同跨度和深度空腔尾流模式流場結(jié)構(gòu)及渦發(fā)展規(guī)律基本一致。

圖 4 不同跨度空腔渦量云圖（T/2）Fig. 4 Vorticity contours of different length cavities（T/2）

圖 5 不同深度空腔渦量云圖（T/2）Fig. 5 Vorticity contours of different depth cavities（T/2）

3.2 空腔振蕩頻率影響分析

3.2.1 空腔跨度對振蕩頻率的影響

Rossiter給出空腔振蕩頻率預(yù)測半經(jīng)驗(yàn)公式：

式中：St為斯特勞哈爾數(shù)；f為振蕩頻率；κ為渦遷移速度與自由流速度之比；α為渦通過與產(chǎn)生壓之間的時(shí)間延遲因子；n為模態(tài)階數(shù)。α=0.25，κ=0.66是常數(shù)[11]。

圖6給出了v=1.5 m/s條件下，模型1（80:20）、模型2（100:20）、模型3（120:20）及模型4（step flow）監(jiān)測點(diǎn)法向速度。從前文流場結(jié)構(gòu)可發(fā)現(xiàn)，空腔振蕩由前緣渦脫落導(dǎo)致，則腔口監(jiān)測點(diǎn)位置頻率即空腔前緣脫渦頻率，也是空腔振蕩頻率。圖中可以觀察到4種不同空腔法向速度周期性均較為明顯，僅起始階段周期性較差，為確保計(jì)算準(zhǔn)確性，選擇空腔振蕩穩(wěn)定后的法向速度，即2～5 s的法向速度，進(jìn)行傅里葉變換得到腔口法向速度頻譜。

圖 6 監(jiān)測點(diǎn)法向速度（v=1.5 m/s）Fig. 6 Normal velocity of monitoring point （v=1.5 m/s）

圖7給出了v=1.5 m/s條件下，模型1（80:20）、模型2（100:20）、模型3（120:20）及模型4（step flow）腔口法向速度頻譜。圖中可以明顯觀察到空腔振蕩存在多階頻率，階數(shù)越高，振蕩幅值越低，其中1階頻率幅值最強(qiáng)，高階頻率幅值依次遞減。監(jiān)測點(diǎn)法向速度頻譜各階頻率，即空腔振蕩各階頻率。

圖8給出了模型1（80:20）、模型2（100：20）、模型3（120:20）及模型4（step flow）不同來流條件下1階和2階振蕩頻率。圖8（a）為1階振蕩頻率，圖中可發(fā)現(xiàn)深度相同、跨度不同的空腔，振蕩頻率隨流速增加而增加，但空腔跨度的改變幾乎不影響振蕩頻率，4種不同跨度空腔的振蕩頻率大小和變化趨勢基本一致；圖8（b）為2階振蕩頻率，可發(fā)現(xiàn)深度相同、跨度不同的空腔結(jié)構(gòu)，其2階振蕩頻率大小和變化趨勢也基本一致。說明尾流模式振蕩頻率與空腔跨度幾乎無關(guān)，而Rossiter公式中空腔振蕩頻率與空腔跨度直接相關(guān)。

圖 7 監(jiān)測點(diǎn)速度頻譜（v=1.5 m/s）Fig. 7 Power spectrum of monitoring point （v=1.5 m/s）

圖 8 空腔振蕩頻率Fig. 8 Frequency of the oscillation cavity

由于腔口自由剪切層由一系列渦組成，這些渦從空腔前緣處脫落，并跟隨主流向下游運(yùn)動，撞擊后緣角產(chǎn)生壓力脈沖，壓力脈沖以當(dāng)?shù)芈曀傧蛏嫌蝹鞑ィ竭_(dá)上游反饋至導(dǎo)邊又會激起新渦產(chǎn)生。在亞聲速和超聲速可壓縮流中，來流運(yùn)動速度和反射壓力脈沖速度量級相當(dāng)，所以反射壓力脈沖對空腔前緣的脫渦頻率有明顯影響，即空腔振蕩頻率與跨度相關(guān)；而不可壓縮流中，反射壓力脈沖速度遠(yuǎn)遠(yuǎn)大于流體運(yùn)動速度，空腔跨度改變對反射壓力脈沖到達(dá)前緣的影響可忽略不計(jì)，所以反射壓力脈沖對空腔前緣的脫渦頻率影響不明顯，即空腔振蕩頻率與跨度幾乎不相關(guān)。

3.2.2 空腔深度對振蕩頻率的影響

圖9給出了v=1.5 m/s條件下，模型1（80:20）和模型5（80:15）腔口監(jiān)測點(diǎn)法向速度。圖中可發(fā)現(xiàn)腔口法向速度有很強(qiáng)的周期性，為保證計(jì)算準(zhǔn)確，選取2～5 s的法向速度進(jìn)行傅里葉變換得到法向速度頻譜。圖10給出了v=1.5 m/s條件下，模型1（80:20）和模型5（80:15）腔口法向速度頻譜?？砂l(fā)現(xiàn)不同深度的空腔，尾流模式也存在多階頻率，階數(shù)越高，振蕩幅值越低。

圖 9 監(jiān)測點(diǎn)法向速度（v=1.5 m/s）Fig. 9 Normal velocity of monitoring point （v=1.5 m/s）

圖 10 監(jiān)測點(diǎn)速度頻譜（80:15, v=1.5 m/s）Fig. 10 Power spectrum of monitoring point （80:15, v=1.5 m/s）

圖11給出了模型1（80:20）和模型5（80:15）不同來流條件下1階和2階振蕩頻率，其中模型5（80:15）4 m/s來流條件時(shí)，空腔振蕩的周期性較差，沒有明顯的主頻，此工況振蕩頻率不作對比。圖8（a）為1階振蕩頻率，圖中可發(fā)現(xiàn)相同跨度、不同深度的空腔結(jié)構(gòu)，其1階振蕩頻率出現(xiàn)了明顯差別，同一工況下模型5（80:15）的振蕩頻率高于模型1（80:20）；圖8（b）為2階振蕩頻率，圖中可發(fā)現(xiàn)相同跨度、不同深度的空腔結(jié)構(gòu)，其2階振蕩頻率也出現(xiàn)明顯差別，同一工況模型5（80:15）的振蕩頻率高于模型1（80:20）的振蕩頻率?？梢娍涨簧疃葘ξ擦髂Ｊ秸袷庮l率有明顯影響，而Rossiter公式中空腔振蕩頻率與空腔的深度不相關(guān)。

前文圖2給出尾流模式一個(gè)周期內(nèi)的渦量云圖，空腔前緣腔口渦在增長過程中尺度接近深度時(shí)，受到腔底限制轉(zhuǎn)而橫向發(fā)展并變得扁平，直至腔口渦脫離導(dǎo)邊。此過程中腔底限制腔口渦發(fā)展是導(dǎo)致腔口渦脫落的關(guān)鍵因素，模型1空腔深度為20 mm，模型5空腔深度為15 mm，模型1深度大于模型5的深度，所以同一工況，模型1腔口渦受腔底限制晚于模型5，導(dǎo)致模型1腔口渦脫離導(dǎo)邊晚于模型5，腔口渦脫離導(dǎo)邊后的發(fā)展兩者幾乎一致，所以模型1一個(gè)周期持續(xù)的時(shí)間長于模型5，即模型1的振蕩頻率小于模型5的振蕩頻率。

圖 11 空腔振蕩頻率Fig. 11 Frequency of the oscillation cavity

4 結(jié) 語

本文采用基于兩方程k-ω-SST模型的IDDES方法，對5類空腔流進(jìn)行數(shù)值模擬，主要對空腔流場結(jié)構(gòu)及振蕩頻率進(jìn)行了詳細(xì)分析，得到以下結(jié)論：

1）IDDES方法可有效模擬尾流模式不可壓縮空腔流，捕捉空腔流的渦結(jié)構(gòu)，反映空腔流的流場特征；

2）5類空腔尾流模式流場結(jié)構(gòu)及渦發(fā)展規(guī)律基本一致，都存在腔口渦和腔角渦，腔口渦由空腔自由剪切層卷起在腔口位置產(chǎn)生，腔角渦出現(xiàn)在腔口渦與空腔前壁底角構(gòu)成的封閉區(qū)域。剪切層卷起形成腔口渦，腔口渦增長并脫離空腔前緣，最終撞擊至后緣角被甩出腔外，腔口渦是產(chǎn)生流激振蕩的主要原因。

3）空腔跨度幾乎不影響振蕩頻率，而深度對振蕩頻率有顯著影響，振蕩頻率隨空腔深度的增長而降低，與可壓縮空腔流Rossiter給出的半經(jīng)驗(yàn)公式有明顯區(qū)別。不同跨度的空腔，不可壓縮流中來自下游的反射壓力脈動速度遠(yuǎn)遠(yuǎn)高于來流流速，空腔跨度改變對反射壓力脈沖到達(dá)前緣的影響可忽略不計(jì)，所以反射壓力脈沖對空腔前緣脫渦頻率影響不明顯，即空腔振蕩頻率與跨度幾乎不相關(guān)。不同深度的空腔，腔底對腔口渦限制是腔口渦脫離導(dǎo)邊的關(guān)鍵因素，空腔深度越大，腔口渦增長尺度越大，一個(gè)周期持續(xù)時(shí)間也越長，導(dǎo)致振蕩頻率越低。

[ 1 ]尹立中, 王本利, 鄒經(jīng)湘. 航天器液體晃動與液固耦合動力學(xué)研究概述[J]. 哈爾濱工業(yè)大學(xué)學(xué)報(bào), 1999, 31(2): 118–122.

[ 2 ]SUPONITSKY V, AVITAL E, GASTER M. On threedimensionality and control of incompressible cavity flow[J].Physics of Fluids, 2005, 17: 104–130.

[ 3 ]DE ROECK W, RUBIO G, BAELMANS R, et al. Toward accurate hybrid prediction techniques for cavity flow noise applications[J]. Int J Numer Meth Fluids, 2009, 61:1363–1387.

[ 4 ]HAIGERMOSER C, VESELY L, NOVARA M, et al. A timeresoved particle image velocimetry investigation of a cavity flow with a thick incoming turbulent boundary layer[J].Submitted to Physics of Fluids, 2008.

[ 5 ]GRACE SM, DEWAR WG, WROBLEWSKI DE.Experimental investigation of the flow characteristics within a shallow wall cavity for both laminar and turbulent upstream boundary layers[J]. Experiments in Fluids, 2004, (36):791–804.

[ 6 ]AHUJA KK, MENDOZA J. Effects of cavity dimensions,boundary layer, and temperature on cavity noise with emphasis on benchmark data to validate computational aeroacoustic codes[J]. NASA contractor report: Final report contract NAS1-19061, task 13.

[ 7 ]ROSSITER J E. Wind-tunnel experiments on the flow over rectangular cavities at subsonic and transonic speeds[J]. Royal Aircraft Establishment Technical Report. 1964, 3438: 8–12.

[ 8 ]GHARIB M, ROSHKO. The effect of flow oscillations on cavity drag[J]. Fluid Mech, 1987, 177: 501–530.

[ 9 ]RIZZETTA D P. Numerical simulation of supersonic flow over a three-dimensional Cavity [J]. AIAA Journal 1988, 26(7):799–807.

[10]SHUR M L, SPALART P R, STRELETS M, TRAVIN A. A hybrid RANS-LES approach with delayed-DES and wallmodelled LES capabilities [J]. International Journal of Heat and Fluid Flow, 2008, 29(6): 1638–1649.

[11]程艷偉. 二維非定?？蓧嚎s湍流的大渦數(shù)值模擬[D]. 長沙:國防科學(xué)技術(shù)大學(xué), 2002.