0.   引言

在基于模型的机器人控制研究中，精确的逆动力学模型是必不可少的。因此需要对实际机器人开展动力学参数辨识工作，得到较为准确的参数，从而得到对应轨迹下的电机输出力矩。

关于机器人动力学参数的获取，目前主要的途径有：直接测量法^[1]、CAD计算法^[2]、实验辨识法^[3]。其中，实验辨识法是指通过辨识实验获得所需参数，该方法被广泛应用^[4-5]。实验流程主要为：线性化动力学建模，优化激励轨迹，采集处理所需数据，参数估计与验证^[6]。

机器人的关节模型是动力学建模中的重要部分。以往的众多研究成果中，在关节建模上方法疏于探索，常采用库伦- 黏滞摩擦模型或Stribeck摩擦模型来代表关节整体的摩擦^[7-9]。但对于协作机器人而言，其关节模组中采用了紧凑的谐波减速器，该减速器能够保证较高的传动精度，并且可以实现更小的体积及更大的传动比。但同时也具有较为复杂的摩擦来源以及关节柔性等问题^[10]。对于采用了谐波减速器的机器人关节，上述研究中的方法并不能准确反映关节内部的摩擦。需要对机器人的关节结构进行细分，以考虑谐波减速器等组件对动力学模型的贡献。

因此，对于以协作机器人为代表的串联机器人而言，需要建立其关节模型，以提高整体模型精度。Kircanski等^[11]对关节摩擦与关节刚度开展了辨识工作。张奇等^[12]对安装谐波减速器的机器人关节摩擦及刚度等开展了辨识研究。但上述研究只探究了关节转速对关节整体摩擦的作用，没有从谐波减速器本身原理出发对其摩擦展开研究。

有研究指出，谐波减速器摩擦的组成为：传动过程中的啮合摩擦、波发生器的滚动摩擦、滚动轴承摩擦、润滑液搅动等^[13]。针对谐波减速器的摩擦现象，Tuttle^[14]、Gandhi^[15]、石崟^[16]等进行了相关研究，并提出相应的摩擦模型。但上述摩擦模型形式较为复杂，参数较多，难以应用在机器人整体的动力学参数辨识研究中。

因此，有必要从谐波减速器的原理出发，得到更为简洁的摩擦力矩理论模型，从而建立更为精准的机器人关节模型，以进一步提高参数辨识的精度。

1.   机器人动力学模型

根据相关原理，可将n自由度机器人动力学方程表示为基于最小惯性参数的线性形式^[17-18]：

其中，τ_l∈为各连杆输出力矩，q_l∈为各关节位置；为最小惯性参数对应观测矩阵，为各连杆最小惯性参数。

但在实际运行中，若要得到关节电机的输出力矩，还需要考虑关节摩擦力矩的作用。

1.1   机器人简化关节模型

在一些参数辨识研究中，通常对机器人关节进行简化，便于开展研究。

图 1为关节i的简化模型，其中，和分别为电机端和连杆端角速度，τ_qⁱ和τ_lⁱ分别为电机力矩和连杆输出力矩，τ_fⁱ为关节摩擦力矩在连杆一侧的等效力矩，η_i为减速比。基于图 1所述的关节模型可得：

图  1  机器人简化关节模型

Fig. 1.  Simplified robot joint model

下载: 全尺寸图片 幻灯片

其中，、、τ_q、τ_l均为各关节对应变量所组成的n×1维的列向量；η为各关节减速比所组成的对角矩阵，即η=diag(η₁，η₂，…，η_n)。对于摩擦力矩τ_f，通常选用可线性化表示的库伦- 黏滞摩擦：

其中，f_c和f_v均为对角矩阵，分别表示机器人n个关节的库伦-黏滞摩擦参数。

综合式(1)~(3)，可得基于图 1所示关节模型的动力学方程：

其中，，为机器人最小惯性参数和摩擦参数；，为对应观测矩阵。

1.2   机器人改进关节模型

协作机器人的关节模组主要包括电机、谐波减速器、位置编码器、力矩传感器等组件。其中谐波减速器结构复杂，包含波发生器、柔轮及刚轮等多个组件。在关节运转时，除电机转子与连杆转动存在相应摩擦外，减速器传动过程中的啮合摩擦以及波发生器滚动摩擦也是关节摩擦的重要来源。

