水陸兩棲飛機全機模型靜水滑行水動性能預報

張杜文, 李新穎, 吳彬, 曹楷

(1.航空工業特種飛行器研究所, 荊門 448035; 2.高速水動力航空科技重點實驗室, 荊門 448035)

21世紀以來,隨著全球化的逐漸深入、海洋資源的不斷開發以及現代海戰向立體化、多層次發展,人類對于海洋的認知和重視程度不斷提升,世界各國紛紛大力發展自身的海洋產業和裝備,以此掌握海上制空權,維護海洋權益,發展海洋經濟,保障領海安全。面對復雜的海洋形勢,中國近年來開始從國家層面上重視海洋權益,通過各種途徑來保證中國海洋國土安全,維護海洋的主權和領土完整,積極開展水面飛行器的研制工作[1]。

水面飛行器水面執行任務時會面臨復雜的水面環境和飛行環境,例如,較高的海況等級時,飛行器可能會有傾覆、側翻、墜毀等風險。因此,在水面飛行器設計階段,水面起飛性能是水面飛行器的基本性能,也是總體技術的核心,涉及多個學科領域。一直以來,水面飛行器性能評估與分析主要依靠物理水池模型試驗實現,隨著現代流體力學的發展和計算機處理超大規模數據能力的提高計算流體力學(computational fluid dynamics,CFD)已經成為和試驗、理論計算并列的重要研究手段[2]。理論解析方法一般只能解決簡單幾何模型并基于各種假設的簡化物理問題,試驗方法往往面臨成本高昂、測量數據有限、誤差控制和尺度效應等難點。基于此,CFD數值計算在實際應用和科學研究中越來越受到重視,與理論解析方法相比,數值計算能夠處理復雜的三維幾何模型和實際的工程問題;與試驗相比,數值計算成本較低、數據全面,其計算結果對于整體和局部細致的理解使實際問題更加清晰透徹,并且可以代替部分危險、昂貴和難以實施的試驗。

水陸兩棲飛機在水面滑行時,由于受風、流、浮力等外力作用使得各種運動之間相互影響,即存在非線性耦合作用。對水陸兩棲飛機的水動力學研究可以借鑒船舶研究的成果,水陸兩棲飛機水面滑行是一個動態平衡的過程,模型的運動必然會導致周圍網格的變形。采用數值模擬的方法研究水陸兩棲飛機的非線性運動,必然涉及較為復雜的包含動邊界的流動問題。水陸兩棲飛機全機模型水動性能研究對比單船身模型水動性能研究,不需要垂向氣動升力卸載,結果更加真實可靠。

針對水陸兩棲飛機水面滑行水動力性能數值預報存在的問題,現結合計算流體力學、船舶理論以及計算軟件數值模擬技術,從自由液面水氣分布、機身底部壓力分布特征、阻力、升沉、姿態角等方面對水陸兩棲飛機全機模型水面滑行進行較為深入的研究。

1 數值計算方法

任何形式的CFD計算都是基于流體力學基本控制方程:連續性方程、動量方程和能量方程,這些方程表述的是物理原理,它們是所有流體力學都必須遵循的三大基本物理定律的數學表述[3]。

1.1 連續性方程(質量守恒方程)

連續性方程可表述為控制體控制面的凈質量流出等于控制體內質量隨時間的減少量,其表達式為

(1)

式(1)中:ρ為流體密度,kg/m3;t為時間,s;u、v、w分別為直角坐標系中x、y、z方向上的速度分量。

矢量形式為

(2)

式(2)中:V為速度矢量。

1.2 動量方程(動量守恒方程)

動量方程可表述為控制體上的合力等于控制體的質量與控制體的加速度的乘積,其表達式為

(3)

式(3)中:fx、fy、fz分別為單位質量流體在x、y、z3個方向上的體積力;p為壓力時均值;μ為動力黏性系數。

1.3 湍流模型

由于RANS方程引入了未知的雷諾應力項,所以方程不是閉合的。根據Boussinesq提出的渦黏假設,需要引入湍流模型建立相對雷諾應力與平均速度梯度的關系。方程為

(4)

