小型固定翼飛機(jī)飛行動(dòng)力學(xué)仿真

陳衛(wèi)東，寧雷，劉淼群，吳限德，汪超亮

(1.哈爾濱工程大學(xué) 航天與建筑工程學(xué)院，黑龍江 哈爾濱150001;2.中國(guó)科學(xué)院 定量遙感信息技術(shù)重點(diǎn)實(shí)驗(yàn)室，北京100094)

0 引言

計(jì)算機(jī)仿真是系統(tǒng)仿真技術(shù)、計(jì)算機(jī)技術(shù)及相關(guān)專(zhuān)業(yè)技術(shù)等多領(lǐng)域結(jié)合的產(chǎn)物，它隨著計(jì)算機(jī)技術(shù)的發(fā)展而發(fā)展起來(lái)，并在諸多領(lǐng)域得到廣泛的應(yīng)用。飛行動(dòng)力學(xué)仿真技術(shù)基于飛機(jī)飛行動(dòng)力學(xué)、空氣動(dòng)力學(xué)以及飛行控制原理等學(xué)科，是飛機(jī)飛行仿真的關(guān)鍵技術(shù)［1］。

早期，由于飛機(jī)設(shè)計(jì)與制造公司無(wú)法提供準(zhǔn)確、充足的飛行仿真所需的設(shè)計(jì)和試飛數(shù)據(jù)，飛機(jī)的飛行建模方法比較粗糙簡(jiǎn)單，模型的逼真度也較低。進(jìn)入20世紀(jì)90年代后，由于航空界對(duì)飛行模擬重視程度的提高，一些大型飛機(jī)制造公司(如波音公司和空中客車(chē)公司)開(kāi)始在飛機(jī)設(shè)計(jì)和試飛中同步生產(chǎn)出飛機(jī)飛行模擬機(jī)所需要的專(zhuān)用數(shù)據(jù)，使得飛行仿真建模越來(lái)越精準(zhǔn)［2］。隨著研究工作的不斷深入，以及對(duì)飛機(jī)模型準(zhǔn)確性要求的不斷提高，其系統(tǒng)龐大、運(yùn)動(dòng)特性復(fù)雜、氣動(dòng)數(shù)據(jù)繁多的特點(diǎn)使得準(zhǔn)確建立飛行動(dòng)力學(xué)模型的難度越來(lái)越高［3］，因此在飛機(jī)建模時(shí)，應(yīng)充分考慮飛機(jī)運(yùn)動(dòng)特性復(fù)雜和模型數(shù)據(jù)量龐大的特點(diǎn)，并對(duì)模型作適當(dāng)簡(jiǎn)化，以求能夠準(zhǔn)確反映飛機(jī)的實(shí)際運(yùn)動(dòng)特性。

本文以Y12固定翼輕型運(yùn)輸機(jī)為研究背景，對(duì)該飛機(jī)進(jìn)行動(dòng)力學(xué)仿真，驗(yàn)證飛機(jī)的飛行特性，同時(shí)，對(duì)其他以飛機(jī)為對(duì)象的研究工作也具有參考價(jià)值。該飛機(jī)一般為低中空、低速(馬赫數(shù)0.3以?xún)?nèi))飛行，并由升降舵、方向舵及副翼控制飛行姿態(tài)。在動(dòng)力學(xué)建模時(shí)，為使仿真模型準(zhǔn)確性更高，首先利用有限元分析軟件，獲得飛機(jī)的各個(gè)氣動(dòng)參數(shù)，并考慮風(fēng)對(duì)飛行的影響，然后使用Matlab語(yǔ)言編寫(xiě)仿真程序，采用標(biāo)準(zhǔn)四階Runge-Kutta法逐步求解飛機(jī)飛行動(dòng)力學(xué)方程組，利用飛機(jī)實(shí)際飛行數(shù)據(jù)與仿真結(jié)果進(jìn)行了對(duì)比。對(duì)比結(jié)果表明，該方法適用于飛機(jī)動(dòng)力學(xué)仿真，獲得的結(jié)果也相對(duì)準(zhǔn)確。

