榮吉利，王 超，趙 瑞，李海旭，嚴(yán) 昊

（1.北京理工大學(xué) 宇航學(xué)院力學(xué)系，北京100081；2.中國船舶工業(yè)系統(tǒng)工程研究院，北京100036）

0 引言

艦載直升機(jī)在執(zhí)行海上補(bǔ)給、搜救、反潛等任務(wù)中表現(xiàn)出巨大的優(yōu)勢，因此越來越多的艦載直升機(jī)被裝備于非航空型艦船上。但與陸基直升機(jī)不同，由于艦船運動以及艉部流場結(jié)構(gòu)復(fù)雜等因素，使直升機(jī)起降受到嚴(yán)重的影響[1-3]。當(dāng)氣流流經(jīng)艦船會形成湍流邊界層，隨后會與飛行甲板前方的機(jī)庫壁面發(fā)生碰撞，由于幾何形狀的突然變化，會導(dǎo)致流動發(fā)生膨脹，此時湍流邊界層在機(jī)庫的分離點處分離，形成自由剪切層，在機(jī)庫后方的飛行甲板處形成巨大的回流區(qū)，該回流區(qū)包含了隨時間變化的渦結(jié)構(gòu)，且具有極強(qiáng)的非定常特性[4-6]。

相關(guān)學(xué)者使用計算流體力學(xué)（CFD）對影響艦船艉流場的因素進(jìn)行了數(shù)值模擬研究。黃斌等人[7]通過建立一個基于NS 方程的艦載直升機(jī)著艦域流場分析方法，探究了機(jī)庫開合狀態(tài)對艦船艉部著艦域的影響，研究結(jié)果表明打開機(jī)庫門更有利于直升機(jī)安全著艦；洪偉宏等人[8]以LHA 兩棲攻擊艦作為參考通過CFD 定常計算研究發(fā)現(xiàn)，艦船的上層建筑是誘發(fā)艦船艦面空氣湍流的因素之一，而島式上層建筑靠前的艦船可以改善后部甲板的流場環(huán)境；但國外相關(guān)的計算流體力學(xué)技術(shù)研究表明，時間精確的非定常計算能夠精確捕捉由艦船船體結(jié)構(gòu)產(chǎn)生的非穩(wěn)態(tài)湍流氣流[9-11]。趙瑞等人[12]通過高精度的DES 方法對SFS 與SFS2 等船型展開非定常數(shù)值模擬，得到了詳細(xì)的船體流場渦脫落和時間精確的流場數(shù)據(jù)，并從數(shù)值上探究了有無上層建筑與機(jī)庫門開合狀態(tài)對艉流場的影響；劉長猛等人[13]使用FLUENT 的非穩(wěn)態(tài)標(biāo)準(zhǔn)k-ε 模型進(jìn)行計算，發(fā)現(xiàn)機(jī)庫的形狀與尺寸會對飛行甲板周圍渦旋的強(qiáng)度和位置產(chǎn)生影響。

由此可見，艦船艉部空氣流場的變化取決于諸多因素，其中類似后臺階類型的艦船艉部會因為機(jī)庫外形尺寸的不同對其流場產(chǎn)生一定的影響。本文使用基于一方程SA（Spalart-Allmaras）模型[14]的脫體渦模擬（Detached-Eddy Simulation，DES）[15]對典型護(hù)衛(wèi)艦船型進(jìn)行數(shù)值模擬，研究了不同機(jī)庫外形對艉流場的影響，并且根據(jù)非定常計算結(jié)果推導(dǎo)了預(yù)測艉流場回流區(qū)范圍的經(jīng)驗公式。

1 計算模型與方法

1.1 計算方法

1.1.1 Spalart-Allmaras 湍流模型

基于一方程SA 模型的DES 方法，其湍流模型求解有關(guān)渦粘性變量的偏微分方程為：

1.1.2 DES 模型

DES 方法是將大渦模擬（LES）方法和雷諾平均（RANS）方法相結(jié)合的湍流模擬方法，通過對Spalart-Allmaras RANS 湍流模型進(jìn)行修正，使在離壁面近處的區(qū)域采用RANS 方法模擬邊界層內(nèi)小尺度渦的運動，而在遠(yuǎn)離壁面的區(qū)域處采用LES 方法模擬，對小尺度渦采用亞格子模型模擬，對大尺度渦進(jìn)行直接模擬。通過結(jié)合LES 方法和RANS 方法的優(yōu)點，DES 方法可以有效處理邊界層內(nèi)小尺度渦的脈動，同時也降低了計算資源和計算時間[16]。