图 2为采用谐波减速器的单个关节结构简图^[12]。其中可按照关节结构将关节模型拆分为3个部分：伺服电机端、谐波减速器、连杆端。其中J_mⁱ为电机转子和波发生器的转动惯量，下文统称为转动惯量；τ_mⁱ和τ_rⁱ分别为电机输出力矩和减速器输出力矩，其中τ_mⁱ为1.1节中的τ_qⁱ；τ_mfⁱ和τ_mfⁱ分别为电机转子和连杆转动时的摩擦力矩；τ_rfⁱ为谐波减速器结构摩擦在靠近连杆一侧的等效力矩；K_fⁱ和K_tⁱ分别为谐波减速器与力矩传感器的刚度系数，统称为Kⁱ。

图  2  改进的机器人关节模型

Fig. 2.  Improved robot joint model

下载: 全尺寸图片 幻灯片

力矩传感器测量力矩为τ_rⁱ，包含了连杆输出力矩τ_lⁱ和连杆转动摩擦力矩τ_lfⁱ。对于机器人整体，改进关节模型可如式(5)所示：

其中，τ_m、τ_mf、τ_lf、τ_r、τ′_r、τ_rf均由各关节相应的变量所组成，均为n×1维列向量；J_m为各关节转动惯量所成的对角矩阵，即J_m=diag(J_m¹，J_m²，…，J_mⁿ)。对于摩擦力矩τ_mf和τ_lf，则可通过式(3)所示的库伦- 黏滞摩擦模型来进行建模。

对于谐波减速器中的摩擦力矩τ_rf，文[19]研究表明，当谐波减速器工作时，其内部摩擦主要来源于减速器工作过程中的啮合摩擦与滚动摩擦。其摩擦力矩τ_rf与输入力矩τ′_r满足：

其中，ξ_gf=diag(ξ_gf¹，ξ_gf²，…，ξ_gfⁿ)为各关节啮合摩擦系数所成对角矩阵；ξ_sf=diag(ξ_sf¹，ξ_sf²，…，ξ_sfⁿ)为各关节滚动摩擦系数所成对角矩阵。对于单个关节i而言，其减速器摩擦系数可表示为

其中，k、m、α分别为当前谐波齿轮的等效变形系数、模数和齿面角；R₀、D、d分别为减速器等速曲线曲率半径、未形变时柔轮直径和波发生器滚子直径；μ和μ_r则分别为齿面的摩擦系数与滚动摩擦系数。

由式(7)可知，当关节i的谐波减速器的结构确定时，其摩擦系数ξ_gfⁱ和ξ_sfⁱ是唯一确定的，因此可以认为减速器的摩擦力矩τ_rf与输入力矩τ′_r之间近似满足关系：

基于式(8)，可通过减速器摩擦系数来估计其摩擦力矩。该方法参数相对更少，表达形式更为简洁，更适用于参数辨识的整体工作中。

2.   参数辨识方案设计

2.1   激励轨迹设计

作为参数辨识工作中的关键部分，激励轨迹的选择至关重要。SWEVERS^[20]提出采用多阶傅里叶级数作为关节的激励轨迹，具有良好的效果。因此，在本次研究选择上述方案作为关节轨迹。对于关节i而言，其关节位置可表示为

其中，a_ik和b_ik分别为正弦、余弦函数的幅值；ω为函数基频；q_i0为初始位置的补偿量。本次研究取i=5，即关节轨迹为5阶傅里叶级数。

为减小辨识理论误差，需优化上述轨迹中的参数。目前广泛以观测矩阵条件数即cond(W)为优化目标，通过减小观测矩阵的条件数来减轻噪声对辨识精度的影响^[21]。关节激励轨迹的优化本质上是一个非线性优化过程，本文通过Matlab中fmincon函数对上述问题开展优化，以得到理论激励轨迹。

2.2   数据采集与处理

通过控制器中的PID控制模块控制机器人关节，以使各关节按照理论轨迹中的位置运动。

在本研究中，为进一步开展辨识，需要对机器人连杆端位置q_l、减速器输出力矩τ_r、电机等效输出力矩τ′_m(τ′_q)等参数进行采集。

由于关节内置力矩传感器，因此力矩τ_r可通过传感器直接测得。而电机等效输出力矩τ′_m则需要借助采集伺服电机的工作电流来得到。电机的工作电流与其输出力矩近似成线性关系，对关节i而言：

其中，η_i为关节减速比，在本研究中机器人各关节减速器的减速比为101；β_i是电机的扭矩常数，其中关节1、关节2为0.124 N·m/A，关节3~关节6为0.096 N·m/A。

受数据采集过程中的噪声影响，实际所得数据存在较多的“毛刺”，需进行相关数据的处理工作。为了去除力矩τ_r和τ′_m中的“毛刺”，可选择五点三次法进行平滑^[22]，以降低噪声对数据的影响。

