欠約束多機協調吊運系統動力學工作空間分析

蘇 程， 葉佳楠， 李 巍， 丁旺才， 趙志剛

(蘭州交通大學 機電工程學院， 蘭州 730070)

用小型吊運設備組合形成一個大負載緊耦合多機聯合吊運系統，從實用性、經濟性等方面考慮都是值得期待的.為了使系統結構簡單，一般情況下該類緊耦合系統為典型的欠約束定位系統.可應用于戰時橋梁的快速架設，以及民用橋梁的分段精確定位架接等.從類型上講，多機協調并聯吊運系統屬于柔索并聯機器人系統，但與傳統意義柔索并聯機器人又有很大的區別，后者將柔索的一端固定做成特定的運動機構，通過調整柔索長度以實現被吊運物的期望運動位姿，其本質是一個單體機器人，故系統結構可重組性、通用性和模塊化顯得不足.多機協調并聯吊運系統是由模塊化串聯吊運機構組成，模塊化串聯吊運機構可以是飛行機器人、工業機器人和常規起重設備等.欠約束多機協調吊運系統通過調整串聯吊運機器人機構末端位置和柔索長度實現被吊運物的位姿，故欠約束多機協調吊運系統具有很大的靈活性、通用性和工作空間.

目前國內外關于欠約束多機協調吊運系統的研究剛剛起步，柔索并聯機器人的研究成果可為該類多機器人系統的研究提供有價值的借鑒.欠約束多機協調吊運系統的一個主要研究方向是利用多架無人直升機來實現三維空間吊運操作.有許多學者探討了多直升機吊運系統的運動學、動力學模型，并分析了系統的空間配置及穩定性[1-5].

另一個研究方向為多地面機器人系統對物體的吊運操作.Zi等[6]建立了3臺起重機吊運一個重物的運動學和動力學模型，并分析了其工作空間，但只考慮了起重機末端的投影點始終保持正三角形的情況，且僅能變化繩長來改變重物姿態，起重機的末端是固定不動的，這種情況下多機系統與單體繩牽引并聯機器人類似.也有學者應用運動旋量斜率對多機器人吊運系統的運動穩定性進行了研究[7].近年來趙志剛等[8-12]建立了多機協調并聯吊運系統的正逆運動學、動力學模型，并探討各種構型下系統的求解方法，同時分析求解了系統的靜態工作空間.

對于欠約束多機協調吊運系統來說，動力學工作空間可以拓寬系統的運動范圍，到達一些特定的位姿，故研究系統的動力學工作空間是非常有必要的.動力學工作空間是指被吊運物滿足動力學方程時，被吊運物質心所在的空間位姿的集合.文獻[13]中利用Farkas和Stiemke引理對不完全約束和完全約束兩類繩牽引并聯機器人工作空間的邊界進行了研究，其在仿真過程中，固定了參數值，這樣只給出了二維的層圖.本文將在其基礎上，利用Farkas和Stiemke引理，求解前人研究較少的欠約束多機協調吊運系統的動力學工作空間，同時利用蒙特·卡羅方法，數值求解得到系統動力學工作空間的整體效果圖.

圖1 系統示意圖Fig.1 Schematic diagram of system

1 系統動力學建模

圖1所示為m臺吊機通過繩索牽引n自由度被吊運物的示意圖.圖中，每個吊機末端通過一根繩索與被吊運物連接.繩索與吊機末端的連接點為Pm，可通過改變連接點Pm的位置和改變繩長來實現被吊運物的位姿變化.繩索與被吊運物的連接點為Bm.

建立全局坐標系O-XYZ，在被吊運物質心上建立局部坐標O′-X′Y′Z′.被吊運物質心O′在全局坐標系O-XYZ的位置為

r=[xyz]T

(1)

設被吊運物的質量為M，重力加速度為g，則重力用零節距的單位旋量在全局坐標系表示為

G=-Mg[ir×i]T,i=[0 0 1]T

(2)

被吊運物在全局坐標系的速度和角速度分別為

被吊運物在局部坐標系中的慣性矩陣為I，則在全局坐標系中的慣性矩陣為

I′=RIRT

(5)

(6)

式中：R為局部坐標系相對于全局坐標系的旋轉變換矩陣；γ、β和α分別為局部坐標系相對于全局坐標系的X、Y和Z軸的旋轉角度；c和s分別為cos和sin.設每根繩的張力分別為T1,T2, …,Tm，對被吊運物利用牛頓-歐拉方程可得

(7)