在DES 方法中，SA 模型中的長度尺度，即流場中任意點到最近物面的距離d，用d~來代替，其中d~表達(dá)為

在當(dāng)前計算中，Δ=max（Δx，Δy，Δz）為網(wǎng)格中心到相鄰單元中心距離中最大的一個，CDES=0.65。其他細(xì)節(jié)可以參考文獻(xiàn)[15]。

1.1.3 計算設(shè)置

本文通過CFD++軟件，采用有限體積法對控制方程進(jìn)行離散，并基于具有空間二階精度的多尺度TVD（Total Variation Diminishing）限制器。時間推進(jìn)采用具有二階精度的雙時間步長（Runge-Kutta）全隱式算法。定常計算采用一方程SA 模型并將其計算的結(jié)果作為非定常湍流計算的初始流場。非定常計算3 000 步，其中時間步長Δt*=1.88×10-2s。

1.2 幾何模型

為了研究不同機(jī)庫外形對艦船艉部流場的影響，并且減少計算資源的消耗，本文使用三種簡化的護(hù)衛(wèi)艦CFD 模型，如圖1所示，分別命名為A 船、B 船與C 船。三種簡化船型在外形和尺寸上均保持一致，僅對艉部機(jī)庫的外形進(jìn)行改動。其中，A 船機(jī)庫橫截面為“口”字型，機(jī)庫高度6.5 m；B 船機(jī)庫橫截面為“凸”字型，最大機(jī)庫高度6.5 m，最小機(jī)庫高度4.5 m；C 船機(jī)庫橫截面為“口”字型，機(jī)庫高度較低，為4.5 m。

圖1 三種簡化船型A 船（左）、B 船（中）與C 船（右）及其尺寸細(xì)節(jié)Fig.1 Dimensions of three simplified frigate ships

1.3 網(wǎng)格劃分

整套網(wǎng)格為圓柱形區(qū)域，如圖2所示，船體模型被放置于圓柱形區(qū)域的中央，半徑r=4.5ls，高d=0.75ls，其中l(wèi)s為船長。本文對三種簡化模型采用混合網(wǎng)格進(jìn)行劃分，該混合網(wǎng)格由兩部分組成：一是整個流場的非結(jié)構(gòu)網(wǎng)格，船體表面利用三角形單元形成的40 層棱柱結(jié)構(gòu)（垂直于壁面的第一層網(wǎng)格高度為5 mm，保證y+=o（10），增長率為1.3）來保證粘性層流動的捕捉；二是艉部飛行甲板上方區(qū)域的結(jié)構(gòu)網(wǎng)格，根據(jù)Spalart 的“重點區(qū)域”理論[17]，為了能充分地發(fā)揮DES 方法的優(yōu)點，網(wǎng)格單元應(yīng)該盡可能地小并且各向同性，故在艉流場重點關(guān)注區(qū)域處使用結(jié)構(gòu)網(wǎng)格能更好地控制網(wǎng)格的精細(xì)程度。網(wǎng)格細(xì)節(jié)圖如圖3所示。A 船、B 船與C 船的網(wǎng)格量分別為1 278 萬、1 240 萬和1 217 萬。

圖3 網(wǎng)格細(xì)節(jié)Fig.3 Details of simplified frigate ships’ grid

2 計算方法驗證

為了準(zhǔn)確地研究艦船艉流場特性并驗證CFD 的計算效果，國際項目TTCP[18]得到了簡化的護(hù)衛(wèi)艦船型SFS2，其模型及尺寸如圖4所示。加拿大國家研究委員會對兩個簡化模型都做過一系列的風(fēng)洞實驗，得到了飛行甲板及前部上層建筑中特定位置的平均速度場及湍流數(shù)據(jù)[19-20]。

為了驗證本文CFD 方法的有效性，以SFS2 模型作為標(biāo)準(zhǔn)算例，針對0°風(fēng)向角下的流場進(jìn)行了計算。圖5給出了艉部甲板上方直線的時均速度分布曲線與風(fēng)洞實測數(shù)據(jù)的對比，該直線與機(jī)庫等高度，位于甲板y-z 平面內(nèi)正中間位置，橫坐標(biāo)為坐標(biāo)值與甲板寬度d 的比值，縱坐標(biāo)為速度分量u、v、w與來流速度U∞的比值。從圖5 中可以看出本文的計算結(jié)果總體上與風(fēng)洞實測數(shù)據(jù)相吻合，這表明本文采用的計算方法是有效的。

