典型超聲速/高超聲速動導(dǎo)數(shù)計算方法研究

秦之軒, 史愛明, 安富強, 裴揚

(西北工業(yè)大學(xué) 翼型葉柵空氣動力學(xué)國家重點實驗室, 陜西 西安 710072)

典型超聲速/高超聲速動導(dǎo)數(shù)計算方法研究

秦之軒, 史愛明, 安富強, 裴揚

(西北工業(yè)大學(xué) 翼型葉柵空氣動力學(xué)國家重點實驗室, 陜西 西安 710072)

摘要:以當(dāng)?shù)亓骰钊碚摓榛A(chǔ)完成超聲速/高超聲速動導(dǎo)數(shù)計算公式推導(dǎo),發(fā)展了一種高效的動導(dǎo)數(shù)計算方法。選取1個超聲速、2個高超聲速動導(dǎo)數(shù)算例標(biāo)模,計算了超聲速/高超聲速下飛行器動導(dǎo)數(shù)變化規(guī)律。研究結(jié)果表明:BFM超聲速流動算例的俯仰組合動導(dǎo)數(shù)絕對值隨馬赫數(shù)的增大而減小,計算重心位置后移亦會導(dǎo)致俯仰組合動導(dǎo)數(shù)絕對值減小,縱向動穩(wěn)定性降低;在0°～20°迎角范圍內(nèi),高超聲速流動中的尖錐和鈍頭旋成體俯仰組合動導(dǎo)數(shù)絕對值隨迎角的增大而增大,縱向動穩(wěn)定性增強;所提方法計算時間約為雙時間動導(dǎo)數(shù)計算方法的1/27。

關(guān)鍵詞:超聲速; 高超聲速; 當(dāng)?shù)亓骰钊碚? 動導(dǎo)數(shù)

0引言

飛行器動導(dǎo)數(shù)是飛機運動中單位速度或角速度變化引起的氣動力或力矩的變化,是飛行器研制工程中必不可少的動態(tài)氣動參數(shù),動導(dǎo)數(shù)的準(zhǔn)確評估對飛機設(shè)計和飛行都具有重要的意義。隨著現(xiàn)代先進戰(zhàn)斗機和導(dǎo)彈飛行包線的擴展,人們對復(fù)雜流場中的飛行器動穩(wěn)定性研究更加重視,對飛行器各種飛行條件下的動導(dǎo)數(shù)都需要評估,工作量大;因此精確、高效地評估飛行器動導(dǎo)數(shù)變得尤為重要。

目前獲取動導(dǎo)數(shù)的主要方法有風(fēng)洞試驗和理論計算。風(fēng)洞試驗是動導(dǎo)數(shù)特性判定的主要依據(jù),同時可以為計算提供基礎(chǔ)和參考;但風(fēng)洞試驗費用較高,且受到風(fēng)洞規(guī)模的限制,試驗雷諾數(shù)較低,無法對設(shè)計者關(guān)心的狀態(tài)進行試驗,需要與理論計算相互補充、互相驗證。隨著CFD的發(fā)展,采用非定常CFD方法計算動導(dǎo)數(shù)已經(jīng)成為目前飛行器動態(tài)特性研究的一種方法[1],但是其計算量大,效率較低。基于CFD的當(dāng)?shù)亓骰钊碚摲椒╗2-3]是一種將高精度的CFD方法和當(dāng)?shù)亓骰钊碚撓嘟Y(jié)合的非定常CFD方法。該方法兼顧了CFD方法流場計算的高精度特點,同時相比于完全的非定常CFD計算耗費的計算時間明顯減少,是一種效率和精度都相對較高的方法。

