米百剛，詹 浩

（西北工業(yè)大學(xué) 航空學(xué)院，西安 710072）

近地、水面時的飛行器動態(tài)穩(wěn)定特性數(shù)值模擬

米百剛，詹 浩

（西北工業(yè)大學(xué) 航空學(xué)院，西安 710072）

基于非結(jié)構(gòu)動態(tài)網(wǎng)格重構(gòu)技術(shù)，采用SA模型以及VOF方法，分別對NACA0012二維翼型自由空間、近地效應(yīng)以及近水面效應(yīng)的非定常運動流場進(jìn)行了數(shù)值計算分析，并根據(jù)建立的基于小幅度強迫振動的組合動導(dǎo)數(shù)辨識方法計算了各個攻角下動導(dǎo)數(shù)。計算結(jié)果分析表明動導(dǎo)數(shù)計算方法準(zhǔn)確可靠；中小攻角下，地面與水面以管道效應(yīng)影響流場，均損失了動態(tài)阻尼，柔性的水面損失更多，動導(dǎo)數(shù)絕對值更小；大攻角時，兩者表現(xiàn)為阻塞效應(yīng)，柔性水面下的動態(tài)阻尼增大，動導(dǎo)數(shù)絕對值也變大；相比地面效應(yīng)，柔性的水面影響范圍更廣。

自由空間；地面效應(yīng)；水面效應(yīng)；動導(dǎo)數(shù)；數(shù)值模擬；計算流體力學(xué)

0 引 言

當(dāng)飛行器貼近地面或者水面飛行時，受到地面或者水面效應(yīng)的影響[1-2]，飛行器的氣動力會發(fā)生不同于無限空間自由流場的變化，充分地利用這些效應(yīng)可以設(shè)計出地效飛行器，此類飛行器載重量大，經(jīng)濟性好，無需專門修建的起降場地，低空飛行隱蔽性好，無論是在民用還是軍用上都有著非常廣闊的前景，因而得到了世界各國的重視。同時，隨著飛行器氣動設(shè)計的快速發(fā)展，非定常動態(tài)特性的研究成為一個熱點，與其相關(guān)的現(xiàn)象如飛行器的顫振，動態(tài)操縱以及動態(tài)穩(wěn)定性等在當(dāng)前的飛行器設(shè)計中表現(xiàn)得越來越重要，無論是針對常規(guī)飛行器還是地效飛行器。

目前對于地效飛行器的氣動特性的研究多集中在常規(guī)的氣動力分析計算上。比如國內(nèi)的盛其虎[3]使用面元法來對風(fēng)浪擾動中的飛行器非定常氣動力進(jìn)行了計算分析；鄒光遠(yuǎn)等[4]發(fā)展了邊界元法來求解近水面飛行時的全機氣動力水波效應(yīng)，給出了水波效應(yīng)隨攻角、飛行高度、入射角等參數(shù)的變化；劉沛清等[5-6]使用基于求解NS方程的方法對飛行器的波浪地面以及波浪水面上的氣動力進(jìn)行了數(shù)值計算，分析了飛行器在不同波長、波高等因素影響下的氣動力變化特性，較為精確地揭示了地面以及水面對飛行器的氣動設(shè)計的影響作用。國外的研究中，Schweikhard[7]采用離散渦模型計算了二維機翼的動態(tài)地面效應(yīng)，NAhait[8]則采用渦格法建立了計算三維地效機翼的數(shù)值方法并考慮了厚度的影響；Kornev等[9]基于渦格法以及非定常線化理論假設(shè)發(fā)展了估算地效下的復(fù)雜構(gòu)型氣動導(dǎo)數(shù)的數(shù)值方法。可以看出，目前的研究對飛行器地、水面效應(yīng)下的非定常動態(tài)運動特性涉及很少，而地效飛行器長時間貼近地面或者水面飛行，其動態(tài)運動時引起的非定常氣動力特性與常規(guī)飛行器有很大差別，因而，分析研究地效飛行器動態(tài)運動氣動特性十分迫切。