圖4 SFS2 模型及尺寸細(xì)節(jié)Fig.4 The geometry of SFS2 and its dimensions

圖5 0°風(fēng)向角計算結(jié)果對比Fig.5 Headwind mean velocities normalized by U∞

3 不同機(jī)庫外形對艦船艉流場的影響

機(jī)庫的幾何形狀是影響艉流場回流區(qū)的重要因素之一，保證船體整體尺寸外形一致，探究不同機(jī)庫縱截面外形對艉流場的影響。

圖6 0°風(fēng)向角下y/b=0 截面處速度云圖與流線圖Fig.6 Streamlines and velocity magnitude contours for a headwind

圖6 為0°風(fēng)向角下y/b=0 截面的速度流線圖與云圖。從圖中可以看出，A 船回流區(qū)縱向長度約為2.05 倍機(jī)庫高度，B 船回流區(qū)縱向長度約為1.95 倍機(jī)庫高度，C 船回流區(qū)縱向長度約為2.27 倍機(jī)庫高度。回流區(qū)范圍差距不大，直升機(jī)著艦域均在回流區(qū)內(nèi)。但從速度云圖可以明顯看出，B 船飛行甲板上方的下洗速度要明顯小于A 船；在飛行甲板長度相同的情況下，A 船與B 船回流區(qū)縱向長度占據(jù)了整個甲板的約74%，而C 船約為57%；從速度云圖也可以看出C 船回流區(qū)更靠近機(jī)庫，直升機(jī)著艦域位于回流區(qū)邊緣處。

圖7 右舷45°風(fēng)向角下橫截面合速度梯度云圖Fig.7 Contours of velocity magnitude for the Green 45° case

圖7給出了右舷45°風(fēng)向角下橫截面x/l=0，25%，50%，75%，100%處的合速度梯度云圖?？梢钥吹剑诩装宓挠覀?cè)出現(xiàn)了y 方向的“回流區(qū)”，A 船速度梯度的變化程度明顯大于B 船，特別是在近機(jī)庫區(qū)域處（x/l=25%與x/l=50%），而相同位置處C 船速度梯度的變化程度要小于A 船。從由高到低的速度梯度變化中可以得到機(jī)庫上方邊界產(chǎn)生的剪切層位置，同位置處速度梯度的大幅變化會對直升機(jī)安全著艦產(chǎn)生一定的影響。

圖8 0°風(fēng)向角下俯視截面z/h=1.15 處合速度梯度云圖（上：時均結(jié)果，下：瞬時結(jié)果）Fig.8 Contours of mean (top) and instantaneous (bottom) velocity magnitude for a headwind,plotted at z/h=1.15 above the deck

圖8 與圖9 分別為0°與右舷45°風(fēng)向角下俯視截面z/h=1.15 處合速度的梯度云圖。由時均結(jié)果可以看出：0°風(fēng)向角下時均流場結(jié)構(gòu)較為整齊劃一，渦結(jié)構(gòu)簡單，流動會在建筑的邊緣產(chǎn)生分離，導(dǎo)致了自由剪切層；B 船機(jī)庫兩側(cè)較低高度的建筑外形減弱了后方艉流場的速度分離，與A 船相比速度梯度變化較為平緩；45°風(fēng)向角下，A 船艉部右舷處的速度梯度變化程度明顯大于B 船，且影響范圍更廣。由瞬時結(jié)果可以看出：0°風(fēng)向角下合速度在背風(fēng)處以及艉部飛行甲板位置處尾跡的不對稱性表明了船體的幾何建筑構(gòu)型會出現(xiàn)相干的渦脫落特性。0°與45°風(fēng)向角下，A 船與C 船飛行甲板上空的流動均較為紊亂，且渦區(qū)的覆蓋范圍更廣，而該位置一般為直升機(jī)在艉部飛行甲板上方高空盤旋的高度，說明等高機(jī)庫的外形使得艉流場具有相對紊亂的湍流結(jié)構(gòu)，會對直升機(jī)在飛行甲板上方安全盤旋造成不利的影響。

圖9 右舷45°風(fēng)向角下俯視截面z/h=1.15處合速度梯度云圖（上：時均結(jié)果，下：瞬時結(jié)果）Fig.9 Contours of mean (top) and instantaneous (bottom) velocity magnitude for Green 45° case,plotted at z/h=1.15 above the deck

