張 磊, 楊鳳田, 周文雅, 黃 俊

(1.沈陽(yáng)飛機(jī)設(shè)計(jì)研究所， 沈陽(yáng) 110135； 2.沈陽(yáng)航空航天大學(xué)遼寧省通用航空重點(diǎn)實(shí)驗(yàn)室， 沈陽(yáng) 110136； 3.遼寧通用航空研究院， 沈陽(yáng) 110136； 4.大連理工大學(xué)航空航天學(xué)院, 大連 116024； 5.北京航空航天大學(xué)航空科學(xué)與工程學(xué)院, 北京 100191)

隨著航空技術(shù)的發(fā)展，無(wú)人機(jī)越來(lái)越多地受到各國(guó)的重視，無(wú)人機(jī)具有使用范圍廣、應(yīng)用限制少、經(jīng)濟(jì)性好等優(yōu)點(diǎn)[1-2]。在無(wú)人機(jī)初步設(shè)計(jì)階段，對(duì)設(shè)計(jì)周期要求嚴(yán)格[3]，需要快速計(jì)算無(wú)人機(jī)的氣動(dòng)特性參數(shù)，來(lái)分析無(wú)人機(jī)的穩(wěn)定性和操縱性，以及無(wú)人機(jī)起飛爬升、降落巡航等飛行性能[4]，并為下一步設(shè)計(jì)無(wú)人機(jī)飛行控制系統(tǒng)提供參考數(shù)據(jù)[5]。在設(shè)計(jì)無(wú)人機(jī)飛行控制系統(tǒng)前，需要將無(wú)人機(jī)升力、阻力、俯仰力矩等氣動(dòng)系數(shù)代入運(yùn)動(dòng)方程，分析穩(wěn)定性和操縱性，進(jìn)而設(shè)計(jì)飛行控制規(guī)律，因此快速計(jì)算無(wú)人機(jī)氣動(dòng)特性參數(shù)并得到線性化方程至關(guān)重要。但計(jì)算流體力學(xué)(computational fluid dynamics，CFD)軟件計(jì)算量較大，風(fēng)洞實(shí)驗(yàn)成本較高，這就需要一款能夠具有較好準(zhǔn)確性、較高計(jì)算速度的氣動(dòng)計(jì)算軟件。DATCOM程序采用了美國(guó)近幾十年來(lái)的氣動(dòng)試驗(yàn)數(shù)據(jù)，適用于大多數(shù)常規(guī)布局飛行器，具有計(jì)算快速、適應(yīng)性強(qiáng)的優(yōu)點(diǎn)[6]。文獻(xiàn)[7-8]將DATCOM軟件計(jì)算得到的氣動(dòng)參數(shù)和風(fēng)洞實(shí)驗(yàn)數(shù)據(jù)進(jìn)行了比較，表明使用DATCOM軟件可以得到較準(zhǔn)確的氣動(dòng)參數(shù)。因此采用美國(guó)氣動(dòng)力工程估算程序DATCOM估算無(wú)人機(jī)氣動(dòng)系數(shù)[9]。

DATCOM計(jì)算程序包含大量美國(guó)飛機(jī)飛行試驗(yàn)數(shù)據(jù)和圖表公式，適用于大多數(shù)常規(guī)布局飛機(jī)，對(duì)于飛翼布局飛機(jī)，也可通過(guò)一定技巧組合計(jì)算。DATCOM程序有四大類(lèi)輸入?yún)?shù)：Ⅰ類(lèi)是飛行條件、Ⅱ類(lèi)是基本構(gòu)型、Ⅲ類(lèi)是附加構(gòu)型、Ⅳ類(lèi)是控制參數(shù)。根據(jù)各個(gè)模塊要求，輸入無(wú)人機(jī)各項(xiàng)外形尺寸等數(shù)據(jù)，可以通過(guò)NACA系列標(biāo)準(zhǔn)翼型代號(hào)定義飛機(jī)翼型，也可以通過(guò)輸入翼型外形坐標(biāo)點(diǎn)來(lái)定義翼型。輸入飛機(jī)參數(shù)后，DATCOM軟件能夠快速計(jì)算出各種氣動(dòng)系數(shù)(升力系數(shù)、阻力系數(shù)、側(cè)力系數(shù)、俯仰力矩系數(shù)、滾轉(zhuǎn)力矩系數(shù)、偏航力矩系數(shù)和舵面效率等)，并以文本格式輸出結(jié)果，其中氣動(dòng)系數(shù)的單位可以選擇公制或英制。DATCOM軟件計(jì)算速度快，使用簡(jiǎn)便。