本文結(jié)合非定常動態(tài)非結(jié)構(gòu)網(wǎng)格重構(gòu)技術(shù)，基于求解非定常NS方程，對二維翼型貼近地面以及水面飛行時的非定常動態(tài)運動流場進(jìn)行數(shù)值模擬研究，對表征非定常動態(tài)穩(wěn)定性的指標(biāo)—動導(dǎo)數(shù)進(jìn)行計算，分析不同近地/近水高度對動導(dǎo)數(shù)的影響。結(jié)果表明，地面/水面的存在劇烈地影響了翼型的動態(tài)穩(wěn)定特性，近地/近水高度不同，影響程度也不一樣，同樣高度下，地面與水面的影響程度也存在差異，說明在地效飛行器的非定?？諝鈩恿υO(shè)計中，地面、水面的干擾應(yīng)該被綜合統(tǒng)籌考慮。

1 控制方程

文中采用二維積分形式的非定常雷諾平均N-S方程，即

其中：V0為任意控制域；W是守恒變量；F為無粘 （對流）通矢量項；Fv為粘性通量；?V0為控制域的邊界；n為控制域邊界單位外法向矢量；Re為計算的雷諾數(shù)。

空間離散采用二階迎風(fēng)格式—量差分分裂（Roe-FDS）格式，使用用LU-SGS雙時間時間推進(jìn)法求解非定常控制方程。湍流模型使用SA一方程模型。

2 VOF模型

本文的計算涉及到水面上的非定常氣動力計算，此時的介質(zhì)包括水和空氣，因此使用VOF模型定義自由水面。VOF模型是一種在固定的歐拉網(wǎng)格下的表面跟蹤方法，通過求解單獨的動量方程和處理穿過區(qū)域的每一流體的體積分?jǐn)?shù)來模擬兩種或三種不能混合的流體。當(dāng)需要得到一種或多種互不相融流體間的交界面時，可以采用這種模型。典型的應(yīng)用包括分層流、射流破碎、流體中的大泡運動、自由表面流動等。

采用VOF方法模擬自由水面，用幾何重構(gòu)方案模擬氣液界面，并對水氣交界區(qū)域網(wǎng)格進(jìn)行了加密。VOF模型中定義了一個體積函數(shù)αq，代表第q種流體占控制體的體積分?jǐn)?shù)，在每個控制體當(dāng)中，所有流體的體積分?jǐn)?shù)之和等于1，即αq=1表示流體當(dāng)中充滿第q種流體，αq=0表示流體中不含有第q種流體，0＜αq＜1表示控制體中含有第q種流體的體積分?jǐn)?shù)為αq。自由界面處滿足體積分?jǐn)?shù)方程：和計算控制條件兩種流體混合的單元采用如下公式計算密度：

動量方程為

3 動穩(wěn)定性導(dǎo)數(shù)計算方法及其驗證

3.1 計算方法

動導(dǎo)數(shù)[10-12]是表征飛行器動態(tài)穩(wěn)定特性的指標(biāo)，它是飛行器六個氣動系數(shù)（即三維直角坐標(biāo)系中的三個力和三個力矩系數(shù)）對飛行器姿態(tài)參數(shù)時間的變化率的導(dǎo)數(shù)，是飛行器導(dǎo)航系統(tǒng)，控制系統(tǒng)以及動態(tài)品質(zhì)分析不可缺少的原始?xì)鈩訁?shù)，直接影響著飛行器的飛行品質(zhì)。

以縱向為例，闡述本文求解動導(dǎo)數(shù)的方法。強迫飛行器做小幅度俯仰振蕩時，剛體飛行器所受非定常氣動力矩可以表示為：

上式實際上是對剛體飛行器小振幅強迫振動所受非定常氣動力矩進(jìn)行了泰勒展開，其中為氣動力矩靜導(dǎo)數(shù)，為氣動力矩對迎角的一階動導(dǎo)數(shù)，為氣動力矩對俯仰角速度的零階動導(dǎo)數(shù)，為氣動力矩對俯仰角速度的一階動導(dǎo)數(shù)，為高階導(dǎo)數(shù)項。

