零耦合度空间2T1R并联机构运动学与刚度建模分析

沈惠平 朱忠颀 孟庆梅 吴广磊 邓嘉鸣

(1.常州大学现代机构学研究中心， 常州 213016； 2.大连理工大学机械工程学院， 大连 116024)

0 引言

三自由度两平移一转动并联机器人机构驱动元件少、结构紧凑、制造方便，在工业生产中具有较高实用价值[1]，其设计需考虑机构的运动学、刚度等因素，国内外学者对此进行了大量研究。GOSSELIN[2]提出了并联机构的刚度映射矩阵,并分析了3-RPR平面并联机构和6-SPS空间并联机构在不同参数下机构的刚度特性；李树军等[3]研究了基于新的守恒协调转换刚度矩阵，对3-RRR平面并联机构进行刚度特性分析；周玉林等[4]基于小变形叠加原理推导得到动平台、球心点位移与外力的关系，通过静力学分析得到机构的整体刚度矩阵；CARBONE等[5]对多个刚度性能评价指标进行了比较；WU等[6-7]对3-PPR并联机械手进行了刚度分析，并用数值法将笛卡尔刚度矩阵解耦为平移刚度矩阵和旋转刚度矩阵，基于虚拟弹簧法对3、4自由度且含有4条相同支链并联机器人的刚度特性进行了分析，指出Ragnar机器人的结构刚度低于以Quattro结构为基础设计的其他机器人；YANG等[8]基于应变能理论建立了并联机器人弹性静态刚度分析模型；项超群等[9]利用最小二乘法建立了单根气动肌肉气压、位移及刚度的关系模型。

并联机构刚度分析主要有以下几种方法：有限元分析法(FEA)[10-11]、螺旋理论分析法[12]、矩阵结构分析法(MSA)[13]、虚拟关节分析法(VJM)[14-15]等。FEA法的优势在于其杆件和关节建模都具有精确的物理模型[12，16]；MSA法综合了FEA法的优点，将杆件和关节看作单元，降低了运算量，但无法直接得到笛卡尔系的刚度矩阵[17-18]；应用较广的VJM法将杆件视为刚体，并虚设柔性关节(为了累积杆件或者关节所有类型的柔性)，较好地反映了机构的实际变形；文献[19]将驱动器柔性作为机构柔性的主要来源，推导出驱动刚度与机构刚度之间的关系，且考虑了杆件柔性[20]，并通过不同的假设和数值方法应用于不同的机构中；文献[21]采用虚拟弹簧法对相同DOF和POC的不同并联机构进行刚度建模分析，并优选出刚度性能较好的机构。

本文根据基于方位特征(POC)方程的并联机构拓扑设计理论[22]，提出一种零耦合度、含1条被动冗余支链的空间2T1R并联机构。首先，进行拓扑分析和位置分析，然后基于虚拟弹簧法建立机构的刚度模型，分析该机构各支链及整体的刚度特性，并讨论冗余支链对整体刚度性能的影响。

1 机构设计与拓扑分析

1.1 机构拓扑设计

设计的空间2T1R机构如图1所示，o′为动平台中心，动平台1由3条支链连接于静平台0，其中，混合支链Ⅰ(HSOC1)为一平面五杆机构(R11-P12-R5-P22-R21)，也可视为由2个相同A、B(RPR)支链并联而成；支链Ⅱ、Ⅲ由3个相互平行的R副串联而成；静平台0上4个转动副R11、R21、R31、R41的轴线平行于Y轴，动平台1上的转动副R6为复合铰链。

图1 零耦合度空间2T1R并联机构Fig.1 Spatial 2T1R PM with zero coupling degree

1.2 拓扑结构分析

1.2.1POC集计算

机构POC方程为

(1)

(2)

式中Mji——第i个运动副的POC集

Mbi——第i条支链末端的POC集

Mpa——机构动平台的 POC集

第Ⅱ、Ⅲ条支链为R‖R‖R支链，记为

SOCi={-Ri1‖Ri2‖Ri3-} (i=2,3)

由式(1)可得，上述3条支链的POC集分别为

由式(2)可得

(3)

(4)

式(3)表明：由支链Ⅰ、Ⅱ组成的子并联机构的动平台，已可实现xoz平面内的两平移和绕R5(或R6)轴线的一维转动(2T1R)的输出运动；因此，此时，并联连接的支链Ⅲ不影响机构的POC。

