馬曙光，于 強(qiáng)，王培棟，王天舒

(清華大學(xué) 航天航空學(xué)院·北京·100084)

0 引 言

在軌服務(wù)[1]是近年來(lái)航天領(lǐng)域的研究熱點(diǎn)之一。通過(guò)在軌服務(wù)可以對(duì)航天器執(zhí)行加注燃料、更換零件、拆解組裝和軌道轉(zhuǎn)移等操作，以幫助在軌航天器保持良好工作狀態(tài)和延長(zhǎng)使用壽命，從而降低航天工程的經(jīng)濟(jì)成本?？臻g機(jī)械臂由于其技術(shù)相對(duì)成熟和通用性較強(qiáng)等優(yōu)點(diǎn)，是最常用的在軌服務(wù)執(zhí)行機(jī)構(gòu)之一。關(guān)于空間機(jī)械臂的工程驗(yàn)證項(xiàng)目包括日本的“工程實(shí)驗(yàn)衛(wèi)星-7”(ETS-VII)[2]和美國(guó)的“軌道快車(chē)”[3]等。與依靠關(guān)節(jié)電機(jī)驅(qū)動(dòng)的機(jī)械臂相比，繩驅(qū)機(jī)械臂[4-5]的尺寸可以設(shè)計(jì)得更小。本文展示一種繩驅(qū)空間柔性機(jī)械臂，該機(jī)械臂由多根臂桿串聯(lián)而成，可靈活改變形狀，適合在狹小空間內(nèi)進(jìn)行精細(xì)在軌服務(wù)操作。

常見(jiàn)的機(jī)械臂動(dòng)力學(xué)建模方法包括虛擬臂法[6]、等效臂法[7]和多體動(dòng)力學(xué)方法[8-10]等。其中，虛擬臂法和等效臂法屬于近似方法，其精確程度不如多體動(dòng)力學(xué)方法。對(duì)于本文要討論的繩驅(qū)空間機(jī)械臂，由于其結(jié)構(gòu)較為復(fù)雜且自由度較高，近似方法無(wú)法精確分析其動(dòng)力學(xué)特性，多體動(dòng)力學(xué)方法例如Lagrange 方程、R-W方法、Kane方程和遞推算法[9-10]可能會(huì)更加合適。在上述動(dòng)力學(xué)方法中，遞推算法由于計(jì)算效率較高(時(shí)間復(fù)雜度可達(dá)O(N)或O(logN))，比較適合自由度較高的多體系統(tǒng)。本文采用遞推算法中的一種即空間算子代數(shù)方法[10]來(lái)對(duì)繩驅(qū)空間機(jī)械臂進(jìn)行動(dòng)力學(xué)建模。

數(shù)值仿真效率往往由單個(gè)時(shí)間步的計(jì)算量和時(shí)間步長(zhǎng)大小來(lái)決定，對(duì)于本文所討論的動(dòng)力學(xué)模型而言，其數(shù)值仿真效率主要取決于時(shí)間步長(zhǎng)大小。在多體動(dòng)力學(xué)模型中，通常將動(dòng)力學(xué)方程寫(xiě)成微分方程的形式，然后再利用微分方程的數(shù)值算法[11]進(jìn)行求解。對(duì)于常見(jiàn)的微分方程數(shù)值算法(例如Runge-Kutta法)，由于存在穩(wěn)定區(qū)域[11]對(duì)時(shí)間步長(zhǎng)的限制，時(shí)間步長(zhǎng)過(guò)大時(shí)仿真結(jié)果會(huì)很快發(fā)散，所以時(shí)間步長(zhǎng)的選取存在臨界值，而該臨界值將限制動(dòng)力學(xué)模型的數(shù)值仿真效率。本文通過(guò)對(duì)“穩(wěn)定區(qū)域”的計(jì)算來(lái)分析動(dòng)力學(xué)模型中影響時(shí)間步長(zhǎng)臨界值的因素，并對(duì)動(dòng)力學(xué)模型做出相應(yīng)改進(jìn)以增大時(shí)間步長(zhǎng)臨界值，從而提高動(dòng)力學(xué)模型的數(shù)值仿真效率。

