李 琴，賈英崎，黃玉峰，李 剛，葉 闖

（1.西南石油大學(xué)機(jī)電工程學(xué)院，四川 成都 610500；2.中國(guó)石油天然氣集團(tuán)有限公司東方地球物理勘探有限責(zé)任公司，河北 涿州 072750）

近幾十年來(lái)，隨著自動(dòng)化技術(shù)的飛速發(fā)展，工業(yè)機(jī)器人在生產(chǎn)線上的應(yīng)用愈加廣泛。軌跡規(guī)劃是機(jī)器人運(yùn)動(dòng)控制的基礎(chǔ)，良好的軌跡規(guī)劃可以提高機(jī)器人的工作效率、降低機(jī)器人的能耗以及延長(zhǎng)機(jī)器人的使用壽命。

目前，針對(duì)機(jī)器人軌跡規(guī)劃優(yōu)化的研究主要集中于采用多種算法相結(jié)合或改進(jìn)原算法來(lái)獲取最優(yōu)的規(guī)劃軌跡。例如：居鶴華等［1］根據(jù)機(jī)械臂的運(yùn)動(dòng)學(xué)約束，采用遺傳算法對(duì)3-5-3多項(xiàng)式進(jìn)行軌跡優(yōu)化，但是3-5-3多項(xiàng)式存在加加速度突變的問(wèn)題，且傳統(tǒng)遺傳算法容易使解陷入局部最優(yōu)。Huang等［2］以5次非均勻有理B樣條（non-uniform rational B-splines，NURBS）曲線作為工業(yè)機(jī)器人的軌跡規(guī)劃曲線，并采用非支配排序遺傳算法-II（non-dominated sorting genetic algorithm-II，NSGA-II）對(duì)運(yùn)動(dòng)時(shí)間和沖擊磨損這2個(gè)目標(biāo)進(jìn)行優(yōu)化，結(jié)果表明基于5次NURBS曲線構(gòu)造的軌跡曲線較平滑，但未解決NSGA-II的不足，最終的結(jié)果仍易陷入局部最優(yōu)。Rout等［3］提出了一種成就標(biāo)量化函數(shù)與NSGA-II相結(jié)合的混合多目標(biāo)算法，并對(duì)六軸工業(yè)機(jī)器人進(jìn)行多目標(biāo)軌跡優(yōu)化，但采用的三次多項(xiàng)式并未實(shí)際解決機(jī)器人加速度、加加速度突變的問(wèn)題。

基于此，筆者針對(duì)文獻(xiàn)［1，3］中機(jī)器人加加速度不連續(xù)的問(wèn)題，以工業(yè)機(jī)器人為研究對(duì)象，采用5次NURBS曲線作為其軌跡規(guī)劃曲線，并結(jié)合布谷鳥(niǎo)搜索（cuckoo search，CS）算法中的Levy飛行策略與NSGA-II，提出了一種用于機(jī)器人多目標(biāo)軌跡優(yōu)化的混合算法（簡(jiǎn)稱為CSNSGA-II），并對(duì)初始種群進(jìn)行Tent混沌優(yōu)化處理，最終得到分布均勻的Pareto最優(yōu)解，旨在提高算法的收斂速度以及解決NSGA-II易陷入早熟的問(wèn)題。

1 工業(yè)機(jī)器人多目標(biāo)軌跡優(yōu)化模型

在工業(yè)生產(chǎn)中，機(jī)器人通常按提前預(yù)設(shè)的軌跡進(jìn)行移動(dòng)，通過(guò)求解逆運(yùn)動(dòng)學(xué)方程得到機(jī)器人各關(guān)節(jié)的軌跡點(diǎn)位置矢量pi=[pi，1pi，2…pi，N]（i=0，1，…，n；3≤N≤6）。其中：i為關(guān)節(jié)的軌跡點(diǎn)，共n+1個(gè)；N為關(guān)節(jié)數(shù)，即自由度數(shù)。各關(guān)節(jié)的軌跡點(diǎn)構(gòu)成了工業(yè)機(jī)器人的運(yùn)行軌跡。