圖10 與圖11 分別為0°與右舷45°風(fēng)向角下艉部飛行甲板上方中心直線的下洗速度分布曲線，該直線的位置選擇基于直升機(jī)在甲板上空懸停時的位置。由結(jié)果可見，0°風(fēng)向角下，A 船最大下洗速度為-0.18U∞，B 船為-0.15U∞，C 船為-0.14U∞，最大下洗速度的位置均約與機(jī)庫同高。45°風(fēng)向角時下洗速度有更為明顯的梯度變化，約在0.36 倍機(jī)庫高度處達(dá)到極值，A 船最大下洗速度為-0.20U∞，而B 船為-0.15U∞，C 船為-0.12U∞。由此可見，在不同風(fēng)向角下B 船的最大下洗速度均小于A 船，且具有較低的機(jī)庫高度的C 船回流區(qū)關(guān)鍵區(qū)域處的下洗速度明顯減弱。

圖10 0°風(fēng)向角下甲板上方中心直線下洗速度分布曲線Fig.10 The downwash velocity component of the flight deck center for a headwind case

圖11 45°風(fēng)向角下甲板上方中心直線下洗速度分布曲線Fig.11 The downwash velocity component of the f light deck center for Green 45° case

0°與右舷45°風(fēng)向角下監(jiān)測點的瞬時速度時域特性如圖12 與圖13所示，監(jiān)測點位于飛行甲板正后方（x/l=0.8，y/b=0，z/H=0.6），該位置約為直升機(jī)進(jìn)艦的途經(jīng)點?？梢钥闯?，瞬時速度分量隨著時間發(fā)生周期性的脈動，表明了流動的非定常特性。0°風(fēng)向角下，A 船垂向速度最大振幅達(dá)到了-0.4U∞～+0.35U∞，B 船為-0.25U∞～+0.2U∞，C 船為-0.05U∞～+0.1U∞，A 船速度震蕩劇烈且幅值要遠(yuǎn)遠(yuǎn)大于B 船與C 船，且約為B 船的2 倍、C 船的4 倍；相同的趨勢也存在于45°風(fēng)向角下，此時A 船垂向速度的最大振幅達(dá)到了-0.6U∞～+0.8U∞，B 船為-0.05U∞～+0.6U∞，C 船為-0.05U∞～+0.5U∞。因此，A 船甲板后方的艉流場速度變化極快，會對直升機(jī)進(jìn)艦的操控性帶來嚴(yán)重的影響。

圖14 為0°風(fēng)向角下某瞬時時刻流場的渦量等值面圖。渦量的判別依據(jù)為λ2準(zhǔn)則，λ2為S2+ Ω2的2階特征值，其中S 為應(yīng)變率幅值，是速度梯度的對稱張量；Ω 為渦量幅值，是速度梯度的反對稱張量。從圖中可以明顯地看出，因艦艏上層建筑而產(chǎn)生的渦結(jié)構(gòu)強(qiáng)度隨著氣流向艉部流動而逐漸減小，但是遇到艉部機(jī)庫后渦結(jié)構(gòu)強(qiáng)度又會加強(qiáng)。A 船與C 船整體的渦結(jié)構(gòu)較為相似，而C 船在甲板處的渦結(jié)構(gòu)強(qiáng)度有所減弱。相比擁有與A 船同高機(jī)庫外形的B 船，B 船因其兩側(cè)較低的機(jī)庫高度使其在甲板處渦結(jié)構(gòu)加強(qiáng)的程度也有所減弱。

圖12 0°風(fēng)向角下監(jiān)測點瞬時速度時域特性Fig.12 Time domain characteristics of spectra point for a headwind case

圖13 45°風(fēng)向角下監(jiān)測點瞬時速度時域特性Fig.13 Time domain characteristics of spectra point for Green 45° case

圖14 0°風(fēng)向角下某瞬時時刻流場渦結(jié)構(gòu)Fig.14 Visualization of time-accurate vortex in the airwake (flooded by the magnitude of downwash flow)

4 艉流場回流區(qū)范圍經(jīng)驗公式探究

Schulman[21]總結(jié)并推導(dǎo)了0°風(fēng)向角下二維后臺階問題中回流區(qū)范圍的經(jīng)驗公式。后臺階的建筑外形會影響流場回流區(qū)的整體結(jié)構(gòu)，建筑外形包括其高度（H），迎風(fēng)寬度（W）以及氣流流經(jīng)臺面的長度（L）。Wilson[22]給出了評價后臺階的長度尺度R 為

