基于計算流體動力學(xué)方法的兩棲飛機規(guī)則波中運動特性研究

胡開業(yè), 石林飛, 周輝, 毛麗君

(哈爾濱工程大學(xué) 船舶工程學(xué)院,黑龍江 哈爾濱 150001)

水陸兩棲飛機具備獨特的環(huán)境適用性能,可在水面和陸上實現(xiàn)起降,廣泛應(yīng)用于消防、海上救援、運輸、巡查、旅游開發(fā)等領(lǐng)域,彌補了救援直升機及船舶航速慢、航程短、載量小等方面的不足。目前部分國家研發(fā)設(shè)計了具備優(yōu)秀性能的水陸兩棲飛機,如俄羅斯別利耶夫設(shè)計局研發(fā)的高航速水上飛機Be-200,日本新明和公司研發(fā)的US-2水上飛機,以及加拿大福恩公司、三叉戟公司研發(fā)的CL-215和CL-415型水上飛機[1-2]。近年來我國研制的AG600水陸兩棲飛機,代表我國水陸兩棲飛機自主設(shè)計制造能力躋身世界先進水平。

水陸兩棲飛機的深V型船底結(jié)構(gòu)與斷階滑行艇相似,李曉文等[3-4]基于Star-CCM+軟件對滑行艇水面直航進行數(shù)值模擬,驗證了該平臺在研究滑行艇波浪中運動水動力性能時具有一定可行性和準(zhǔn)確性。黃淼等[5-7]對水陸兩棲飛機抗浪能力做了大量試驗研究,分析船底結(jié)構(gòu)及波浪要素對飛機滑行運動及機身載荷的影響,得到了模型在波浪上運動響應(yīng)曲線隨遭遇頻次、波浪、速度的變化規(guī)律,以及飛機在波浪水面發(fā)生劇烈運動響應(yīng)時對應(yīng)的波浪范圍。蔣榮等[8]采用船舶耐波性分析中常用的層次分析法對水上飛機的抗浪性能進行了評估。Chambliss等[9]研究了V型剖面斜升角對運動性能的影響。段旭鵬等[10]則提出了一種基于OpenFoam開源平臺的水面飛行器波浪中數(shù)值模擬技術(shù)。對于水陸兩棲飛機著水性能的研究,往往是基于數(shù)值計算方法模擬初始條件以及波浪因素對水陸兩棲飛機著水性能的影響[11-15]。水陸兩棲飛機在水上滑行起降時,飛機升沉、俯仰、尾流現(xiàn)象、噴濺特性、起降距離和機身底部壓力變化特性,以及數(shù)值模擬技術(shù)對于處理水陸兩棲飛機波浪中運動問題的適應(yīng)性,也是以往研究的熱點之一[16-19]。

水陸兩棲飛機波浪中加速航行時縱搖和升沉的耦合運動具有強非線性,受到波浪因素和速度影響較大,導(dǎo)致高海況非定常運動全過程運動的研究難度大。本文基于Star-CCM+平臺對水陸兩棲飛機波浪中滑行起飛過程進行數(shù)值仿真,以波長和波高為變量,探討加速滑行過程中飛機升沉、縱搖運動響應(yīng)以及機身重心處垂向加速度幅值

1 數(shù)值計算模型

1.1 計算模型

本文研究對象為某型水陸兩棲飛機,該飛機具有V型船身結(jié)構(gòu)、大展弦比懸臂上單翼、抑波槽以及船底中部的斷階結(jié)構(gòu),數(shù)值計算幾何模型按1∶ 8.5倍等比例縮小,尺寸如表1所示。水陸兩棲飛機在波浪中滑行時受到水動升力、水動阻力、氣動升力、氣動阻力、浮力、重力、以及發(fā)動機推力共同作用,航行穩(wěn)定性影響因素多,數(shù)值計算中模型為節(jié)省網(wǎng)格和提高計算效率,省略螺旋槳結(jié)構(gòu),幾何模型如圖1所示。

表1 模型飛機主要參數(shù)Table 1 Main parameters of model aircraft

1.2 數(shù)學(xué)模型

流體流動滿足質(zhì)量守恒、動量守恒和能量守恒三大基本守恒定律,引入湍流模型求解RAND方程具備較好的湍流中物體高速運動模擬精度,本文選取SSTk-ω模型,應(yīng)用有限體積法劃分計算域控制體,SIMPLE算法求解速度和壓力的離散方程,對流項采用二階迎風(fēng)格式,應(yīng)用壁面函數(shù)法求解近壁面流場,流體域體積方法捕捉水氣兩相流自由交界面。