通常情況下，工業(yè)機(jī)器人在作業(yè)過(guò)程中大多會(huì)重復(fù)同一種動(dòng)作。因此，為提高機(jī)器人的工作效率、降低機(jī)器人的能耗以及延長(zhǎng)機(jī)器人的使用壽命，基于其運(yùn)動(dòng)學(xué)約束（速度、加速度和加加速度約束），定義目標(biāo)函數(shù)［4］，可表示為：

其中，約束條件為：

式中：S1(T)為機(jī)器人的總運(yùn)動(dòng)時(shí)間；S2(A)為機(jī)器人各關(guān)節(jié)加速度的平均值，用于衡量關(guān)節(jié)的能耗損失；S3(J)為機(jī)器人各關(guān)節(jié)加加速度的平均值，用于衡量關(guān)節(jié)的沖擊損耗；ti為經(jīng)過(guò)第i個(gè)軌跡點(diǎn)時(shí)所對(duì)應(yīng)的時(shí)間；am、jm分別為機(jī)器人第m個(gè)關(guān)節(jié)的加速度和加加速度；vi，m、ai，m和ji，m分別為第i個(gè)軌跡點(diǎn)處機(jī)器人第m個(gè)關(guān)節(jié)的速度、加速度和加加速度；vmax、amax和jmax分別為機(jī)器人關(guān)節(jié)的速度、加速度和加加速度的最大值。

2 基于NURBS曲線的工業(yè)機(jī)器人軌跡規(guī)劃

2.1 NURBS曲線的數(shù)學(xué)表示

根據(jù)工業(yè)機(jī)器人某一關(guān)節(jié)經(jīng)過(guò)的軌跡點(diǎn)位置pi（i=0，1，…，n），采用k次NURBS曲線來(lái)規(guī)劃其軌跡曲線。k次NURBS曲線可表示為分段有理多項(xiàng)式函數(shù)形式：

式中：P(u)為分段有理多項(xiàng)式函數(shù)；ωi為NURBS曲線的權(quán)因子；di為NURBS曲線的控制頂點(diǎn)位置；Ni，k(u)為k次NURBS曲線的B樣條有理基函數(shù)。

Ni，k(u)是由節(jié)點(diǎn)矢量U=[u0u1…un+2k]按Cox De Boor遞推公式推導(dǎo)得到的（規(guī)定0/0=0），可表示為：

設(shè)給定的軌跡點(diǎn)共有n+1個(gè)，為使機(jī)器人各關(guān)節(jié)經(jīng)過(guò)給定的軌跡點(diǎn)，使k次NURBS曲線各分段的連接點(diǎn)與軌跡點(diǎn)一一對(duì)應(yīng)，則需反求NURBS曲線的控制頂點(diǎn)。k次NURBS曲線由n+k個(gè)控制頂點(diǎn)與n+2k+1個(gè)節(jié)點(diǎn)組成。對(duì)于節(jié)點(diǎn)矢量U=[u0u1…un+2k]，本文采用規(guī)范節(jié)點(diǎn)矢量，即令u0=u1=…=uk=0，un+k=un+k+1=…=un+2k=1，其他節(jié)點(diǎn)均采用積累弦長(zhǎng)參數(shù)化法對(duì)時(shí)間間隔矢量進(jìn)行歸一化求解獲得：

式中：Δtj為軌跡點(diǎn)間的時(shí)間間隔。

2.2 5次NURBS曲線的控制頂點(diǎn)反算

由于3次NURBS曲線存在加加速度不連續(xù)變化的情況，本文采用5次NURBS曲線作為工業(yè)機(jī)器人的軌跡規(guī)劃曲線。但是，5次NURBS曲線的插值數(shù)據(jù)量龐大，因此本文主要推導(dǎo)NURBS曲線的矩陣表達(dá)式，以簡(jiǎn)化計(jì)算過(guò)程。

基于上文推導(dǎo)的節(jié)點(diǎn)矢量U=[u0u1…un+2k]和軌跡點(diǎn)位置pi(i=0，1，…，n)，根據(jù)節(jié)點(diǎn)矢量、軌跡點(diǎn)與控制頂點(diǎn)的關(guān)系，反求5次NURBS曲線的控制頂點(diǎn)。由B樣條曲線的局部支撐性可知，5次NURBS曲線在節(jié)點(diǎn)區(qū)間ui≤u≤ui+1內(nèi)可表示為：