隨著無(wú)人機(jī)飛行環(huán)境愈加復(fù)雜，各種外界干擾越來(lái)越多，更需要關(guān)注無(wú)人機(jī)的穩(wěn)定性[10-11]。穩(wěn)定性包括靜穩(wěn)定性和動(dòng)穩(wěn)定性，其運(yùn)動(dòng)分為縱向運(yùn)動(dòng)和橫側(cè)向運(yùn)動(dòng)。靜穩(wěn)定性是指飛機(jī)在平衡狀態(tài)受到擾動(dòng)后，飛機(jī)操縱舵面不動(dòng)，飛機(jī)本身具有自動(dòng)恢復(fù)到原平衡狀態(tài)的趨勢(shì)；動(dòng)穩(wěn)定性是指飛機(jī)受到擾動(dòng)后，能逐漸恢復(fù)到擾動(dòng)之前的運(yùn)動(dòng)狀態(tài)。

雖然電動(dòng)無(wú)人機(jī)動(dòng)力系統(tǒng)不消耗燃油，只消耗動(dòng)力電池電能[12]，動(dòng)力電池質(zhì)量不會(huì)發(fā)生變化，其質(zhì)量和重心在載荷質(zhì)量固定的情況下不隨飛行時(shí)間改變，擁有穩(wěn)定性和操縱性不隨飛行時(shí)間改變的特點(diǎn)。但現(xiàn)分析的承載電動(dòng)無(wú)人機(jī)，任務(wù)載荷質(zhì)量可達(dá)到空機(jī)質(zhì)量的5倍，其穩(wěn)定性和操縱性在承載和空載兩種情況下有較大變化，因此需著重對(duì)比分析這兩種條件下的飛機(jī)穩(wěn)定性和操縱性。

使用CATIA軟件計(jì)算得到電動(dòng)無(wú)人機(jī)的轉(zhuǎn)動(dòng)慣量和質(zhì)量重心，使用DATCOM計(jì)算安裝載荷和空載兩種條件下的電動(dòng)無(wú)人機(jī)升力系數(shù)、阻力系數(shù)、俯仰力矩系數(shù)、升降舵面效率等參數(shù)，使用MATLAB軟件將無(wú)人機(jī)在水平直線飛行狀態(tài)下配平，并計(jì)算得到線性化運(yùn)動(dòng)方程，采用SIMULINK軟件仿真其縱向穩(wěn)定性和操縱性，分析加裝載荷后，電動(dòng)無(wú)人機(jī)縱向穩(wěn)定性和操縱性的變化情況。

1 運(yùn)動(dòng)方程的線性化

采用小擾動(dòng)法將非線性運(yùn)動(dòng)方程簡(jiǎn)化為線性運(yùn)動(dòng)方程進(jìn)行分析[13-14]。小擾動(dòng)法是將無(wú)人機(jī)運(yùn)動(dòng)分為基準(zhǔn)運(yùn)動(dòng)和擾動(dòng)運(yùn)動(dòng)?；鶞?zhǔn)運(yùn)動(dòng)一般選擇定常直線平飛或定常盤(pán)旋等，擾動(dòng)運(yùn)動(dòng)是指在基準(zhǔn)運(yùn)動(dòng)上施加較小的外界干擾。

