陳宣合,張向強,崔燕香

(1.中國科學院空天信息創(chuàng)新研究院,北京 100094;2.中國科學院大學,北京 100049)

0 引言

對于浮空器的氣動穩(wěn)定性,無動力的系留氣球是靜穩(wěn)定的,中低空動力飛艇是靜不穩(wěn)定的。臨近空間飛艇并非中低空飛艇的比例放大版,其顯著特點在于偏航控制機構不同。中低空飛艇的偏航控制通過舵面提供控制力矩,該力矩與空速的平方成正比;而臨近空間的大氣密度僅為中低空的1/10,低速時舵面效率更低;此外,低空飛艇舵面機構的重量、低溫低氣壓特性均不適合臨近空間飛艇。因此,臨近空間飛艇的偏航控制采用螺旋槳提供力矩。由于電機功率的限制,偏航或差動螺旋槳的最大控制力矩基本固定,導致臨近空間飛艇空速超出設計值后,面臨控制力矩不足的局面,此時,若臨近空間飛艇動穩(wěn)定裕度較低,可能出現偏航通道發(fā)散,國內近兩年的平流層飛行試驗中已發(fā)生多次。高抗風需求下的橫側向穩(wěn)定性問題需加強研究。

近年來,國外臨近空間飛艇飛行次數有所減少,臨近空間飛艇操控與氣動穩(wěn)定性方面公開發(fā)表的研究成果很少。本文針對兩種已有的飛艇構型,通過CFD仿真得出靜導數與動導數,進行穩(wěn)定性分析,為分析偏航通道發(fā)散問題的機理提供依據,同時為平流層飛艇氣動設計提供參考。

1 模型與計算處理方法

1.1 飛艇模型

本文利用已有的HALE-D和HiSentinel兩種飛艇模型進行數值仿真,模型示例如圖1所示,模型數據參看參考文獻[1]。使用ICEM CFD生成網格,附面層層數為30層,模型表面Y plus<1,空間六面體網格總數約800萬。

圖1 模型示例

坐標系的定義不同于飛機等常規(guī)飛行器,艇身機體坐標系原點為飛艇體心,方向定義:X向后沿機身縱軸;Z向左垂直于模型對稱面;Y向上垂直于XZ平面。計算輸出的氣流坐標系原點為飛艇重心,方向定義:X向后沿遠場氣流方向;Y向上垂直X軸并位于飛艇對稱面內;Z向左垂直于XY平面。

飛行速度相對于機體坐標軸系的角度定義如下:側滑角為飛行速度與飛艇對稱面的夾角,當飛行速度在橫軸分量為正時側滑角為正;迎角為飛行速度在飛艇對稱面的投影與縱軸的夾角,當飛行速度沿豎軸分量為正時迎角為正。

側向力、所有力矩結果均為體軸系,升阻力、升阻比均為氣流坐標系,動導數結果均為體軸系。

機體坐標系向氣流坐標系轉換的公式為:

Oxayaza=By·Bz·Oxbybzb

(1)

其中:Oxayaza為氣流坐標系坐標;Oxbybzb為機體坐標系坐標;By為側滑角轉換矩陣;Bz為迎角轉換矩陣;α為迎角;β為側滑角。

1.2 穩(wěn)定性參數識別方法

本文通過Fluent軟件自編UDF程序,采用模型強迫振蕩方法來獲得動導數[2-4]。定常流動的NS方程采用SIMPLEC算法求解,非定常流動的NS方程采用Coupled算法求解,計算模型采用SST-Kω模型。

(2)

其中:j為l或m,分別指升力和力矩系數;V為來流速度。俯仰受迫振動有以下運動學關系:

(3)

其中:αA為受迫周期性振蕩振幅;w為運動頻率。整理可得:

(4)

(5)

則

(6)

(7)

公式(6)為靜導數,公式(7)為組合導數,其它振動方向處理方式與俯仰類似。

1.3 動網格方法

模型強迫振蕩方法需要結合動網格方法實現,在CFD仿真軟件中,實現動網格的常用方法有:剛性網格[5]、滑移網格[6]和重構網格[7]。三種動網格生成方式均會在計算域存在網格重疊,插值誤差會隨著計算時間推移而逐漸累積[8]。結合模型結構網格劃分特點,參考飛機動導數及穩(wěn)定性相關文獻,選用剛性網格結果,誤差相比另外兩種小[9]。

2 數值計算及分析

2.1 算例驗證

為保證數值模擬結果可靠,首先,選用Finner[10-11]導彈模型驗證穩(wěn)定性參數求解過程與精度。Finner導彈模型如圖2所示。

圖2 Finner導彈模型

計算條件與參考文獻[12]一致,馬赫數Ma1.58,初始迎角α0=0°,減縮頻率k=0.015 822 6,振幅αm=1°,振蕩規(guī)律可寫為:

α=α0+αmsin(ωt)=1°sin(17t)

(8)

驗證過程從初始迎角α0出發(fā)做定常計算,待殘差收斂及力矩系數穩(wěn)定在某值附近時,開始做振蕩動態(tài)計算。迎角及俯仰力矩系數隨時間變化圖如圖3所示,圖中顯示出了氣動力矩的遲滯效應。

圖3 迎角及俯仰力矩系數隨時間變化圖

根據積分法得到俯仰組合動導數,對照參考文獻[12]的風洞試驗結果,動導數驗證結果如表1所示。

表1 動導數驗證結果

結合參考文獻[12-13],本文驗證相對誤差在3%以下,說明數值模擬結果比較準確,方法可行。

2.2 縱向靜態(tài)氣動力

圖4至圖8為迎角在-24°至24°范圍內,風軸升力系數Cl、風軸阻力系數Cd、升阻比K、體軸俯仰力矩系數Cmz和體軸法向力系數Cy的靜態(tài)參數變化圖。如圖4所示,HiSentinel(以下簡稱H1)的升力系數變化圖形似三次函數,在0°迎角附近存在一段緩慢變化區(qū)域,且迎角繼續(xù)增大時,升力系數也在均勻較快增長;HALE-D(以下簡稱H2)的升力系數變化圖在該范圍內近似直線,升力系數增長變化明顯。

圖4 風軸升力系數Cl靜態(tài)參數變化圖

在圖5所示的阻力系數變化圖中,H1與H2變化圖近似二次函數,隨著迎角的增大,H2的阻力系數增長相較于H1更快。

圖5 風軸阻力系數Cd靜態(tài)參數變化圖

圖6所示中,H1與H2升阻比變化圖形態(tài)基本相似,但H2的最大升阻比所對應的迎角絕對值較H1略小。

圖6 升阻比K靜態(tài)參數變化圖

圖7所示中兩種模型俯仰力矩變化圖形似中心對稱正弦變化圖,在±24°內都不穩(wěn)定,拐點約為-16°和16°。相較于H1模型而言,H2模型的力矩系數變化范圍較小。

圖7 體軸俯仰力矩系數Cmz靜態(tài)參數變化圖

由圖8可以看出,法向力系數變化圖趨勢與升力系數變化圖有較大相似性,相比于升力系數變化圖,法向力系數隨迎角變化更加明顯。

圖8 體軸法向力系數Cy靜態(tài)參數變化圖

2.3 橫向靜態(tài)氣動力

圖9至圖11為側滑角在0°至24°范圍內,風軸阻力系數Cd、體軸側向力系數Cz、體軸偏航力矩系數Cmy的靜態(tài)參數變化圖。

圖9 風軸阻力系數Cd靜態(tài)參數變化圖

圖10 體軸側向力系數Cz靜態(tài)參數變化圖

圖11 體軸偏航力矩系數Cmy靜態(tài)參數變化圖

結合縱向與橫向的靜態(tài)穩(wěn)定性參數變化圖來看,H2模型相較于H1模型,隨相應方向角度的增大,力的系數增長變化較為顯著,而在力矩系數上,兩種模型的變化趨勢大致相同,但H1模型的變化幅度更大。

2.4 俯仰動穩(wěn)定性參數

在得到上述靜態(tài)穩(wěn)定性參數的基礎上,強迫飛艇振動,計算動態(tài)穩(wěn)定性參數。兩種模型的強迫振動頻率為0.1 Hz,振幅為3°,縱向俯仰力矩系數在0°迎角遲滯環(huán)曲線如圖12所示。可以觀測到,遲滯曲線以相應的靜態(tài)穩(wěn)定性參數值為中心,呈橢圓形偏移。

圖12 縱向俯仰力矩系數在0°迎角遲滯環(huán)曲線

由遲滯環(huán)提取動導數,縱向振蕩俯仰力矩組合導數變化圖及縱向振蕩升力組合導數變化圖如圖13和圖14所示。整體上,模型俯仰振蕩俯仰力矩組合導數在整個迎角范圍內為負值,模型為阻尼狀態(tài)。模型俯仰振蕩升力組合導數在整個迎角范圍內為負值。

圖13 縱向振蕩俯仰力矩組合導數變化圖

圖14 縱向振蕩升力組合導數變化圖

2.5 偏航動穩(wěn)定性參數[14]

偏航通道頻率與振幅和俯仰的相同,偏航力矩系數在0°迎角的遲滯環(huán)曲線如圖15所示。橫向振蕩偏航力矩組合導數變化圖和橫向振蕩偏航側向力組合導數變化圖分別如圖16、17所示。在縱向與橫向上,俯仰力矩組合導數及偏航力矩組合導數在對應角度范圍內均小于0,飛艇在對應中心角度范圍附近為振幅收斂的振蕩,是動穩(wěn)定狀態(tài)。

