考慮彈性變形的月面表取機械臂精度補償方法

陳 鋼，孫豐磊，李 彤，邵 文

（北京郵電大學自動化學院，北京100876）

月球是距離地球最近的天體，是人類進行深空 探測的前沿陣地，對月球的研究將促進人類對月球、地外行星乃至太陽系的認識。月壤采樣返回是中國探月工程“繞、落、回”三部曲的最終章[1]。空間機械臂具有技術成熟、靈活性好等優點，適合代替人類執行采樣任務[2]。表取機械臂是空間機械臂的一個重要分支，在月球表面采樣工作中，表取機械臂的任務是多次抓取月球表面的月壤并傾瀉到采樣封裝裝置中，采樣結束后，提取采樣封裝裝置，將其轉移到返回艙。機械臂末端絕對定位精度決定了表取機械臂能否順利完成采樣任務，是衡量機械臂操作性能的重要指標，提高機械臂絕對定位精度對于確保采樣任務完成至關重要[3?4]。

由于受到加工裝配誤差、關節背隙、編碼器誤差等因素影響，機械臂實際運動參數與名義運動參數之間存在誤差，導致機械臂末端絕對定位精度下降[5]。此外，對于表取機械臂而言，兩桿件為細長桿，跨度較大，桿件上的微小變形映射到機械臂末端會產生較大的末端位置誤差，機械臂自重、負載及外力等引起的彈性變形是導致機械臂末端位置誤差的主要原因。

分析彈性變形對機械臂末端位置誤差影響，需要建立彈性變形量對于末端位置誤差映射關系。針對該問題，國內外學者進行了大量研究。Khalil等基于逆向動力學模型對關節、桿件進行受力分析，建立坐標系間彈性傳遞矩陣，并推導構造幾何參數與彈性參數雅可比矩陣[6]；Meggiolaro 等以多項式近似的方法建立了彈性誤差模型[7]；Zhou 等通過分析關節附加力矩與末端負載映射關系建立有關關節變形的彈性誤差模型[8]；劉志等對機械臂桿件彎曲變形進行分析并提出一種考慮結構變形的運動學補償模型[9]；伍小凱通過建立幾何誤差與彈性誤差綜合運動學誤差模型對機器人進行了運動學辨識[10]；陳宵燕采用線性扭簧理論構建機械臂柔度誤差模型，并基于此建立了機器人剛柔耦合模型[11]。

調研結果表明：現階段彈性建模方法并未綜合考慮靜外力/力矩對機械臂桿件、關節產生的拉伸、彎曲、扭轉變形的影響，造成彈性變形量與機械臂末端誤差間映射關系表征不準確，從而導致后續幾何參數標定精度不高。針對表取機械臂多種彈性變形耦合引起的末端絕對定位精度問題，本文建立了受彈性變形影響的機械臂彈性誤差模型；在此基礎上，結合幾何誤差模型建立機械臂綜合誤差模型，并通過幾何參數標定的方法對末端絕對定位誤差進行補償；通過實驗驗證了文中提出精度補償方法的有效性。

1 機械臂運動學建模

有效的運動學模型必須滿足完備性、比例性和連續性的要求[12]。機械臂運動學建模常常使用Denavit?Hartenberg (DH)法，表取機械臂由一 個偏航關節和3 個俯仰關節組成[13?14]。由于3 個俯仰關節均為平行關節，而平行關節坐標系之間運動學參數對兩軸之間的夾角誤差非常敏感[15]，當兩平行關節軸間存在微小夾角時，會引起運動學參數突變，不滿足有效的運動學模型要求。因此，本文采用改進的DH（Modified DH, MDH)方法建立運動學模型[16]，MDH 方 法 采 用αi、ai、θi、di、βi這5 個 參 數 描述桿件間的關系。其中：αi稱為桿件i 扭角，即關節i 和關節i+1 軸線的夾角；ai為桿件i 長度，即關節i 和關節i+1 軸線公垂線長度；θi為關節i 的轉角，即桿件i 和桿件i+1 軸線的夾角；di為關節i 的偏置，即桿件i 和桿件i+1 軸線公垂線偏置；βi與αi類似，同樣用于描述兩關節間軸線夾角。

“嫦娥5”表取機械臂MDH 坐標系及MDH 參數分別如圖1 及表1 所示。