為方便計(jì)算，令第i（i≠0）段曲線為：

由Cox De Boor遞推公式推導(dǎo)得到Ni中第1行的元素，分別為：

規(guī)定：

根據(jù)5次NURBS曲線的插值要求，可得：

其中，對(duì)于Qi(t)，當(dāng)u=ui時(shí)，t=0；u=ui+1時(shí)，t=1。聯(lián)立式（5）和式（6）可得：

令 ：ai+2=n11ωi，bi+2=n12ωi+1，ci+2=n13ωi+2；ei+2=n14ωi+3，fi+2=n15ωi+4；qi=(ai+2+bi+2+ci+2+ei+2+fi+2)pi，則式（7）可以簡(jiǎn)化為：

由式（8）構(gòu)成的矩陣可以得出n+1個(gè)線性方程組，而所求的控制頂點(diǎn)為n+5個(gè)，根據(jù)線性方程組有唯一解的條件，還須根據(jù)邊界條件確定4個(gè)附加方程。指定啟、停速度為vs、ve以及啟、停加速度為as、ae，并以啟、停位置處5次NURBS曲線的一階和二階導(dǎo)數(shù)作為邊界條件構(gòu)建方程。

k次NURBS曲線的求導(dǎo)公式為：

對(duì)5次NURBS曲線進(jìn)行一階、二階求導(dǎo)，可得啟、停時(shí)刻的速度vs、ve和加速度as、ae：

由于直接求導(dǎo)法的計(jì)算量較大且過(guò)程復(fù)雜，本文通過(guò)德布爾遞推算法計(jì)算各階導(dǎo)數(shù)的控制點(diǎn)。NURBS曲線上某一點(diǎn)的r階導(dǎo)數(shù)P(r)(u)為：

結(jié)合式（10）和式（11），對(duì)上文的n+5個(gè)線性方程進(jìn)行求解，由此得到5次NURBS曲線控制頂點(diǎn)的位置矢量：

3 CSNSGA-II構(gòu)建

NSGA-II是目前在多目標(biāo)優(yōu)化問(wèn)題中應(yīng)用最為廣泛的算法之一，其是在NSGA的基礎(chǔ)上引入快速非支配排序方法，并結(jié)合精英策略，利用其對(duì)某些測(cè)試函數(shù)進(jìn)行優(yōu)化可收斂到真正的最優(yōu)Pareto前沿。但隨著多目標(biāo)優(yōu)化問(wèn)題的復(fù)雜性增加，傳統(tǒng)的NSGA-II出現(xiàn)了不足之處，主要體現(xiàn)在以下2個(gè)方面［5］：1）易出現(xiàn)局部早熟的情況；2）搜索后期會(huì)出現(xiàn)局部搜索能力不足等。因此，須對(duì)NSGA-II進(jìn)行改進(jìn)，以獲取更加真實(shí)且均勻的Pareto前沿。

CS算法是一種新型的啟發(fā)式算法，由Yang和Deb于2009年提出［6］。研究表明，CS算法具有局部和全局搜索性能，基于其獲得的最優(yōu)解大多優(yōu)于遺傳算法和粒子群優(yōu)化算法［7］，且已有較多的學(xué)者將CS算法應(yīng)用于機(jī)器人軌跡優(yōu)化［8-10］。基于此，提出了一種改進(jìn)CS算法與NSGA-II相結(jié)合的CSNSGA-II，其中CS算法的Levy飛行策略可以幫助NSGA-II跳出局部最優(yōu)，找到真實(shí)的Pareto前沿。

3.1 基于Tent混沌映射的種群初始化

為提高NSGA-II對(duì)工業(yè)機(jī)器人軌跡的優(yōu)化性能，采用混沌算法對(duì)種群進(jìn)行初始化。文獻(xiàn)［4］采用一維logistic映射對(duì)種群進(jìn)行初始化，但是一維logistic映射并不具有最好的遍歷性，因此本文引入Tent混沌映射來(lái)產(chǎn)生初始種群。Tent混沌映射相比于logistic映射具有更好的均勻性與隨機(jī)性［11］，可以提高算法的全局尋優(yōu)能力。Tent混沌映射的表達(dá)式如下：