圖15 偏航力矩系數在0°迎角的遲滯環(huán)曲線

圖16 橫向振蕩偏航力矩組合導數變化圖

圖17 橫向振蕩偏航側向力組合導數變化圖

3 飛艇橫向穩(wěn)定性分析

3.1 飛艇橫向穩(wěn)定性與控制理論[1]

飛艇的橫向穩(wěn)定性是由飛艇的橫向運動和定向運動之間的相互作用決定的:橫向運動是由飛艇的側滑角產生的側向力造成的,當飛艇橫向平移時,側滑角在平移的相反方向上產生一個側力,側滑角的小擾動會隨著時間的推移而減少;定向運動是由圍繞偏航軸旋轉產生的偏航力矩引起的。通常情況下,飛艇理論上是不穩(wěn)定的,因為空氣動力中心在參考點的前面,側滑角的小擾動所產生的偏航力矩將傾向于增加該側滑角。側向和定向運動都受到飛艇側滑角和偏航率的影響,為了確定側向系統(tǒng)的整體穩(wěn)定性,必須同時考慮這兩個特性。

(9)

其中,V∞為遠場速度。

側滑角速率簡化為體軸側向加速度除以真實空速。側向加速度等于側向力除以質量加上由于體軸相對于慣性軸的旋轉而產生的科里奧利力。側向力用靜態(tài)側向力系數和動態(tài)側向力偏航速率導數表示。

(10)

(11)

3.2 飛艇橫向穩(wěn)定性結果分析

穩(wěn)定性分析變化圖如圖18所示,Ⅰ區(qū)與Ⅱ區(qū)、Ⅱ區(qū)與Ⅲ區(qū)交界的兩條直線分別劃出H1和H2模型,理論上實現穩(wěn)定轉彎的穩(wěn)定區(qū)域和不穩(wěn)定區(qū)域。其中,Ⅰ區(qū)為H1的非穩(wěn)定轉彎區(qū)域,I區(qū)與Ⅱ區(qū)為H2模型轉彎的不穩(wěn)定區(qū)域。由于穩(wěn)定性和操縱性的孿生矛盾,H1模型即便存在不穩(wěn)定區(qū)域,相應地也具備更好的操縱性;而H2模型即便距分界線較遠,穩(wěn)定裕度較大,相應地飛行操縱性則大大降低。橫向偏航穩(wěn)定性需結合日后實際飛行試驗進一步加深研究。

圖18 穩(wěn)定性分析變化圖

4 結論

本文結合實際飛行需要,通過數值模擬方法得出兩種飛艇模型的靜態(tài)參數、動穩(wěn)定參數,并結合所得數據對飛艇偏航的穩(wěn)定性進行分析,得出如下結論:

(1)靜態(tài)參數上,H1模型與H2模型的變化趨勢基本相同。在氣動力系數上,H2模型隨迎角或側滑角的變化更加顯著;在氣動力矩系數上,H1和H2在計算角度內均靜不穩(wěn)定,H1模型變化幅值較H2模型的更大。此外,相同迎角(側滑角)時,H2的阻力較大。

(2)動穩(wěn)定性參數上,H1模型與H2模型變化范圍均較小,俯仰振蕩、偏航振蕩組合導數在整個角度變化范圍內均為負值,模型處于阻尼狀態(tài)。

(3)橫向穩(wěn)定性分析表明,H2模型穩(wěn)定區(qū)間大于H1模型。

按照對常規(guī)飛行器運動穩(wěn)定性的理解,氣動靜不穩(wěn)定即意味著運動不穩(wěn)定,只有依靠控制系統(tǒng)才可保持運動穩(wěn)定[15-16],但從已有的飛艇氣動資料發(fā)現,許多飛艇是靜不穩(wěn)定的[17]。相關文獻[18]指出,由于飛艇獨特的氣動外形,動導數對飛艇的氣動影響是常規(guī)飛行器的數百倍,對飛艇穩(wěn)定性的影響極大。國內平流層飛艇項目整體尚在小規(guī)模實驗階段,技術方法有待完善。HALE-D和HiSentinel系列是國外飛艇經過飛行實驗的典型構型,也是倒“Y”型尾翼和“X”型尾翼的典型代表。對于飛艇橫向穩(wěn)定性的研究,是以尾翼形狀為主體的,艇體外形的影響微乎其微。

由于穩(wěn)定性與操控性的孿生矛盾,飛艇的設計需根據實際需求,綜合考慮阻力、穩(wěn)定裕度和操縱性能之間的平衡。本文所研究的兩種不同尾翼形狀,可以根據不同的任務需求,應用于飛艇控制系統(tǒng)設計中。