當(dāng)剛體飛行器以低頻ω做小幅度振蕩時，其模型運動方程可以簡化為：

展開（2）式并略去高次諧波分量，同時將（4）式代入略去高階分量Δ贊，整理同類項可得：

當(dāng)非定常問題計算足夠長時，令ωt=2nπ，就可以抹去初始效應(yīng)的影響，氣動力矩達(dá)到一個周期性穩(wěn)態(tài)值，于是（5）式可以寫為：

至此，基于傅立葉及泰勒展開方法推導(dǎo)出了縱向組合導(dǎo)數(shù)的表達(dá)式（7）。

根據(jù)國家標(biāo)準(zhǔn)，減縮頻率取為k=ωl/2V*，其中l(wèi)為參考長度，V*為遠(yuǎn)場自由來流速度，代入可得無量綱化動導(dǎo)數(shù)計算公式（8）

基于此方法，通過CFD方法得到穩(wěn)態(tài)氣動力以及非定常周期氣動力即可數(shù)值求解無限空間自由流場、地面效應(yīng)以及水面效應(yīng)影響下的動導(dǎo)數(shù)。

3.2 標(biāo)準(zhǔn)算例驗證

動導(dǎo)數(shù)的計算主要在于準(zhǔn)確地求解非定常氣動力及力矩系數(shù)。本節(jié)對前面開發(fā)的非結(jié)構(gòu)網(wǎng)格的翼型動態(tài)特性的數(shù)值模擬方法進(jìn)行驗證。驗證算例為AGARDCT3[12]，NACA0012翼型的運動狀態(tài)可以寫為

其中：α0=4.86°，αm=2.44°，減縮頻率 k=ωc/2V=0.080 8。

計算馬赫數(shù)為Ma=0.6，基于弦長c的雷諾數(shù)為480萬，翼型繞1/4弦長處做俯仰運動，力矩系數(shù)參考位置也選擇此點。圖1分別為使用本文方法計算得到的非定常法向力系數(shù)及俯仰力矩系數(shù)隨迎角的變化曲線，可以看出，計算值與試驗值符合得很好，本文的非定常強迫振動方法能夠很好地預(yù)測迎角與力/力矩系數(shù)之間的遲滯效應(yīng)。

圖1 非定常氣動力/力矩系數(shù)遲滯環(huán)Fig.1 Hysteresis loops of unsteady aerodynamic force and moment coefficient

根據(jù)本文的動穩(wěn)定性導(dǎo)數(shù)計算方法可以計算得到該狀態(tài)下的動導(dǎo)數(shù)為-3.342 3，與文獻(xiàn)值的-3.41吻合得較好。

4 計算模型、網(wǎng)格以及邊界條件

本文的計算針對三種情況進(jìn)行對比：一是無限空間自由流場下的的動導(dǎo)數(shù)計算；二是考慮地面效應(yīng)的動導(dǎo)數(shù)計算；三是考慮柔性水面效應(yīng)的動導(dǎo)數(shù)計算。

這三種情況的計算模型、網(wǎng)格及邊界條件各不相同，分別論述如下：

第一種情況是常規(guī)的飛行器在無限空間自由流場運動的情形，此時的計算模型及網(wǎng)格分別見圖2（a）、（b），邊界條件包括遠(yuǎn)場（farfield）以及無滑移壁面邊界條件。

第二種情況是考慮地面效應(yīng)影響的情形，此時翼型下部的邊界被地面所代替，所采用的邊界條件為滑移壁面條件，用來模擬地面相對飛行器的運動。翼型表面為無滑移壁面條件，其余邊界為遠(yuǎn)場條件。計算模型和網(wǎng)格分別見圖3（a）、（b），圖中V為自由來流速度，h為翼型尾緣點距離地面（水面）的距離，使用h/c來表示這一高度的相對值。

