蘇華勇，李海艷

（廣東工業(yè)大學(xué)機(jī)電工程學(xué)院，廣州 510006）

0 引言

機(jī)械臂具有控制范圍廣及其操作靈活性高等優(yōu)點(diǎn)，在工業(yè)各領(lǐng)域中得到廣泛應(yīng)用，并且機(jī)械臂系統(tǒng)在當(dāng)下的發(fā)展需求材料更輕質(zhì)、操作范圍更廣及運(yùn)行精度更高[1]。因此，需考慮機(jī)械臂在運(yùn)動(dòng)過(guò)程中各節(jié)臂發(fā)生彈性變形對(duì)位置控制的影響，這就是典型的柔性多體動(dòng)力學(xué)問(wèn)題，其主要特點(diǎn)是機(jī)械臂系統(tǒng)在運(yùn)動(dòng)的過(guò)程中不但經(jīng)歷著大的剛性轉(zhuǎn)動(dòng)，也伴隨著變形運(yùn)動(dòng)，并且這兩種運(yùn)動(dòng)是高度耦合的[2]。機(jī)械臂的動(dòng)力學(xué)建模是機(jī)械臂運(yùn)動(dòng)分析的基礎(chǔ)，典型的動(dòng)力學(xué)分析方法主要有拉格朗日法、牛頓-歐拉法和Kane 方法等等，目前較常用的是拉格朗日法，它在分析過(guò)程中不必求機(jī)械臂的內(nèi)作用力而只需求其對(duì)應(yīng)的速度[3]。由于彈性變形有無(wú)限多自由度，動(dòng)力學(xué)模型的精確解是無(wú)法得到的，通常將彈性變形離散成有限自由度來(lái)近似[4-6]。典型的離散方法有模態(tài)分析法[7]和瑞利-里茲法[8]等等，目前常用的是通過(guò)瑞利-里茲法離散彈性變形。為了方便模型的求解，常取彈性變形離散形式前兩階近似，但沒(méi)有分析不同階數(shù)對(duì)模型求解準(zhǔn)確度的影響。此外，動(dòng)力學(xué)模型屬于高度復(fù)雜的非線性方程組，從而忽略了模型中一些小變形項(xiàng)以簡(jiǎn)化分析的難度，但這同樣影響了動(dòng)力學(xué)模型求解的準(zhǔn)確度。

本文以三節(jié)臂機(jī)械臂為研究對(duì)象，如圖1 所示。通過(guò)拉格朗日法和瑞利-里茲法推導(dǎo)出準(zhǔn)確的動(dòng)力學(xué)模型，分析彈性變形對(duì)機(jī)械臂運(yùn)動(dòng)的影響及彈性變形離散不同階數(shù)對(duì)模型求解準(zhǔn)確度的影響。

圖1 三節(jié)臂機(jī)械臂

1 機(jī)械臂動(dòng)力學(xué)建模

在拉格朗日法中，只要求出機(jī)械臂系統(tǒng)的總動(dòng)能Ek和總勢(shì)能Ep，得到它們的差值L＝Ek-Ep，最后建立拉格朗日動(dòng)力學(xué)方程：

式中：Qi為系統(tǒng)的廣義坐標(biāo)；Q˙i為對(duì)應(yīng)的廣義速度；Fi為系統(tǒng)的廣義力；n為完整約束方程的個(gè)數(shù)。

1.1 參考坐標(biāo)系

如圖2 所示，分析機(jī)械臂系統(tǒng)首先需要建立一系列的坐標(biāo)系來(lái)表示各機(jī)械臂的位置，其中xoy 為整體坐標(biāo)系，x1o1y1、x2o2y2、x3o3y3為局部坐標(biāo)系，局部坐標(biāo)系分別以各臂桿產(chǎn)生彈性變形的首末端的連線作為坐標(biāo)系的橫軸，這樣建立坐標(biāo)系的好處是各臂桿自身產(chǎn)生的彈性變形不會(huì)影響到其他臂桿的位置描述，從而降低建立動(dòng)力模型的復(fù)雜度。

圖2 機(jī)械臂坐標(biāo)系示意圖

在機(jī)械臂系統(tǒng)中，把各臂桿均考慮為均質(zhì)桿，驅(qū)動(dòng)電機(jī)及機(jī)械手夾持裝置均考慮為集中質(zhì)量，并設(shè)臂桿i 的質(zhì)量為mi，長(zhǎng)度為li，電機(jī)i的質(zhì)量為mi-1，機(jī)械手夾持裝置的質(zhì)量為M3，臂桿i 在整體坐標(biāo)系中的偏轉(zhuǎn)角為θi，將臂桿的位置向量從局部坐標(biāo)系過(guò)渡到整體坐標(biāo)系需要借助旋轉(zhuǎn)矩陣Ai：

其中i＝1，2，3。

1.2 機(jī)械臂的總動(dòng)能