1 飛機(jī)模型建立

本文以Y12飛機(jī)為背景，建立表征飛機(jī)運(yùn)動(dòng)規(guī)律的數(shù)學(xué)模型。在建模時(shí)假設(shè)飛機(jī)的運(yùn)動(dòng)為六自由度的剛體運(yùn)動(dòng)，將地面坐標(biāo)系視為慣性參考系，并假設(shè)地面為平面;同時(shí)考慮非平靜大氣的作用。

將飛機(jī)質(zhì)心動(dòng)力學(xué)矢量方程投影到航跡坐標(biāo)系中，可得到如下形式:

將飛機(jī)繞質(zhì)心轉(zhuǎn)動(dòng)的動(dòng)力學(xué)矢量方程投影到機(jī)體坐標(biāo)系上。則方程有如下形式:

式中:Jx，Jy，Jz分別為飛機(jī)對(duì)機(jī)體坐標(biāo)系各軸的轉(zhuǎn)動(dòng)慣量;Jxy為飛機(jī)對(duì)機(jī)體坐標(biāo)系x軸和y軸的慣性積，忽略誤差影響，假設(shè)Oxy為飛機(jī)對(duì)稱(chēng)面，則有慣性積Jyz=Jzx=0。

根據(jù)速度與位置關(guān)系，在地面坐標(biāo)系中，飛機(jī)質(zhì)心的運(yùn)動(dòng)學(xué)方程為:

確定飛機(jī)在空中的姿態(tài)，需要建立描述飛機(jī)相對(duì)地面坐標(biāo)系姿態(tài)變化的運(yùn)動(dòng)學(xué)方程如下:

非平靜大氣中的仿真模型，需要考慮風(fēng)對(duì)飛機(jī)飛行的影響。飛機(jī)相對(duì)地面的速度V、相對(duì)氣流的速度VU與風(fēng)速VW有如下關(guān)系:

其中，速度VU可表示為:

在描述飛機(jī)運(yùn)動(dòng)的模型中，還需有描述飛機(jī)質(zhì)量變化的微分方程。其表達(dá)式為:

式中:Q為燃料消耗率，由油門(mén)開(kāi)度決定［5］。

下面建立飛機(jī)仿真的大氣模型。飛機(jī)工作于對(duì)流層，因此可將大氣視為一種因地球引力場(chǎng)作用而呈固定分布的、服從玻意耳定律關(guān)系式的氣體，并假設(shè)地球引力為常數(shù)。隨著海拔高度的升高，大氣密度從海平面處的1.225 kg/m3開(kāi)始遞減，利用理想氣體狀態(tài)方程，密度可表示為:

式中:大氣壓力p和溫度T的具體表達(dá)式見(jiàn)文獻(xiàn)［6］;M0為平均空氣分子量(28.964 4 kg/kmol);RS為大氣常數(shù)(8 314.32 J/(kmol˙K))。

聲速c依賴(lài)于大氣溫度，計(jì)算公式為:

式中:γ為理想雙原子氣體的比熱容(γ=1.4)。

以上便建立了飛機(jī)的飛行動(dòng)力學(xué)模型，結(jié)合相關(guān)參數(shù)后，便可以對(duì)其進(jìn)行求解。

2 作用在飛機(jī)上的氣動(dòng)力和力矩

作用在飛機(jī)上的空氣動(dòng)力沿氣流坐標(biāo)系分解為三個(gè)分量:阻力D、升力L和側(cè)向力C。實(shí)驗(yàn)分析表明:空氣動(dòng)力的大小與來(lái)流的動(dòng)壓q和飛機(jī)的特征面積S有如下關(guān)系:

式中:CD，CL和CC分別為阻力系數(shù)、升力系數(shù)和側(cè)向力系數(shù)。

作用在飛機(jī)上的氣動(dòng)力矩沿機(jī)體坐標(biāo)系分解為三個(gè)分量:滾轉(zhuǎn)力矩Mx、俯仰力矩My和偏航力矩Mz，可寫(xiě)成如下形式:

式中:Lc和Lz分別為側(cè)向和縱向特征長(zhǎng)度;Cl，Cn和Cm分別為滾轉(zhuǎn)力矩系數(shù)、偏航力矩系數(shù)和俯仰力矩系數(shù)［7］。

本文采用計(jì)算流體力學(xué)有限元分析軟件進(jìn)行固定翼飛機(jī)氣動(dòng)特性計(jì)算。在計(jì)算時(shí)假設(shè)機(jī)翼的安裝角為0°;考慮飛行的真實(shí)過(guò)程以及各種運(yùn)動(dòng)學(xué)參數(shù)對(duì)氣動(dòng)特性的影響，選取主要影響參數(shù)(如迎角、馬赫數(shù)、舵偏角等)對(duì)固定翼飛機(jī)氣動(dòng)特性進(jìn)行計(jì)算，得到諸如升力系數(shù)、阻力系數(shù)等氣動(dòng)參數(shù)，并將各氣動(dòng)參數(shù)表示成多個(gè)在不同影響因素(如側(cè)滑角β、升降舵偏角δe)下關(guān)于馬赫數(shù)Ma和迎角α的二維表格形式［8］。表1和表2分別為β與 δe均為0°時(shí)飛機(jī)的升力系數(shù)和阻力系數(shù)，其他工況由于篇幅所限未全列出。

由氣動(dòng)數(shù)據(jù)表可知，氣動(dòng)參數(shù)是多變量函數(shù)，形式比較復(fù)雜，并且在飛行仿真時(shí)需要根據(jù)采樣點(diǎn)處的函數(shù)值，計(jì)算非采樣點(diǎn)處的數(shù)據(jù)，因而需要進(jìn)行多維函數(shù)插值計(jì)算。

表1 β=0°，δe=0°時(shí)的升力系數(shù)Table 1 Lift coefficient with β =0°，δe=0°

表2 β=0°，δe=0°時(shí)的阻力系數(shù)Table 2 Drag coefficient with β =0°，δe=0°

3 仿真和模型驗(yàn)證

本文所仿真的飛機(jī)飛行軌跡描述如下:在高度3 430 m左右定高飛行，直線(xiàn)平飛一段距離后，完成一個(gè)右轉(zhuǎn)彎動(dòng)作，并沿著平行于第一段飛行直線(xiàn)的軌跡返回，再完成一個(gè)左轉(zhuǎn)彎，回到接近飛行初始點(diǎn)的位置(見(jiàn)圖1)，并重復(fù)上述動(dòng)作。

仿真與數(shù)據(jù)處理程序用Matlab語(yǔ)言編寫(xiě)，首先設(shè)定仿真初始參數(shù)，然后采用標(biāo)準(zhǔn)四階Runge-Kutta法對(duì)方程進(jìn)行求解。求解完成得到結(jié)果后，與飛機(jī)飛行過(guò)程中實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比，下面分別敘述各參數(shù)飛行試驗(yàn)實(shí)測(cè)數(shù)據(jù)與仿真曲線(xiàn)的對(duì)比結(jié)果。

圖1為飛機(jī)水平軌跡對(duì)比圖。由圖可見(jiàn)，仿真飛行軌跡與實(shí)測(cè)軌跡吻合度較高，但是，實(shí)際飛行過(guò)程中的飛行環(huán)境隨機(jī)因素和人為因素會(huì)帶來(lái)兩條曲線(xiàn)間的差異:飛行速度的變化使得實(shí)測(cè)曲線(xiàn)平均轉(zhuǎn)彎半徑相對(duì)較小，兩個(gè)轉(zhuǎn)彎動(dòng)作也不可能完全相同;仿真模型對(duì)兩個(gè)轉(zhuǎn)彎部分設(shè)定相同的參數(shù)，曲線(xiàn)呈對(duì)稱(chēng)形狀，而并不會(huì)考慮上述因素帶來(lái)的誤差。故可認(rèn)為仿真模型能夠完成飛機(jī)實(shí)際的飛行動(dòng)作。