第三種情況是考慮水面影響的情形，此時翼型下部為具有一定深度的液態(tài)水，水面上的空氣部分邊界條件為速度入口，壓力出口，上部以及右部的出口與大氣相通；水面下的部分同樣是速度入口，壓力出口，此時的壓力由使用P=ρghw計算得到不同水深處的壓力，底部的壁面為滑移壁面。為了更為清晰地捕捉水面，交界附近的網(wǎng)格進(jìn)行了加密。圖4（a）、（b）分別顯示了水面影響時的計算模型以及網(wǎng)格。

圖2 無限空間自由流場Fig.2 Free flow field

圖3 地面效應(yīng)Fig.3 Ground effect

圖4 水面效應(yīng)Fig.4 Water surface effect

5 計算結(jié)果及分析

自由空間、地面效應(yīng)以及水面效應(yīng)的計算狀態(tài)如下：

馬赫數(shù)Ma=0.2

計算攻角：0、2、4、6、8 及 10 度

近地/水高度：h/c=1，2，5

減縮頻率為：

非定常運動的運動形式為：

其中：α0為初始攻角。

利用文中建立的基于非結(jié)構(gòu)網(wǎng)格重構(gòu)的方法實現(xiàn)非定常運動流場的求解，得到的非定常氣動力/力矩呈現(xiàn)周期性的變化，并且相比瞬時迎角的變化呈現(xiàn)明顯的遲滯效應(yīng)，這是非定常運動的一個特點，圖5（a）、（b）分別顯示了0°以及8°攻角時的俯仰力矩系數(shù)隨瞬時攻角變化的遲滯環(huán)曲線，可以看出，由于近水面效應(yīng)的復(fù)雜性，非定常運動曲線的光滑程度不如近地以及無限空間自由流場狀態(tài)，同時，當(dāng)h/c=5時，各個攻角下的近地效應(yīng)俯仰力矩系數(shù)遲滯環(huán)曲線與無限空間的結(jié)果基本一致，但是水面影響的曲線依然相差很大，這更加清晰地說明了水面的存在影響范圍更廣。

圖6-9為初始攻角為10°時翼型上仰到瞬時攻角為10.79°的三種計算模型的物面壓力系數(shù)以及流場圖。由于有了地面和水面的影響，翼型的上下翼面的壓力系數(shù)均發(fā)生了變化，而相比之下，地面和水面影響下的上翼面壓力系數(shù)基本一致，但是顯然近水面時的下翼面的壓力系數(shù)更高，從而導(dǎo)致整體的瞬時升力較近地面效應(yīng)大。

圖5 俯仰力矩系數(shù)遲滯環(huán)Fig.5 Hysteresis loops of pitching moment coefficient

圖6 三種模型物面壓力系數(shù) Fig.6 Pressure coefficients of three models

圖7 自由空間壓力云圖Fig.7 Pressure contour of free flow field

圖8 地面影響的壓力云圖Fig.8 Pressure contour of ground effect

圖9 水面影響的壓力云圖Fig.9 Pressure contour of water surface effect

應(yīng)用本文的動導(dǎo)數(shù)計算方法可以計算得到這三種情形下各個攻角的縱向俯仰力矩系數(shù)組合動導(dǎo)數(shù)，結(jié)果見圖10，從圖中可以看出近地/水高度相同的情況下，中小攻角時，受地面/水面的影響，動導(dǎo)數(shù)絕對值減小，動態(tài)穩(wěn)定特性損失，且水面效應(yīng)影響下的動導(dǎo)數(shù)絕對值減小得更為厲害，而在較大攻角時，這兩種效應(yīng)均能夠增加動態(tài)阻尼，導(dǎo)致動導(dǎo)數(shù)絕對值的增大。究其原因，中小攻角時，受到翼型自身形狀的影響，下翼面與地面/水面之間形成的通道表現(xiàn)出文丘里管效應(yīng)，這種效應(yīng)對通過管道的氣流進(jìn)行了加速，疊加上自身的非定常運動，致使翼型重心前后段的壓力分布變化不再穩(wěn)定，損失了整體的動態(tài)阻尼，地面由于呈現(xiàn)剛性，形成的管道效應(yīng)一側(cè)為平面，水面由于柔性作用易發(fā)生變形，形成的管道更為復(fù)雜，因而動導(dǎo)數(shù)絕對值也較地面影響時的小。大攻角時，這一情況剛好相反，下翼面與地面/水面形成的管道則以阻塞效應(yīng)為主，非定常運動因為氣流阻塞的影響而趨于減弱，也即系統(tǒng)的阻尼增大，因而動導(dǎo)數(shù)絕對值增大。

