基于多體系統(tǒng)傳遞矩陣法的四旋翼飛行器振動(dòng)建模

范肖肖，賀嘉璠，戚國(guó)慶，李銀伢，盛安冬

（南京理工大學(xué)自動(dòng)化學(xué)院，江蘇南京210094）

在四旋翼飛行器的控制問(wèn)題研究中，建模方法通常是將四旋翼飛行器簡(jiǎn)化為剛體，從整機(jī)受力出發(fā)，基于牛頓第二運(yùn)動(dòng)定律及歐拉方程建立動(dòng)力學(xué)模型[1-3]。然而依據(jù)實(shí)際情況，此類模型未考慮系統(tǒng)的固有頻率，在飛行中四旋翼飛行器受到固有頻率附近的周期激振力時(shí)，容易產(chǎn)生共振現(xiàn)象。目前對(duì)四旋翼飛行器振動(dòng)的研究相對(duì)較少，共振的研究主要集中在航天器、直升機(jī)這類飛行器中[4-6]。文獻(xiàn)[4]針對(duì)共軸鉸接式旋翼直升機(jī)，提出時(shí)—頻分析法研究上下旋翼擺振和機(jī)體運(yùn)動(dòng)間的耦合，揭示共軸式直升機(jī)地面共振的物理機(jī)理。文獻(xiàn)[5]采用傳遞矩陣方法對(duì)直升機(jī)旋翼試驗(yàn)臺(tái)傳動(dòng)系統(tǒng)的振動(dòng)特性進(jìn)行分析。文獻(xiàn)[6]針對(duì)含有大型柔性附件的航天器，用有限元分析軟件與MATLAB聯(lián)合分析系統(tǒng)振動(dòng)模態(tài)。在飛行過(guò)程中四旋翼飛行器也會(huì)受到外界的氣流干擾，當(dāng)擾動(dòng)頻率接近飛行器的固有頻率時(shí)，容易引起共振，從而影響飛行器姿態(tài)的穩(wěn)定性。如果不能對(duì)其振動(dòng)特性尤其是固有頻率進(jìn)行準(zhǔn)確求解，往往難以設(shè)計(jì)出動(dòng)態(tài)性能良好的控制器，因此對(duì)四旋翼飛行器的振動(dòng)特性進(jìn)行準(zhǔn)確表征有著重要的意義。

經(jīng)典的有限元分析法雖然具備對(duì)旋翼類飛行器進(jìn)行模態(tài)分析的能力，但旋翼類飛行器結(jié)構(gòu)復(fù)雜，采用有限元法存在計(jì)算規(guī)模大、運(yùn)算耗時(shí)長(zhǎng)等問(wèn)題。多體系統(tǒng)傳遞矩陣法在自行火炮發(fā)射動(dòng)力學(xué)、多管火箭發(fā)射動(dòng)力學(xué)與艦炮發(fā)射動(dòng)力學(xué)等工程領(lǐng)域的成功應(yīng)用驗(yàn)證了該理論的有效性[7-9]，其中線性多體系統(tǒng)傳遞矩陣法[10]特別適用于復(fù)雜多剛?cè)狍w系統(tǒng)的振動(dòng)特性求解，該方法由于涉及的矩陣階次低從而避免了復(fù)雜多剛?cè)狍w系統(tǒng)振動(dòng)特性計(jì)算病態(tài)，保證了系統(tǒng)特征值的快速求解[11]。

文中基于線性多體系統(tǒng)傳遞矩陣法，建立四旋翼系統(tǒng)每個(gè)元件的傳遞矩陣以及系統(tǒng)總傳遞矩陣。在MATLAB上求解特征方程，獲得系統(tǒng)的固有頻率，并且通過(guò)實(shí)驗(yàn)驗(yàn)證固有頻率求解的準(zhǔn)確性。在此基礎(chǔ)上，修正了四旋翼飛行器未考慮系統(tǒng)諧振模態(tài)的常規(guī)模型，提出一種帶有諧振模態(tài)的四旋翼飛行器動(dòng)力學(xué)模型，對(duì)后續(xù)控制器設(shè)計(jì)中抑制共振現(xiàn)象有實(shí)際意義。

1 振動(dòng)特性求解

1.1 結(jié)構(gòu)分析