对于连杆端、实际输出角度，可通过关节内部连杆端的编码器采集。考虑到机器人实际跟踪精度可能存在的不足，因此选择对关节角度先差分计算后低通滤波降噪的方式获取连杆实际角速度与角加速度。

采集激励轨迹中各时刻对应的关节运动状态，进而可计算对应轨迹下的观测矩阵^[18]，即。

在得到关节连杆端运动参数q_l、后，可基于式(5)，由τ_l和对减速器输入端即电机端的运动参数进行表示：

2.3   机器人动力学参数辨识流程

2.3.1   基于WLS法与简化关节模型的辨识方法

以文[7-9]中的研究为例，在许多机器人参数辨识研究中，均采用了相似的辨识思路。即基于1.1节中的模型来对关节进行建模，并通过加权最小二乘法来辨识所需参数。为了评估本文所提出的改进方法的效果，可按照上述传统思路开展辨识，得到相应动力学参数作为对照。

采用统一辨识的方法对机器人整体的参数开展辨识。在运行激励轨迹中，共在N个时间点上采集各个关节电机输出力矩和关节端实际输出角度，由式(4)可知：

其中，T_q为各时刻对应的电机等效输出力矩，即T_q=[τ′_q(t₁)，τ′_q(t₂)，…，τ′_q(t_N)]^T；为各时刻对应的观测矩阵所组成的组合矩阵，以下可简称为W。通过加权最小二乘法，可基于式(13)得到动力学参数：

其中，为测量力矩T_q测量方差所成的对角矩阵^[23]，称其为权重矩阵。辨识流程如图 3所示。

图  3  基于简化关节模型的辨识流程

Fig. 3.  Identification process based on the simplified joint model

下载: 全尺寸图片 幻灯片

如图 3所示，控制器将控制命令下发给伺服驱动，进而使各关节根据理论轨迹来运动。在完成数据采集及数据处理工作后，可得到激励轨迹下各个关节的力矩T_q与对应矩阵W，进而对动力学参数进行辨识。

2.3.2   基于WLS法与改进关节模型的辨识方法

在对机器人运行数据进行采集和处理后，可得到所对应的q_l、、τ_r、τ′_m。从而进行参数辨识实验，辨识步骤为：

1) 连杆参数辨识

由图 2所示关节结构可知，关节力矩传感器所测量力矩τ_r包含了τ_l和τ_f，由式(1)可知，τ_l为连杆输出力矩。τ_lf为连杆所受摩擦力矩，可通过式(3) 所示的库伦- 黏滞摩擦模型来建模。因此可采用与2.3.1节中相同的方法对改进模型连杆端的动力学参数进行辨识：

其中，Τ_r为各时刻力矩传感器测量数值，Τ_r=[τ_r(t₁)，τ_r(t₂)，…，τ_r(t_N)]^T；p_r为所需辨识的参数，主要包括连杆的最小惯性参数p_b和及相应的摩擦参数f_lc和f_lv，即p_r=[p_b，f_lc，f_lv]^T。

通过加权最小二乘法(WLS)辨识连杆参数p_r：

其中，为反映测量力矩T_r误差的对角矩阵。

2) 关节模型参数辨识

为完成机器人完整模型的参数辨识，在辨识度得到连杆参数p_r后，需进一步对关节模型的参数开展辨识工作。由式(5)可知，采集电机输出力矩τ_m和力矩传感器力矩τ_r，通过二者之差可得：

由式(6)可知，减速器摩擦力矩τ_rf与其输入力矩τ′_r为线性关系。考虑到无法对τ′_r直接测量，由式(5) 可知，忽略J_m和τ_mf，用电机输出力矩τ′_m作为替代。综合式(3)、式(8)、式(11)后可得：

其中，f_mc和f_vc分别为各关节库伦- 黏滞摩擦参数组合而成的对角矩阵，fmc=diag(f_mc¹，f_mc²，…，f_mcⁿ)，f_vc=diag(f_vc¹，f_vc²，…，f_vcⁿ)；τ′_m、可通过采集处理获得，其余则为待辨识参数。

为简化上述模型的辨识过程，提高参数辨识的效率，可对式(17)中的参数进行合并整理。合并后的各参数分别为

经过对参数的整理，可将式(17)的形式转换为

基于式(19)，在轨迹运行中的N个时刻采集数据后可得：