對比圖10中的不同近地/水高度下的動導(dǎo)數(shù)變化，可知，隨著h/c的增大，地面以及水面的動導(dǎo)數(shù)曲線都有著向自由空間結(jié)果靠攏的趨勢，這表明兩種效應(yīng)的影響在減弱，而地面效應(yīng)衰減得更快，在h/c=5時，其動導(dǎo)數(shù)隨攻角變化曲線已經(jīng)很接近自由空間結(jié)果，此時水面的影響依然很強，這也與之前的遲滯環(huán)變化一致。由此說明，在進(jìn)行現(xiàn)代地效飛行器動態(tài)氣動設(shè)計中，關(guān)于地面、水面以及自由空間的氣動特性需要詳細(xì)研究，區(qū)別對待。

圖10 縱向組合動導(dǎo)數(shù)Fig.10 Longitudinal combined dynamic derivatives

6 結(jié) 論

根據(jù)對自由空間流場、地面以及水面影響下的流場的動態(tài)特性計算結(jié)果進(jìn)行研究和分析，可以得出以下結(jié)論：

（1）本文建立的非定常方法能夠較為準(zhǔn)確地計算翼型運動時的非定常氣動力及其遲滯效應(yīng)，基于小幅度諧和振動的動導(dǎo)數(shù)辨識方法也能夠精確地計算組合動導(dǎo)數(shù)。

（2）地面與水面的存在對于翼型的動態(tài)流場有著很大的影響，進(jìn)而影響了翼型的動態(tài)阻尼，導(dǎo)致動導(dǎo)數(shù)絕對值的差異。但是地面的影響與水面不同。

（3）中小攻角下，地面與水面的影響均以管道效應(yīng)為主，柔性的水面形成的效應(yīng)使得翼型損失了更多的動態(tài)阻尼，動導(dǎo)數(shù)絕對值相比地面效應(yīng)要小。

（4）大攻角時，影響動態(tài)特性的為阻塞效應(yīng)，柔性水面的局部變形形成更為顯著的恢復(fù)力矩以減弱非定常運動，因而動導(dǎo)數(shù)絕對值較地面的影響大。

（5）由于地面、水面影響的差異，在對近地、水面飛行器的動態(tài)氣動特性進(jìn)行研究時，要分別對待，細(xì)致分析不同條件下的流場特性。

[1]楊 美,楊 韡,楊志剛.地效翼地面粘性效應(yīng)風(fēng)洞試驗研究[J].空氣動力學(xué)學(xué)報,2015,33(1):82-86.Yang Mei,Yang Wei,Yang Zhigang.Wind tunnel test of ground viscous effect on wing aerodynamics[J].Acta Aerodynamica Sinica,2015,33(1):82-86.

[2]韓 龍,陳陽陵,王福新,等.高海況機翼波浪地面效應(yīng)數(shù)值模擬與分析[J].上海交通大學(xué)學(xué)報,2014,48(8):1127-1133.Han Long,Chen Yangling,Wang Fuxin,et al.Numerical simulation and aerodynamic analysis of airfoil flying over wave ground at high grade sea condition[J].Journal of shanghai Jiaotong University,2014,48(8):1127-1133.

[3]盛其虎,吳德銘,張 亮.信天翁近海綿飛行時的氣動力研究[J].應(yīng)用數(shù)學(xué)和力學(xué),2005,26(9):1114-1120.Sheng Qihu,Wu Deming,Zhang Liang.Aerodynamic forces acting on an albatross flying above sea-waves[J].Applied Mathematics and Mechanics,2005,26(9):1114-1120.