1.2.2机构自由度(DOF)

并联机构自由度计算式为[22]

(5)

(6)

v=m-n+1

式中F——机构自由度m——运动副数

fi——第i个运动副的自由度

v——独立回路数n——构件数

ξLj——第j个独立回路的独立位移方程数

第1个回路为平面五杆机构(R11-P12-R5-P22-R21)本身，显然，其独立位移方程数为ξL1=3，则自由度为

由上述第1个子并联机构、转动副R5及支链Ⅱ组成第2个回路；由式(6)、(5)得，其独立位移方程数和自由度为

由上述第2个子并联机构和支链Ⅲ组成第3个回路，由式(6)得，其独立位移方程数为

由式(5)得

该机构自由度为3，当取子并联机构上的P12、P22以及静平台0上的R31作为驱动副时，动平台1可实现两平移一转动的运动输出。因此，支链Ⅲ也不影响机构自由度。

由于支链Ⅲ不影响机构的POC、DOF，因此，为冗余支链。

1.2.3机构耦合度

由基于有序单开链(SOC)的机构组成原理可知，机构可分为若干个最小子运动链(Sub kinematics chain，SKC)，每一个SKC仅含一个基本运动链，又可分解为约束度为正、零、负的单开链(SOC)，第j个SOCj的约束度定义为

(7)

式中mj——第j个SOCj的运动副数

Ij——SOCj第j个驱动副

(8)

1.2.2节中已求得3个回路的独立位移方程数，即ξL1=3，ξL2=3，ξL3=3，其约束度分别由式(7)计算得

由式(8)得，机构的耦合度为

因此，该机构包含3个SKC，它们耦合度均为零，即κ1=0,κ2=0,κ3=0。

2 机构位置分析

由基于拓扑特征[23]的机构运动学建模与求解原理[24]可知，该机构的位置正解，可转换为其包含的3个SKC的位置正解求解；对本机构而言，因每个SKC的约束度为零，无需设虚拟变量，其运动位置可由其本身独立求解。

2.1 坐标系建立和参数标注

为便于理解，各运动副的标注如图2所示，设机构静平台0上的4个转动副为A11、A21、A31、A41，位于一边长为2l1的正方形的4个顶点处。

图2 机构运动学建模Fig.2 Kinematic modeling of PM

不失一般性，在静平台0上建立oxyz坐标系，o为静平台重心，x轴与A31A41连线平行，y轴与A11A41连线平行；在动平台1上建立o′x′y′z′坐标系，o′为正三角形动平台重心，x′轴与C33C43连线平行，y′轴在R5轴线上，而z轴和z′轴由右手笛卡尔坐标系法则确定。

设矢量A11C23、A31B32与x轴负方向的夹角分别为θ1、θ2，动平台1绕y′轴逆时针转动时的姿态角为α。该机构构件尺寸参数为

lA11A21=lA21A31=lA31A41=lA11A41=2l1
lA11C23=l2lA21C23=l3lC23C33=lC23C43=lC33C43=2l4
lB32C33=lB42C43=l5lA31B32=lA41B42=l6

2.2 位置正解求解

已知l2、l3、θ2，求解动平台1的位置。

易知各点的坐标为A11=(l1,l1,0),A21=(-l1,l1,0)，A31=(-l1，-l1，0),A41=(l1,-l1,0)，B32=(-l1-l6cosθ2,-l1,l6sinθ2)。

(1)SKC1(R11-P12-R5-P22-R21)位置求解

即在A11-C23-A21中，易求得

C23=(l1-l2cosθ1,l1,l2sinθ1)
C33=(l1-l2cosθ1-l4cosα,-l1,l2sinθ1+l4sinα)

从而求得

(9)

(2)SKC2(R6-R33-R32-R31)位置求解

即在C23-C33-B32-A31中，由几何约束建立位置方程，并整理得

Asinα+Bcosα+C=0

(10)

其中

A=2l4(l2sinθ1-l6sinθ2)

即可完成动平台的位置正解求解。

2.3 位置逆解求解及其验证

已知：动平台1坐标(X1，Y1，Z1)和姿态角α，求输入l2、l3、θ2。

B32、C23、C33的坐标在位置正解中已表示，由杆长约束条件lA11C23=l2及lA21C23=l3，求得