圖1 表取機械臂MDH 坐標系Fig.1 MDH coordinate systems of the manipulator

表1 表取機械臂MDH 參數Table 1 MDH parameters of lunar sampling manipulator

機械臂桿件坐標系間變換可以表示為

式中：i=1,2,3,4；Trot(X,αi)表示繞X 軸旋轉αi變換矩陣；Ttrans(X,ai)表示沿X 軸平移ai變換矩陣；Trot(Z,θi)表示繞Z 軸旋轉θi變換矩陣；Ttrans(Z,di)表示沿Z 軸平移di變換矩陣；Trot(Y,βi)表示繞Y軸旋轉βi變換矩陣。

從慣性坐標系到末端坐標系的轉換關系為

2 機械臂彈性誤差建模

分析桿件、關節因靜外力導致的彈性變形量，基于此推導桿件彈性模型與關節彈性模型，綜合兩者建立機械臂綜合彈性誤差模型。表取機械臂的實際工作環境為月球重力環境，為保證地面模擬實驗精度滿足機械臂實際工作精度要求，需要模擬月球的重力環境，往往采用施加靜態吊點力的方式來補償地球重力環境，在進行受力分析時需要考慮吊點力對關節、桿件彈性變形的影響。

2.1 桿件彈性變形分析

桿件A 與桿件B 跨度較大，在機械臂自身重力、機械臂上攜帶的設備等載荷的重力、操作物的重力等共同作用下發生彎曲、拉伸、扭轉變形。

2.1.1 桿件B 末端彈性變形誤差量

將桿件B 視為首端固定、末端自由的懸臂梁進行受力分析[17]，桿件B 主要受力為自身重力G4、遠攝像機重力Gc、吊點力F1；所受力矩為遠攝像機重力產生力矩Mc、吊點力產生力矩MF1。

如圖2 所示，將桿件4 的重力和機械臂末端負載簡化到桿件B 末端（關節4 軸線與桿件B 軸線的交點），得主矢G4和主矩M4，主矩M4分解為M4P與M4V，則有

式中：M4P數值為M4P=m4g(a5-r4z)；M4V數值為M4V=m4gr4xcosφ4；m4為 桿 件4 質 量；g 為 重 力 加速度；r4x、r4z分別為桿件4 質心在∑E系下的位置矢量沿XE和ZE軸分量。

受力分析如圖3 所示。桿件B 與水平方向的夾角φ3=θ2+θ3，將G4可分解為G4τ和G4n，G4τ使桿件拉伸，G4n使桿件向+y 方向彎曲，M4P可分解為M4Pτ和M4Pn，M4Pτ使 桿 件 沿+x 方 向 扭 轉，M4Pn使桿件向+z 方向彎曲，M4V使桿件向+y 方向彎曲。桿件自身的重力GB形成均勻載荷qB=mBg/lB（mB為桿件B 質量，lB為桿件B 長度）分解為qBτ和qBn，qBτ使 桿 件 拉 伸，qBn使 桿 件 沿+y 方 向彎曲。遠攝像機的重力簡化到桿件軸線上，得主矢Gc=mcg 和主矩Mc=-mcgrcz（mc為遠攝像機質量，rcz為遠攝像機質心在∑E系下的位置矢量沿ZE軸 分 量 的 模），Gc、Mc引 起 的 桿 件 變 形同G4、M4。

圖2 桿件4 重力向桿件3 末端等效示意圖Fig.2 Equivalent diagram of gravity of link 4 toward the end of link 3

圖3 桿件B 受力分析Fig.3 Schematic drawing of force analysis of link B

吊點在桿件B 上加一個豎直向上的力F1和力矩MF1，其中：MF1=FdFsinφ3，dF為吊點力作用點到桿件B 軸線的垂直距離。F1可分解為F1τ和F1n，F1τ使 桿 件 拉 伸，F1n使 桿 件 向-y 方 向 彎 曲；MF1使桿件向-y 方向彎曲。

桿件B 的變形示意圖如圖4 所示，桿件B 的末端由B1點移動到B2點，則∑3e系的微分移動和轉動DBe為

圖4 桿件B 的變形示意圖Fig.4 Schematic drawing of deformation of link B

