劉 洋，彭世國，馬宏志，廖維新

（1. 廣東工業(yè)大學(xué) 自動化學(xué)院, 廣東 廣州 510006；2. 香港中文大學(xué) 機械與自動化系, 香港 999077）

人口老年化、交通意外和腦中風(fēng)導(dǎo)致越來越多的人出現(xiàn)下肢運動功能障礙。下肢外骨骼機器人(Lower Limb Exoskeleton Robot, LLER) 能夠幫助下肢傷殘、偏癱、有運動功能障礙的患者和身體虛弱的老年人進行康復(fù)訓(xùn)練和助力行走。因此，備受業(yè)內(nèi)專家學(xué)者和健康醫(yī)療機構(gòu)關(guān)注[1-2]。剛出臺的《“十四五”機器人產(chǎn)業(yè)發(fā)展規(guī)劃》和《“十四五”醫(yī)療裝備產(chǎn)業(yè)發(fā)展規(guī)劃》對于外骨骼機器人相關(guān)技術(shù)的研究給予高度重視和大力支持。下肢外骨骼機器人是一個多變量、強耦合的非線性系統(tǒng)，為了使其動作更協(xié)調(diào)柔順、步態(tài)更自然平穩(wěn)、軌跡更平滑精準，建立人機協(xié)同系統(tǒng)精確的動力學(xué)模型和參數(shù)辨識是首要和關(guān)鍵步驟[3-4]。目前，國內(nèi)外學(xué)者提出的機器人參數(shù)辨識方法有：(1) 解體測量計算法[5]：無需運動數(shù)據(jù)采集和算法設(shè)計，但不易精確測量質(zhì)心位置，且拆解后無法考慮關(guān)節(jié)摩擦力的影響。(2) CAD輔助法[6]：慣性參數(shù)由軟件或理論公式自動計算得到的理想值，但在實際加工及裝配過程中會產(chǎn)生誤差。(3) 線性化整體辨識法[7-8]：先對動力學(xué)模型進行線性化處理，得到線性回歸矩陣和慣性參數(shù)矩陣，再驅(qū)動機器人運動，一次性辨識出所有參數(shù)，但計算量較大，存在重復(fù)辨識和誤差累積等問題。(4) 多次辨識法[9]：對整體辨識法加以改進，將機器人按關(guān)節(jié)分組后進行多次驅(qū)動和參數(shù)辨識，但忽略了關(guān)節(jié)摩擦的影響。王靖森等[10]和席萬強等[11]考慮摩擦，對六自由度機器人進行建模和辨識，提高了模型精度。(5) 內(nèi)部傳感器法：采用編碼器測量關(guān)節(jié)角度，根據(jù)驅(qū)動器電流計算轉(zhuǎn)矩，但轉(zhuǎn)矩常數(shù)隨著溫度和壽命而變化存在誤差。耿令波[12]采用測驅(qū)動電流和關(guān)節(jié)力矩的綜合辨識法，實現(xiàn)慣性參數(shù)和摩擦解耦。(6) 外部傳感器法[13]：通過高精度視覺系統(tǒng)捕獲運動數(shù)據(jù)，通過安裝在基座的力傳感器測力矩，但并不能得到所有的慣性參數(shù)，且機器人基座需要多次定位，位置坐標(biāo)變換計算復(fù)雜。(7) 負載辨識法：陳友東等[14]僅驅(qū)動機器人的部分關(guān)節(jié)減少關(guān)節(jié)耦合誤差，運行激勵軌跡，在負載變化時進行動力學(xué)參數(shù)辨識。陳柏等[15]通過優(yōu)化激勵軌跡對六自由度工業(yè)機器人負載動力學(xué)模型進行辨識，再補償已知機器人本體參數(shù)模型，得到本體和負載結(jié)合的動力學(xué)模型。