由几何约束条件lB32C33=l5，有

A2sinθ2+B2cosθ2+C=0

令tan(θ2/2)=k2，则有

(11)

其中

A2=-2l6(Z1+l4sinα)

综上可知，输入角θ2有2组解，因此，该机构有2种构型。

由式(9)、(10)，运用Matlab计算得机构的位置正解为(X,Y,Z,α)=(0,-66.67 mm，346.41 mm，-71.8°)；将此组数据代入逆解式(11)，得到的值与输入值一致，表明正逆解推导正确。

2.4 奇异位置分析

2.4.1雅可比矩阵

将由3个杆长约束条件(lA11C23=l2，lA21C23=l3，lB32C33=lB42C43=l5)建立的位置约束方程对时间t求导，可得到此机构末端执行器输入速度v和主动关节的输入速度ω的关系为

Jpν=Jqω

2.4.2奇异性分析

依据矩阵Jp、Jq是否奇异，将机构的奇异位形分为如下3类：

(1)输入奇异

机构发生输入奇异，意味着每条支链靠近驱动杆的2根杆处于折叠或完全展开状态，动平台的自由度减少，此时，det(Jq)=0，方程解集合A为

A={(u11=0)∪(u22=0)∪(u33=0)}

由于在实际应用中u11=l2，u22=l3皆不能为零，所以只存在唯一解，即

u33=l6sinθ2(X1-l4cosα+l1)+l6cosθ2(Z1+l4sinα)=0

满足u33方程解的三维构型，如图3所示。

图3 输入奇异位形Fig.3 Input singular configuration

(2)输出奇异

当机构发生输出奇异，意味着每条支链靠近动平台的杆处于折叠或完全展开的状态，此时，动平台自由度数增多，即使锁住输入，动平台也可能存在自由度输出。

设[fi1fi2fi3]=ei(i=1,2,3)，若det(Jp)=0，则向量e1、e2、e3有2种情况：

①存在2个向量线性相关

若e1=ke2，即满足

[f11f12f13]=k[f21f22f23]

由于e1=[X1-l1Z10]，e2=[X1-l1Z10]，e1≠ke2，所以，此类情况不存在。

若e1=ke3，即满足[f11f12f13]=k[f31f32f33]，此时机构到达空间左边界处，其奇异位形如图4a所示。

若e2=ke3，即满足[f21f22f23]=k[f31f32f33]，此时，机构到达空间右边界处，其奇异位形如图4b所示。

图4 输出奇异位形Fig.4 Output singular configuration

②存在3个向量线性相关

设e2=k1e1+k2e3(k1k2≠0)，此时有

[f21f22f23]=k1[f11f12f13]+k2[f31f32f33]

通过Matlab计算表明，该种情况下k1、k2无解，因此，此种情况不存在。

(3)构型奇异

当det(Jp)=det(Jq)=0时，也就是输入奇异与输出奇异同时发生；此时，机构的驱动关节和末端执行器都存在着瞬时互不影响的非零输入和输出，对应的位姿即构型奇异，处于该类奇异时，机构将失去自由度，在机构实际阶段应予以避免。

样机调试或机构工作过程中，一旦上述奇异位置发生，启动冗余支链Ⅲ产生动作，从而避免奇异位置发生，这对样机调试时的轨迹规划与运动控制具有参考价值。

3 刚度模型建立

3.1 单杆刚度矩阵

在基于虚拟弹簧的刚度模型中，杆件被视作梁单元分析其末端变形并求解其刚度矩阵，当杆件受到力/力矩时，由材料力学中的梁理论可得到杆件的挠曲线方程。在弯曲变形很小且材料服从胡克定律的情况下，挠曲线方程是线性的，考虑杆件所受力和力矩的耦合情况，采用叠加法计算杆件在力和力矩作用下的柔度矩阵(刚度矩阵的逆矩阵)，从而求得空间中杆件的刚度矩阵为[25]

(12)

式中Gi——杨氏模量Ei——弹性模量

Ai——截面积Ixi、Iyi、Izi——惯性矩

li——杆长

对于移动副作为驱动的杆件，在杆的刚度建模中，还需考虑导轨的变形，由文献[7]可得，该2T1R机构移动导轨刚度矩阵为

(13)

式中Kdg——导轨刚度矩阵