式中：I3為3×3單位矩陣；T=[T1T2…Tm]為m根繩的拉力組成的m維列矩陣.

對式(7)進行變形，則被吊運物的動力學方程可表示為

AT=C

(8)

式中：C為n維列矩陣，表示作用在被吊運物上的所有外力旋量之和(包括被吊運物的重力、慣性力等)，

(9)

滿足

A=[e1e2…em]∈Rn×m

(10)

矩陣A是由m個旋量ei構成，

(11)

li=Pi-RBi-r

(12)

i=1,2,…,m

式中：Pi為機器人位置Pi在全局坐標系下的位置矢量；Bi為繩與被吊運物的連接點Bi在全局坐標系下的位置矢量.

2 動力學工作空間條件分析

首先對m≥1+n時的過約束、完全約束系統和m=n時的不完全約束系統的動力學工作空間做簡要分析.隨后主要討論m

(1) 當m≥1+n時，可利用矢量封閉原理的求解方法，結合蒙特·卡羅算法隨機產生單個機器人工作范圍內的位置坐標或者利用掃描法逐行產生單個機器人的位置坐標.當m≥1+n時滿足矢量封閉原理的條件，在結構矩陣A不發生降秩的條件下，AT=C可以轉換成如下形式:

T=A+C+(Im-A+A)k

(13)

式中：A+為結構矩陣A的M-P偽逆；k為m維的任一向量；Im為m×m的單位矩陣.要保證每根繩的拉力Ti都大于零，根據矩陣的相關知識，只要式(13)中的右邊第2項足夠大時，就能使T中的每個元素都為正.因此，(I-A+A)中每個元素都為正時，能保證T中的每個元素都為正.

(2) 當m=n時，可以將重力看成一根特殊的虛擬繩索，該繩的方向始終為重力方向，拉力大小為重力大小，這樣相當于有m+1根繩，使其滿足力封閉條件，這時動力學方程可變換成：

(14)

式中：j=[0 0 -1 -yx0]T.這樣也可以利用力封閉求解算法，并結合蒙特·卡羅算法隨機產生單個機器人工作范圍內的位置坐標或者利用掃描法逐行產生單個機器人的位置坐標.

(3) 當m

由于動力學方程AT=C的形式滿足動力學空間的標準是每根繩的拉力都大于零(即T>0)，在研究過程中可以考慮T=0的時候作為工作空間的邊界，所以動力學工作空間的條件可以表述為

AT=C，T≥0

(15)

本文基于Farkas和Stiemke引理分析欠約束多機器人繩牽引系統的動力學工作空間條件：

(1) 基于Farkas引理分析.設A∈Rn×m，C∈Rn，則下面2組關系式：

有且僅有一組有解.

Farkas引理的幾何解釋是考慮q作為一個分離超平面的法向量，這個超平面分離了2個集合，一個是由A中各個列向量組合而成的凸錐，另外一個是向量C.這個凸錐是A中所有單位旋量ei組成的空間.如果向量C不屬于這個凸錐的空間，則說明這個超平面存在，這樣式(17)有解，也就說明式(16)無解.因此，要使動力學工作空間的條件成立，即式(16)有解，就必須不能形成分離的超平面.

這其中有一種特殊情況，C位于超平面上，這樣CTq=0，顯然式(17)無解，這樣式(16)也有解.由于q是法向量，是一些向量的叉積而成，不可能為零向量，對于CTq=0，會有2種情況：①C≠0，這樣會導致部分的繩拉力大小為零，可以作為工作空間的邊界.②C=0，零向量必共面，如果C=0，則式(15)變成了如下形式：

AT=0，T≥0

(18)

從式(18)可以看出，總存在一組解T=0使式(18)成立.如果T=0即每根繩的張力都為0，但從實際的角度，不可能每根繩的張力都為0卻能對被吊運物實施牽引運動.因此，式(18)將變成如下形式：

AT=0，T>0

(19)

對于式(19)可以利用Stiemke引理進行分析.

(2) 基于Stiemke引理分析.設A∈Rn×m，則下面2組關系式：

有且僅有一組有解.

Stiemke引理的幾何解釋是考慮q作為一個支撐超平面的法向量，使凸錐A位于這個支撐超平面的一邊.如果這個支撐超平面存在，這樣式(21)有解，也就說明式(20)無解.因此，要使動力學工作空間的條件成立，即式(20)有解，就必須不能形成支撐超平面.

總之，通過對Farkas和Stiemke引理的幾何解釋可以看出，要使達到動力學工作空間的條件即式(16)或式(20)有解，就不能形成分離超平面或支撐超平面.