SSTk-ω(shear stress transportk-ω)模型稱為剪應力輸運k-ω模型,此模型對湍流運動中的剪切應力做了相應的考慮,對模型中的常數項進行了修正,增加了剪應力輸運項,為可以模擬近壁面流動的低雷諾數模型,因此本文湍流模型采用SSTk-ω結合壁函數進行處理。

目前的湍流數值模擬方法根據簡化方法以及簡化程度不同可分為直接數值模擬法(direct nume-rical simulation,DNS)、大渦模擬法(large eddy simulation,LES)和雷諾平均法(Reynolds average navier-stokes,RANS)[4]。RANS方法將流場中的瞬時變量分為時均量和脈動量,然后對N-S方程進行時間平均得到了RANS方程。在許多工程實踐中,該方法可以在數值精度和計算成本之間實現更好的平衡,因此RANS方法是當前使用最廣泛的湍流數值模擬方法。水面飛行器的數值模擬均可使用該方法進行。

VOF方法(volume of fluid)具有容易實現、計算量小和精度高等優點并且可以處理自由面的大變形和自由面的拓撲結構發生變化等復雜三維情況,因此自由液面捕捉采用VOF方法[5]。網格采用非結構動網格技術進行相關的處理,近壁面處理采用壁面函數法,數值離散方法采用有限體積法,離散格式采用二階迎風格式。

2 水陸兩棲飛機全機模型水動性能數值計算

2.1 幾何模型

以某水陸兩棲飛機簡化全機模型為研究對象,略去抑波槽、襟翼等結構,該飛機外形三維示意圖如圖1所示。

圖1 某水陸兩棲飛機三維示意圖Fig.1 A three-dimensional diagram of an amphibious hydroplane

在進行水動力性能計算時,由于存在機翼升力,不需要在船體重心位置施加一個垂向力。對于模型兩相流運動,利用DFBI模型(dynamic fluid body interaction)根據流動引起的作用力來模擬飛機模型的運動。在DFBI的運動設置中選擇Y軸旋轉,Z軸平移,流體沿X軸正方向以一定的速度流過飛機模型。

2.2 網格劃分

由于全機模型為左右對稱結構,為節省計算資源,僅對半模進行仿真計算,將模型導入后,建立流場計算域。由于水面飛行器在滑行過程中姿態變化較大,機身會出現較大的縱傾,因此采用重疊網格的方式進行網格劃分[6]。建立一個小的幾何體將模型包圍作為重疊域,并為其建立局部笛卡爾坐標系,計算過程中僅重疊域發生姿態變化。整個流場區域由小型重疊域及流場大域組成。

重疊域與流場大域的網格進行獨立劃分。網格劃分時選取的計算域大小對計算精度有一定的影響,足夠大的計算域能避免因邊界過小而造成計算精度問題。計算域選定之后,對網格進行劃分,應遵循如下幾點要求:①流場大域的最大網格單元尺寸設定為L/10(L為機身長度);②重疊域網格基本尺寸等于L/300;③保證模型長度方向上網格數量不小于200個,邊界層按Y+=60設定;④對自由液面進行加密,單層網格高度按L/1 000設定,長細比為4∶1。

按照以上要求,劃分的得到的計算區域網格如圖2所示,模型表面網格如圖3所示。

圖2 計算區域網格示意圖Fig.2 The calculate regional grid diagram

圖3 模型表面網格示意圖Fig.3 Model surface grid diagram

2.3 時間步長的選擇

庫朗數(Courant number,CFL)定義如式(5)所示,為來流速度U與網格尺度Δx的比值,通常為了計算的穩定性,庫朗數一般要求小于或等于1。然而,通常在隱式非定常的計算中,時間步長一般由流體特性決定而不僅僅是庫朗數。為了在合理的運行時間內得到準確的預報結果,通常時間步長根據ITTC(2011)推薦設置,本文研究中的時間步長均參照式(6)設置。

(5)

Δt=0.005～0.02L/U

(6)

2.4 邊界條件

在仿真計算中,邊界條件的設置非常重要,必須確保邊界不會影響到黏性流場的求解。在計算前需要根據經驗對飛機的運動姿態進行預估[7]。在設置邊界條件時,由于模型幾何對稱,采用半機身模型進行仿真模擬,因此對稱面采用對稱邊界條件。其余邊界的邊界條件:來流入口、底部、頂部和側面采用速度入口,尾流出口采用壓力出口,在機身表面定義無滑移的壁面邊界條件[5],計算域的尺寸和邊界條件設置如表1和圖4所示。

