脈沖推力作用下上樓梯雙足機器人行走的建模與動力學分析

陳嘉睿，凌 琳，蔣貴榮

(桂林電子科技大學 數學與計算科學學院, 廣西 桂林 541004)

1990年，McGeer[1]首次提出“被動動態行走”概念。與傳統雙足機器人不同，被動行走機器人能夠在沒有任何輸入情況下，以循環步態在平緩的下坡上行走?？紤]機器人行走的基本特征，學者們建立機器人行走的模型，進而研究其復雜的動力學行為[2-4]。

與平路雙足機器人相比，上樓梯雙足機器人的研究更復雜且更具挑戰性，這是由于樓梯具有復雜幾何形狀，且雙足機器人在爬升樓梯過程中抬腿高度增加、能量損失增多，使得雙足機器人的行走穩定性降低。為解決這些難題，人們設計了帶膝關節或帶伸縮腿的上樓梯雙足機器人，并給出步態規劃[13-14]。Shih等[15]設計一個具有最簡單雙足結構的7-DOF 雙足機器人，該機器人有2條可變長度的腿，用于模擬人類行走時的膝關節；文獻[16]采用運動力控制方案結合遺傳算法(GA)，討論機器人斜坡上行走的最佳步行速度；文獻[17]開發了機器人步態的基于轉矩的比例積分-微分(PID)控制器。

綜上，上樓梯雙足機器人行走的步態規劃和控制算法方面取得了很多成果，但行走動力學的理論研究有待深入。Kuo[11]和HAO等[18]認為擺動階段中在關節處施加力矩，可以保證系統極限環的存在，從而提高機器人行走的穩定性。由于存在擺動腿與樓梯的碰撞切換，上樓梯雙足機器人行走也是一種碰撞模型[19]，具有復雜的動力學行為。

本文考慮伸縮腿結構以解決上樓梯的觸碰臺階等問題，采用幅值限制和與支撐腿角速度相關的脈沖推力作為上樓梯雙足機器人的動力源。擺動階段在關節處施加力矩，提高機器人行走的穩定性。通過理論分析和數值模擬，研究上樓梯雙足機器人行走的復雜動力學行為，比較可調節脈沖推力和固定脈沖推力的性能，驗證可調節脈沖推力的有效性，討論脈沖參數與系統自身結構參數對系統穩定性的影響。

1 上樓梯雙足機器人行走模型

圖1 上樓梯雙足機器人模型示意Fig. 1 Schematic representation of the model of ascending stair biped robot

雙足機器人左右對稱、前后對稱，以支撐腿的質點為原點O建立二維角坐標系，A為髖關節質點，B為擺動腿踝關節質點。機器人初始姿態為兩只腳分別落在階梯中間 (支撐腿位于擺動腿的前方)，髖關節質點A投影在樓梯邊緣上。在運動過程中，雙足機器人的支撐腿的腳部始終貼緊平面且不發生滑動，腳踝以下部分不發生轉動。

如圖2所示，雙足機器人上樓梯過程有4個階段：支撐腿伸長前擺、支撐腿伸長瞬間、支撐腿伸長后擺以及碰撞切換瞬間。本文提出的上樓梯過程描述如下：(I-II) 此階段雙足機器人雙腿等長，即Lst=Lsw=l0(l0為常數)，做類似于最簡雙足機器人擺動階段的運動，支撐腿繞踝關節向前擺動直至與臺階平面垂直；(II-III)為支撐腿伸長階段，此階段發生于θst=0, 支撐腿的長度由l0伸長至L0,且該伸長過程是瞬時完成；(III-IV) 支撐腿繼續繞踝關節向前擺動，直至擺動腿運動到下一級階梯平面；(IV)為碰撞階段, 在原支撐腿踝關節處施加沿其指向髖關節方向的脈沖推力p,在p的作用下支撐腿離地，完成雙腿角色轉換。另外，為了避免新擺動腿與樓梯發生擦碰，碰撞發生的同時，新擺動腿由L0伸長至l0，且踝關節以下部分被收回。階段IV完成后，新擺動腿與新支撐腿處于下一步行走的初始狀態。

圖2 上樓梯雙足機器人步態示意Fig. 2 Gait diagram of ascending stair biped robot。

雙足機器人在上樓梯過程中單腿支撐的擺動階段分為支撐腿伸長前擺動階段(I-II)和支撐腿伸長后擺動階段(III-IV)，這2個階段都可以根據拉格朗日方程推出機械系統的擺動階段動力方程，機器人的動態方程由連續微分方程表示。在支撐腿伸長階段和碰撞階段，機器人的狀態變化可以根據動量矩守恒通過代數映射函數來表示。