3 欠約束系統動力學工作空間求解

3.1 求解方法

動力學工作空間的求解方法是先通過動力學方程組判斷屬于哪類超平面，再利用單位旋量ei通過不同的組合形成待測定的超平面，然后再測試剩余的單位旋量ei和C來判斷待測定的超平面是否是分離超平面或支撐超平面，在測試的過程中需要測試所有可能的組合形式，來排除形成分離超平面或支撐超平面的可能性.具體測試步驟如下：

(1)A是由m個單位旋量ei構成，每個單位旋量是一個n維的列向量.由于超平面是n-1維子空間的平移，因此，待測定的超平面由n-1個線性無關的單位旋量構成，超平面的法向量是n-1個單位旋量的叉積，可以表示為

(22)

式中：kn為標準的基向量；ei,j為第i個單位旋量的第j個元素.

(2) 根據2個引理測試待測定的超平面是不是分離超平面或支撐超平面：

3.2 特殊情況分析

在測試的過程中存在sgn=0的特殊情況，對于Farkas引理，如上述提及的當C位于超平面上的情況，這樣sgn(CTq)=0.由于法向量q是n-1個單位旋量的叉積構成，同時A中的單位旋量是線性無關的，所以q不可能為零向量.那么存在sgn=0的特殊情況的原因除了上述提到的C=0時的情況，另一種解釋是C≠0，但位于超平面上，剩余沒有形成q的單位旋量必位于超平面的一邊.由于AT=C，是由單位旋量通過非負的組合形成C，同時單位旋量不可能為零向量，但是位于超平面一邊的剩余單位旋量非負組合形成一個向量時不可能與超平面共面，也就不能形成超平面上的C，因此只可能是剩余的單位旋量對應的T中部分拉力為零.對于部分拉力為零的情況，可以作為工作空間的邊界.同理，對于Stiemke引理時，一些剩余單位旋量位于超平面上，可考慮部分繩拉力為零，這樣雖不滿足式(19)，但能使式(15)成立.

在特殊情況測試中，除了能形成法向量q，還需要計算剩余的單位旋量sgn，因此，對于Farkas引理，應滿足m-n+1≥1，即m≥n.對于Stiemke引理，應滿足m-n+1≥2，即m≥n+1.因此，當m≥n時，可以直接利用Farkas引理或者Stiemke引理.但是本文研究的欠約束時即m

m=1或m=2時很難實現對被吊運物的控制，故主要研究的是m≥3且3

對于Stiemke引理的條件m≥n+1，部分情況下無法使用；同時也可以看出，當C=0時，使用的范圍也有局限，適用于被吊運物的重力、慣性力等被其它外力平衡掉或者為零.

4 動力學工作空間計算分析

4.1 計算步驟

在仿真過程中，對于Farkas引理，可以給定被吊運物的位移加速度、角速度、角加速度的范圍，再利用蒙特·卡羅方法隨機產生被吊運物的位移加速度、角速度、角加速度具體值，這樣C就已知了.動力學工作空間的仿真步驟如下：

(1) 先通過已知條件判斷是利用Farkas引理還是Stiemke引理，將動力學方程拆分成2部分.

(2) 確定每個機器人的空間范圍即Pm點的范圍.

(3) 通過機器人的空間范圍確定被吊運物的位置r范圍和姿態范圍.對于Farkas引理，給定被吊運物的位移加速度、角速度、角加速度的范圍.

(4) 利用蒙特·卡羅方法隨機循環N次，每次循環隨機產生每個機器人的位置Pm、被吊運物的位姿和被吊運物的位移加速度、角速度、角加速度.

(5) 在每次循環周期內，計算出每根繩長，并分別計算出拆分后的2組結構矩陣A，判斷A是否發生降秩.如果發生降秩，本循環停止，進入下一循環.如果不發生降秩，則對試待測定的超平面進行判斷，是否為分離超平面或者支撐超平面.如果是，則本循環停止，進入下一循環.如果不是，則重復測試，直到測試完該循環下的全部組合形式都不是分離超平面或者支撐超平面，同時拆分后的2組方程組應同時滿足要求，否則本循環停止，進入下一循環.

(6) 每次循環周期內，如果全部組合形式下都不是分離超平面或者支撐超平面，則此次循環內，被吊運物的質心位置屬于動力學工作空間的點，記錄下此時的空間坐標.

(7) 重復步驟(5)和(6)，直至N次循環完成，最后輸出動力學工作空間的圖形.