式中：β為混沌系數(shù)，取β=0.5；γ=1，2，…，v，v為決策變量數(shù)。

將混沌種群映射到對(duì)應(yīng)的決策變量空間（tγmin，tγmax）中，得到初始種群的變量：

3.2 約束不可行度處理

在NSGA-II處理約束優(yōu)化問(wèn)題時(shí)，常用到對(duì)不可行度的處理，以保留與可行解接近的解，舍棄與可行解相差較大的解。為了使不等式約束轉(zhuǎn)化為等式約束［12］，定義不可行度Cx（x=1，2，…，G，G為當(dāng)前種群中解的數(shù)量）：

式中：gy、lz分別為不等式約束和等式約束；L、I分別為多目標(biāo)優(yōu)化模型中不等式和等式約束的數(shù)量。

不可行度Cx表示超過(guò)約束條件的距離，當(dāng)個(gè)體在約束條件內(nèi)時(shí)，Cx=0。

同時(shí)定義約束違規(guī)閾值Cxmean：

超過(guò)約束違規(guī)閾值的解定義為不可行解，滿足約束違規(guī)閾值的解定義為可行解。

3.3 自適應(yīng)Levy飛行策略

傳統(tǒng)CS算法中的Levy飛行策略存在種群多樣性不易調(diào)節(jié)、搜索方式不能合理控制等問(wèn)題，導(dǎo)致該算法的精度不高及收斂速度較慢，因此本文提出了一種自適應(yīng)Levy飛行策略，其表達(dá)式如下：

式中：Xs，δ為第s只布谷鳥(niǎo)的第δ代鳥(niǎo)巢的位置；τ為動(dòng)態(tài)慣性權(quán)重；α為變步長(zhǎng)因子；⊕表示矩陣之間的各元素對(duì)應(yīng)相乘，Levy(λ)為L(zhǎng)evy隨機(jī)搜索策略。

Levy隨機(jī)搜索策略可表示為：

式中：U、V服從標(biāo)準(zhǔn)正態(tài)分布，μ=1.5。

動(dòng)態(tài)慣性權(quán)重τ越大，CS算法的全局搜索能力越好；相反，τ越小全局搜索能力越差。因此采用自動(dòng)調(diào)整慣性權(quán)重的方法，以使算法初期即具有較好的搜索能力。τ的表達(dá)式為：

式中：τmax為慣性權(quán)重的最大值；τmin為慣性權(quán)重的最小值；N()為服從(0，1)上均勻分布的隨機(jī)數(shù)。

變步長(zhǎng)因子α可表示為：

其中：

式中：Xs，best為隨機(jī)選取的當(dāng)前Pareto前沿的最優(yōu)解；δmax為最大迭代次數(shù)。

3.4 CSNSGA-II流程

本文所提出的CSNSGA-II的流程如圖1所示，其中遺傳操作中的交叉算子采用獨(dú)立點(diǎn)交叉算子，變異算子采用多項(xiàng)式變異算子。

圖1 CSNSGA-II流程Fig.1 Flow chart of CSNSGA-II

由圖1可知，CSNSGA-II的具體步驟為：

步驟1 設(shè)定交叉、變異系數(shù)。采用Tent混沌映射生成初始化時(shí)間種群，并計(jì)算種群的適應(yīng)度函數(shù)，在進(jìn)行非支配排序后對(duì)排序后的解集進(jìn)行擁擠距離計(jì)算，并按照擁擠距離進(jìn)行排序。

步驟2 對(duì)排序后的所有解進(jìn)行不可行度計(jì)算，得到可行解集W1和不可行解集V1。

步驟3 對(duì)可行解集W1進(jìn)行二聯(lián)錦標(biāo)賽選擇，對(duì)篩選的解進(jìn)行下一代的交叉、變異操作，最終得到可行解集W2。

步驟4 采用自適應(yīng)Levy飛行策略對(duì)不可行解集V1進(jìn)行更新，得到新的解集V2。