1.1 支撐腿伸長前的擺動階段

在階段I-II，機器人雙腿的長度均為l0, 支撐腿與地面無彈起和相對滑動，擺動腿向前擺動。本階段采用拉格朗日方程得到雙足機器人擺動運動時的動力學模型。以支撐腿踝關節質點為原點，建立直角坐標系，定義拉格朗日方程為

式中K和V分別為系統的動能和勢能。根據哈密頓原理，運動方程可以從以下拉格朗日方程得到

(1)

式中Fi為施加到質點上的廣義力。分別計算系統的動能和勢能，可得拉格朗日表達式為

(2)

將式(2)代入式(1)，整理可得

(3)

1.2 支撐腿伸長階段

在階段III中，支撐腿運動到與臺階垂直位置的瞬間，支撐腿的伸縮裝置開始作用，該階段的條件可表示為

HS1={θ∈R2:T(θ)=θst=0}。

除重力外，雙足機器人整體只在支撐腿踝關節受到外力作用，機器人整體關于O點的角動量守恒，如式(4)。擺動腿僅受到髖關節A處的外力作用，擺動腿關于髖關節A的角動量守恒，如式(5)。

(4)

(5)

支撐腿伸長前后，雙足機器人關于支撐腿的角動量以及擺動腿關于髖關節的角動量都是守恒的，即滿足

1.3 支撐腿伸長后的擺動階段

在階段IV，機器人支撐腿長度為L0,擺動腿長度為l0，擺動腿向前擺動。本階段采用拉格朗日方程得到雙足機器人擺動運動時的動力學模型。

分別計算系統的動能和勢能，可得拉格朗日表達式為

(6)

將式(6)代入式(1)，整理可得

1.4 擺動腿與樓梯碰撞切換階段

在上述假設條件下，當滿足以下3個條件時，即為沖擊階段。

① 擺動腿到達下一級樓梯表面；

② 擺動腿向下移動；

③ 擺動腿運動到支撐腿前方,且兩腳質點水平相距2S。

這3種沖擊條件可以分別由以下集合HS2轉換

(7)

由假設可知，雙足機器人發生碰撞時，擺動腿與樓梯表面的接觸是剛性的，即非彈性的。雙足機器人在行走周期結束時發生剛性碰撞，碰撞是瞬間發生的，階梯對于擺動腿的作用力被建模為脈沖矢量，擺動腿與行走下一級階梯瞬時碰撞會導致兩條腿和站姿發生轉換，擺動腿變成站姿腿，反之亦然。在此過程中，雙足機器人的配置變量不變，但廣義速度會發生跳躍。瞬間碰撞切換模型如圖3所示(V-和V+分別為施加脈沖推力前髖關節的運動速度和地面沖擊后髖關節的運動速度)。

圖3 雙足機器人的瞬時碰撞模型Fig. 3 Instantaneous collision switching model for bipedal robots

由圖3的幾何條件，很容易求出碰撞前后的各自由度變換關系以及雙腿長度間的關系如下

(8)

(9)

(10)

由擺動腿關于髖關節的動量矩守恒可得如下關系：

(11)

結合式(11)，雙足機器人碰撞前后的狀態變化關系可以用代數方程描述，寫成矩陣的形式如下

1.5 上樓梯雙足機器人步態的脈沖混雜模型

帶伸縮裝置的上樓梯雙足機器人的周期步態包括4個階段：支撐腿未伸長的擺動階段、支撐腿伸長階段、支撐腿伸長后的擺動階段以及擺動腿與階梯碰撞切換階段。利用拉格朗日動力學原理和角動量守恒原理建立上樓梯雙足機器人的動力學方程，其中支撐腿伸長前后的擺動階段都可以由連續型非線性微分方程描述；支撐腿伸長過程和擺動腿與階梯碰撞過程都可以由離散型代數切換映射描述，因此上樓梯雙足機器人的動力學模型可表示為由多個連續型非線性微分方程和多個離散型代數映射構成的脈沖混雜系統。

(12)

式中：Ω1?R4,Ω2?R4;HS1,HS2分別定義如下：

HS1={x∈R2:T(x)=θst=0},

HS2={x∈R4:Υ1(x)=H,Υ2(x)<0,Υ3(x)=2S},

fi(x)、Δi(x)以及g(x)分別為

2 龐加萊映射

2.1 脈沖推力

(13)

即pk滿足

2.2 龐加萊映射