本文采用基于CFD的當(dāng)?shù)亓骰钊碚摲椒ㄟM行超聲速/高超聲速飛行器動導(dǎo)數(shù)計算。通過對3種典型的超聲速/高超聲速動導(dǎo)數(shù)模型的計算,研究其動導(dǎo)數(shù)的變化規(guī)律。動導(dǎo)數(shù)計算模型中,Basic Finner Missile(BFM)模型是超聲速導(dǎo)彈的基本外形,大多數(shù)空空導(dǎo)彈屬于此類;高超聲速尖錐和Hyper Ballistic Shape(HBS)模型為高超聲速外形,火箭、再入式彈道導(dǎo)彈、飛船返回艙等屬于這類外形。

1非定常氣動力計算方法

根據(jù)當(dāng)?shù)亓骰钊碚?物面壓力可表示為[4-5]:

(1)

式中:Pl為定常壓力;W為物面變形法向速度大小和振動法向速度的模之和;n為物面變形后的外法線單位矢量;Vb為物面振動速度;n0為物面變形前的外法線單位矢量。

設(shè)飛行器俯仰振動的模態(tài)為zx=ai+bj+ck,對應(yīng)的廣義坐標(biāo)為ξ(t),當(dāng)ξ(t)=1時飛行器發(fā)生單位變形后的物面外法線單位矢量為n0,則有:

從而有:

則物面的壓力表達式為:

P-Pl=

定常壓力項Pl并不反映計算對象的運動特性,故將其略去,并在整個飛行器表面做積分,得到非定常廣義氣動力計算公式:

(2)

其中:

A={}dS

B=

Dp=Pl/P∞,Dρ=ρl/ρ∞

當(dāng)ξ(t)為迎角時,對應(yīng)的非定常廣義氣動力Q為俯仰力矩。

2動導(dǎo)數(shù)計算方法

對于做如下簡諧振動的飛行器[6]:

(3)

式中:αm為振幅;ω為強迫振動角頻率。

將式(3)代入式(2),則可得俯仰力矩系數(shù)為:

其中:

式中:S和b分別為參考面積和參考長度。

對俯仰力矩在一個周期的內(nèi)變化進行積分,得到俯仰力矩系數(shù)動導(dǎo)數(shù)的計算公式:

(4)

式中:k=ωb/(2V∞)為減縮頻率。從式(4)可以看出,最終推導(dǎo)結(jié)果中減縮頻率k不再出現(xiàn),因此該方法不計入減縮頻率對動導(dǎo)數(shù)的影響。

3BFM動導(dǎo)數(shù)計算分析

“十”字翼導(dǎo)彈BFM為一個尖錐形頭部、圓柱形彈身并帶有4個矩形尾翼的外形,尾翼為十字布局,其彈身截面直徑為d,彈體總長為10d,尾翼弦長為0.94d,尾翼厚度為0.08d,尾翼半翼展為d。計算網(wǎng)格采用非結(jié)構(gòu)網(wǎng)格,網(wǎng)格單元總數(shù)1 690 219,導(dǎo)彈表面網(wǎng)格如圖1所示。

圖1　BFM外形表面網(wǎng)格Fig.1　Surface grid of BFM

取參考長度為彈身截面直徑d=1.0 m,參考面積為彈身截面面積,分別計算了0°迎角下重心在xcg=5.0d,6.1d時,俯仰力矩系數(shù)動導(dǎo)數(shù)隨馬赫數(shù)的變化,并與文獻[7]中的試驗值和文獻[3]參考值作對比,對比分析結(jié)果如圖2所示。

圖2　BFM俯仰力矩系數(shù)動導(dǎo)數(shù)Fig.2　Pitching dynamic derivative of BFM

由圖2可以看出,本文動導(dǎo)數(shù)計算結(jié)果與試驗值比較接近,變化趨勢吻合良好,計算得到的動導(dǎo)數(shù)數(shù)值略高于試驗值。從圖中曲線的變化趨勢可以看出,隨著馬赫數(shù)的增加,俯仰力矩系數(shù)動導(dǎo)數(shù)絕對值逐漸減小,縱向氣動阻尼穩(wěn)定性降低;而隨著俯仰轉(zhuǎn)軸位置后移,動導(dǎo)數(shù)絕對值也在減小,即俯仰轉(zhuǎn)軸位置后移也使縱向氣動阻尼穩(wěn)定性降低。