本文研究的四旋翼飛行器系統(tǒng)為常見(jiàn)的十字形對(duì)稱布局，簡(jiǎn)化的三維模型如圖1所示。根據(jù)其結(jié)構(gòu)特點(diǎn)，四旋翼飛行器可視為四端輸入一端輸出的空間振動(dòng)多剛?cè)狍w[12]。應(yīng)用線性多體系統(tǒng)傳遞矩陣法的思想，可推導(dǎo)出該模型輸入與輸出之間的關(guān)系，即得到傳遞方程和傳遞矩陣。

圖1 四旋翼飛行器簡(jiǎn)化模型

根據(jù)四旋翼飛行器的結(jié)構(gòu)特征，將其劃分為若干部件并依據(jù)各個(gè)部件的特性，進(jìn)行簡(jiǎn)化處理。機(jī)身簡(jiǎn)化為剛體元件，機(jī)身的彈性效應(yīng)等效為空間彈性鉸元件。機(jī)臂等效成Euler-Bernoulli彈性梁元件[12]。驅(qū)動(dòng)電機(jī)簡(jiǎn)化為剛體元件，驅(qū)動(dòng)電機(jī)與機(jī)臂之間的聯(lián)接簡(jiǎn)化為空間彈性鉸元件。根據(jù)多體系統(tǒng)傳遞矩陣法“體”元件和“鉸”元件統(tǒng)一編號(hào)的原則，將各元件分別編號(hào)，系統(tǒng)邊界均編號(hào)為0，元件的簡(jiǎn)化方式及編號(hào)見(jiàn)表1。

根據(jù)線性多體系統(tǒng)傳遞矩陣法的原理[12]，該多剛?cè)狍w系統(tǒng)有4個(gè)輸入端，分別為槳葉22、槳葉23、槳葉24、槳葉25；機(jī)身1的中心為系統(tǒng)的輸出端，沿輸入端到輸出端的方向即為系統(tǒng)的傳遞方向。沿著傳遞方向，進(jìn)入元件的點(diǎn)稱為輸入點(diǎn)，用I表示，傳遞完成后，離開(kāi)元件的點(diǎn)稱為輸出點(diǎn)，用O表示。在每條傳遞路線上元件之間的聯(lián)接點(diǎn)不僅是前一元件的輸出點(diǎn)，也是后一元件的輸入點(diǎn)，用Pi，j表示該節(jié)點(diǎn)，其中i表示體元件的編號(hào)，j表示鉸元件的編號(hào)。

表1 元件簡(jiǎn)化方式及編號(hào)

1.2 狀態(tài)矢量

圖1中的四旋翼系統(tǒng)含有24個(gè)聯(lián)接點(diǎn)和5個(gè)邊界點(diǎn)。其中聯(lián)接點(diǎn)依次為P1,i(i=2～5)?、Pi,i-4(i=6～9,14～17,22～25)?、Pi,i+4(i=6～9,14～17)?，5 個(gè)邊界點(diǎn)依次為P0,22,P0,23,P0,24,P0,25,P1,0。定義各點(diǎn)的狀態(tài)矢量如下

其中各狀態(tài)變量分別用來(lái)表示各點(diǎn)在模態(tài)坐標(biāo)下的線位移、角位移、內(nèi)力矩和內(nèi)力[12]。其中狀態(tài)矢量Z0,22、Ζ22,18、Z14,18、Z14,10、Z6,10、Z6,2、Z1,2、Z1,3、Z1,4、Z1,5、Z1,0定義在I1x1y1z1坐標(biāo)系內(nèi)，Z0,23、Z23,19、Z15,19、Z15,11、Z7,11、Z7,3定義在I2x2y2z2坐標(biāo)系內(nèi)，Z0,24、Z24,20、Z16,20、Z16,12、Z8,12、Z8,4定義在I3x3y3z3坐標(biāo)系內(nèi)，Z0,25、Z25,21、Z17,21、Z17,13、Z9,13、Z9,5定義在I4x4y4z4坐標(biāo)系內(nèi)。

1.3 傳遞方程

1.3.1 空間振動(dòng)剛體

槳葉22～25和驅(qū)動(dòng)電機(jī)14～17均為單端輸入單端輸出的空間振動(dòng)剛體，傳遞方程為：