KψY、KψZ——Y、Z方向的扭转刚度矩阵

KlY、KlZ——Y、Z方向的线性刚度矩阵

3.2 支链运动方程

一般支链通常由驱动器，主、从动臂以及被动副组成，其虚拟弹簧模型如图5所示。

图5 一般支链的虚拟弹簧模型Fig.5 Virtual spring model of general branch chain

其中，1-dof的虚拟弹簧表示驱动关节为转动副R的伺服刚度，其变形可表示为Δθα；6-dof的虚拟弹簧表示对应连杆在笛卡尔坐标系中3自由度旋转变形特性和3自由度拉伸变形特性，主动臂和从动臂上弹簧的变形量可分别表示为(Δθ1,Δθ2，…,Δθ6)和(Δθ7,Δθ8，…,Δθ12)。

由图5可得，支链中弹簧变形和被动关节变形到末端变形之间运动方程的一般形式为

(14)

其中

式中 Δt——笛卡尔坐标系中机构末端的变形，由沿坐标轴方向的3个旋转变形和3个拉伸变形组成

Δθi——末端变形Δt的映射

将式(14)的支链运动方程表示成螺旋形式，即

(15)

式中 $——机构末端参考点变形的旋量，即支链末端参考点相对于虚拟弹簧和被动副的螺旋运动

3.3 支链的静力学方程与笛卡尔刚度矩阵

为得到支链的静力学方程，设fi为支链i所受的外力/外力矩，wi为虚拟弹簧所受的力/力矩，δθi和δψi分别表示虚拟弹簧和被动关节在受力后产生的微小变形，则末端的变形为

利用虚功原理，外力所做的虚功之和为

(16)

(17)

被动关节受力后会发生被动运动，所以在静平衡状态下，被动运动不做功，只有弹簧受力做功，即

(18)

且δθi与δψi为相互独立的变量，联立式(16)和式(17)，消去δθi和δψi，可以得到

(19)

由式(16)～(19)可得，支链的静力平衡方程为

(20)

由式(20)可得，静力平衡方程的矩阵形式为

(21)

其中

由文献[26]可得

(22)

式中Ki——支链i的笛卡尔刚度矩阵

若将一个具有n条支链的机构末端点的受力f分解到每条支链上，则有

因此，机构整体刚度矩阵为

(23)

3.4 2T1R并联机构的刚度建模

对该2T1R并联机构进行刚度计算，设机构各连杆的参数如表1所示，为提高机构整体的强度并减轻机构的质量，各连杆均选用碳纤维材料。

机构各个输入采用相同的驱动电机，其驱动刚度由实验测得，选取其刚度为5×104N/m。

3.4.1混合支链Ⅰ的刚度建模

机构中第Ⅰ条支链为混合支链(HSOC1)，它由两个相同的A、B支链并联而成。

由式(23)可知，KHSOC1=KA+KB，求解混合支链的刚度矩阵，只需求两条子支链A、B的刚度矩阵。

为了便于运动方程的建立与分析，先在a1、b1两点处建立了两个坐标系(图6)，oxyz为静平台坐标系，设o′为动平台中心点。由图5建立支链A的虚拟弹簧模型，如图7所示，其中，1-dof虚拟弹簧表示驱动关节的变形；6-dof虚拟弹簧表示杆a1c1和杆b1d1在笛卡尔坐标系中的旋转和线性变形；4-dof表示导轨在y、z方向上的旋转和线性变形；R表示转动副，Ac为驱动副。

由式(22)得，支链A的静力方程为

式中KA——支链A的刚度矩阵

O——零矩阵

Krodi——单杆刚度矩阵，i=1，2

Kdg——移动导轨的刚度矩阵

fA——末端点o′所受外力在支链A上的分量

利用Matlab软件，由机构的位置逆解，可得支链A在工作空间的刚度分布(图8)。

图8 支链A的刚度分布Fig.8 Stiffness distribution of branch chain A

由于支链B与支链A相同，可得支链B的静力方程为

式中KB——支链B的刚度矩阵

fB——末端点o′所受外力在支链B上的分量

利用Matlab软件，由机构的位置逆解，可得支链B在工作空间的刚度分布(图9)。

图9 支链B的刚度分布Fig.9 Stiffness distribution of branch chain B

3.4.2R-R-R支链的刚度建模