質(zhì)量守恒方程又稱連續(xù)性方程,若速度得3個方向分量表示為uνw,則有不可壓縮流體滿足:

(1)

流體域體積模型是建立在歐拉網(wǎng)格下的表面跟蹤方法,廣泛適用于不相容的多項流交界面捕捉,適用于研究具有噴濺現(xiàn)象發(fā)生的強非線性問題。在STAR-CCM+中采用流體域體積多相模型來捕捉流體交界面,并可預(yù)報不相容相交界面的分布和移動。不相容相的位置和分布由各個相的體積分?jǐn)?shù)q來描述。相1的體積分?jǐn)?shù)定義為:

q1=V1/V

(2)

式中:q1為相1的體積分?jǐn)?shù);V1為網(wǎng)格單元中相1的體積;V為網(wǎng)格單元的體積。

一個網(wǎng)格單元所包含的所有相的體積分?jǐn)?shù)的總和一定為1。

針對本文采用的水相和空氣相組成的多相流,若水相為相1,則當(dāng)q1=1時,表示此網(wǎng)格單元內(nèi)只有水相;當(dāng)q1=0時,表示此網(wǎng)格單元內(nèi)只有空氣相;當(dāng)0

1.3 流體域邊界條件

本文數(shù)值模擬涉及水陸兩棲飛機滑行起飛過程中的升沉和縱搖運動,Star-CCM+平臺提供動態(tài)流體相互作用模塊(DFBI)可對飛機3個自由度運動進行模擬[20]。如圖2所示,因本文需要模擬兩棲飛機迎浪起降過程,故流體域設(shè)置速度入口距離飛機中心長度64 m,壓力出口距離飛機中心16 m,消波方法采用阻尼消波,阻尼消波區(qū)不小于一個波長。計算域頂部距離液面6 m,底部距離液面8 m,兩側(cè)距離飛機中心6.5 m,邊界條件均設(shè)置為速度入口。

圖2 流體域邊界示意Fig.2 Schematic diagram of fluid domain boundary

1.4 網(wǎng)格劃分

受航速和波浪影響,水陸兩棲飛機縱搖和升沉運動運動變化劇烈,重疊網(wǎng)格可以巧妙解決飛機大傾角運動問題,其中重疊區(qū)域隨剛體做大幅運動,而自由液面始終處于流體域網(wǎng)格加密區(qū),自由液面捕捉精度較高。

計算域網(wǎng)格采用具有較高精度和計算效率的結(jié)構(gòu)性網(wǎng)格,機體近壁面采用12層增長率為1.2的棱柱層網(wǎng)格,對機體周圍及液面區(qū)做網(wǎng)格加密處理。計算域采用半?？梢欢ǔ潭葴p少網(wǎng)格數(shù)量,但無法精確計算船底龍骨線附近湍流壓力變化,因此本文研究采用全模計算域,最終得到網(wǎng)格數(shù)量為1 578萬,計算采用64核雙路至強工作站,計算域網(wǎng)格劃分方案如圖3所示。

圖3 流場網(wǎng)格劃分方案Fig.3 Grid division scheme of flow field

2 數(shù)值計算驗證

2.1 數(shù)值造波

由于流體粘性和網(wǎng)格因素的影響,長域規(guī)則波傳遞過程中出現(xiàn)波浪衰減現(xiàn)象[21-22],為模擬實際尺度海況,將數(shù)值水池域尺寸設(shè)置為1 500 m×200 m×200 m,使用邊界造波法模擬數(shù)值規(guī)則波,壓力出口處設(shè)置阻尼消波區(qū),阻尼消波區(qū)長度為1.5倍波長,沿水池中縱剖面每間隔100 m設(shè)置了一個波浪監(jiān)測點,自由液面波形監(jiān)測點位置分布如圖4所示。

圖4 流體域網(wǎng)格及波形監(jiān)測點位置分布Fig.4 The distribution of fluid domain grid and waveform monitoring points

設(shè)置數(shù)值波浪波長λ為74 m(約2倍實機長度),波高h分別為1.5 m,波陡為0.010 1。本文對水線面附近網(wǎng)格進行3層加密以減少傳遞過程中的波浪衰減,其中最里層加密區(qū)波高方向網(wǎng)格不小于20個,波長方向網(wǎng)格不小于80個,時間步長設(shè)置兼顧庫朗數(shù)和波浪周期要求,庫朗數(shù)值小于0.2,周期值不小于T/500[23],通過對比時間步長為0.002 5、0.002、0.001 5、0.001 s時波浪衰減情況,最終選擇時間步長取值為0.002 s。取監(jiān)測點6處數(shù)值模擬波高和周期如圖5所示。