針對下肢外骨骼機器人人機協(xié)同系統(tǒng)，本文借鑒上述工業(yè)機器人的參數(shù)辨識方法，并進行綜合與改進，提出一種將外骨骼參數(shù)辨識實驗方法與人體下肢參數(shù)計算方法相結(jié)合的動力學(xué)參數(shù)辨識方案。首先，為了考慮摩擦力對參數(shù)辨識的影響，在文獻[16]建立的二連桿結(jié)構(gòu)機器人動力學(xué)模型中加入摩擦力項，得到轉(zhuǎn)矩表達式和一組待辨識的最小慣性參數(shù)。其次，提出一種靜態(tài)與動態(tài)結(jié)合的實驗辨識方法：先多次驅(qū)動單個關(guān)節(jié)，使大腿和小腿分別運動到不同位置后保持靜止，使轉(zhuǎn)矩表達式中的角速度項和角加速度項變?yōu)榱?，計算更簡單；再以最簡單的正弦函?shù)作為最優(yōu)激勵軌跡驅(qū)動下肢連續(xù)擺動；并利用MATLAB進行在線數(shù)據(jù)采集、計算和曲線繪制，依據(jù)曲線特征和最小二乘法完成外骨骼參數(shù)辨識。再次，選取余弦函數(shù)作為激勵軌跡進行參數(shù)辨識結(jié)果驗證。然后，依據(jù)人體生理特征對下肢參數(shù)開展近似計算，并綜合實驗辨識結(jié)果建立人機協(xié)同系統(tǒng)精確的動力學(xué)模型。最后，采用基于模型上界的滑?？刂?，并引入低通濾波器開展步態(tài)軌跡跟蹤研究。

1 動力學(xué)模型

筆者課題組研制了一款如圖1(a～b) 所示的下肢外骨骼服[16]。當(dāng)使用者穿戴下肢外骨骼服，單腿支撐且另一條腿擺動行走時，為了便于分析和計算，將腳和小腿作為一個整體，下肢的擺動近似簡化為二連桿在矢狀面XOY內(nèi)的運動，如圖1(c) 所示。以擺動腿的髖關(guān)節(jié)O為坐標(biāo)原點，大腿繞其轉(zhuǎn)動，OA長度為L1，質(zhì)量為m1， 質(zhì)心C1距 離O點l1，與豎直方向的夾角為θ1，θ1∈[?30?,120?]。小腿繞膝關(guān)節(jié)A轉(zhuǎn)動，長度為L2，質(zhì)量為m2， 質(zhì)心C2距離A點l2，與大腿延長線方向的夾角為θ2，θ2∈[?120?,0?]。大腿繞通過髖關(guān)節(jié)O并垂直于XOY平面的Z軸轉(zhuǎn)動的轉(zhuǎn)動慣量為I1，小腿繞通過膝關(guān)節(jié)A并平行于Z軸轉(zhuǎn)動的轉(zhuǎn)動慣量為I2。

圖1 下肢外骨骼機器人及二連桿結(jié)構(gòu)Fig.1 The diagram of LLER and two links structure

一般，一個關(guān)節(jié)機器人的動態(tài)性能可由如式(1)二階非線性微分方程描述[17]。

式中：θ ，θ˙，θ¨ ∈Rn×1分別為關(guān)節(jié)的角位置、角速度和角加速度，M(θ)∈Rn×n為機器人的正定慣性矩陣，C(θ,θ˙)∈Rn×1為 離心力和哥氏力矩陣，G(θ)∈Rn×1為重力項，F(xiàn)(θ˙)∈Rn×1為摩擦力項，若只考慮黏滯摩擦和庫侖摩擦，F(xiàn)(θ˙)=fvθ˙+fcsgn(θ˙)，fv∈Rn×n為黏滯摩擦系數(shù)矩陣，fc∈Rn×n為庫倫摩擦系數(shù)矩陣，T∈Rn×1為關(guān)節(jié)驅(qū)動力矩。重力加速度g取9.8 m/s2。

因此，二連桿結(jié)構(gòu)的機器人動力學(xué)方程為

式 (3) 中，共有12個未知參數(shù)，為了計算簡便，重新定義5個新參數(shù)作為一組動力學(xué)最小慣性參數(shù)：

則，式(2) 可以寫成

如圖1(b) 所示，安裝在關(guān)節(jié)處的成套減速電機(包括編碼器、直流伺服電機和行星減速齒輪)較重，且足部踏板質(zhì)量不可忽略，導(dǎo)致各連桿的質(zhì)心不在幾何中心位置，即式(3) 中的l1和l2難以準確描述。同時即使已知各連桿的質(zhì)量，而I1和I2也未知，且摩擦力項F(θ˙)無法通過傳感器檢測。顯然，僅通過測量與計算無法建立精確的系統(tǒng)模型。因此，進行系統(tǒng)動力學(xué)參數(shù)辨識十分必要。

2 參數(shù)辨識實驗

2.1 實驗方案及準備

