王明海,張威,劉香辰
(1.沈陽航空航天大學(xué)機電工程學(xué)院,遼寧沈陽 110136;2.沈陽航空航天大學(xué)航空制造工藝數(shù)字化國防重點學(xué)科實驗室,遼寧沈陽 110136;3.沈陽航空航天大學(xué)航空宇航學(xué)院,遼寧沈陽 110136)
近年來,機械臂制孔在航空工業(yè)領(lǐng)域的應(yīng)用越來越普遍,隨著機器人技術(shù)的不斷進步,對其加工精度也提出了更高的要求。然而工業(yè)機器人固有的弱剛性卻嚴(yán)重阻礙了行業(yè)的發(fā)展,故提升其作業(yè)剛度、保證工件加工質(zhì)量很必要。針對此問題,國內(nèi)外學(xué)者已經(jīng)在工業(yè)機器人剛性方面做過相關(guān)探索,CVITANIC等[1]利用靜態(tài)和動態(tài)剛度模型對不同切割場景下的機器人姿態(tài)優(yōu)化進行了比較研究,結(jié)果表明優(yōu)化可減少執(zhí)行器偏轉(zhuǎn),但不能使其回落到可重復(fù)性范圍內(nèi)。華中科技大學(xué)的楊靖等人[2]提出一種對銑削機器人主剛度定向的方法,深入研究機械臂的冗余度,改進加工位姿以減少銑削過程的顫振。XIONG等[3]提出機械臂冗余度,采用離散搜索方法將五軸刀具數(shù)控路徑轉(zhuǎn)換為六軸工業(yè)機器人軌跡進行優(yōu)化。陳欽韜等[4]將銑削力橢圓平面的各向同性度作為評價指標(biāo),運用遺傳算法進行位姿優(yōu)化,使被加工件銑削平面度提升45%。諸多學(xué)者對工業(yè)機器人奇異問題、剛度辨識、誤差補償?shù)榷歼M行了相關(guān)理論研究,但針對鉆鉚機器人實際工況下剛度優(yōu)化探索較少。本文作者以位姿的調(diào)整為切入點,結(jié)合實際,利用改進后的智能優(yōu)化算法,對串聯(lián)機器人重載下的整體剛度進行深入研究,并證實了其有效性。
工業(yè)機器人運動學(xué)模型可以清晰表示工業(yè)機器人加工時所處位姿,反映機械臂連桿間的幾何關(guān)系。本文作者以KUKA KR600 R2830型工業(yè)機器人為例進行研究,該型號機器人本體質(zhì)量約2 650 kg,額定負(fù)載為6 000 N,負(fù)載高、工作范圍大,姿態(tài)重復(fù)精確度可達±0.08 mm,使其適用于高精密、大負(fù)載條件下制孔加工任務(wù)。機器人整體結(jié)構(gòu)由6個旋轉(zhuǎn)關(guān)節(jié)串聯(lián)而成,稱為串聯(lián)工業(yè)六自由度機器人,可以采用改進的D-H法[5]來建立機器人運動學(xué)正解,對機器人參數(shù)化。圖1所示為機器人在進行板類部件鉆鉚作業(yè)時的場景。
圖1 工業(yè)機器人鉆鉚
依照工業(yè)機器人作業(yè)位姿在其關(guān)節(jié)連接部位創(chuàng)建笛卡爾坐標(biāo)系,根據(jù)連桿數(shù)據(jù)建立運動學(xué)數(shù)學(xué)模型,如圖2所示。參數(shù)見表1。
表1 機器人D-H參數(shù)
圖2 運動學(xué)模型
模型中的4個關(guān)鍵參數(shù)為ai-1、αi-1、di和θi,分別是連桿長度、連桿扭轉(zhuǎn)角、連桿偏距和關(guān)節(jié)角。坐標(biāo)系{i}可以由坐標(biāo)系{i-1}經(jīng)過4次旋轉(zhuǎn)或平移運動變換得到,將4步轉(zhuǎn)換矩陣連乘可以獲得兩相鄰連桿間的變換關(guān)系,可寫為
i-1iT=R(X,αi-1)T(X,αi-1)R(Z,θi)
T(Z,di)=
(1)
雅可比矩陣的數(shù)學(xué)意義是表示從機器人關(guān)節(jié)速度到機器人末端操作速度的一種映射關(guān)系[6],需要注意的是該廣義傳動比非定值。雅可比矩陣的求解方法通常有矢量積法[7]和微分變換法,本文作者通過對比兩種方式編程及結(jié)果的復(fù)雜度優(yōu)先選用矢量積法進行計算。假設(shè)機器人末端線速度和角速度分別用v和ω表示,則關(guān)節(jié)速度與末端速度的關(guān)系如式(2)所示:
(2)
(3)
綜上得到KUKA機器人的雅可比矩陣公式為
(4)
由公式(4)可以在MATLAB中編程求解R2830型機器人雅可比矩陣。
剛度反映的是機械臂在外力作用下抵抗變形的能力,關(guān)節(jié)剛度、連桿剛度以及所處位姿都影響機器人末端操作剛度[8]。在大負(fù)載條件下,由連桿變形所造成的末端變形遠(yuǎn)小于關(guān)節(jié)所造成的末端變形,為簡化計算可將機器人連桿視為剛性,在此條件下機械臂末端操作剛度僅與各關(guān)節(jié)剛度和所處位姿相關(guān)。在這一假設(shè)下,討論機器人關(guān)節(jié)剛度與末端操作剛度之間的關(guān)系,討論機器人末端變形與末端受力之間的關(guān)系,建立數(shù)學(xué)模型,從而得到機器人的末端操作剛度模型為
F=J-TKθJ-1X
(5)
其中:F為末端力矢量;J為雅可比矩陣;Kθ為關(guān)節(jié)剛度矩陣;X為末端變形矢量。將式(5)變形處理可得:
(6)
整理后得:
(7)
針對工業(yè)機器人關(guān)節(jié)剛度獲取,主要有以下兩種方法:(1)通過對關(guān)節(jié)處伺服電機、減速器等部件的結(jié)構(gòu)、材料、規(guī)格等數(shù)據(jù)進行計算[9],求其等效剛度值。(2)設(shè)計關(guān)節(jié)剛度辨識實驗[10],處理實驗數(shù)據(jù)計算剛度值。其中方法一過于理想化,機器人關(guān)節(jié)部件的各項參數(shù)很難獲得且計算繁瑣,故選擇方法二求解工業(yè)機器人關(guān)節(jié)剛度。
式(6)為關(guān)節(jié)剛度識別實驗的原理,加工時由扭矩造成的旋轉(zhuǎn)位移較小,可忽略不計,故辨識實驗只取變形的前3項即末端平移位移[dx,dy,dz]T進行計算[11]。辨識實驗所搭建的實驗平臺如圖3所示,其中激光跟蹤儀型號為Radian Core(50),測力儀為瑞士KISTLER及配套Dynoware軟件。實驗中激光跟蹤儀和測力儀所測得的數(shù)據(jù)均是在自身坐標(biāo)系原點處,故先將激光跟蹤儀坐標(biāo)系轉(zhuǎn)換至機器人基坐標(biāo)系下,測力儀還需二次轉(zhuǎn)換至機器人末端法蘭中心處進行標(biāo)定。
圖3 關(guān)節(jié)剛度辨識實驗平臺
通過測量工業(yè)機器人在多個位姿時的末端變形及力矢量代入式(7),利用最小二乘法[12]在MATLAB中進行編程可以求得機器人關(guān)節(jié)剛度矩陣Kθ。求解出的工業(yè)機器人各關(guān)節(jié)剛度如表2所示。
表2 工業(yè)機器人關(guān)節(jié)剛度 單位:N·m·rad-1
機器人末端操作剛度Kx為六階矩陣,是一張量測度[13],故需一準(zhǔn)確值來表述機器人末端的剛度特性。工業(yè)機器人剛度性能評價指標(biāo)[14]有3種,分別是矩陣瑞利熵、變形橢球(剛度橢球)和力橢球(柔度橢球)[14]。在使用力橢球進行計算時,無需計算雅可比矩陣的逆,計算較為簡便且易編程,故選擇力橢球作為剛度性能評價指標(biāo)。將式(5)展開并經(jīng)變形處理可得:
(8)
由于機器人末端在外力下的扭轉(zhuǎn)變形量遠(yuǎn)小于位移變形量且無法精確測得,故省略扭轉(zhuǎn)變形[15],將δ設(shè)為0,由式(8)可得:
d=cfd·f
(9)
其中:
cfd=
(10)
假設(shè)機械臂法蘭端受到外力作用時造成單位變形,即d=1,可得:
(11)
通過分析可知,式(11)表示一橢球體,其方程可以寫為
(12)
借助蒙特卡洛法可以求解工業(yè)機器人的工作空間及空間內(nèi)的剛度性能[16],原理是利用機器人末端所能達到所有點的集合來表征機器人可達空間,如圖4所示。
圖4 機器人工作空間云圖
結(jié)合柔度矩陣及剛度性能評價指標(biāo),在MATLAB中編程仿真可以得到機器人工作空間剛度性能云圖如圖5所示,圖中藍色較深區(qū)域工業(yè)機器人剛度性能指標(biāo)較優(yōu),故所處位姿末端剛度較高,故機器人在加工時應(yīng)盡量選擇末端處于藍色區(qū)域的位姿進行加工。
圖5 剛度性能云圖(xOz向)
尋求最優(yōu)解即求函數(shù)極值的問題,重點在于3個關(guān)鍵條件的建立,即數(shù)學(xué)模型的目標(biāo)函數(shù)、設(shè)計變量及約束條件。下面以工業(yè)機器人作業(yè)時剛度性能最大化為目標(biāo)建立數(shù)學(xué)模型,解決此優(yōu)化問題的基本思想是通過整個可達空間中各個位姿之間的協(xié)作和信息共享來尋找最優(yōu)位姿。選用粒子群算法的原因在于其思想簡單,求解過程迅速。建立目標(biāo)函數(shù)時,選擇柔度矩陣cfd的最大特征值λmax作為剛度性能評價指標(biāo)k,λmax的值越小,機器人末端產(chǎn)生單位變形所需最小外力越大,機器人剛度性能越好,反之剛度性能越差。
機器人所處位姿與關(guān)節(jié)角度相關(guān),故變量為關(guān)節(jié)角θ。其中θ1及θ6對工業(yè)機器人末端操作剛度無影響[17],θ2及θ3主要影響機械臂末端位置,θ4及θ5主要影響末端方向。在自動化鉆鉚系統(tǒng)中,工業(yè)機器人多配備滑軌進行作業(yè),θ4對末端位置及方向的影響可以通過機械臂在滑軌上的平移替代。由此粒子群優(yōu)化算法中自變量的個數(shù)可以減少為3個,可加快算法求解速度。
將關(guān)節(jié)角上下限作為界限約束,規(guī)定變量取值范圍。實際工況下機械臂末端刀具頂點C需與工件作業(yè)點P重合,令點C在基坐標(biāo)系{B}中的表達為[Cx,Cy,Cz]T,點P在基坐標(biāo)系{B}中的表達為[Px,Py,Pz]T,則需保證:
ΔH1=(Cx-Px)2+(Cy-Py)2+(Cz-Pz)2
(13)
其次刀具軸線與目標(biāo)點法方向需一致,故刀具坐標(biāo)系{C}相對于基坐標(biāo)系{B}與目標(biāo)點坐標(biāo)系{P}相對于基坐標(biāo)系{B}的位姿偏差[Δax,Δay,Δaz]T為0,如下:
ΔH2=(Δax)2+(Δay)2+(Δaz)2
(14)
故可將工作中末端姿態(tài)偏差H(θ)作為性態(tài)約束,可得機器人數(shù)學(xué)模型為
(15)
粒子群算法(PSO)操作簡單、收斂速度快,在函數(shù)優(yōu)化中全局最優(yōu)解搜索能力較強,故選擇粒子群算法進行計算。算法具體流程如圖6所示。
圖6 粒子群算法流程
優(yōu)化參數(shù)中迭代次數(shù)設(shè)為1 000,自變量個數(shù)n=3,設(shè)置群體粒子個數(shù)N=100并隨機初始化每個個體的速度和位置。定義適應(yīng)度函數(shù),個體極值是第i個粒子搜尋到的最優(yōu)位置,全部個體極值的最優(yōu)解為此次全局最優(yōu)解,結(jié)果與歷史全局最優(yōu)解進行比對,然后更新。此次計算ω取0.8,C1和C2取2。粒子速度、位置迭代的公式為
v(i+1)i=ωvid+C1rand(0,1)(Pid-xid)+C2rand(0,1)(Pgd-xid)
(16)
x(i+1)d=xid+vid
(17)
設(shè)計好算法相關(guān)參數(shù)和流程后,可以在MATLAB中編寫程序求解,最終得到最優(yōu)位姿。為了驗證優(yōu)化結(jié)果的正確性,參照工業(yè)機器人鉆鉚作業(yè)實況,隨機選取一組機器人關(guān)節(jié)角結(jié)合上述優(yōu)化算法進行求解,得到優(yōu)化后的工業(yè)機器人加工位姿(關(guān)節(jié)角配置),迭代結(jié)果如圖7所示。
圖7 位姿優(yōu)化迭代過程
根據(jù)鉆鉚工業(yè)機器人實際工況隨機選擇一組關(guān)節(jié)角(見表3)與算法優(yōu)化后的關(guān)節(jié)角進行末端理論變形量的對比,將兩組關(guān)節(jié)角分別代入式(9)并假設(shè)在末端施加z軸方向的力F=100~1 100 N,圖8為優(yōu)化前后末端變形量的對比??梢园l(fā)現(xiàn):優(yōu)化后的末端變形量更小,效果明顯,優(yōu)化程度在60%以上且精度可以控制在0.03 mm以內(nèi),驗證了改變加工位姿對提升工業(yè)機器人末端剛度的有效性,且可使加工孔位精度小于0.5 mm,法向精度小于0.5°,滿足了飛機對制孔位置的要求。
表3 優(yōu)化前后關(guān)節(jié)角及剛度參數(shù)
圖8 位姿優(yōu)化前、后變形量
由表3可知KUKA KR600 R2830型串聯(lián)機械臂剛度最大時關(guān)節(jié)角的準(zhǔn)確值,可為機器人作業(yè)時加工位姿的選擇提供參考。在SolidWorks中繪制工業(yè)機器人模型并導(dǎo)出到MATLAB中,可清晰表達工業(yè)機器人剛度最優(yōu)的位姿,如圖9所示。
通過對六自由度機械手進行運動學(xué)分析、工作空間分析及剛度分析,采用參數(shù)優(yōu)化后的粒子群算法來尋找機器人剛度最大的加工位姿,減小機械臂末端變形,優(yōu)化機器人的加工質(zhì)量,保證工業(yè)機器人在高精度、大負(fù)載條件下的加工剛度,而且考慮實際工況改良的優(yōu)化方法適用于同類型各種型號的串聯(lián)機器人。在航空航天部件鉆鉚加工中,面對配備大負(fù)載末端執(zhí)行器或加工合金或復(fù)合材料工況時,通過調(diào)整工業(yè)機器人的位姿來提高剛度,滿足加工精度要求,既節(jié)省資源又有良好的優(yōu)化作用。其中的分析過程可為后續(xù)研究提供基礎(chǔ)。