2.1.2 桿件A 末端彈性變形誤差量

桿件A 主要受力為桿件3、桿件4 重力G34以及自身重力均勻載荷qA，對桿件A 按照桿件B 的分析思路進行處理，得到桿件A 受力分析圖如圖5 所示。桿件A 末端坐標系∑2e微分移動和轉動DAe為

式中：qAτ、qAn為qA分力，分別使桿件拉伸和沿+y方向彎曲；G34τ、G34n為G34分力，使桿件拉伸變形和沿+y 方向彎曲；M34P、M34V為桿件3、桿件4 重力向桿件A 末端轉換時的等效力矩分力矩，M34Pτ、M34Pn則 為M34P分 力 矩，M34Pτ使 桿 件 沿+x 方 向 扭 轉，M34Pn使 桿 件 向+z 方 向 彎 曲，M34V使 桿 件 向+y 方向 彎 曲；EAA為 桿 件A 抗 拉 剛 度，EIA為 桿 件A 彎曲 剛 度，GIpA為 桿 件A 扭 轉 剛 度，lA為 桿 件A 長度。微分變換為

圖5 桿件A 受力分析圖Fig.5 Schematic drawing of force analysis of link A

2.2 關節彈性變形分析

關節彈性變形主要為繞關節轉動方向的扭轉變形，拉伸變形與彎曲變形可忽略不計。關節受力主要為末端負載重力、吊點力、其余關節以及桿件重力。

2.2.1 關節4 坐標系彈性變形誤差量

如圖6 所示，關節4 主要受到桿件4 重力以及末端負載影響。

圖6 桿件4 重力向關節4 簡化示意圖Fig.6 Equivalent diagram of gravity of link 4 toward joint 4

將桿件4 的重力和機械臂末端負載簡化到關節4 末端（關節4 軸線與桿件4 軸線的交點），得主矢G4=m4g 和 主 矩M4=M4P+M4V。 式 中：M4V=m4gr4xcosφ4，φ4=θ2+θ3+θ4，M4V使關節4 沿+z 方 向 扭 轉；m4為 桿 件4 質 量，r4x、r4z分 別 為桿件4 質心在∑4系下的位置矢量沿X4和Z4軸分量。

關節4 受力分析圖如圖7 所示。則∑4系相對于∑3e系的微分移動和轉動為

2.2.2 其余關節坐標系彈性變形誤差量

關節3 在桿件B、關節4、桿件4 的重力和機械臂末端負載合力G34以及吊點力F1作用下產生扭轉變形，其受力分析如圖8 所示。∑3系相對于∑2e系的微分移動和轉動為

圖7 關節4 受力分析圖Fig.7 Schematic drawing of force analysis of joint 4

圖8 關節3 受力分析圖Fig.8 Schematic drawing of force analysis of joint 3

關節2 扭轉變形為桿件A、關節3、桿件B、關節4、桿件4 的重力和機械臂末端負載合力G234以及吊點力共同作用的結果，其受力分析圖如圖9所示。

圖9 關節2 受力分析圖Fig.9 Schematic drawing of force analysis of joint 2

∑2系相對于∑1系的微分移動和轉動為

G234、吊點力向關節2 末端（關節2 軸線與桿件A 軸線的交點）轉換等效力矩分力矩，GIpJ2為關節2 扭轉剛度。

關節1 受力分析圖如圖10 所示。

圖10 關節1 受力分析圖Fig.10 Schematic drawing of force analysis of joint 1

由受力分析圖可知，關節1 無扭轉變形，則∑1相對于∑0的微分移動和轉動為

2.3 彈性誤差建模

機械臂彈性變形包括關節彈性變形與桿件彈性變形兩部分，關節、桿件彈性變形將導致末端坐標系發生偏移，與運動學參數造成的末端位置誤差發生耦合，導致幾何參數無法準確辨識出來。去除彈性誤差影響，對于提高幾何參數標定結果精度至關重要。

2.3.1 關節彈性模型

設T 和T′分別為慣性系到末端坐標系的名義變換矩陣和實際變換矩陣，則

關節彈性變形后機械臂運動學模型為

由式(12)和式(13)可得關節彈性變形導致的末端變換矩陣微分dT 為

進一步推導可得到關節彈性變形導致的末端位置誤差Df_J為