機器人動力學(xué)參數(shù)辨識是一個復(fù)雜的系統(tǒng)工程，包括動力學(xué)建模、激勵軌跡的選取、實驗設(shè)計、數(shù)據(jù)采集與處理、參數(shù)估計、模型驗證等環(huán)節(jié)[12，19]。本文動力學(xué)參數(shù)辨識與實驗驗證方案如圖2所示。其中，X表示外骨骼動力學(xué)最小慣性參數(shù)集，Xload表示人體下肢動力學(xué)最小慣性參數(shù)集，T表示采用滑?？刂戚敵龅年P(guān)節(jié)驅(qū)動轉(zhuǎn)矩，T'表示經(jīng)過低通濾波后輸出的關(guān)節(jié)驅(qū)動轉(zhuǎn)矩。e,e˙,e¨表示關(guān)節(jié)角度跟蹤誤差及誤差一階和二階導(dǎo)數(shù)。

圖2 動力學(xué)參數(shù)辨識方案設(shè)計Fig.2 The design of dynamic parameter identification of LEER

如圖3所示，實驗前用繩子將外骨骼懸空吊起，用支架將其固定，使其軀干與雙腿都處于豎直靜止?fàn)顟B(tài)?？刂朴彝润y關(guān)節(jié)和膝關(guān)節(jié)處的電機來驅(qū)動大腿和小腿(含足部) 運動，模擬外骨骼左腿豎直支撐且右腿擺動行走的過程。

圖3 LLER參數(shù)辨識實驗過程Fig.3 The process of LLER parameter identification

2.2 靜態(tài)實驗

對于髖關(guān)節(jié)和膝關(guān)節(jié)，定義關(guān)節(jié)離開初始位置為反向運動，返回初始位置為正向運動。輪流讓其中一個關(guān)節(jié)在初始位置固定不動，另一個關(guān)節(jié)運動到不同位姿后保持靜止，此刻角速度和角加速度都為零，轉(zhuǎn)矩表達式中僅剩重力項和摩擦力項，易于參數(shù)辨識。

圖4所示為關(guān)節(jié)轉(zhuǎn)矩與角度的關(guān)系。

第1步：讓右腿豎直且膝關(guān)節(jié)保持靜止， θ2=0?，髖關(guān)節(jié)初始角度 θ10=0?。由上位機給電機控制器發(fā)送髖關(guān)節(jié)角度期望信號θ1d0=30?，控制器在簡單比例控制作用下輸出控制電壓u1=kp1(θ1d0?θ10)驅(qū)動右大腿擺動到 30?位置后停止，記錄該位置的關(guān)節(jié)驅(qū)動轉(zhuǎn)矩T10。 每次將θ1di增 加2?， 直到θ1d30=90?為止，記錄此刻關(guān)節(jié)處驅(qū)動轉(zhuǎn)矩T130， 得到θ1反向增加時髖關(guān)節(jié)驅(qū)動轉(zhuǎn)矩曲線T1?

i，如圖4(a) 所示。

第4步：由上位機給電機控制器發(fā)送膝關(guān)節(jié)角度期望信號 θ2d31=?68?，電機驅(qū)動右小腿反向擺動到?68?位 置后停止，記錄此刻關(guān)節(jié)處驅(qū)動轉(zhuǎn)矩T231。依次將 θ2di增 加2?， 直到θ2d60=?10?，記錄此刻關(guān)節(jié)處驅(qū)動轉(zhuǎn)矩T260， 得到θ2正向增加時膝關(guān)節(jié)驅(qū)動轉(zhuǎn)矩T2+i曲線，如圖4(b) 所示。

根據(jù)式(16) ，結(jié)合實驗測得的60組數(shù)據(jù)，求均值得到fc2=?1.521。根據(jù)式(17) ，結(jié)合實驗測得的60組數(shù)據(jù)，求均值得到X5=0.513。根據(jù)式(12) ，結(jié)合實驗測得的60組數(shù)據(jù)，求均值得到X4=1.871。已知外骨骼大腿的長度L1=0.5， 由式(4) 得X3=X5L1，計算得出X3=0.257。

2.3 動態(tài)實驗

對于髖關(guān)節(jié)和膝關(guān)節(jié)，輪流讓其中一個關(guān)節(jié)在初始位置固定不動，另一個關(guān)節(jié)先勻速運動使其角加速度為零，繪制摩擦曲線可求黏滯摩擦系數(shù)。再通過正弦函數(shù)激勵軌跡讓關(guān)節(jié)周期性擺動求轉(zhuǎn)動慣量參數(shù)。

第1步：讓右側(cè)小腿與大腿保持直線，即膝關(guān)節(jié)不運動，θ2=0 ，θ ˙2=0 ，θ ¨2=0。由上位機給電機控制器發(fā)送髖關(guān)節(jié)角度期望信號θ1d(t)=kt，通過簡單比例控制算法使得電機驅(qū)動大腿帶動小腿和足部一起擺動。髖關(guān)節(jié)勻速運動，θ˙1=k∈[?25,25]， θ ¨1=0。根據(jù)式(6) 得，