傳遞矩陣為：

其中m為剛體的質(zhì)量；在以輸入端I為坐標(biāo)原點(diǎn)的連體坐標(biāo)系中，JI為剛體相對(duì)輸入點(diǎn)I的轉(zhuǎn)動(dòng)慣量矩陣，lIO為剛體輸出端相對(duì)I點(diǎn)的位移，lIC為剛體的質(zhì)心相對(duì)I點(diǎn)的位移，分別為lIO和lIC的叉乘矩陣。

機(jī)體1為四端輸入一端輸出空間振動(dòng)剛體，令P1,r+1(r=1～4)為第r輸入點(diǎn)Ir，則從輸入點(diǎn)P1,2～5到輸出點(diǎn)P1,0的傳遞方程為

在以輸入端I1為原點(diǎn)的連體坐標(biāo)系中，為剛體相對(duì)I1點(diǎn)的傳動(dòng)慣量矩陣，lIrO為剛體輸出端相對(duì)其它輸入端Ir的位移，是其叉乘矩陣。Hα、Hβ、H?分別為坐標(biāo)系I2x2y2z2、I3x3y3z3、I4x4y4z4到I1x1y1z1的坐標(biāo)變換矩陣[12]。

1.3.2 空間彈性鉸

元件 2～5，10～13，18～21為空間彈性鉸，傳遞方程為：

傳遞矩陣為：

式中

Kr,x、Kr,y、Kr,z分別為彈簧在x、y、z3個(gè)方向的彈性剛度，分別為扭簧在x、y、z3個(gè)方向上的扭轉(zhuǎn)剛度。

1.3.3 空間彈性梁

機(jī)臂6，7，8，9為空間振動(dòng)彈性梁，考慮x軸、y軸、z軸的橫向空間振動(dòng)，其傳遞方程為：

式中，Ur∈?12×12，除了以下說(shuō)明元素，Ur中的其余元素均為零：

式中：0≤xi≤li，ρ為梁的密度，為梁的線質(zhì)量密度，Iy,i、Iz,i為截面的慣性積，E為彈性模量，li為梁的長(zhǎng)度。S(x)、T(x)、U(x)、V(x)

為 Крылов 函數(shù)[12]。

1.3.4 系統(tǒng)總傳遞矩陣

根據(jù)文獻(xiàn)[13]中多體系統(tǒng)傳遞矩陣法總傳遞方程的推導(dǎo)定理，可推導(dǎo)出四旋翼飛行器系統(tǒng)的總傳遞方程為：

分別為四旋翼飛行器系統(tǒng)邊界點(diǎn)組成的狀態(tài)矢量和系統(tǒng)總傳遞矩陣，并有：

1.4 固有頻率求解

利用多體系統(tǒng)傳遞矩陣法，拼接各元件的傳遞矩陣得到系統(tǒng)總的傳遞矩陣，代入邊界條件得到系統(tǒng)的特征方程，對(duì)其求解即可得到系統(tǒng)的固有頻率。

不加任何外部激勵(lì)的情況下，對(duì)四旋翼飛行器系統(tǒng)進(jìn)行自由振動(dòng)分析，槳葉和機(jī)身輸出端均為自由端，邊界條件表示如下：

將邊界條件帶入總傳遞方程，刪除狀態(tài)矢量Zall中的零元素記為，刪除狀態(tài)矢量Uall中對(duì)應(yīng)零元素的行和列記為，改寫(xiě)式（13）為如下形式：

從以上過(guò)程可知，和Uall都只與系統(tǒng)自身的結(jié)構(gòu)參數(shù)及固有頻率ωk(k=1,2,3…)有關(guān)。對(duì)于實(shí)際的四旋翼飛行器系統(tǒng)，方程（19）必有非零解，則系統(tǒng)的固有振動(dòng)頻率ωk(k=1,2,3…)對(duì)應(yīng)的矩陣Uˉall需滿足下面的條件

式（20）即為系統(tǒng)的特征方程。求解特征方程得到四旋翼飛行器系統(tǒng)的固有振動(dòng)頻率[14]。

2 數(shù)值仿真與實(shí)驗(yàn)

2.1 固有頻率仿真與分析