表1 邊界條件設置Table 1 Boundary condition setting

圖4 計算區域劃分示意圖Fig.4 Calculate area division diagram

對全機模型進行數值計算時,基本參數設置如下:①模型重量、重心和慣量;②采用隱式非定常;③湍流模型為SSTk-ω湍流模型;④多相模型組合中選擇流體域體積(VOF);⑤最大非線性迭代次數為10。

2.5 數值計算結果

對各個速度下全機模型的靜水面滑行進行數值模擬,得到的數值計算結果如表2所示,選取長度傅汝德數為橫坐標繪制阻力、升沉和姿態曲線圖如圖5所示。可以看出,模型升沉值隨著速度的增加而增大,而阻力和縱傾角整體上隨速度的增加呈現出先增大后減小的趨勢,阻力峰出現在約長度傅汝德數Fr=1.294的速度狀態下,此時R/W=0.207,R為阻力,W為排水量。

表2 模型數值計算結果Table 2 Numerical results of the model

圖5 數值計算結果曲線圖Fig.5 Numerical calculation results curve

第一阻力峰出現在Fr=0.728時,這是因為低速段的縱傾和升沉并未產生大幅變化,導致斷階后方并未產生氣穴區,因此,Fr=0.728時的阻力較Fr=0.647時有了較大增幅。而隨著速度的增大,水動升力隨之增大,下機身底部斷階部位出現空穴區,使得相鄰速度的阻力出現減小的趨勢。而越峰后,下機身底部空穴區域面積隨著速度的增大而不斷增大,對飛機摩擦阻力的降低起到了積極的作用。因此,在Fr=1.294后,飛機水阻力出現了逐漸降低的趨勢。

3 計算結果的試驗驗證

3.1 全機模型水池拖曳試驗

為了驗證數值計算方法的可行性,按Fr相似,縮尺比為1∶10進行水陸兩棲飛機全機模型的制作和物理水池拖曳試驗。

全機模型經過了陽模-陰模-蒙皮-合模-拆模-噴漆整個生產制作過程。該全機模型采用金屬框架與航空層板框架相配合的方式,并采用碳管進行對接穿插,蒙皮采用玻璃鋼材料,經檢驗,模型制作完成后外形尺寸合格,工藝良好,滿足技術設計要求。按照該水陸兩棲飛機的相關參數對模型重量、重心、慣量進行調試,調試完成后將模型安裝到拖車上進行試驗。試驗過程中保持模型航向、俯仰和垂蕩3個自由度,利用阻力傳感器、傾角傳感器等監測模型航行過程中的阻力、姿態。試驗安裝示意圖如圖6所示,模型水池試驗示意圖如圖7所示。

WL表示水線位置

圖7 水陸兩棲飛機模型水池試驗Fig.7 Amphibious hydroplane model tank test

3.2 數值計算結果驗證

為了驗證數值計算結果的準確性,將數值計算結果與水池試驗結果進行對比如圖8所示。

圖8 計算與試驗結果對比圖Fig.8 Comparison of calculation and test results

由圖8可知,數值計算結果與試驗結果較為吻合,阻力值誤差在10%以內,縱傾角誤差在0.8°以內,升沉趨勢較為吻合,誤差在15 mm以內,數值計算方法可行,數值計算結果真實可信。

4 全機模型數值計算流場分析

4.1 自由液面波形云圖

由圖9～圖17所示的不同速度下自由液面分布云圖,可以看出,隨著速度的提高,飛機兩側的興波以及尾部的波形發生了不同的變化。在排水航行階段,各速度下的行波基本相似,其開爾文角未發生明顯變化,而隨著速度的增大,模型縱傾角及升沉均有較大程度的提高,導致浸濕面積減小,尾部開爾文波的寬度變窄,開爾文角逐漸變小。同時可以看出,隨著速度的提升,“雞尾流”出現的相對位置逐漸后移[8],這與實際的流動現象是一致的。

圖9 Fr=0.485自由液面波形云圖Fig.9 Free surface waveform cloud diagram at Fr=0.485