圖5 監(jiān)測點6處數(shù)值模擬波形與理論波形時歷曲線對比Fig.5 Comparison of the time history curve of the numerical simulation waveform and the theoretical waveform at the monitoring position 6

根據(jù)統(tǒng)計對比結(jié)果得出,監(jiān)測位置越靠近速度入口端,波浪數(shù)值模擬精度越高,各監(jiān)測點模擬波高相對于理論波高誤差在4.7%以內(nèi),模擬周期相對于理論周期誤差在2.6%以內(nèi),該數(shù)值造波方法精度較高。

2.2 網(wǎng)格無關(guān)性驗證

水陸兩棲飛機在起降滑行過程中速度變化極快,計算精度受近壁面棱柱層厚度影響較大,本文設(shè)置棱柱層增長率為1.2,總層數(shù)為12層,壁面率y+分別取30、100、200、300生成4套不同的網(wǎng)格。設(shè)置規(guī)則波波高h=0.06 m,波長λ=8 m,模型航速v=5 m/s。得到水陸兩棲飛機升沉、縱搖曲線與中試驗數(shù)據(jù)對比如圖6所示。各方案與試驗偏差如表2所示。

圖6 水陸兩棲飛機升沉、縱搖曲線與試驗數(shù)據(jù)對比Fig.6 Comparison of heave and pitch curves of amphibious aircraft with test data

表2 網(wǎng)格方案與試驗對比及偏差Table 2 Comparison and deviation between grid scheme and experiment

由表2可得方案B飛機滑行過程縱搖和升沉與試驗值偏差相對最小,符合工程實驗要求,后續(xù)數(shù)值仿真實驗中y+值取100,空氣和水的運動粘度相差為10倍,該數(shù)值同樣適用于機翼氣動升力的y+取值范圍。

3 數(shù)值計算分析

3.1 加速滑行起飛過程運動狀態(tài)分析

飛機水上滑行過程運動狀態(tài)可分為4個階段,分別為排水航行階段、過渡階段、滑行階段、起飛階段[24]。取工況波長為10 m,波高為0.165 m,飛機初始2 s加速度值勻速增加到2 m/s2,然后以恒定加速度2 m/s2運動,給定來流速度為5 m/s,設(shè)置初始流場穩(wěn)定時間10 s,采用前述CFD仿真方法,計算得到該工況水陸兩棲飛機波浪中加速起飛過程運動時歷曲線如圖7所示。

圖7 水陸兩棲飛機波浪中滑行時歷曲線Fig.7 Taxiing time history curve of amphibious aircraft in waves