線性化方法的假設(shè)有：無(wú)人機(jī)油門(mén)固定不變，大氣是平靜無(wú)風(fēng)狀態(tài)；基準(zhǔn)運(yùn)動(dòng)是處于配平條件下的定常直線平飛，無(wú)人機(jī)受到的擾動(dòng)是小量；無(wú)人機(jī)的縱向運(yùn)動(dòng)和橫側(cè)向運(yùn)動(dòng)互不影響。

基于以上假設(shè)，結(jié)合無(wú)人機(jī)縱向運(yùn)動(dòng)方程，經(jīng)過(guò)簡(jiǎn)化可得到其縱向小擾動(dòng)方程為

(1)

(2)

(3)

(4)

式中：V為速度；α為迎角；θ為俯仰角；q為俯仰角速率；δe為升降舵偏角；g為重力加速度；μe為航跡傾斜角；A為狀態(tài)矩陣；B為控制矩陣[15]；E為雅克比矩陣。式中氣流坐標(biāo)軸系內(nèi)力和力矩的量綱導(dǎo)數(shù)定義如下：

2 穩(wěn)定性理論

(5)

則有

(6)

當(dāng)Kn>0時(shí)，即重心在氣動(dòng)焦點(diǎn)之前，Cmα<0，則在迎角α增大時(shí)，負(fù)的Cmα將產(chǎn)生低頭力矩，使無(wú)人機(jī)低頭，進(jìn)而減少α，因此無(wú)人機(jī)是縱向靜穩(wěn)定的。若Kn<0，則Cmα>0，則在α增大時(shí)，產(chǎn)生抬頭力矩，進(jìn)而增加α，則無(wú)人機(jī)是縱向靜不穩(wěn)定的。若Kn=0，則無(wú)人機(jī)是縱向中立穩(wěn)定的。

將無(wú)人機(jī)縱向小擾動(dòng)方程轉(zhuǎn)化為特征方程，通過(guò)求解特征方程的根來(lái)分析無(wú)人機(jī)動(dòng)穩(wěn)定性。特征方程形式如下：

D(S)=S4+a1S3+a2S2+a3S+a4=0

(7)

式(7)中：S為縱向運(yùn)動(dòng)特征方程的復(fù)變量；D(S)為縱向運(yùn)動(dòng)的特征多項(xiàng)式；D(S)=0為縱向運(yùn)動(dòng)的特征方程。

無(wú)人機(jī)受到擾動(dòng)后，其各運(yùn)動(dòng)參數(shù)隨時(shí)間的變化，通常由一些典型模態(tài)所對(duì)應(yīng)的簡(jiǎn)單運(yùn)動(dòng)疊加而成，通過(guò)分析各個(gè)典型運(yùn)動(dòng)模態(tài)特征參數(shù)，研究無(wú)人機(jī)的動(dòng)穩(wěn)定性。飛機(jī)的模態(tài)參數(shù)是根據(jù)小擾動(dòng)運(yùn)動(dòng)方程的特征根進(jìn)行計(jì)算的。常規(guī)布局飛機(jī)的特征方程一般有一大一小兩對(duì)共軛復(fù)根，大的對(duì)應(yīng)短周期模態(tài)，小的對(duì)應(yīng)長(zhǎng)周期模態(tài)。動(dòng)穩(wěn)定性用飛機(jī)運(yùn)動(dòng)的特征方程根是否在S平面的左半平面來(lái)度量。

3 仿真分析

3.1 無(wú)人機(jī)基本參數(shù)

電動(dòng)無(wú)人機(jī)空載時(shí)整機(jī)質(zhì)量為3.2 kg，在無(wú)人機(jī)重心處可均勻放置15 kg質(zhì)量的任務(wù)載荷，承載后氣動(dòng)外形保持不變，重心橫向位置保持不變。無(wú)人機(jī)翼展4.5 m，機(jī)身長(zhǎng)度2.1 m，飛行高度200 m，由電機(jī)及螺旋槳驅(qū)動(dòng)。外形圖如圖1所示。

圖1 無(wú)人機(jī)外形圖Fig.1 The outline drawing of UVA