假設(shè)機(jī)械臂系統(tǒng)各臂桿在相對(duì)坐標(biāo)系上一點(diǎn)的位置向量為pi[ xi, 0 ]T，則產(chǎn)生彈性變形wi( xi, t )后的位置為pi′＝[xi,wi(xi,t)]T，末端位置表示為 pli＝[li,0]T，借助旋轉(zhuǎn)矩陣可將在局部坐標(biāo)系上的位置轉(zhuǎn)換成在整體坐標(biāo)系上的位置，轉(zhuǎn)換坐標(biāo)系后的位置向量分別表示如下：

通過(guò)對(duì)位置向量關(guān)于時(shí)間t 求導(dǎo)后便可得到對(duì)應(yīng)的速度。因此，臂桿的總動(dòng)能表示為：

式（6）的前一項(xiàng)表示臂桿i 在運(yùn)動(dòng)時(shí)產(chǎn)生的動(dòng)能，后一項(xiàng)表示表示相應(yīng)臂桿末端負(fù)載的動(dòng)能。

1.3 機(jī)械臂的總勢(shì)能

機(jī)械臂在運(yùn)動(dòng)過(guò)程中，總勢(shì)能的變化除了各臂桿運(yùn)動(dòng)時(shí)產(chǎn)生的重力勢(shì)能Eg1，還有臂桿產(chǎn)生彈性變形產(chǎn)生的重力勢(shì)能Eg2和變形產(chǎn)生的應(yīng)變能Es，各勢(shì)能可以表示為：

式中：EIi為臂桿i的剛度；E為彈性模量；I為慣性矩。

因此，機(jī)械臂總勢(shì)能可以表示為：

1.4 瑞利-里茲法

由于機(jī)械臂臂桿的彈性變形是與距離和時(shí)間有關(guān)的量，通常需要將彈性變形離散成有限自由度來(lái)近似，典型的離散方法有模態(tài)分析法和瑞利-里茲法。目前常用的是瑞利-里茲法，它也是彈性連續(xù)體力學(xué)中求近似解的最基本方法。該方法是將彈性變形項(xiàng)離散成里茲基函數(shù)與對(duì)應(yīng)時(shí)間系數(shù)乘積的級(jí)數(shù)展開(kāi)式，以往的方法一般都是取級(jí)數(shù)展開(kāi)式的前2階進(jìn)行近似。為了進(jìn)行對(duì)比驗(yàn)證，此次選取了級(jí)數(shù)展開(kāi)式前4階來(lái)近似彈性變形項(xiàng)，各臂桿離散后的彈性變形可表示如下：

1.5 建立動(dòng)力學(xué)方程

將彈性變形的離散形式分別代入到動(dòng)能和勢(shì)能的求解式中，最后得到的總動(dòng)能和總勢(shì)能表達(dá)式都是僅與時(shí)間t有關(guān)的多項(xiàng)式。

拉格朗日函數(shù)可表示為L(zhǎng)＝Ek-Ep，方程的廣義坐標(biāo)表示為Q，機(jī)械臂系統(tǒng)的廣義坐標(biāo)除了臂桿的偏轉(zhuǎn)角外，還有考慮各臂桿發(fā)生彈性變形造成的變化，它也是受時(shí)間t影響的量，在彈性變形的離散形式中以各階基函數(shù)對(duì)應(yīng)的時(shí)間系數(shù)表示。方程的廣義力表示為F，在旋轉(zhuǎn)運(yùn)動(dòng)中的廣義力為驅(qū)動(dòng)轉(zhuǎn)矩，并設(shè)電機(jī)的轉(zhuǎn)矩為T(mén)，臂桿1 受到電機(jī)1 的驅(qū)動(dòng)轉(zhuǎn)矩和電機(jī)2反作用的驅(qū)動(dòng)轉(zhuǎn)矩，臂桿2受到電機(jī)2的驅(qū)動(dòng)轉(zhuǎn)矩和電機(jī)3的反作用驅(qū)動(dòng)轉(zhuǎn)矩，臂桿3 受到電機(jī)3 的驅(qū)動(dòng)轉(zhuǎn)矩。因此臂桿i受到的驅(qū)動(dòng)轉(zhuǎn)矩τi為：

最后，機(jī)械臂的拉格朗日動(dòng)力學(xué)方程表示如下：