式中：BS為較小的迎風(fēng)高度H，BL為較大的迎風(fēng)寬度W。

當(dāng)L＞0.9R 時，流動會在建筑邊緣分離后形成自由剪切層，隨著流動的發(fā)展會在某一位置與后方壁面發(fā)生再附，此時回流區(qū)的最大高度為HR=H。Fackrell[23]給出了再附點與建筑背風(fēng)面之間的距離，即回流區(qū)縱向范圍LR為

用機(jī)庫高度H 對其進(jìn)行無量綱化，得LR，H為

由A 船與C 船的機(jī)庫外形尺寸：

HA=6.5 m，HC=4.5 m，WA=WC=16 m，LA=LC=18 m

可得A 船與C 船的長度尺度RA=8.76 m 和RC=6.88 m，由

LA=18 m＞0.9RA=7.88 m，LC=18 m＞0.9Rc=6.19 m

可知當(dāng)氣流流經(jīng)A 船與C 船的機(jī)庫后會與甲板產(chǎn)生再附，其位置為

這與前文所述結(jié)論完全吻合。

Schulman[21]給出了二維后臺階問題中回流區(qū)垂向范圍變化公式：

對于艦船艉流場回流區(qū)問題，由公式（3）可以確定不同機(jī)庫外形的長度尺度，即機(jī)庫外形的幾何因子R。本文根據(jù)艦船特殊的建筑構(gòu)型與公式（5），并對非定常數(shù)值仿真計算結(jié)果進(jìn)行擬合，得到了確定等高機(jī)庫艉部不同橫向位置處回流區(qū)縱向長度的經(jīng)驗公式LR* 和不同縱向位置處回流區(qū)垂向高度的經(jīng)驗公式zH：

式中：下標(biāo)H 是用機(jī)庫高度無量綱化的值。對于公式（8），當(dāng)y＜0 時取其絕對值。

圖15 A 船經(jīng)驗公式驗證Fig.15 Verification of empirical equation for ship A

通過經(jīng)驗公式LR*和zH可以對等高機(jī)庫外形的艉流場回流區(qū)范圍進(jìn)行預(yù)估。圖15 與圖16 分別是對A 船與C 船艉流場回流區(qū)范圍經(jīng)驗公式的驗證，由圖可以看出本文所推導(dǎo)的經(jīng)驗公式能夠準(zhǔn)確預(yù)測回流區(qū)的范圍。

圖16 C 船經(jīng)驗公式驗證Fig.16 Verification of empirical equation for ship C

5 結(jié)論

本文采用基于SA 湍流模型的DES 方法對艦船艉流場進(jìn)行了非定常數(shù)值模擬，詳細(xì)分析了0°與右舷45°風(fēng)向角下不同機(jī)庫外形尺寸對艉部流場的影響，得到的主要結(jié)論如下：

（1）雖然具有等高機(jī)庫外形的A 船與兩側(cè)具有較低機(jī)庫高度外形的B 船在艉部產(chǎn)生的回流區(qū)范圍大體一致，但是和B 船相比，A 船更不利于直升機(jī)在艉部飛行甲板上方安全起降與懸停，具體表現(xiàn)為：A 船艉流場速度梯度變化程度大，飛行甲板上空流動更為紊亂，渦區(qū)的覆蓋范圍更廣且機(jī)庫后方渦結(jié)構(gòu)強(qiáng)度更強(qiáng)，飛行甲板中心直線上的下洗速度比B 船分別大3%（0°風(fēng)向角）與5%（右舷45°風(fēng)向角），飛行甲板艉部上方的垂向速度振幅過大，相當(dāng)于B 船的2 倍；

（2）降低機(jī)庫高度會使回流區(qū)位置相對前移，因此通過合理設(shè)計機(jī)庫高度可以使回流區(qū)有效避開直升機(jī)的著艦域。較低的機(jī)庫高度使得艉部速度梯度變化趨于平緩，且大幅度降低了關(guān)鍵位置處的下洗速度和艉部渦結(jié)構(gòu)的強(qiáng)度，同時提升了直升機(jī)進(jìn)艦的穩(wěn)定性；

（3）通過經(jīng)驗公式可以快速預(yù)測艉流場回流區(qū)的范圍與位置，為直升機(jī)著艦提供相應(yīng)的參考信息，增加其著艦的安全性。