式中：T1、X4、X5、fc1均 為已知，計算并繪制F(θ˙1)如圖5(a) 所示。根據(jù)黏滯摩擦和庫倫摩擦特性，將圖5(a) 中實測數(shù)據(jù)擬合成黑色折線，利用最小二乘法進行估計得到fv1=?0.062， 擾動f1=?1.796。

第2步：讓右側(cè)大腿重新保持豎直和靜止，即髖關(guān)節(jié)不運動，θ1=0 ，θ ˙1=0 ，θ ¨1=0。由上位機給電機控制器發(fā)送膝關(guān)節(jié)角度期望信號 θ2d(t)=kt，通過簡單比例控制算法使得電機驅(qū)動小腿帶動足部一起擺動。膝關(guān)節(jié)勻速運動，θ ˙2=k∈[?25,25]， θ¨2=0。

根據(jù)式(7) 得，

式中：T2、X5、fc2均 為已知，計算并繪制F(θ˙2)如圖5(b)所示。根據(jù)黏滯摩擦和庫倫摩擦特性，將圖5(b) 中實測數(shù)據(jù)擬合成黑色折線，利用最小二乘法進行估計得到fv2=?0.503。

圖5 關(guān)節(jié)摩擦力矩F (θ˙)曲線Fig.5 The curves of joint friction torque

第3步：讓右側(cè)大腿重新保持豎直和靜止，即髖關(guān)節(jié)不運動 θ1=0 ，θ ˙1=0 ，θ ¨1=0。膝關(guān)節(jié)激勵軌跡選?。害?d(t)=?30?+20?sin(2πft)，通過簡單比例控制算法使得電機驅(qū)動小腿和足部一起擺動。其中f為0.2 Hz，以頻率100 Hz采集T2和 θ2，計算出每個采樣點的角速度θ ˙2和角加速度θ ¨2。

根據(jù)式(7) 得，

式中：T1、X3、X4、X5、fc1、fv1均為已知，利用最小二乘法辨識可得X1=9.506。

綜上所述，參數(shù)辨識結(jié)果如表1所示。

表1 參數(shù)辨識結(jié)果Table 1 The results of parameter indentification

2.4 辨識結(jié)果驗證

根據(jù)髖關(guān)節(jié)和膝關(guān)節(jié)實際角度的變化范圍：θ1∈[?30?,120?]，θ2∈[?120?,0?]，選取驗證激勵軌跡：

式中，f=1 Hz。將式(26) 與參數(shù)辨識結(jié)果代入式(6～7) 計算，并繪制曲線。在式(26) 激勵軌跡驅(qū)動下，外骨骼髖、膝關(guān)節(jié)同時運動，根據(jù)電機驅(qū)動器電流和已知的轉(zhuǎn)矩系數(shù)，計算得到關(guān)節(jié)驅(qū)動轉(zhuǎn)矩實際值，并繪制曲線。如圖6所示，將轉(zhuǎn)矩計算值與實際測量值進行對比分析發(fā)現(xiàn)，關(guān)節(jié)轉(zhuǎn)矩計算值(紅色曲線)都比實際測量值(藍色曲線)略小，原因在于實際系統(tǒng)運行時，存在外部擾動轉(zhuǎn)矩，而動力學(xué)模型對此忽略不計。在一定誤差范圍內(nèi)，二者波形基本一致，表明外骨骼參數(shù)辨識結(jié)果正確。

圖6 關(guān)節(jié)轉(zhuǎn)矩實際值與計算值Fig.6 The tested and calculated values of joint torque

3 步態(tài)跟蹤滑?？刂?/h2>

3.1 人機協(xié)同系統(tǒng)模型

在本項目中[16,20]，外骨骼服的質(zhì)量m0=20 kg，穿戴者身高168 cm，與外骨骼服尺寸完全相符合，L1=0.5 m，L2=0.38 m ，體重mload=70 kg，根據(jù)成年人各部分肢體長度與身高的比例關(guān)系和各部分肢體質(zhì)量與體重的關(guān)系來確定人體下肢的幾何參數(shù)及慣性參數(shù)[21]。大腿的質(zhì)量為

小腿與足部的總質(zhì)量為

由于人體肌肉組織分布比較均勻，為了計算簡便，將大腿和小腿(含足部)的質(zhì)心都取在幾何中心位置，則有

3.2 基于模型上界的滑?？刂?/h3>

滑模控制(Sliding mode control, SMC)具有響應(yīng)快、對參數(shù)變化和擾動不靈敏、魯棒性強等特點，但其不連續(xù)的開關(guān)特性會引起系統(tǒng)高頻抖振。