對(duì)四旋翼飛行器振動(dòng)分析能夠指導(dǎo)后續(xù)控制系統(tǒng)的設(shè)計(jì)，更好地抑制共振，從而提高四旋翼飛行器姿態(tài)控制系統(tǒng)的性能。

基于實(shí)驗(yàn)室現(xiàn)有的四旋翼飛行器實(shí)驗(yàn)平臺(tái)，得到四旋翼飛行器各部件的主要參數(shù)如表2所示。在線性多體系統(tǒng)傳遞矩陣法對(duì)四旋翼飛行器進(jìn)行模態(tài)分析的基礎(chǔ)上，結(jié)合MATLAB計(jì)算得到系統(tǒng)的各階固有頻率，如表3所示?？梢钥闯?，前10階振動(dòng)模態(tài)分布在12.3～5131.3 Hz之間，一階固有頻率最低，固有頻率隨著模態(tài)階數(shù)的升高呈遞增趨勢(shì)。固有頻率越低，越容易激發(fā)出相應(yīng)的模態(tài)發(fā)生共振，因此工程上著重考慮低階固有頻率。

表2 四旋翼飛行器的主要參數(shù)

表3 四旋翼飛行器振動(dòng)各階固有頻率

2.2 四旋翼飛行器振動(dòng)實(shí)驗(yàn)與分析

為驗(yàn)證固有頻率求解的正確性，對(duì)四旋翼飛行器進(jìn)行實(shí)驗(yàn)研究，實(shí)驗(yàn)平臺(tái)如圖2所示。采用F450四軸飛行器，配有ArduPilot Mega自動(dòng)駕駛儀，通過(guò)地面控制站Mission Planner讀取MPU6000加速度傳感器，得到飛行過(guò)程中的姿態(tài)信息。

圖2 四旋翼飛行器實(shí)驗(yàn)平臺(tái)

實(shí)驗(yàn)過(guò)程中調(diào)節(jié)油門(mén)的大小來(lái)改變電機(jī)的轉(zhuǎn)速，從而使四旋翼飛行器在不同的頻率下工作。為了降低外界環(huán)境對(duì)實(shí)驗(yàn)的干擾，在室內(nèi)環(huán)境下，緩慢增大油門(mén)使旋翼轉(zhuǎn)速逐漸升高，得到四旋翼飛行器在機(jī)體坐標(biāo)系下X、Y、Z3個(gè)方向的加速度時(shí)域曲線，如圖3、5、7所示，其中采樣頻率為50 Hz。

由X、Y、Z3個(gè)方向的時(shí)域曲線可知，在旋翼轉(zhuǎn)速上升初始階段四旋翼飛行器姿態(tài)保持穩(wěn)定，當(dāng)四旋翼飛行器工作在一階固有頻率附近時(shí)，X、Y、Z方向的振動(dòng)范圍均超過(guò)±2 m/s2，四旋翼飛行器發(fā)生共振現(xiàn)象。

傅里葉變換能夠準(zhǔn)確辨識(shí)出信號(hào)的諧振情況，因此對(duì)共振時(shí)間內(nèi)X、Y、Z 3個(gè)方向的加速度時(shí)域信號(hào)進(jìn)行FFT頻譜分析，結(jié)果如圖4、6、8所示，可以看出圖中主頻均在10 Hz左右，對(duì)應(yīng)四旋翼飛行器的一階固有頻率，次頻均為23 Hz左右，對(duì)應(yīng)四旋翼飛行器的二階固有頻率。實(shí)驗(yàn)結(jié)果與本文表3理論計(jì)算的結(jié)果基本一致，證明本文所用的多體系統(tǒng)傳遞矩陣法在四旋翼飛行器振動(dòng)建模方面具有較好的精度。

圖3 X方向振動(dòng)時(shí)域曲線

圖4 X方向75～85s共振頻域曲線

3 四旋翼飛行器模型修正

在研究四旋翼飛行器的控制問(wèn)題中，文獻(xiàn)[15]給出了一類常見(jiàn)的建模方法，根據(jù)牛頓-歐拉方程得出飛行器機(jī)體的非線性動(dòng)力學(xué)模型之后，通過(guò)引入4個(gè)控制量，如式（21）所示，將復(fù)雜的非線性耦合系統(tǒng)解耦為4個(gè)獨(dú)立的控制通道，在懸停狀態(tài)下根據(jù)小擾動(dòng)線性化原理將非線性模型進(jìn)行線性化處理，從而得到每個(gè)通道的線性化模型。