3.2 縱向穩(wěn)定性

通過(guò)DATCOM軟件，計(jì)算得到空載時(shí)電動(dòng)無(wú)人機(jī)靜穩(wěn)定導(dǎo)數(shù)Cmα為-0.997 7，由于Cmα<0，說(shuō)明該無(wú)人機(jī)是靜穩(wěn)定的。

將DATCOM計(jì)算得到的升力系數(shù)、阻力系數(shù)、俯仰力矩系數(shù)、升降舵效參數(shù)代入電動(dòng)無(wú)人機(jī)非線性運(yùn)動(dòng)方程，利用MATLAB軟件配平線性化得到線性方程，并計(jì)算得到運(yùn)動(dòng)特征方程。在給定油門(mén)狀態(tài)下，空載時(shí)電動(dòng)無(wú)人機(jī)的縱向運(yùn)動(dòng)特征方程為

D(S)=S4+20.567 1S3+112.116S2+

39.951 9S+79.563 7=0

(8)

求解縱向運(yùn)動(dòng)特征方程，其根為-10.16±1.8i，-0.116±0.86i，由于其特征根均在S平面左半部分，說(shuō)明其是穩(wěn)定的。

在給定油門(mén)狀態(tài)下，承載時(shí)電動(dòng)無(wú)人機(jī)的縱向運(yùn)動(dòng)特征方程為

D(S)=S4+16.645 2S3+121.071S2+

15.623 9S+54.774 0=0

(9)

求解縱向運(yùn)動(dòng)特征方程，其根為-8.289±7.1268i，-0.0336±0.6762i，由于其特征根均在S平面左半部分，說(shuō)明其是穩(wěn)定的。

給空載和承載兩種條件下的電動(dòng)無(wú)人機(jī)施加小擾動(dòng)，用Simulink軟件仿真得到V、α、θ和q受擾動(dòng)后的變化仿真曲線，如圖2所示。

圖2 受擾動(dòng)后V、α、θ、q的變化仿真曲線Fig.2 The simulation curve of V, α, θ and q after disturbance

由圖2可以看出，空載時(shí)電動(dòng)無(wú)人機(jī)的速度和俯仰角變化量經(jīng)過(guò)40 s振蕩后趨于零。迎角變化量經(jīng)過(guò)0.3 s急速變小，再經(jīng)過(guò)15 s小幅度振蕩后趨于零。俯仰角速度經(jīng)過(guò)0.5 s急速變小，再經(jīng)15 s小幅度振蕩后趨于零?？蛰d電動(dòng)無(wú)人機(jī)受擾動(dòng)后有恢復(fù)到原飛行狀態(tài)的能力。

相同擾動(dòng)量承載后的電動(dòng)無(wú)人機(jī)速度，俯仰角需要經(jīng)過(guò)100 s振蕩后趨于零。迎角變化量經(jīng)過(guò)0.25 s急速變小，再經(jīng)過(guò)40 s小幅度振蕩后趨于零。俯仰角速度經(jīng)過(guò)0.25 s急速變小，再經(jīng)過(guò)40 s小幅度振蕩后趨于零。承載荷后電動(dòng)無(wú)人機(jī)受擾動(dòng)后也有恢復(fù)到原飛行狀態(tài)的能力。

承載后電動(dòng)無(wú)人機(jī)速度和俯仰角在擾動(dòng)后恢復(fù)時(shí)間為空載時(shí)的2.5倍，承載后電動(dòng)無(wú)人機(jī)的迎角和俯仰角速度快速下降時(shí)間略有減少，小幅度振蕩的恢復(fù)時(shí)間為空載時(shí)的2.7倍。