单开链SOC1(记为支链Ⅱ)、SOC2(记为支链Ⅲ)均为R-R-R支链，且支链Ⅲ为冗余支链。

对支链Ⅱ进行运动方程的分析与建立，o′为动平台中点，其虚拟弹簧模型如图10所示。

图10 支链Ⅱ虚拟弹簧模型Fig.10 Virtual spring model of branch chainⅡ

因此，由式(22)得支链Ⅱ的静力方程为

式中K3——支链Ⅱ的刚度矩阵

f3——末端点o′所受外力在支链Ⅱ上的分量

利用Matlab软件,由机构的位置逆解，可得支链Ⅱ在工作空间的的刚度分布(图11)。

图11 支链Ⅱ的刚度分布Fig.11 Stiffness distribution of branch chain Ⅱ

3.4.3冗余支链的刚度建模

由于冗余支链在正常工作状态下处于随动状态，所以在刚度建模时，无需考虑驱动副，3个R副均为被动副，其虚拟弹簧模型如图12所示。

图12 支链Ⅲ虚拟弹簧模型Fig.12 Virtual spring model of branch chain Ⅲ

因此，由式(22)得支链Ⅲ的静力方程为

式中K4——支链Ⅲ的刚度矩阵

f4——末端点o′所受外力在支链Ⅲ上的分量

利用Matlab软件，由机构的位置逆解，求得支链Ⅲ在工作空间的刚度分布(图13)。

图13 支链Ⅲ的刚度分布Fig.13 Stiffness distribution of branch chain Ⅲ

4 刚度分析

4.1 不含和含冗余支链时机构总刚度比较

利用式(23)，通过Matlab可以求出不含和含冗余支链时，机构的总体刚度，如图14所示。

图14 不含和含冗余支链时机构的刚度分布Fig.14 Stiffness distributions of PM without/with redundant branch chain

由图14可知，当机构存在冗余支链时，机构在工作空间内的总刚度增加，约提升了22%。

4.2 不含和含冗余支链时机构扭转刚度、线性刚度比较

由式(23)得到机构整体的笛卡尔空间刚度矩阵K6×6；其主对角线前3项分别表示x、y、z方向的扭转刚度，而后3项分别对应各方向的线性刚度，单均为N/m。

利用Matlab软件，计算该2T1R并联机构在不含和含冗余支链时，机构的扭转刚度和线性刚度，如图15、16所示。

图15 不含和含冗余支链时机构的扭转刚度分布Fig.15 Torsional stiffness distributions of PM without/with redundant branch chain

图16 不含和含冗余支链时机构的线性刚度分布Fig.16 Linear stiffness distributions of PM without/with redundant branch chain

4.3 数值验证

若取末端支链的坐标为P(0,0.067 m,0.3 m),其所在位置不含和含冗余支链两种情况下的扭转与线性刚度由Matlab导出，如表2所示。由于线性刚度远大于扭转刚度，可不再考虑冗余机构对机构的扭转刚度的影响，即冗余支链提供该机构更多的线性刚度，且由图16b可知，x、z轴方向的坐标值越大，刚度越大，表明机构向上运动时，刚度性能更好，动平台更稳定。

由表2可知，冗余支链Ⅲ使得机构的线性刚度增加约23.7%，与4.1节中冗余支链Ⅲ对机构总刚度提升22%相吻合。

表2 机构的扭转刚度和线性刚度Tab.2 Torsional and linear stiffness of PM

5 结论

(1)根据基于方位特征(POC)方程的并联机构拓扑设计理论和方法，提出一种空间两平移一转动并联机构，其具有的优势包括：零耦合度使机构具有符号式位置正解；3个SKC且驱动副分布在2个不同的SKC中，使机构具有部分运动解耦特性；全部由低副(R、P)构成，使机构制造容易；被动冗余支链能避免奇异位置，改善刚度。

(2)基于机构的符号式运动学，在工作空间中反映了笛卡尔空间刚度矩阵的分布，并通过Matlab得出在不含和含冗余支链两种情况下机构的刚度分布图；由比较不同截面刚度性能可知，在x、z方向偏移越大，整体刚度性能越好。

(3)通过对比不含和含冗余支链机构的整体刚度特性可知，加入冗余支链后机构整体刚度特性提升23.7%，由于该冗余支链为3R支链，所以，机构的线性刚度增加较为明显，提升约22%。