圖1 水平軌跡對(duì)比Fig.1 Horizontal trajectory comparison

圖2 為飛行高度的變化曲線(xiàn)。受到飛行環(huán)境中隨機(jī)因素與人為因素的干擾，飛行試驗(yàn)的高度變化在一定范圍內(nèi)產(chǎn)生波動(dòng)，但這種波動(dòng)依然是在盡量保持高度3 430 m飛行的條件下產(chǎn)生的，仿真模型忽略了上述因素的干擾，對(duì)高度的保持控制較為理想，波動(dòng)較小。從總體上看，仿真很好地體現(xiàn)了定高飛行的要求。

圖3為偏航角對(duì)比圖，在完成兩次直線(xiàn)和轉(zhuǎn)彎飛行的時(shí)間上，實(shí)測(cè)數(shù)據(jù)與仿真結(jié)果略有差別。這同樣受到不可控因素的影響，駕駛員在駕駛飛機(jī)時(shí)速度并不能保證恒定，由圖1中也可見(jiàn)，實(shí)測(cè)曲線(xiàn)在轉(zhuǎn)彎處的飛行距離要小于仿真數(shù)據(jù)，因此在圖3的后半段，仿真曲線(xiàn)的變化也會(huì)相應(yīng)延后，但這并不影響整體飛行軌跡的完成和仿真對(duì)飛行特性的體現(xiàn)。在完成轉(zhuǎn)彎過(guò)程中，偏航角發(fā)生變化，由0°變化到180°的過(guò)程中，有一個(gè)先減小再增大的過(guò)程，與圖1中的軌跡正好吻合。因此，偏航角仿真能夠很好地反映飛行試驗(yàn)中的真實(shí)情況。

圖2 飛行高度對(duì)比Fig.2 Flight altitude comparison

圖3 偏航角對(duì)比Fig.3 Yaw angle comparison

圖4 為滾轉(zhuǎn)角和俯仰角對(duì)比圖。由圖可見(jiàn)，實(shí)測(cè)與仿真曲線(xiàn)吻合度較高，仿真結(jié)果能夠極好地模擬直線(xiàn)平飛與轉(zhuǎn)彎兩個(gè)過(guò)程中偏航角和滾轉(zhuǎn)角的對(duì)應(yīng)變化情況。在平飛階段，滾轉(zhuǎn)角為0°，進(jìn)入轉(zhuǎn)彎時(shí)，滾轉(zhuǎn)角發(fā)生變化，隨著轉(zhuǎn)彎的方向，飛機(jī)機(jī)身發(fā)生相應(yīng)傾斜，圖中實(shí)測(cè)與仿真曲線(xiàn)都體現(xiàn)了這種變化。與偏航角相同，在飛行仿真后段，仿真模型的飛行距離大于飛機(jī)飛行距離，轉(zhuǎn)彎時(shí)滾轉(zhuǎn)角的變化也相應(yīng)延后。俯仰角方面，進(jìn)入轉(zhuǎn)彎時(shí)飛機(jī)姿態(tài)發(fā)生變化，相應(yīng)的俯仰角也做出調(diào)整，俯仰角在飛機(jī)進(jìn)入轉(zhuǎn)彎和完成轉(zhuǎn)彎時(shí)均發(fā)生變化，與飛行實(shí)測(cè)數(shù)據(jù)中的姿態(tài)略有差異，這是由駕駛員選擇轉(zhuǎn)彎的方式在仿真中被忽略的原因造成，但是對(duì)飛行軌跡和高度造成的影響較小，對(duì)仿真反映真實(shí)飛行規(guī)律和特點(diǎn)的影響也較小。

圖4 滾轉(zhuǎn)角和俯仰角對(duì)比Fig.4 Roll angle and pitch angle comparison