由于數(shù)值模擬實驗飛機襟翼、升降舵角度固定,以及受模型尺度和波浪因素影響,數(shù)值模擬過程起飛離水速度Vga為20 m/s左右。流場穩(wěn)定時間內(nèi)無剛體運動,飛機排水航行階段(Vx≤0.25Vga,Vx表示飛機前進方向的速度)不被監(jiān)測。過渡階段(0.25Vga

3.2 波長對飛機滑行起飛過程運動特性的影響

研究波長對水陸兩棲飛機滑行起飛過程運動狀態(tài)的影響,根據(jù)機身長度和中國近海海域自然海況條件,選取波長范圍為2～20 m,每隔2 m選取為一組,波高選取0.147、0.182 m。

本文研究飛機加速度恒定為2 m/s2,滑行起飛過程中升沉、縱搖運動時歷曲線呈不穩(wěn)定非周期性。圖8、9給出了飛機升沉運動最值曲線和飛機縱搖運動最值曲線。由圖8～9可得,隨著波長增加,飛機升沉運動最大值、雙幅值最大值先增加,在波長為10 m附近出現(xiàn)小幅峰值,而后變化緩慢趨于穩(wěn)定。縱搖運動最大值、雙幅值最大值呈現(xiàn)先增加后緩慢減小趨勢,在波長為10 m附近出現(xiàn)峰值,縱搖最大值超過11°,雙幅值最大值超過13°。圖10、圖11給出了水陸兩棲飛機波浪上滑行起飛過程中垂向加速度和縱搖角加速度的最大值,隨著波長增加,飛機重心處垂向加速度、縱搖角加速度呈現(xiàn)出先增大后緩慢減小趨勢,在波長為8～10 m時出現(xiàn)峰值。由圖12可知,飛機垂向加速度值超過垂向加速度g的工況波長超過機身長度,在波長為10 m附近頻次最大。

圖8 飛機升沉運動最值曲線Fig.8 Maximum curve of aircraft heave motion

圖9 飛機縱搖運動最值曲線Fig.9 The curve of the maximum value of aircraft pitch motion

圖10 飛機重心處垂向加速度最值Fig.10 The maximum vertical acceleration at the center of gravity of the aircraft

圖11 縱搖角加速度最值Fig.11 Maximum pitch angular acceleration

圖12 垂向加速度大于g頻次分布Fig.12 Frequency distribution of vertical acceleration greater than g

由此可見,波長為2～3倍機身長度時,飛機波浪中加速起飛過程運動受波浪影響較大,對機身結(jié)構(gòu)和機組人員安全不利,應(yīng)盡量避免在波長為2～3倍機身長度海況中進行滑行起降工作。

3.3 波高對飛機滑行起飛過程運動特性的影響

研究波高對水陸兩棲飛機滑行起飛過程運動狀態(tài)的影響,根據(jù)機身長度和中國近海海域條件,選取波長為4、12 m,波高選取0.147、0.165、0.182、0.2、0.218、0.235 m。

水陸兩棲飛機升沉運動最大值,隨著波浪高度增加,12 m長波中飛機升沉最大值增加,在波高為0.182～0.218 m區(qū)間漲幅最大,而短波中飛機升沉運動最大值隨波高增加緩慢,如圖13所示。飛機縱搖運動最大值隨著波高增加而增加,如圖14所示,波長為12 m的波浪波陡較小,單個波浪周期內(nèi)飛機有較長加速距離,所以同一波高內(nèi)長波中飛機升沉、縱搖運動最大值較短波中更大。飛機重心垂向加速度最大值集中在高速滑行階段,長波中飛機躍出波峰時具有更大的垂向速度,隨波高增大垂向速度增大,故飛機升沉幅值隨波高增加而增加,船身接水時受到水動升力更大,重心垂向加速度更大,如圖15所示,重心垂向加速度超過重力加速度g的頻次隨著波高增大而增大。

圖13 飛機升沉運動最值曲線Fig.13 The maximum value curve of aircraft heave motion

圖14 飛機縱搖運動最值曲線Fig.14 The curve of the maximum value of aircraft pitch motion

圖15 飛機重心垂向加速度曲線Fig.15 Vertical acceleration curve of the aircraft′s center of gravity

由此可見。隨著波浪高度增加,飛機滑行起飛過程穩(wěn)定性下降,升沉、縱搖運動響應(yīng)加劇,對飛機結(jié)構(gòu)和機組人員安全影響不利。

3.4 機翼氣動升力計算分析

取工況波長為2 m,波高為0.147 m,得到數(shù)值模擬主翼氣動升力隨速度變化曲線,圖16給出機翼升力系數(shù)C隨機翼攻角變化曲線,代入氣動升力計算公式得到理論機翼升力變化曲線。氣動升力計算公式為:

圖16 飛機機翼升力系數(shù)曲線Fig.16 The lift coefficient curve of the aircraft wing

Y=ρCSv2/2

式中:ρ為空氣密度;C為升力系數(shù);S為機翼面積;ν為飛機速度。

由于波浪中機翼隨著機體做俯仰、升沉運動,升力數(shù)值模擬曲線為上升波浪線,由圖17觀察到,理論計算得到的機翼升力曲線與數(shù)值模擬值一致,最大偏差小于10%,同時驗證了本文數(shù)值計算結(jié)果的準(zhǔn)確度較高,本文研究的數(shù)值模擬方案可用于水陸兩棲飛機波浪中加速起飛過程機翼氣動特性計算。

圖17 飛機氣動升力速度曲線Fig.17 Aerodynamic lift speed curve of aircraft

4 結(jié)論

1)隨著波長增大,飛機升沉運動最大值、最大雙幅值隨之增大,在波長為2～3倍機身長度時出現(xiàn)小幅峰值,隨后緩慢變化,在波長為2～3倍機身長度時,飛機縱搖最大值和雙幅值出現(xiàn)峰值,隨后緩慢減小。飛機重心處垂向加速度、縱搖角加速度變化規(guī)律與縱搖運動一致,均在在波長為2～3倍機身長度時出現(xiàn)峰值,且垂向加速度超過g的頻次在該工況附近較大。

2)隨著波高增大,飛機升沉、縱搖運動及重心垂向加速度最大值隨之增大,且長波中的變化值較短波中更大。

3)機翼氣動升力數(shù)值模擬結(jié)果與理論計算值一致,同時驗證了本文方案的可行性,本文數(shù)值模擬方案可用于水陸兩棲飛機波浪中加速起飛過程機翼氣動力特性研究。