1 動(dòng)力學(xué)模型

1.1 機(jī)械結(jié)構(gòu)

如圖1所示，繩驅(qū)空間機(jī)械臂的主體機(jī)械結(jié)構(gòu)分為根部、中段和端部三部分。其中，根部用于連接機(jī)械臂和漂浮基座；中段主要包含五組臂桿，每組臂桿由六根相同的臂段串聯(lián)而成，全部臂段通過(guò)萬(wàn)向節(jié)進(jìn)行相連；端部用于放置機(jī)械臂的載荷。除30根臂桿之外，機(jī)械臂還包括15根驅(qū)動(dòng)繩索和60根聯(lián)動(dòng)繩索。如圖2所示，驅(qū)動(dòng)繩索穿過(guò)臂桿的過(guò)孔并帶動(dòng)臂桿進(jìn)行運(yùn)動(dòng)。如圖3所示，對(duì)于同一組臂段而言，其組內(nèi)的萬(wàn)向節(jié)由聯(lián)動(dòng)繩索來(lái)保證轉(zhuǎn)角的同步性，該同步性使得機(jī)械臂在運(yùn)動(dòng)時(shí)具有呈分段等曲率的幾何特征，從而增強(qiáng)機(jī)械臂的操控性能。

圖1 空間機(jī)械臂結(jié)構(gòu)示意圖Fig.1 Structure of the space manipulator

圖2 驅(qū)動(dòng)繩索示意圖Fig.2 Driving ropes of the space manipulator

圖3 聯(lián)動(dòng)繩索示意圖Fig.3 Connecting ropes of the space manipulator

1.2 動(dòng)力學(xué)方程

如圖4所示，機(jī)械臂可視為鏈?zhǔn)蕉鄤傮w系統(tǒng)，其動(dòng)力學(xué)方程可表示為

(1)

式(1)中的廣義坐標(biāo)q包含漂浮基座的6個(gè)自由度和機(jī)械臂的60個(gè)轉(zhuǎn)動(dòng)自由度，因此q的維數(shù)為66。

式(1)中的廣義力Q由三部分構(gòu)成：

Q=Qc+Qd+Qe

(2)

式(2)中的Qe為漂浮基座和空間機(jī)械臂受到的外力，可通過(guò)給定的邊界條件計(jì)算得出；Qc為聯(lián)動(dòng)繩索的等效作用力，Qd為驅(qū)動(dòng)繩索的等效作用力，二者需要建立相應(yīng)的繩索模型來(lái)計(jì)算得出。

圖4 空間機(jī)械臂動(dòng)力學(xué)模型示意圖Fig.4 The model of space manipulator as a rigid multibody system

對(duì)于聯(lián)動(dòng)繩索，可將其中的一對(duì)聯(lián)動(dòng)繩索等效為扭簧：

(3)

式(3)中的ΔT為扭簧的等效扭矩，kc為扭簧的等效彈性系數(shù)，cc為扭簧的等效阻尼系數(shù)，Δθ為轉(zhuǎn)角之差，式(3)的推導(dǎo)可通過(guò)將聯(lián)動(dòng)繩索類(lèi)比皮帶傳動(dòng)得到，具體的推導(dǎo)過(guò)程在此略去。

將式(3)同運(yùn)動(dòng)學(xué)方程結(jié)合即可計(jì)算出聯(lián)動(dòng)繩索的等效作用力Qc。

對(duì)于驅(qū)動(dòng)繩索，可將其中某一根繩索視為由臂桿過(guò)孔分割的多段彈性繩連接而成，同一根繩索相鄰兩段的張力大小關(guān)系為

Tj+1=Tjesjθj,sj=-μjsgn(vj)

(4)

式(4)中的μj為臂桿過(guò)孔處的摩擦系數(shù)，θj為兩根繩索之間的角度差，vj為繩索相對(duì)過(guò)孔的運(yùn)動(dòng)速度大小。