龐加萊映射是一種研究周期運動穩定性及其隨參數變化的分岔特征的幾何方法。為了建立上樓梯雙足機器人動態行走步態對應的龐加萊映射，選取擺動腿與階梯碰撞后瞬間機器人的狀態所在的超平面作為龐加萊截面，記為Σ; 定義龐加萊映射為Φ:Σ→Σ，其中

在龐加萊截面上取一點x0∈Σ作為系統的初始狀態，記xk為系統第k次返回截面∑的狀態變量，則龐加萊映射Φ:Σ→Σ可表示為

xk+1=Φ(xk)。

若存在一點xd∈Σ，滿足

xd=Φ(xd),

則稱xd為龐加萊映射的不動點，生成的龐加萊映射的不動點可用于分析行走模式和構建模型分岔圖。

2.3 穩定性分析

運用龐加萊映射方法，通過分析系統不動點的穩定性來分析行走穩定性。

假設系統不動點為xd，對其施加微小擾動Δxd，擾動點經龐加萊映射作用后為x*, 采用泰勒級數展開， 則有

x*=Φ(xd+Δxd)=Φ(xd)+JΔxd,

式中J為龐加萊映射函數相對于系統狀態的梯度，稱為雅可比矩陣，其表達式為

行走中的誤差通過雅可比矩陣J傳遞，J決定了不動點的穩定性[9]，若J的最大特征值λm的模小于1或所有特征值均在復平面單位圓里，則誤差會逐漸消除，行走是穩定的；反之，行走的誤差會隨著步數增加而被放大，行走不穩定。特別地，當J存在一個特征值λ0=-1時，則有一個超臨界的Flip分岔，此時一個周期-2解從周期-1分岔出來。

3 仿真實驗及分析

3.1 周期步態的初始條件及參數設置

假設一周期步行時間為T，雙足機器人上樓梯的過程被分為4個階段：支撐腿伸長前擺動階段(0～Ts)、支撐腿伸長階段(Ts)、支撐腿伸長后擺動階段(Ts～T)以及碰撞階段(T)。

(14)

周期步態參數設置的具體流程如圖4所示。

圖4 周期步態參數設置的流程Fig. 4 Flowchart for setting parameters of the cycle gait

3.2 上樓梯雙足機器人模型的穩定周期步態

如圖5(a)所示，將不動點xd作為初始狀態，上樓梯雙足機器人的狀態在相空間中形成一個穩定的極限環，其中點R1(K1)對應支撐腿(擺動腿)運動的初始點；點R2(K2)對應支撐腿(擺動腿)在圖2中(I-II)和(II-III)狀態之間的切換點；點R3(K3)對應支撐腿(擺動腿)在圖2中(II-III)和(III-IV)狀態之間的切換點；點R4(K4)對應支撐腿(擺動腿)在擺動足發生瞬間的碰撞切換點。點Ri和點Ki(i=1, 2, 3, 4)分別表示支撐腿和擺動腿在一個周期步態內4個階段的切換點。圖5(b)為雙足機器人穩定的周期-1步態下雙腿的角度和角位移的時間序列，可以看出在周期-1步態中支撐腿的角位移由正值逐漸減小到零后變為負值再逐漸增大，擺動腿角位移則是由負值向正方向運動直至最大角位移處，隨后擺動腿做“回擺”運動，此時擺動角位移逐漸減小。圖6所表示的是擾動解的相圖與時間序列，由圖6(a)可知，系統從出發的軌線趨向于對應不動點的周期-1解。圖6(b)表示以擾動點為初始點的雙足機器人雙腿的角度和角位移隨時間變化逐漸穩定。圖7為雙足機器人上樓梯的棍棒示意。結合圖5～7可以看出該雙足機器人在這組參數下具有穩定的周期步態。

圖5 系統(12)在b=3.12,d=3.12,c=2.9 時的周期-1解Fig. 5 Period -1 gait for system (12) with b=3.12,d=3.12,c=2.9

圖6 系統(12)在b=3.12,d=3.12,c=2.9，初始點為[0.694 5,-0.627 6,-0.840 0,-0.170 0]T時的擾動解Fig. 6 For system(12) with b=3.12,d=3.12,c=2.9, a gait of perturbation solutions starting from the point [0.694 5,-0.627 6,-0.840 ,-0.170 0]T

圖7 以[0.694 5,-0.627 6,-0.904 0,-0.185 0]T為初始點,b=3.12,d=3.12,c=2.9 的雙足機器人上樓梯的棍棒圖Fig. 7 Stick diagram of ascending stair biped robot with b=3.12,d=3.12,c=2.9,and the initial point [0.694 5,-0.627 6,-0.904 0,-0.185 0]T