4尖錐外形動導(dǎo)數(shù)計算分析

高超聲速尖錐外形是國際上計算動導(dǎo)數(shù)的標(biāo)準(zhǔn)模型,該模型是半錐角為10°的尖錐,錐體總長L=141 mm,錐底面直徑為50 mm,計算網(wǎng)格為非結(jié)構(gòu)網(wǎng)格,網(wǎng)格單元數(shù)481 630,如圖3所示。

圖3　尖錐表面網(wǎng)格Fig.3　Surface grid of cone

取錐體參考長度為總長L=141 mm,參考面積為錐體底部面積S=1 963.5 mm2,重心位置為0.7L,計算了馬赫數(shù)為6.85,迎角為0～18°時各個狀態(tài)下的俯仰力矩動導(dǎo)數(shù),并與文獻[8]中的試驗值、牛頓內(nèi)伏理論值及文獻[3]的參考值進行比較,結(jié)果如圖4所示。可以看出,在迎角小于10°時,本文計算結(jié)果與試驗值有較大差距,且誤差大于其他兩個參考計算結(jié)果;在迎角大于10°時,本文計算結(jié)果較好,優(yōu)于另外兩個計算結(jié)果。計算結(jié)果的總趨勢與參考計算結(jié)果一致,與試驗結(jié)果的對比中,最大誤差不超過50%。

圖4　尖錐俯仰力矩系數(shù)動導(dǎo)數(shù)Fig.4　Pitching dynamic derivative of cone

5HBS動導(dǎo)數(shù)計算分析

彈道外形(Hyperballistic Shape,HBS)為一個半球鈍頭、帶有兩段擴張裙部的外形,彈道前部半球直徑為d,第一個裙部擴張角為5°,第二個裙部擴張角為15°,兩個裙部長度均為1.5d,彈道外形總長4.5d。計算網(wǎng)格采用非結(jié)構(gòu)網(wǎng)格,網(wǎng)格單元數(shù)609 052,物面網(wǎng)格如圖5所示。

圖5　HBS表面網(wǎng)格Fig.5　Surface grid of HBS

針對該外形,取直徑d=1.0 m,計算了Ma=6.85,初始迎角α分別為0°,5°,10°,15°,20°時的俯仰力矩動導(dǎo)數(shù),并與文獻[8]中的試驗值、牛頓內(nèi)伏理論值及文獻[6]的參考值進行比較,結(jié)果如圖6所示。

圖6　HBS俯仰力矩系數(shù)動導(dǎo)數(shù)Fig.6　Pitching dynamic derivative of HBS

可以看出,本文計算結(jié)果與試驗值及其他兩個參考計算結(jié)果都較吻合,在迎角小于15°時,動導(dǎo)數(shù)變化幅度較小,絕對值隨迎角增大而緩慢增大;15°以后，動導(dǎo)數(shù)絕對值隨迎角增大而增大的速度變快,本文計算結(jié)果略高于試驗值。

6計算效率對比

以HBS模型算例為例,在CPU為Inter core i7的計算機上計算,單核運行,以迭代時間步1 000步為限。使用本文的方法計算一個狀態(tài)下的動導(dǎo)數(shù)用時約40 min,使用非定常的雙時間方法計算(每個周期分為16個時間步,共計算兩個周期,每個時間步進行300次子迭代),在相同計算條件下計算一個狀態(tài)下俯仰阻尼動導(dǎo)數(shù)約需要6 h,計算俯仰、偏航、滾轉(zhuǎn)阻尼動導(dǎo)數(shù)則需要18 h。本文方法計算用時僅為雙時間方法的1/27,由此可見該方法能顯著提高計算效率。

7結(jié)束語