步驟5 將解集W2與解集V2合并，并對(duì)新解集進(jìn)行適應(yīng)度計(jì)算、非支配排序和擁擠距離計(jì)算。

步驟6 重復(fù)步驟2、3和4，直到滿足結(jié)束條件。

4 6R勃朗特機(jī)器人多目標(biāo)軌跡優(yōu)化仿真分析

4.1 仿真參數(shù)設(shè)置

在MATLAB軟件中對(duì)6R勃朗特機(jī)器人（見(jiàn)圖2）進(jìn)行運(yùn)動(dòng)學(xué)建模，并采用本文提出的CSNSGA-II對(duì)其軌跡進(jìn)行優(yōu)化。設(shè)置6R勃朗特機(jī)器人啟、停時(shí)刻的位置（笛卡爾坐標(biāo)系中），并采用5次NURBS曲線構(gòu)建其關(guān)節(jié)的軌跡，同時(shí)定義啟、停時(shí)刻的速度、加速度均為0，指定5次NURBS曲線的各項(xiàng)權(quán)因子均為1，初始化種群數(shù)量為200，最大迭代次數(shù)為100次。根據(jù)經(jīng)驗(yàn)，遺傳算法的交叉概率設(shè)為0.6，變異概率設(shè)為0.1。6R勃朗特機(jī)器人各關(guān)節(jié)軌跡點(diǎn)的位置如表1所示，基于MATLAB軟件可得其啟、停時(shí)刻位置分別如圖3和圖4所示，各關(guān)節(jié)的運(yùn)動(dòng)學(xué)約束如表2所示。

表1 6R勃朗特機(jī)器人各關(guān)節(jié)軌跡點(diǎn)的位置Table 1 Position of each joint trajectory point of 6R Bronte robot 單位：（°）

表2 6R勃朗特機(jī)器人各關(guān)節(jié)的運(yùn)動(dòng)學(xué)約束Table 2 Kinematic constraints of each joint of 6R Bronte robot

圖2 6R勃朗特機(jī)器人實(shí)物Fig.2 Physical object of 6R Bronte robot

圖3 6R勃朗特機(jī)器人啟動(dòng)時(shí)刻的位置Fig.3 Position of 6R Bronte robot at the moment of start

圖4 6R勃朗特機(jī)器人停止時(shí)刻的位置Fig.4 Position of 6R Bronte robot at the moment of stop

4.2 仿真結(jié)果分析

圖5所示為分別采用CSNSGA-II、NSGA-II和多目標(biāo)粒子群優(yōu)化（multi-objective particle swarm optimization，MOPSO）算法對(duì)6R勃朗特機(jī)器人進(jìn)行軌跡優(yōu)化后得到的Pareto前沿分布結(jié)果。從圖中可以看出：相比于NSGA-II、MOPSO算法，基于CSNSGA-II得到的Pareto前沿分布得更加均勻，且更加靠近真實(shí)的Pareto前沿。

圖5 基于不同算法的6R勃朗特機(jī)器人Pareto前沿分布對(duì)比Fig.5 Comparison of Pareto front distribution of 6R Bronte robot based on different algorithms

為了進(jìn)一步驗(yàn)證CSNSGA-II的優(yōu)勢(shì)，在基于3種算法優(yōu)化后獲得的6R勃朗特機(jī)器人Pareto前沿中分別取3組解，結(jié)果如表3所示。其中：A1、B1、C1為CSNSGA-II優(yōu)化后得到的3組解；A2、B2、C2為NSGA-II優(yōu)化后得到的3組解；A3、B3、C3為MOPSO算法優(yōu)化后得到的3組解。通過(guò)對(duì)比分析可得出，解A1支配解A2、A3，解B1支配解B2、B3，解C1、C2、C3互為非支配解。結(jié)果表明，相比于NSGA-II、MOPSO算法，基于CSNSGA-II能得到更優(yōu)的Pareto前沿。

表3 基于不同算法的6R勃朗特機(jī)器人軌跡優(yōu)化結(jié)果對(duì)比Table 3 Comparison of trajectory optimization results of 6R Bronte robot based on different algorithms

為了進(jìn)一步分析CSNSGA-II的收斂性能，給定一個(gè)優(yōu)化目標(biāo)綜合函數(shù)S：