將式(4)同運(yùn)動(dòng)學(xué)方程結(jié)合即可計(jì)算得到所有驅(qū)動(dòng)繩索的張力分布規(guī)律，進(jìn)而得到驅(qū)動(dòng)繩索的等效作用力Qd。

當(dāng)Qc、Qd和Qe都得到之后，可將式(2)代入式(1)可以得出空間機(jī)械臂的動(dòng)力學(xué)方程

(5)

式(5)是空間機(jī)械臂動(dòng)力學(xué)分析的基礎(chǔ)。需要說(shuō)明的是，在求解式(5)時(shí)，無(wú)需計(jì)算出質(zhì)量矩陣M后求逆，而是利用空間算子代數(shù)計(jì)算中間變量后利用遞推方法求解。關(guān)于空間算子代數(shù)的算法細(xì)節(jié)可參考相關(guān)文獻(xiàn)[10]。

式(5)在形式上可以改寫(xiě)為

(6)

式(6)代表一組常微分方程，可用相關(guān)數(shù)值算法例如Runge-Kutta法進(jìn)行求解。本文采用四階Runge-Kutta法來(lái)求解式中的微分方程，圖5給出了相應(yīng)的數(shù)值仿真流程。

圖5 數(shù)值仿真流程示意圖Fig.5 The scheme of the numerical simulation

2 數(shù)值仿真效率分析

動(dòng)力學(xué)模型的數(shù)值仿真效率由單個(gè)時(shí)間步的計(jì)算量和時(shí)間步長(zhǎng)大小來(lái)共同決定。在圖5所示的仿真流程中，每一時(shí)間步需要重新計(jì)算各驅(qū)動(dòng)繩索的張力分布規(guī)律，該部分計(jì)算量占據(jù)單個(gè)時(shí)間步總計(jì)算量的2/3左右，而且這些計(jì)算量很難通過(guò)模型簡(jiǎn)化來(lái)避免。對(duì)于本文所討論的動(dòng)力學(xué)模型而言，其數(shù)值仿真效率主要取決于時(shí)間步長(zhǎng)大小。

圖6 各階Runge-Kutta法穩(wěn)定區(qū)域Fig.6 Stability regions for the explicit Runge-Kutta methods

2.1 限制時(shí)間步長(zhǎng)的因素

基于式(5)和式(6)，可利用攝動(dòng)法來(lái)推導(dǎo)機(jī)械臂在平衡狀態(tài)下做微振動(dòng)的動(dòng)力學(xué)方程

(7)

式(7)中，Mv是質(zhì)量矩陣，Kv是剛度矩陣，T是系統(tǒng)總動(dòng)能，E是繩索總的彈性勢(shì)能。

式(7)中各系數(shù)的值可由動(dòng)力學(xué)模型給出，其具體推導(dǎo)過(guò)程在這里略去?？捎墒?7)計(jì)算得出微振動(dòng)的最高階頻率為1217Hz。式(7)可以改寫(xiě)為

(8)

式(8)是一階線性微分方程組，其系數(shù)矩陣對(duì)應(yīng)的特征根滿足|λI-A|=0。

由圖6可知，對(duì)于四階Runge-Kutta法，系數(shù)λ和時(shí)間步長(zhǎng)h需滿足λh<2.8i，可得時(shí)間步長(zhǎng)h與最高階頻率fmax的關(guān)系為

(9)

將fmax=1217代入到式(9)中，可計(jì)算得出時(shí)間步長(zhǎng)臨界值的理論值約為3.8×10-4s，這與通過(guò)數(shù)值仿真實(shí)驗(yàn)得到的經(jīng)驗(yàn)值4×10-4s很接近。

