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

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

0 引言

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

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

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

1 飛機模型建立

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

將飛機質(zhì)心動力學矢量方程投影到航跡坐標系中，可得到如下形式:

將飛機繞質(zhì)心轉(zhuǎn)動的動力學矢量方程投影到機體坐標系上。則方程有如下形式:

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

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

確定飛機在空中的姿態(tài)，需要建立描述飛機相對地面坐標系姿態(tài)變化的運動學方程如下:

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

其中，速度VU可表示為:

在描述飛機運動的模型中，還需有描述飛機質(zhì)量變化的微分方程。其表達式為:

式中:Q為燃料消耗率，由油門開度決定［5］。

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

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

聲速c依賴于大氣溫度，計算公式為:

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

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

2 作用在飛機上的氣動力和力矩

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

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

作用在飛機上的氣動力矩沿機體坐標系分解為三個分量:滾轉(zhuǎn)力矩Mx、俯仰力矩My和偏航力矩Mz，可寫成如下形式:

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

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

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

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

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

3 仿真和模型驗證

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

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

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

圖1 水平軌跡對比Fig.1 Horizontal trajectory comparison

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

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

圖2 飛行高度對比Fig.2 Flight altitude comparison

圖3 偏航角對比Fig.3 Yaw angle comparison

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

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

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

圖5 飛行速度對比Fig.5 Flight velocity comparison

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

4 結(jié)束語

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

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

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

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

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

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

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

［7］ 李新國，方群.有翼導彈飛行動力學［M］.西安:西北工業(yè)大學出版社，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.