由機器人系統(tǒng)動力學(xué)特性可知[17]：存在一個線性回歸矩陣 Φ(θ,θ˙,θ˙r,θ¨r)∈R2×3和依賴于機器人質(zhì)量特性的慣性參數(shù)向量p∈R3×1，使得機器人系統(tǒng)動力學(xué)模型式(1)滿足式(34)的線性關(guān)系：

下肢外骨骼機器人在單腿支撐且另一條腿擺動過程中，髖關(guān)節(jié)和膝關(guān)節(jié)實際角度θ =[ θ1θ2]T，期望角度θd=[ θ1dθ2d]T， 誤差e=θd?θ。Λ 為一個二階正對角矩陣，定義參考角速度 θ ˙r= θ˙d+Λe，則參考角加速度θ ¨r= θ¨d+Λe˙，取滑模面的切換函數(shù)為

3.3 穩(wěn)定性分析

3.4 仿真結(jié)果及分析

圖7表明，采用滑模控制的外骨骼人機協(xié)同系統(tǒng)，除了角速度存在抖動外，步態(tài)軌跡跟蹤效果良好。在零時刻，髖關(guān)節(jié)期望初始位置 θ1d(0)=?30?，膝關(guān)節(jié)期望初始位置 θ2d(0)=0?，而外骨骼雙腿豎直向下，髖關(guān)節(jié)實際初始位置θ1(0)=0?，膝關(guān)節(jié)實際初始位置 θ2(0)=0?。系統(tǒng)起動后，外骨骼膝關(guān)節(jié)立刻就能完全跟隨期望軌跡變化；而髖關(guān)節(jié)在約0.15 s內(nèi)先反向擺動微小角度再立刻正向擺回初始位置，直到θ1d(0)=0?時，髖關(guān)節(jié)才會完全跟隨期望軌跡變化。在實際過程中，因存在慣性，髖關(guān)節(jié)電機產(chǎn)生瞬時微弱震動。

如圖8(a) 所示，滑模控制器輸出轉(zhuǎn)矩具有較大的抖振。在滑?？刂破鞔?lián)一階低通濾波器后，轉(zhuǎn)矩抖振大幅度減小,如圖8(b) 所示。表明低通濾波器能夠有效地抑制或減小滑?？刂埔鸬亩墩?。

圖8 滑?？刂频屯V波前后輸出轉(zhuǎn)矩對比Fig.8 The comparison of output torque before and after SMC low-pass filtering

4 結(jié)論與展望

本文先對外骨骼開展參數(shù)辨識實驗，再對人體下肢進行慣性參數(shù)計算，最后綜合得到二連桿結(jié)構(gòu)的人機協(xié)同系統(tǒng)精確的動力學(xué)模型?；谀Ｐ蜕辖绲幕？刂颇軐崿F(xiàn)步態(tài)軌跡精準跟蹤，低通濾波器能有效減小滑模控制引起的高頻抖振。這為下肢外骨骼動力學(xué)參數(shù)辨識提供了一種解決方案，為基于模型的控制方法提供了一種參考模型,也為下肢外骨骼人機協(xié)同系統(tǒng)的步態(tài)軌跡精準跟蹤控制提供了一種參考方法。

將下肢外骨骼近似簡化為二連桿結(jié)構(gòu)，模型更簡單，待辨識參數(shù)更少。在靜態(tài)與動態(tài)結(jié)合的實驗辨識過程中，利用MATLAB在線采集電機驅(qū)動器的電流和重新標(biāo)定后的轉(zhuǎn)矩系數(shù)來實時計算驅(qū)動轉(zhuǎn)矩，減小誤差，也無需在關(guān)節(jié)內(nèi)安裝力矩傳感器，元件更少，結(jié)構(gòu)更簡單，結(jié)果更精確。與通常采用的求線性化回歸矩陣的參數(shù)辨識方法相比，不僅待辨識參數(shù)更少，而且實現(xiàn)了慣性參數(shù)和摩擦力解耦，還避免了對回歸矩陣求逆。但人體下肢慣性參數(shù)近似計算可能引起模型參數(shù)誤差?？砂惭b人機交互力傳感器構(gòu)成反饋環(huán)節(jié)，通過人工智能算法對模型參數(shù)進行自動調(diào)整，使物理模型與參數(shù)模型的高精度一致。引入一階低通濾波器減小滑?？刂埔鸬亩墩瘢椒ê唵斡行?。還可通過改進滑模控制律或者人工智能算法來減小抖振。另外，外骨骼機器人助力行走的步態(tài)規(guī)劃問題值得進一步研究。