其中，ΔT=[Δτ(t₁)，Δτ(t₁)，…，Δτ(t_N)]^T为各时刻力矩差Δτ所组成向量；T_m=[τ′_m(t₁)，τ′_m(t₁)，…，τ′_m(t_N)]^T为各时刻τ′_m所组成的向量；为各时刻所组成的向量；表示各时刻所成向量；A_i=diag(a_i，…，a_i)∈ 为总的关节系数矩阵。

基于上述得到的式(20)，可对各个关节模型的参数开展辨识工作，具体的流程如图 4所示。

图  4  关节模型辨识流程

Fig. 4.  Identification process of the joint model

下载: 全尺寸图片 幻灯片

3.   参数验证方案设计

在完成参数辨识工作后，需要对参数的精度进行评估验证。值得注意的是，参数辨识本身的意义是在给定任意轨迹下，都可基于所辨识参数得到电机输出力矩τ′_m预测值，进而得到关节电机的控制电流。

1) 基于WLS法与简化关节模型的参数验证

对于基于简化模型辨识的参数而言，可根据理论运动参数得到验证轨迹对应的理论观测矩阵，从而基于式(12)计算关节预测力矩：

在计算电机输出力矩预测向量T_pre后，可通过计算预测值与实测力矩T_q间的均方差(RMS)对其进行评估：

2) 基WLS法与改进关节模型的参数验证

在辨识连杆参数时，其方法步骤与简化关节模型中的相同，参数验证方法也相同。因此，可由式(14)计算连杆端预测力矩向量T_rp：

为了获得改进模型的电机输出力矩预测值T_mp，需要预测力矩损耗ΔT_p与T_rp相加。但ΔT_p的计算却需要事先获得电机输出力矩T_m。为此，可以采用基于简化关节模型的辨识结果对T_m进行预测，即用T_pre来替代T_m：

虽然该方法会引入T_pre的误差，但相对于简化模型的辨识结果，精度依然有所提升。图 5为所得参数的验证过程。

图  5  参数验证流程

Fig. 5.  Parameter verification process

下载: 全尺寸图片 幻灯片

4.   结果分析与验证

图 6为实验所需的6轴协作型串联机器人，其D-H坐标如图 7所示。

图  6  6轴串联机器人

Fig. 6.  Six-axis serial robot

下载: 全尺寸图片 幻灯片

图  7  机器人D-H坐标系

Fig. 7.  D-H coordinate system of robot

下载: 全尺寸图片 幻灯片

该机器人的关节结构与1.2节中所述模型相同，主要为位置编码器、谐波减速器、驱动器、伺服电机等。控制系统与机器人之间通过Ether CAT总线进行通信，控制信号发送间隔为1 ms，信号采集间隔为10 ms。如表 1所示为D-H坐标参数，表中参数θ_i对应如图 7所示的机器人俯卧状态。

表  1  机器人D-H参数

Tab. 1.  D-H parameters of robot

i	αi-1/(°)	ai-1/mm	di-1/mm	θi/(°)
1	0	0	122.3	0
2	90	0	0	0
3	0	-270	0	0
4	0	-253	123.3	0
5	90	0	107.1	0
6	-90	0	99.1	0

下载: 导出CSV 

| 显示表格

令激励轨迹的基频为0.1π，轨迹周期为20 s，以cond(W)为优化指标，对傅里叶级数中的若干参数开展优化。可采用Matlab中用于求解非线性函数最值的fmincon函数进行优化，得到相应参数。所得各关节激励轨迹如图 8所示。

图  8  机器人激励轨迹

Fig. 8.  Robot excitation trajectory

下载: 全尺寸图片 幻灯片

图 8中的关节位置对应了D-H坐标系中的角度θ_i，其中关节2、关节4的始末位置为-90°，此时机器人运动过程中的始末状态对应了图 6中的直立姿态。

控制器通过其PID控制模块，生成相应控制命令，驱动机器人关节运动。同时采集所需数据并进行处理，其中电机端输出力矩和连杆传感器力矩信号进行5点3次平滑处理，取单个关节力矩的处理次数为20次。通过Matlab中的Butterworth滤波函数，对差分后所得的角速度与角加速度低通滤波处理，以减少差分处理造成的噪声影响。

在对数据进行采集和处理后，按照图 3和图 4的步骤分别对参数、p_r和{a₁，a₂，a₃，a₄，a₅}进行辨识。其中关节模型的参数如表 2所示。

表  2  关节模型参数

Tab. 2.  Model parameters of the joints

i	a1i	a2i	a3i	a4i	a5i
1	0.402	1.578	2.937	0.515	0.014
2	0.804	1.139	7.870	0.477	0.009
3	0.464	0.822	3.047	0.638	0.029
4	0.266	0.925	2.721	0.534	0.033
5	0.603	1.916	3.690	0.123	-0.034
6	0.282	0.795	2.689	0.614	-0.076