空載時(shí)電動(dòng)無(wú)人機(jī)XV為-0.209 1,承載后的XV為-0.087 6，空載時(shí)是承載后的2.4倍。由于該項(xiàng)表明速度擾動(dòng)量為正值時(shí)，產(chǎn)生負(fù)的速度增量，抑制速度增加?？蛰d時(shí)小特征根的實(shí)部為-0.116，虛部為±0.86，承載后為-0.033 6，虛部為±0.676，空載時(shí)實(shí)部是承載后的3.45倍，虛部相差不大。承載后的電動(dòng)無(wú)人機(jī)在發(fā)生速度擾動(dòng)后，速度振蕩恢復(fù)趨于零的時(shí)間為空載時(shí)的2.5倍，近似與空載時(shí)的XV和承載時(shí)的XV的比值一致，即速度擾動(dòng)后恢復(fù)時(shí)間與XV正相關(guān)。

空載時(shí)Zα/V為-13.788 5，承載時(shí)Zα/V為-5.777 3，空載時(shí)是承載時(shí)的2.4倍。該項(xiàng)表明迎角擾動(dòng)量為正值時(shí)，產(chǎn)生負(fù)的迎角增量，抑制迎角增加?？蛰d時(shí)大特征根的實(shí)部為-10.16，虛部為±1.8，承載后為-8.289，虛部為±7.12，承載時(shí)虛部是空載后的4倍，空載時(shí)實(shí)部是承載后的1.4倍。承載后的電動(dòng)無(wú)人機(jī)在發(fā)生迎角擾動(dòng)后，迎角振蕩恢復(fù)趨于零的時(shí)間為空載時(shí)的2.7倍,近似與空載時(shí)的和承載時(shí)的Zα/V的比值一致。即迎角擾動(dòng)后恢復(fù)時(shí)間與Zα/V正相關(guān)。

3.3 縱向操縱性

施加1°升降舵舵偏角后，V，α，θ和q的變化仿真曲線如圖3所示。

圖3 施加1°升降舵的V、α、θ和q的變化仿真曲線Fig.3 The simulation curve of V, α, θ and q after 1° elevator

由圖3可以看出，空載時(shí)電動(dòng)無(wú)人機(jī)升降舵偏轉(zhuǎn)1°，速度最終增加0.06 m/s，迎角減小0.14°，俯仰角減小0.14°，俯仰角速度經(jīng)過(guò)40 s趨于零。

承載后電動(dòng)無(wú)人機(jī)升降舵偏轉(zhuǎn)1°，速度增加0.15 m/s，迎角減小0.14°，俯仰角減小0.14°，俯仰角速度經(jīng)過(guò)100 s趨于零。

承載后無(wú)人機(jī)施加1°升降舵偏角后，迎角、俯仰角變化量與空載無(wú)人機(jī)相同，承載后速度增加量是空載時(shí)的2.5倍，承載后的俯仰角速度恢復(fù)為空載的2.5倍。

4 結(jié)論

使用DATCOM軟件分別計(jì)算承載和空載兩種狀態(tài)下電動(dòng)無(wú)人機(jī)氣動(dòng)系數(shù)，根據(jù)小擾動(dòng)理論，將其非線性運(yùn)動(dòng)方程轉(zhuǎn)化為線性方程，使用MATLAB軟件進(jìn)行仿真，得到以下結(jié)論。

(1)空載和承載電動(dòng)無(wú)人機(jī)均具有縱向靜穩(wěn)定性和動(dòng)穩(wěn)定性，均具有良好操縱性。

(2)在電動(dòng)無(wú)人機(jī)受擾動(dòng)時(shí)，承載后速度的恢復(fù)時(shí)間為空載時(shí)的2.5倍，與XV正相關(guān)；迎角的恢復(fù)時(shí)間與Zα/V正相關(guān)。電動(dòng)無(wú)人機(jī)升降舵偏轉(zhuǎn)相同角度時(shí)，承載后速度增量是空載時(shí)的2.5倍，迎角變化量基本相同。

(3)提出的承載前后的電動(dòng)無(wú)人機(jī)的縱向穩(wěn)定性和操縱性的分析方法，對(duì)各種尺寸和重量構(gòu)型的無(wú)人機(jī)飛行控制系統(tǒng)初步設(shè)計(jì)具有重要的參考價(jià)值。