圖5 Y方向振動(dòng)時(shí)域曲線

圖6 Y方向75～85s共振頻域曲線

圖7 Z方向振動(dòng)時(shí)域曲線

圖8 Z方向75～85s共振頻域曲線

其中b為升力系數(shù)，d為反扭矩系數(shù)，Ωi為第i個(gè)旋翼的轉(zhuǎn)速。

根據(jù)文獻(xiàn)[15]的小擾動(dòng)線性化理論以及本文中實(shí)際的四旋翼飛行器參數(shù)，俯仰通道開(kāi)環(huán)傳遞函數(shù)可以表示為:

實(shí)際系統(tǒng)中，四旋翼機(jī)架的柔性、連接件的彈性等非剛體因素引起的諧振效應(yīng)不可忽略，然而式（22）表示的俯仰角傳遞函數(shù)模型未對(duì)系統(tǒng)的諧振模態(tài)進(jìn)行建模。后續(xù)通過(guò)設(shè)計(jì)復(fù)雜的控制器來(lái)補(bǔ)償這種建模誤差，雖然反饋控制元件能夠反饋系統(tǒng)的響應(yīng)，但動(dòng)力學(xué)模型的過(guò)度簡(jiǎn)化常常導(dǎo)致難以預(yù)測(cè)的運(yùn)動(dòng)規(guī)律和多余的系統(tǒng)高頻行為，從而嚴(yán)重影響控制精度和速度，因此引入系統(tǒng)的諧振模態(tài)對(duì)此類模型進(jìn)行修正非常有必要。由于實(shí)際系統(tǒng)中高階模態(tài)較難激發(fā)，因此工程上主要考慮系統(tǒng)的一階模態(tài)。結(jié)合表3得到的四旋翼系統(tǒng)的固有頻率，通過(guò)引入諧振環(huán)節(jié)，修正式（22）的傳遞函數(shù)模型，修正后的俯仰通道開(kāi)環(huán)傳遞函數(shù)表示如下:

其中ω1為四旋翼系統(tǒng)的一階固有頻率，單位rad/s。

修正后的四旋翼動(dòng)力學(xué)模型描述了系統(tǒng)的振動(dòng)特性。圖9實(shí)線部分表示修正后俯仰通道的開(kāi)環(huán)頻率特性，虛線表示加入PID控制器后系統(tǒng)的開(kāi)環(huán)頻率特性。對(duì)比發(fā)現(xiàn)在ω=ω1頻率點(diǎn)附近系統(tǒng)增益均很大，即一階模態(tài)被激發(fā)后系統(tǒng)可能出現(xiàn)較強(qiáng)的震蕩，單純通過(guò)經(jīng)典的PID控制器難以對(duì)機(jī)械諧振抑制產(chǎn)生顯著效果[16]，很難消除系統(tǒng)的震蕩傾向，后續(xù)必須引入其他補(bǔ)償措施濾除此諧振頻率從而抑制共振[17]現(xiàn)象。

圖9 修正后四旋翼俯仰通道的頻率特性

4 結(jié) 論

文中針對(duì)四旋翼飛行器在飛行過(guò)程中的共振現(xiàn)象，以線性多體系統(tǒng)傳遞矩陣?yán)碚摓榛A(chǔ)，建立四旋翼系統(tǒng)的總傳遞矩陣模型，根據(jù)各元件的基本結(jié)構(gòu)參數(shù)獲得四旋翼飛行器的前十階固有頻率[18]，并通過(guò)實(shí)驗(yàn)驗(yàn)證了固有頻率求解的準(zhǔn)確性。在此基礎(chǔ)上，修正了四旋翼飛行器未考慮系統(tǒng)諧振模態(tài)的常規(guī)動(dòng)力學(xué)模型，提出了帶有一個(gè)諧振峰的四旋翼飛行器模型。此模型不僅較為精確地刻畫(huà)了四旋翼飛行器的振動(dòng)特性，也為引入其他補(bǔ)償措施抑制共振現(xiàn)象提供了理論基礎(chǔ)，對(duì)四旋翼飛行器的建模和控制工作具有一定的指導(dǎo)意義。