圖10 Fr=0.647自由液面波形云圖Fig.10 free surface waveform cloud diagram at Fr=0.647

圖11 Fr=0.728自由液面波形云圖Fig.11 Free surface waveform cloud diagram at Fr=0.728

圖12 Fr=0.890自由液面波形云圖Fig.12 Free surface waveform cloud diagram at Fr=0.890

圖13 Fr=1.051自由液面波形云圖Fig.13 Free surface waveform cloud diagram at Fr=1.051

圖14 Fr=1.294自由液面波形云圖Fig.14 free surface waveform cloud diagram at Fr=1.294

圖15 Fr=1.617自由液面波形云圖Fig.15 Free surface waveform cloud diagram at Fr=1.617

圖17 Fr=2.426自由液面波形云圖Fig.17 Free surface waveform cloud diagram at Fr=2.426

對全機模型中縱剖面的水氣分布進行分析,可以看出,在Fr=0.485～0.728低速狀態下,下機身底部并未發生較為明顯的水氣分離,因此,導致第一阻力峰出現在了Fr=0.728速度處。速度自Fr=0.809起,下機身底部開始出現較為明顯的氣穴,且隨著速度的提升氣穴長度也在不斷增加,這也是模型水阻力在Fr=1.294速度后未再隨著速度的增大而增加的主要原因。

4.2 下機身流體體積分布云圖

圖18所示的是各個速度下的下機身水體積分布云圖。在較低的速度(Fr=0.485～0.647)下,由于機身的縱傾角及升沉的變化較小,斷階并未起到明顯的阻斷水流的作用,后體的下機身表面主要以水相為主。而隨著速度的增大,飛機的縱傾角增大,斷階處產生了流動分離,導致斷階后形成的氣穴長度和氣穴區面積不斷增加,水氣分離明顯,增大了滑行面的展弦比,降低了下機身浸濕面積,這也是水面飛行器斷階的減阻機理所在[9]。下機身水相分布形狀與壓力駐線較為相似,均呈“三角”狀,這是因為飛機穩定滑行后,下機身與水面的交線(即壓力駐線)至尾部均為水相。

圖18 下機身流體體積分布云圖Fig.18 Lower fuselage fluid volume distribution cloud

4.3 下機身底部壓力云圖

圖19所示的是各個速度下的下機身底部壓力分布云圖。可以看出,下機身底部壓力的分布主要呈現出“三角”的形狀,這是因高速水面航行體的流動特性所致。在全機模型穩定滑行的過程中,下機身底部與水面的交線形成兩道壓力等值線,即壓力駐線[10],下機身水下部分的前緣點為壓力駐點,即“三角”駐線的尖點。通常情況下,下機身底部壓力駐線與舭線的交點處存在高壓,這是因為滑行過程中水流產生的噴濺會沿壓力駐線向機身兩側飛濺,水流的沖擊作用導致了高壓的產生。

圖19 下機身底部壓力云圖Fig.19 Lower fuselage pressure cloud view

5 結論

采用數值計算與模型試驗結合的方法對水陸兩棲飛機全機模型的水動力性能進行了預報與深入研究。首先結合水陸兩棲飛機水面高速滑行特點,利用RANS數值方法和重疊網格技術對水陸兩棲飛機全機模型開展數值仿真模擬,分析了自由液面水氣分布、機身底部壓力分布特征,并將阻力、姿態和升沉與水池模型試驗結果進行對比,可以得出以下結論。

(1)計算水陸兩棲飛機全機模型水動力性能時,模型繞流場區域采取分區混合網格劃分方案。網格的質量要求較高,保證模型長度方向上網格數量不小于200個,對自由液面進行加密,單層網格高度按L/1 000設定,長細比為4∶1效果較好。

(2)基于仿真計算結果,從自由液面、流體體積分布、壓力分布等方面對水陸兩棲飛機全機模型的水動力特性進行了分析,較好地解釋了該模型的水動力性能規律。

(3)對水陸兩棲飛機全機模型進行了試驗驗證,通過對比計算結果與試驗結果發現:阻力值誤差在10%以內,縱傾角誤差在0.8°以內,升沉變化趨勢較為吻合,誤差在15 mm以內,數值計算方法可行,數值計算結果真實可信。