定義關節i 彈性變形引起的位置誤差與末端位置誤差之間的傳遞矩陣為

則式(17)可以寫為

式中：Jf_J=[JJ1JJ2JJ3JJ4]表示關節彈性誤差與末端位置誤差之間的誤差傳遞矩陣，df_J=[dJ1dJ2dJ3dJ4]T為關節彈性誤差向量矩陣，各個誤差向量如式(8)至式(11)所示。

2.3.2 桿件彈性模型

由式(4)和式(6)可以得到∑2e系和∑3e系分別有位置誤差DAe、DBe。考慮桿件變形的機械臂坐標系如圖11 所示。對于∑2系與∑3系的變換矩陣T3，設∑2系與∑2e系的變換矩陣為T31，∑2e系與∑3系的變換矩陣為T32；對于∑3系與∑4系的變換矩陣T4，設∑3系與∑3e系的變換矩陣為T41，∑3e系與∑4系的變化矩陣為T42。則

圖11 考慮桿件變形的機械臂坐標系Fig.11 Coordinate systems considering deformation of links

由于桿件變形，使得矩陣T31和T41分別有微分變換Δk2e和Δk3e，即矩陣T31和T41的實際值為

桿件變形后機械臂末端坐標系相對于慣性系的變換矩陣為

忽略高階無窮小，結合式(12)可得

定義dT =TΔL，則

可得由桿件變形產生的末端位置誤差Df_L在慣性系下的表示為

2.3.3 綜合彈性誤差模型

定義Jf=[Jf_LJf_J]，表示彈性誤差與末端位置誤差之間的誤差傳遞矩陣，df=[df_Ldf_J]T為彈性誤差向量矩陣，則綜合桿件彈性誤差模型與關節彈性誤差模型可得到彈性變形引起的末端位置誤差向量為

3 機械臂末端綜合精度補償

機械臂末端位置誤差包括由幾何參數誤差引起的位置誤差和彈性變形引起的位置誤差兩部分。上文已完成彈性誤差建模工作，下面主要對幾何參數誤差建模，然后對機械臂末端位置誤差進行耦合性分析，分離由關節、桿件變形導致的末端位置誤差，最終對幾何參數進行標定。

3.1 幾何參數誤差建模

對于桿件坐標系，設Ti和Ti′ 分別為{i+1}系相對于{i}系的名義變換矩陣和實際變換矩陣，對變換矩陣微分進行一階近似，可得

由式(28)可得相鄰桿件坐標系相對于慣性系的位置誤差矩陣為

對變換矩陣微分進行一階近似，可由式（30）計算[18]

計算得

設T 和T′分別為慣性系到末端坐標系的名義變換矩陣和實際變換矩陣，則

由式(28)和式(32)可得

令Ui=TiTi+1…Tn(i=1,2,3,4)，Un+1=I4，式(33)可化為

則Δ 可寫為

定義Dg為由幾何參數誤差引起的末端位置誤差向量

由式(35)和式(36)可推導出

令Jgi表示幾何參數誤差與末端位置誤差之間的誤差傳遞矩陣； ei=[ ΔαiΔaiΔθiΔdiΔβi]T（i=1，2，3，4）表示幾何參數誤差。式(37)可寫為

3.2 模型參數耦合性分析

綜合幾何誤差模型與彈性誤差模型，可得到機械臂末端位置誤差向量為

式中：DEnd=PEmeas-PEnom，PEmeas為末端位置測量值，PEnom為末端位置理論值。

由式(39)表示的機械臂末端位置誤差向量中，末端位置理論值PEnom、由彈性誤差引起的末端位置誤差Df、由幾何參數誤差引起的末端位置誤差Dg、彈性誤差與末端位置誤差之間的傳遞矩陣Jf、幾何參數誤差與末端位置誤差之間的傳遞矩陣Jg，均為幾何參數的函數；此外，由于關節和桿件所受力矩陣FLi和FJi均為動力學參數的函數，因此由彈性變形引起的末端位置誤差向量df為幾何參數和彈性參數的函數。

以E 表示全部幾何參數，S 表示全部彈性參數，則有

式中：E′=E+e，S′=S+Δs。式(40)中幾何參數誤差與彈性參數誤差存在耦合，忽略高階無窮小量，可將式(40)簡化為