4.2 仿真結果分析

為了實驗的方便性和安全性，以3臺直角坐標機器人協同吊運6自由度被吊運物為實體實驗平臺，單個直角坐標機器人末端的工作空間范圍為長方體.設直角坐標系機器人X、Y、Z方向的長度分別為 0.8、0.5 和 0.5 m.在地面上建立全局坐標系，每個機器人Z軸方向高度為 1.5 m. 3臺機器人在全局坐標系中XOY平面擺放示意圖如圖2所示.

圖2 3臺機器人XOY平面示意圖Fig.2 Schematic diagram of the XOY plane of the three robots

設被吊運物為一正三角形物體，三角形的頂點分別為B1、B2和B3，也是繩索與被吊運物的連接點，三角形中心到各頂點的距離d=0.1 m，被吊運物的質量M=10 kg.每個機器人末端點Pm的空間范圍和初設的被吊運物質心O′的空間范圍在全局坐標系中的表達式：

(23)

(24)

(25)

(26)

式中：被吊運物的姿態角范圍為[0, 1] rad；被吊運物的X、Y、Z軸方向位移加速度范圍都為[-1, 1] m/s2；被吊運物分別繞X、Y、Z軸轉動的角加速度范圍都為[0, 1] rad/s2；被吊運物分別繞X、Y、Z軸轉動的角速度范圍都為[0, 1] rad/s.

圖3 被吊運物的動力學工作空間Fig.3 Dynamic workspace of the payload

按照第4.1節中的動力學工作空間的計算步驟，共循環1×105次，共得到 10 159 個滿足要求的被吊運物的位置點.被吊運物動力學工作空間整體效果圖如圖3所示.該圖是保證繩的拉力T≥0的工作空間，因此是包括部分繩的拉力為零作為邊界的工作空間.動力學工作空間在XOY平面、YOZ平面以及XOZ平面的投影圖如圖4～6所示.

從圖4的XOY平面可以看出，工作空間整體大致關于y=x的平面對稱，且離原點越遠越密集.從圖5和圖6可以看出，工作空間的下端比上端要密集些.

圖4 XOY平面投影圖Fig.4 XOY plane projection

圖5 YOZ平面投影圖Fig.5 YOZ plane projection

圖6 XOZ平面投影圖Fig.6 XOZ plane projection

在仿真中，確定了空間參數以后，還限定了角速度、角加速度和線性加速度的范圍，但上述動力學工作空間的分析中并沒有限制繩索的拉力范圍，只是需要判斷拉力的正負.故可推斷：在確定的空間和結構參數條件下，被吊運物角加速度和線性加速度的范圍的變化對動力學工作空間影響不大.

為了驗證以上推論再做一組仿真，將上述的被吊運物的X、Y和Z軸方向位移加速度范圍都變為[-5, 5] m/s2，被吊運物分別繞X、Y和Z軸轉動的角加速度范圍都變為[0, 2] rad/s2，被吊運物分別繞X、Y和Z軸轉動的角速度范圍都變為[0, 2] rad/s，其他參數都不變，這樣C的范圍發生了變化.同樣，按照第4.1節中的動力學工作空間的仿真步驟，共循環1×105次，共得到 10 244 個滿足要求的被吊運物的位置點.

分別對比2組仿真結果圖，發現形狀基本一樣，滿足要求的點的個數也相差不多.由于篇幅的限制，第2組仿真結果就不再圖示.說明C的范圍對動力學工作空間影響不大，這是因為動力學工作空間重點考慮的是保證繩的拉力T≥0，至于繩的拉力具體多大未予考慮，同時從測試步驟中sgn的計算也可以看出，重點只考慮了符號的正負和是否形成分離超平面或者支撐超平面，大小則無關.因此，C的范圍只是影響繩的拉力值大小，而無法影響繩拉力的正負.

5 結論

通過對欠約束多機協調吊運系統動力學工作空間的分析，得出以下結論：

(1) 在動力學工作空間中，通過比較動力學方程與Farkas和Stiemke引理的關系，得出要使動力學工作空間的條件成立，就必須不能形成分離超平面或支撐超平面.

(2) 在動力學工作空間中，當系統為欠約束時，無法直接使用Farkas和Stiemke引理來求解動力學工作空間，但對m≥3時，可以通過拆分動力學方程進行調整，最后能夠得到動力學工作空間.

(3) 通過對比不同被吊運物的加速度范圍時的動力學工作空間仿真圖，發現被吊運物的加速度對動力學工作空間影響不大.