[4]鄒光遠(yuǎn),譚 華,張立柱.解極近水面飛行時全機氣動力中水波效應(yīng)的邊界元法[J].水動力學(xué)研究與進(jìn)展A輯,2003,18(4):379-396.Zou Guangyuan,Tan Hua,Zhang Lizhu.The boundary element method(BEM)solving the water wave effect on the aerodynamic of a vehicle flying close to water surface[J].Journal of Hydrodynamics Ser.A,2003,18(4):379-396.

[5]秦旭國,劉沛清,屈秋林,等.三維多段機翼地面效應(yīng)數(shù)值模擬[J].航空學(xué)報,2011,32(2):257-264.Qin Xuguo,Liu Peiqing,Qu Qiulin,et al.Numerical simulation on 3D multi-element wings in ground effect[J].Acta Aeronautica et Astronautica Sinica,2011,32(2):257-264.

[6]秦旭國,劉沛清,屈秋林.翼型波浪水面巡航地面效應(yīng)數(shù)值模擬[J].北京航空航天大學(xué)學(xué)報,2011,37(3):295-304.Qin Xuguo,Liu Peiqing,Qu Qiulin.Numerical simulation on aerodynamics of airfoil flying over wavy water surface[J].Journal of Beijing University of Aeronautics and Astronautics,2011,37(3):295-304.

[7]Chen Y S,Schweikhard W G.Dynamic ground effect on a two-dimensional flat plate[J].Journal of Aircraft,1985,22(7):638-640.

[8]Konstadinopoulos P,Thrasher D F,Mook D T,et al.A vortex-lattice method for general,unsteady aerodynamics[J].Journal of Aircraft,1985,22(1):43-49.

[9]Kornev N,Matveev K.Complex numerical modeling of dynamics crashes of wing-in-ground vehicles[J].AIAA 2003-600,2003.

[10]葉正寅,楊永年,陳迎春.有地面干擾的動導(dǎo)數(shù)計算方法[J].飛行力學(xué),1998,16(3):68-72.Ye Zhengyin,Yang Yongnian,Chen Yingchun.Numerical calculation of dynamic stability derivatives in ground effect[J].Flight Dynamics,1998,16(3):68-72.

[11]Ronch A D,Vallespin D.Computation of dynamic derivatives using CFD[R].AIAA 2010-4817,2010.

[12]Jean-Francois L,Stephane M,Dominique F.Static and dynamic derivatives on generic UCAV without and with leading edge control[R].AIAA 2014-2391,2014.

Numerical simulation of aircraft dynamic stability characteristics flying over ground and water surface

MI Bai-gang,ZHAN Hao
(School of Aeronautics,Northwestern Polytechnical University,Xi’an 710072,China)

Based on unstructured dynamic mesh technique,SA turbulence model and volume of fluid(VOF)method,the 2-dimensional airfoil(NACA0012)unsteady aerodynamics flying over free flow field,ground and water surface were numerically studied,also the method to identify dynamic stability derivative was developed.It is shown that the method to calculate dynamic derivative is reliable.The absolute value of combined longitudinal dynamic derivative with water effect is the least at small angle of attack(AoA)as dynamic damping mostly decreases due to pipeline effect,which also happens on ground effect model,however,as the blockage effect at large AOA,both the dynamic stability of airfoil flying over ground and water surface are improved,besides,compared to ground,water affects a wider range to unsteady aerodynamics of airfoil.

free flow field;ground effect;water surface effect;dynamic derivative;numerical simulation;computational fluid dynamics(CFD)

V211.3

A

10.3969/j.issn.1007-7294.2017.11.004

1007-7294（2017）11-1348-08

2017-07-26

航空科學(xué)基金資助項目（2015ZA53013）；中央高?；究蒲袠I(yè)務(wù)費專項資金資助項目（3102016HKJJ004

米百剛（1989-），男，博士研究生，E-mail：mibaigang@163.com；

詹 浩（1972-），男，教授，博士生導(dǎo)師。