3.3 綜合精度補償

經過上述分析可知，末端位置誤差是幾何參數和彈性參數共同作用的結果，幾何參數的標定會受彈性參數誤差的影響，因此需先對彈性參數進行標定，然后再對幾何參數誤差進行標定，精度補償流程如圖12 所示。

由于幾何參數誤差與末端位置誤差之間的關系是線性的，因此采用最小二乘法進行標定。對于n 自由度串聯機器人，利用MDH 建立的運動學模型中含有5n 個幾何參數，取m 組末端位置構建超定方程組，如式(42)所示。為保證順利求解，m 應取大于5n/6 的整數。

本文只關注機械臂末端絕對定位誤差，即式(42)中Dg1…Dgm以及Jg1…Jgm均只有前三行有效，后三行中對應的末端姿態誤差均為零。因此，為了順利求解全部的幾何參數，要求m 應取大于5n/3的整數。

圖12 精度補償算法流程圖Fig.12 Flowchart of algorithm for accuracy compensation

利用最小二乘法求解式(43)，可得

為保證得到的幾何參數具有足夠的精度，需要利用最小二乘法進行有限次迭代，以Ej=[αjajθjdj]T( j=1,…,N )表示第j 次標定過程輸入的幾何參數，則有

式中：當j=1 時（即第一次標定過程），取e1=0，E1為初始幾何參數，E′j=Ej+ej。式(45)的收斂條件為ej+1=0。

4 實驗與結果分析

為了對本文所提出的考慮彈性變形的精度補償方法的正確性和有效性進行驗證，選取28 組關鍵構型，利用激光跟蹤儀測得28 組構型絕對末端位置，將測得的末端位置作為輸入，利用本文提出的精度補償方法求解坐標系運動學參數誤差，進而得到機械臂坐標系運動學參數。如表2 所示為補償前后機械臂MDH 參數值以及誤差值。

表2 表取機械臂補償前后的MDH 參數Table 2 MDH parameters before and after compensation

對標定結果進行驗證，選取50 組驗證構型，對比精度補償前后機械臂末端位置誤差。同時，忽略彈性變形對機械臂末端位置誤差的影響，僅對運動學參數進行辨識，并基于辨識獲得的坐標系運動學參數計算50 組驗證構型末端位置誤差。補償前、考慮彈性變形進行補償、忽略彈性變形僅進行運動學參數辨識3 種情況末端位置誤差沿X、Y、Z 坐標軸的分量以及三軸合成誤差如圖13 所示。

對實驗數據進行處理，可以得到補償前后末端位置誤差在各個坐標軸上的最大值、最小值和平均值如表3 所示。

分析表3 中數據可以看出，補償后末端位置誤差沿各個坐標軸的平均值比標定前有明顯下降，誤差沿X、Y、Z 三軸分量的平均值分別下降88.90%、92.81%、93.13%，三軸合成的位置誤差平均值下降92.16%。補償后末端位置誤差沿各個坐標軸的最大值比標定前分別下降87.67%、84.87%、89.01%，三軸合成的位置誤差最大值下降88.63%。而且，綜合精度補償后末端位置誤差較幾何標定末端位置誤差更小、效果更優。通過分析可知，利用本文提出考慮彈性變形的精度補償方法得到的機械臂運動學參數精度較高，能夠大幅提高表取機械臂末端絕對定位精度。

圖13 精度補償前后末端位置誤差Fig.13 Positioning error before and after compensation

表3 精度補償前后末端位置誤差對比Table 3 Comparison of positioning error in different sit?uations

5 結 論

針對大跨度、低剛度的表取機械臂存在的末端絕對定位精度問題，本文首先對機械臂關節、桿件進行了受力分析和形變分析，并推導了機械臂彈性誤差和末端位置誤差間的映射關系，從而建立了機械臂彈性誤差模型；然后，基于機械臂的MDH 運動學模型建立了機械臂幾何誤差模型，表征機械臂運動學參數誤差和末端位置誤差間的映射關系；進而，建立彈性機械臂綜合誤差模型并且設計了可適用于表取機械臂的精度補償方法，并進行仿真驗證。結果表明，該方法能夠獲得準確的機械臂運動學參數，大幅提高機械臂末端絕對定位精度。