下载: 导出CSV 

| 显示表格

表 2中，关节5、关节6的参数a₅⁵、a₅⁶(即K^-1)为负。这是由于关节5、关节6的负载小，关节发生的弹性形变很小，可以基本视作刚体，因此a₅⁵、a₅⁶数值较小。当采用最小二乘拟合时，为了最小化整体误差，会造成上述较小的参数存在负数的情况。因此本实验中，关节5、关节6的等效刚度难以准确辨识。但本文目的是对电机输出力矩进行预测，该参数对力矩的作用很小，因此可忽略其影响。

分别基于式(21)、式(23)、式(20)计算激励轨迹下预测力矩向量T_pre、T_rp、ΔT_p，进而得到T_mp。基于式(22)，可计算基于简化模型的预测力矩T_pre和基于改进模型的预测力矩T_mp与实测力矩T_m(即T_q)间的均方差(RMS)，分别为1.685 9 N·m和0.681 5 N·m。对于激励轨迹下的关节力矩，基于改进模型的预测精度相对于简化模型提高了59.58%。

对T_pre、T_mp、T_m进行整理，可得各关节对应力矩，如图 9所示。其中传统法预测力矩为采用2.3.1节中辨识方法所得参数对应的关节预测力矩；而改进法预测力矩则对应了2.3.2节中辨识方法所得参数对应的预测力矩。

图  9  不同辨识方法预测力矩

Fig. 9.  Predicted torque of identification methods

下载: 全尺寸图片 幻灯片

将两种辨识方法所得的预测力矩与实测关节力矩进行比较可知，对于基于改进关节模型辨识所得参数，其对应的预测力矩与实测关节力矩更为接近。这也间接说明改进的关节模型具有更高的精度。

同时，由图 9可知，对于关节2和关节3而言，改进关节模型的精度提升效果更为明显。这可能是由于关节2和关节3所承受转矩较大，因此在传动过程中，减速器的啮合与滚动摩擦则更加明显。并且较大的转矩载荷会造成关节内润滑介质的挤压，进而造成内部润滑状态的转变。内部的润滑状态由完全润滑转变为混合润滑状态及边界润滑状态。此时关节的摩擦并不能简单地通过库伦- 黏滞摩擦来表述。

为了进一步对辨识所得参数的效果进行评估，需要在与辨识轨迹不同的轨迹下对比效果。在实际情况下，往往需要对给定目标轨迹下的关节力矩进行预测，此时的T_mp未知。因此无法通过式(20)来获得，而需基于式(24)来对其进行近似估计。生成的验证轨迹如图 10所示。

图  10  机器人验证轨迹

Fig. 10.  Robot verification trajectory

下载: 全尺寸图片 幻灯片

对所得力矩向量进行整理，可得各关节在验证轨迹下的关节力矩，如图 11所示。

图  11  验证轨迹预测力矩

Fig. 11.  Predicted torque of verification trajectory

下载: 全尺寸图片 幻灯片

图 11中结果与图 9类似，对于验证轨迹而言，通过改进法所得的关节预测力矩更接近实测数值。这也进一步证明改进后的关节模型具有更高的精度，并验证了相应参数辨识方法的效果。

通过控制器驱动机器人关节按照验证轨迹进行跟踪，并按图 5所示的验证流程分别基于式(21)、式(23)、式(24)计算此时验证轨迹下预测力矩向量T_pre、T_rp、ΔT_p。计算T_pre和T_mp与实测力矩T_m的均方差(RMS)，分别为1.704 5 N·m和1.576 2 N·m。由结果可知，对于验证轨迹对应力矩，基于改进模型的辨识精度相对于简化模型提高了7.52%。

5.   结论

本文研究了谐波减速器的摩擦来源，提出了一种改进的机器人关节模型。并提出了相应的参数辨识方案。通过辨识与验证实验，对该关节模型与相应辨识方法的效果进行了评估。实验结果说明，改进的关节模型具有良好的效果，相应的参数辨识方法具有一定的工程意义。

但在对验证轨迹下的关节力矩进行预测时，需要借助传统辨识方法所得参数对谐波减速器的输入力矩进行近似计算，因此难以摆脱传统方法中误差的影响，在一定程度上限制了精度提升的效果。针对上述问题，可在后续研究中，通过力矩传感器等设备对减速器的输入力矩进行测量。同时可针对谐波减速器摩擦来源及其估计方法开展更深入的研究，以进一步提高关节模型的精度与动力学参数辨识的效果。