式中：ψ1、ψ2、ψ3為加權(quán)因子，在實(shí)際生產(chǎn)應(yīng)用中，根據(jù)用戶的需求可以確定各加權(quán)因子的系數(shù)。在本文中，為了對(duì)比算法收斂的速度和收斂性能，取ψ1=ψ2=ψ3=1，η1=1，η2=100，η3=100。

利用式（21）計(jì)算得到基于不同算法的6R勃朗特機(jī)器人軌跡優(yōu)化目標(biāo)的歸一化加權(quán)迭代最優(yōu)值，對(duì)比結(jié)果如圖6所示。圖6結(jié)果顯示，基于CSNSGA-II得到的機(jī)器人軌跡優(yōu)化目標(biāo)的最優(yōu)值優(yōu)于基于NSGA-II和MOPSO算法的，由此說(shuō)明CSNSGA-II有助于跳出局部最優(yōu)，且CSNSGA-II達(dá)到穩(wěn)定收斂時(shí)的迭代次數(shù)比NSGA-II和MOPSO算法分別減少了63.16%和59.62%，有效提高了優(yōu)化效率。

圖6 基于不同算法的6R勃朗特機(jī)器人軌跡優(yōu)化目標(biāo)歸一化加權(quán)迭代最優(yōu)值對(duì)比Fig.6 Comparison of normalized weighted iterative optimal values of trajectory optimization objectives of 6R Bronte robot based on different algorithms

取解B1對(duì)應(yīng)的時(shí)間間隔矢量Δt=［0 5.059 7 3.231 5 1.609 9 1.309 8 2.495 2 1.988 1 3.686 2］，其對(duì)應(yīng)5次NURBS曲線的節(jié)點(diǎn)矢量U=［0000 0 0 0.261 0 0.427 8 0.510 9 0.578 4 0.707 2 0.809 7 1 1 1 1 1 1］，其對(duì)應(yīng)的各關(guān)節(jié)的位置、速度、加速度和加加速度隨時(shí)間的變化曲線如圖7至圖10所示。從圖中可以看出，利用CSNSGA-II優(yōu)化后，6R勃朗特機(jī)器人各關(guān)節(jié)的運(yùn)動(dòng)軌跡平滑且無(wú)沖擊，相比于文獻(xiàn)［13-15］所得結(jié)果更具有平穩(wěn)性。

圖7 6R勃朗特機(jī)器人各關(guān)節(jié)的位置—時(shí)間曲線Fig.7 Position-time curve of each joint of 6R Bronte robot

圖8 6R勃朗特機(jī)器人各關(guān)節(jié)的速度—時(shí)間曲線Fig.8 Velocity-time curve of each joint of 6R Bronte robot

圖9 6R勃朗特機(jī)器人各關(guān)節(jié)的加速度—時(shí)間曲線Fig.9 Acceleration-time curve of each joint of 6R Bronte robot

圖10 6R勃朗特機(jī)器人各關(guān)節(jié)的加加速度—時(shí)間曲線Fig.10 Jerk-time curve of each joint of 6R Bronte robot

5 總 結(jié)

本文采用5次NURBS曲線作為工業(yè)機(jī)器人的軌跡規(guī)劃曲線，以運(yùn)動(dòng)時(shí)間、能耗、沖擊磨損這3個(gè)指標(biāo)作為優(yōu)化目標(biāo)，并在運(yùn)動(dòng)學(xué)約束下進(jìn)行優(yōu)化求解。提出了一種基于改進(jìn)CS算法與NSGA-II結(jié)合的混合算法——CSNSGA-II。通過(guò)對(duì)6R勃朗特機(jī)器人的軌跡進(jìn)行優(yōu)化后可知，CSNSGA-II能成功地跳出局部最優(yōu)，解決了NSGA-II易陷入“早熟”的問(wèn)題，同時(shí)基于該算法得到的Pareto前沿分布得更加均勻，有效提高了收斂速度。結(jié)果表明，所提出的CSNSGA-II具有較高的實(shí)用價(jià)值，可為工業(yè)機(jī)器人高效、持久、可靠地作業(yè)提供指導(dǎo)。