圖5 為飛行速度對(duì)比圖。飛機(jī)在完成上述飛行過(guò)程時(shí)，駕駛員的駕駛習(xí)慣及其他環(huán)境干擾因素使得速度的波動(dòng)較大，但是依然能夠保持在60 m/s至70 m/s的范圍內(nèi)，在中后段的飛行中，飛機(jī)速度隨著波動(dòng)逐漸變大，因此造成了一系列如滾轉(zhuǎn)角、偏航角等仿真結(jié)果較飛行實(shí)測(cè)數(shù)據(jù)延遲的效果，這些誤差也與圖1中轉(zhuǎn)彎處飛機(jī)平均轉(zhuǎn)彎半徑小于仿真結(jié)果的情況相符。仿真拋開(kāi)了隨機(jī)因素和人為因素的影響，設(shè)定飛行速度為63 m/s，且仿真全過(guò)程中保持較好，僅在轉(zhuǎn)彎時(shí)稍有波動(dòng)，總體速度與實(shí)測(cè)的平均速度相當(dāng)，很好地完成了速度仿真。

圖5 飛行速度對(duì)比Fig.5 Flight velocity comparison

由以上對(duì)比結(jié)果可知，各個(gè)仿真參數(shù)曲線(xiàn)與飛機(jī)飛行試驗(yàn)數(shù)據(jù)所反映的飛行狀態(tài)吻合度較高，仿真能夠較好地模擬固定翼飛機(jī)飛行的基本規(guī)律與真實(shí)情況，從而可為飛機(jī)飛行性能分析與評(píng)估提供有益的參考。

4 結(jié)束語(yǔ)

本文研究了固定翼飛機(jī)飛行動(dòng)力學(xué)仿真問(wèn)題。在保證模型精確度和一定的假設(shè)與簡(jiǎn)化基礎(chǔ)上，建立了非平靜大氣中固定翼飛機(jī)的飛行動(dòng)力學(xué)模型;采用有限元分析軟件計(jì)算固定翼飛機(jī)的氣動(dòng)參數(shù)，從而相對(duì)客觀地描述了飛機(jī)的氣動(dòng)力學(xué)環(huán)境。采用Matlab程序?qū)ι鲜鲲w機(jī)模型進(jìn)行仿真，并將仿真結(jié)果與飛行試驗(yàn)的實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比，結(jié)果顯示仿真能夠體現(xiàn)飛機(jī)的飛行特性，從而驗(yàn)證了模型與仿真的有效性和正確性，為相關(guān)領(lǐng)域的研究提供了有益的參考。

［1］ 劉春，曹碩，王志超.一種固定翼飛機(jī)飛行仿真動(dòng)力學(xué)組件的開(kāi)發(fā)［J］.計(jì)算機(jī)仿真，2013，30(4):49-53.

［2］ 汪泓，謝東來(lái).運(yùn)輸機(jī)飛行仿真技術(shù)及應(yīng)用［M］.北京:清華大學(xué)出版社，2013:114.

［3］ 楊新，王小虎，申功璋，等.飛機(jī)六自由度模型及仿真研究［J］.系統(tǒng)仿真學(xué)報(bào)，2000，12(3):210-213.

［4］ 肖業(yè)倫.飛行器運(yùn)動(dòng)方程［M］.北京:航空工業(yè)出版社，1987:14-17.

［5］ 肖業(yè)倫.航空航天器運(yùn)動(dòng)的建?！w行動(dòng)力學(xué)的理論基礎(chǔ)［M］.北京:北京航空航天大學(xué)出版社，2003:30.

［6］ Allerton D.飛行仿真原理［M］.劉興科，譯.北京:電子工業(yè)出版社，2013:84-87.

［7］ 李新國(guó)，方群.有翼導(dǎo)彈飛行動(dòng)力學(xué)［M］.西安:西北工業(yè)大學(xué)出版社，2005:9.

［8］ Gyllhem D，Mohseni K，Lawrence D，et al.Numerical simulation of flow around the Colorado micro aerial vehicle［R］.AIAA-2005-4757，2005.