本文詳細(xì)推導(dǎo)了基于CFD的當(dāng)?shù)亓骰钊碚摰膭訉?dǎo)數(shù)計算公式,并通過1個超聲速、2個高超聲速算例,驗證了本方法在計算高超聲速動導(dǎo)數(shù)時的精度和效率。從最終結(jié)果來看,該方法能夠較好地揭示超聲速/高超聲速流動下動導(dǎo)數(shù)的變化規(guī)律,相比于完全非定常的CFD時域推進方法可以節(jié)約大量計算時間,適用于工程計算。

參考文獻：

[1]史愛明,楊永年,葉正寅.結(jié)合CFD技術(shù)的跨音速動導(dǎo)數(shù)計算方法研究[J].西北工業(yè)大學(xué)學(xué)報,2008,26(1):11-14.

[2]盧學(xué)成,葉正寅,張偉偉.超聲速/高超聲速飛行器動導(dǎo)數(shù)的高效計算方法[J].航空計算技術(shù),2008,38(3):28-31.

[3]劉溢浪,張偉偉,田八林,等.一種超聲速高超聲速動導(dǎo)數(shù)的高效計算方法[J].西北工業(yè)大學(xué)學(xué)報,2013,31(5):825-828.

[4]韓漢橋,張陳安,王發(fā)民.一種高空高超聲速非定常氣動力近似模型[J].力學(xué)學(xué)報,2013, 45(5):690-698.

[5]張陳安,張偉偉,葉正寅.基于當(dāng)?shù)亓骰钊碚摰臍鈩訌椥苑€(wěn)定性分析方法研究[J].工程力學(xué),2007,24(2):22-27.

[6]范晶晶,閻超,李躍軍.飛行器大迎角下俯仰靜、動導(dǎo)數(shù)的數(shù)值計算[J].航空學(xué)報,2009,30(10):1846-1850.

[7]Tong B G,Hui W H.Unsteady embedded Newton-Busemann flow theory[J].Journal of Spacecraft and Rockets,1986,23(3):129-135.

[8]East R A,Hutt G R.Comparison of predictions and experimental data for hypersonic pitching motion stability[J].Journal of Spacecraft and Rockets,1988,25(3):225-233.

(編輯：李怡)

Computing dynamic derivatives for supersonic and hypersonic models based on local piston theory

QIN Zhi-xuan, SHI Ai-ming, AN Fu-qiang, PEI Yang

(National Key Laboratory of Aerodynamic Design and Research, NWPU, Xi’an 710072, China)

Abstract:Based on the local piston theory, the formulas derivation of dynamic derivatives for supersonic and hypersonic models was completed and an effective dynamic derivatives method was developed. A supersonic model and two hypersonic models were selected and their dynamic derivatives change laws were calculated. Research results show that, the absolute value of pitching dynamic derivatives of BFM supersonic flow modulus are decreased with the increasing of Mach number and the center of gravity moved backward, the longitudinal dynamic stability is reduced. The absolute value of pitching dynamic derivatives of the cone and HBS hypersonic flow modulus of pitching derivative are increased with the increasing of AOA when the AOA are located at 0°～ 20°, the longitudinal dynamic stability is enhanced. The method takes only 1/27 calculation time of dual time domain dynamic derivatives computing method.

Key words:supersonic; hypersonic; local piston theory; dynamic derivative

收稿日期:2015-08-18;

修訂日期:2016-01-25; 網(wǎng)絡(luò)出版時間:2016-02-29 16:37

基金項目:國家自然科學(xué)基金資助(10602046)；航天CAST-BISEE基金資助(2015MC0175)；西北工業(yè)大學(xué)畢業(yè)設(shè)計重點項目資助(W002223)

作者簡介:秦之軒(1988-)，男，內(nèi)蒙古五原人，博士研究生，研究方向為計算流體力學(xué)和飛行器設(shè)計； 史愛明(1977-)，男，江蘇金壇人，副教授，博士，研究方向為空氣動力學(xué)和流固耦合力學(xué)。

中圖分類號：V211.3

文獻標(biāo)識碼：A

文章編號：1002-0853(2016)03-0017-04