從上述分析中，可以看出時(shí)間步長(zhǎng)臨界值主要由繩驅(qū)空間機(jī)械臂微振動(dòng)的最高階頻率所決定，而該最高階頻率與振動(dòng)方程中的質(zhì)量矩陣和剛度矩陣有關(guān)。對(duì)于繩驅(qū)空間機(jī)械臂而言，質(zhì)量矩陣主要由質(zhì)量分布所決定，而剛度矩陣則由繩索彈性所決定。由于質(zhì)量分布不太受動(dòng)力學(xué)模型的影響，因此時(shí)間步長(zhǎng)臨界值主要取決于剛度矩陣，即影響仿真效率的主要因素為動(dòng)力學(xué)模型中的繩索模型。

2.2 動(dòng)力學(xué)模型的改進(jìn)

通過(guò)2.1的分析可以得知時(shí)間步長(zhǎng)的臨界值主要由繩索彈性所決定。在1.2中給出的動(dòng)力學(xué)模型中，聯(lián)動(dòng)繩索被等效成了扭簧，由于扭簧的等效剛度較高，導(dǎo)致最高階振動(dòng)頻率較高，時(shí)間步長(zhǎng)臨界值較小。為增大時(shí)間步長(zhǎng)的臨界值，可對(duì)聯(lián)動(dòng)繩索的模型進(jìn)行改進(jìn)。

將聯(lián)動(dòng)繩索視為幾何約束，此時(shí)系統(tǒng)的動(dòng)力學(xué)方程可改寫(xiě)為

(10)

式(10)中的R(q)=0為聯(lián)動(dòng)繩索對(duì)應(yīng)的幾何約束。

基于式(10)可利用攝動(dòng)法來(lái)推導(dǎo)空間機(jī)械臂在平衡狀態(tài)下做微振動(dòng)的動(dòng)力學(xué)方程。計(jì)算結(jié)果給出式(10)的最高階頻率值為223Hz，對(duì)應(yīng)的時(shí)間步長(zhǎng)臨界值的理論值為2.0×10-3s，這與對(duì)應(yīng)的經(jīng)驗(yàn)值2×10-3s相符。對(duì)于式(10)所代表的動(dòng)力學(xué)模型，時(shí)間步長(zhǎng)臨界值增大了5倍，即仿真效率最多可以提高5倍。

為了進(jìn)一步提高仿真效率，還可以對(duì)驅(qū)動(dòng)繩索的模型進(jìn)行改進(jìn)。將驅(qū)動(dòng)繩索視為剛性繩，此時(shí)動(dòng)力學(xué)模型中不存在彈性項(xiàng)，理論上將不存在時(shí)間步長(zhǎng)臨界值，時(shí)間步長(zhǎng)的選取將由計(jì)算精度等因素來(lái)決定。

如果將聯(lián)動(dòng)繩索視為幾何約束，同時(shí)將驅(qū)動(dòng)繩索視為剛性繩，數(shù)值仿真實(shí)驗(yàn)顯示時(shí)間步長(zhǎng)可以取至0.1s且仿真結(jié)果未出現(xiàn)發(fā)散現(xiàn)象。與1.2中的動(dòng)力學(xué)模型相比，改進(jìn)后的動(dòng)力學(xué)模型仿真效率提高了200倍以上。特別的，當(dāng)時(shí)間步長(zhǎng)取0.1s時(shí)，動(dòng)力學(xué)模型可以實(shí)現(xiàn)實(shí)時(shí)仿真。

3 結(jié) 論

本文建立了一種繩驅(qū)空間機(jī)械臂的動(dòng)力學(xué)模型，并在此基礎(chǔ)上研究了影響數(shù)值仿真效率的因素。通過(guò)分析微分方程數(shù)值算法的穩(wěn)定區(qū)域，發(fā)現(xiàn)仿真步長(zhǎng)的臨界值主要取決于動(dòng)力學(xué)模型中的繩索彈性。基于上述結(jié)論對(duì)動(dòng)力學(xué)模型進(jìn)行改進(jìn)，仿真算例表明改進(jìn)后的動(dòng)力學(xué)模型相比最初的動(dòng)力學(xué)模型仿真效率可提高200倍以上，并可以實(shí)現(xiàn)實(shí)時(shí)仿真，從而增加了該動(dòng)力學(xué)模型的工程實(shí)際應(yīng)用價(jià)值。