3.3 可調節脈沖力和固定脈沖力的性能比較

依然選取xd=[0.694 5,-0.627 6,-0.904 0,-0.185 0]T為初始狀態向量，其他初始參數不變，比較輸入不同種類的脈沖力對系統的影響。

圖8 系統(12) 分別在可調節脈沖(b=3.12,d=3.12,c=2.9)和固定脈沖pd=7.251作用下的效果Fig. 8 Effects of system (12) under the action of adjustable impluse (b=3.12,d=3.12,c=2.9) and fixed impluse pd=7.251, respectively

圖9 系統(12)分別在可調節脈沖(b=3.12,d=3.12,c=2.9)和固定脈沖pd=7.251作用下碰撞后支撐腿的角速度變化Fig. 9 Variation diagram of angular velocity of the standing leg just after the collision for system (12) under the action of adjustable impluse (b=3.12,d=3.12,c=2.9) and fixed impluse pd=7.251, respectively

3.4 周期步態的分岔

在環境參數或自身結構參數發生改變的情況下，當參數變化到一個臨界值時，系統變為不穩定，雙足機器人動態行走呈現倍周期步態或混沌步態。

3.4.1 特征值分布

基于2.3節中的理論，采用數值方法求解J的特征值，并分析其在復平面空間的分布情況。

圖10和圖11分別為穩定周期-1解和Flip分岔的特征值在復平面的空間分布圖。如圖10 所示，對3.2節中的不動點xd=[0.694 5,-0.627 6,-0.904 0,-0.185 0]T施加微小擾動，其雅可比矩陣J的所有特征值均在復平面單位圓內，即最大特征值的模小于1，所以不動點是穩定的，機器人具有穩定周期-1步態。圖11反映了Flip分岔現象的特征值分布，采用相同初始狀態施加微小擾動，其他參數固定不變，輸入脈沖參數為d=-3,c=5,b=- 0.99時系統達到臨界值。由圖11(b)可知，當b=-0.99 時，雅可比矩陣存在一個λ*≈-1，當b=-1時，有λ*<-1,此時，一個周期-2解從周期-1分岔出來。

圖10 系統(12)在b=3.12,d=3.12,c=2.9時周期軌道的所有特征值分布Fig. 10 Distribution diagram of all eigenvalues of the periodic orbit for system (12) with b=3.12,d=3.12,c=2.9

圖11 系統(12)在d=-3,c=5，b分別取-0.99和-1時的特征值分布Fig. 11 Eigenvalue distribution for system (12) with d=3.12,c=2.9, and b is taken as -0.99 and -1, respectively

3.4.2 脈沖參數分析

圖12 系統(12)在d=-3,c=5時的周期解的分岔Fig. 12 Bifurcation diagram of system (12) with d=-3,c=5

圖13 系統(12)在b=-1,c=5時周期解的分岔Fig. 13 Bifurcation diagram of system (12) with b=-1,c=5

3.4.3 髖關節質點參數分析

影響雙足機器人穩定步態的因素還有自身結構參數，接下來分析髖關節質點對系統穩定性的影響。將髖關節質量mH視為變量，固定其他參數(b=5，d=-7,c=5,m=0.5), 如圖14所示，隨著mH由10.8減小至約為10.02時，系統出現分岔。

圖14 系統(12)在b=5，d=-7,c=5,m=0.5時的周期解分岔，mH∈(9.6,10.8)Fig. 14 Bifurcation diagram of system (12) with b=5，d=-7,c=5,m=0.5, mH∈(9.6,10.8)

4 結語

本文研究脈沖推力作用下帶伸縮腿上樓梯雙足機器人行走的復雜動力學。引入一種與支撐腿角速度相關且存在幅值限制的脈沖推力，運用拉格朗日方程和角動量守恒定律建立上樓梯雙足機器人的動力學模型。利用龐加萊映射分析帶伸縮腿上樓梯雙足機器人行走的穩定性,討論脈沖參數b和d以及自身結構參數mH對周期-1步態存在、穩定以及分岔的影響。

理論分析和仿真實驗結論表明：當選擇合適的參數時，系統存在穩定的周期-1步態，即雙足機器人能穩定地上樓梯；與支撐腿角速度相關且存在幅值限制的脈沖推力優于固定的常值脈沖推力，能讓上樓梯雙足機器人快速進入穩定的行走狀態。這些結論為上樓梯雙足機器人實際行走的穩定與控制提供了理論依據。