式中：Q＝[θ,q]T， θ＝[θ1,θ2,θ3]， q ＝ [q11,q12,q13,q14,q21,，。

從式（13）可以看到，機(jī)械臂的動(dòng)力學(xué)模型包含了15個(gè)方程，并且該模型屬于復(fù)雜非線性強(qiáng)耦合的方程組。它既包含了機(jī)械臂各臂桿間剛性運(yùn)動(dòng)的耦合，還包含了剛性運(yùn)動(dòng)與臂桿彈性變形的耦合。前者表示機(jī)械臂的大范圍剛性運(yùn)動(dòng)，而后者會(huì)使機(jī)械臂系統(tǒng)運(yùn)動(dòng)產(chǎn)生擾動(dòng)。該模型考慮了彈性變形的影響，因此對(duì)機(jī)械臂系統(tǒng)描述更加全面，但是也很大程度地增加模型的求解難度。

圖3 動(dòng)力學(xué)模型求解流程圖

1.6 動(dòng)力學(xué)方程求解

以往的方法對(duì)方程組的求解過(guò)程中，為了方便求解就一定程度地簡(jiǎn)化了模型，把方程中的小變形項(xiàng)都省略了，但是這無(wú)疑對(duì)方程的準(zhǔn)確度有一定影響。因此，動(dòng)力學(xué)方程的求解不應(yīng)該舍棄小變形項(xiàng)對(duì)模型的貢獻(xiàn)。

通過(guò)Mathematica對(duì)公式推導(dǎo)和整理，拉格朗日方程組最后可以整合成My¨＝z 的形式，其中M15×15是個(gè)滿秩的對(duì)稱矩陣，是廣義坐標(biāo)的二階導(dǎo)數(shù)項(xiàng)，將其他剩余的項(xiàng)整合到z15×1里面，且方程組存在唯一解。

由于該方程組屬于非線性的常微分方程組，可以借助MATLAB 軟件的ode 函數(shù)求解，但是方程一般需要滿足格式y(tǒng)˙＝f (t,y) 的形式，因此要對(duì)動(dòng)力學(xué)方程作進(jìn)一步處理并求解，流程如圖3所示。

2 數(shù)值仿真

現(xiàn)在通過(guò)具體的參數(shù)進(jìn)行數(shù)值仿真，三節(jié)臂機(jī)械臂的長(zhǎng)度分別為 l1＝l2＝l3＝1.5 m，各臂桿的質(zhì)量分別為 m1＝m2＝3 kg，m3＝2.5 kg，電機(jī)的質(zhì)量分別為M0＝M1＝M2＝5.5 kg，末端夾持裝置的質(zhì)量 M3＝4.5 kg，重力加速度 g＝9.8 m/s2，剛度 EI1＝EI2＝EI3＝406 N·m2，電機(jī)驅(qū)動(dòng)力矩T1＝1 150 sin(2πt)3+850(N·m)；T2＝452sin(2πt)3+450(N·m)；T3＝155 sin(2πt)3+150(N·m)。仿真時(shí)間為1.5 s，并設(shè)初始條件Q ＝Q˙＝ 0 。

為了分析產(chǎn)生變形前后對(duì)位置的影響，考慮將機(jī)械臂的末端位置進(jìn)行對(duì)比，可以表示為：

如圖4 所示，機(jī)械臂末端位置在變形前后的曲線變化存在一定的差異，在變形前的曲線是很光滑的，機(jī)械臂彈性變形后的曲線在局部放大處更是可以看到很明顯的振動(dòng)，這正是彈性變形產(chǎn)生的抖動(dòng)現(xiàn)象導(dǎo)致的。另外，機(jī)械臂各臂桿的偏轉(zhuǎn)角在變形前變化后角度也產(chǎn)生了一定的偏差。在變形前的曲線也是很光滑的，在機(jī)械臂發(fā)生彈性變形后在局部放大處也是可以看到很明顯的鋸齒狀。這些現(xiàn)象都可以說(shuō)明機(jī)械臂系統(tǒng)發(fā)生的彈性變形對(duì)位置控制有一定程度影響。

圖4 機(jī)械臂變形前后對(duì)比

由于彈性變形有無(wú)限多自由度，動(dòng)力學(xué)模型的精確解是無(wú)法得到的，通常將彈性變形離散成有限自由度來(lái)近似。為了分析不同階數(shù)近似的彈性變形對(duì)求解的影響，用彈性變形近似的前2階和前4階進(jìn)行對(duì)比，其誤差分析如圖5所示。

可以看出，選取前2階和前4階近似彈性變形對(duì)模型求解準(zhǔn)確度的影響并不大。但是，選取展開(kāi)式近似的階數(shù)越大，動(dòng)力學(xué)模型規(guī)模就會(huì)更大更復(fù)雜，求解的難度也會(huì)陡增，因此選取階數(shù)多少需要根據(jù)實(shí)際情況進(jìn)行綜合權(quán)衡。

圖5 彈性變形2階與4階近似誤差

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

本文通過(guò)以三節(jié)臂的機(jī)械臂作為研究對(duì)象，考慮彈性變形對(duì)機(jī)械臂系統(tǒng)的影響，建立動(dòng)力學(xué)模型并進(jìn)行數(shù)值仿真求解與分析，可以得到如下的結(jié)論。

（1）機(jī)械臂系統(tǒng)在運(yùn)動(dòng)過(guò)程中產(chǎn)生的彈性變形對(duì)臂桿的位置控制產(chǎn)生的擾動(dòng)較大。因此，機(jī)械臂系統(tǒng)的控制有必要考慮彈性變形的影響，彈性變形的擾動(dòng)特性需要進(jìn)一步分析。

（2）機(jī)械臂的彈性變形離散近似取前2階和前4階對(duì)模型的求解準(zhǔn)確度影響不大，并且階數(shù)增加會(huì)導(dǎo)致模型更加復(fù)雜，階數(shù)的選取也需要根據(jù)實(shí)際情況綜合權(quán)衡。