面向工业机器人的数字孪生建模精度补偿方法*

康瑞浩，胡俊山，田 威，张嘉伟，马创业

（1.南京航空航天大学，南京 210016；2.中国空空导弹研究院，洛阳 471000）

随着机器人技术的快速发展，以工业机器人为载体的智能装备在航空航天智能生产线中的应用愈加广泛，如飞机翼面钻铆、型面铣削、部件装配等场合[1]。同时，对机器人作业过程的动态监测需求也变得更加迫切。数字孪生技术 （DT）能够对装备、产线乃至全工厂进行实时、有效的监测，而数字孪生能否准确地映射物理空间状态取决于数字孪生建模的准确度。现有机器人数字孪生建模方法大多依赖理论参数来构建模型，在传感器数据驱动下，往往导致数字孪生模型的运动轨迹与物理空间的实际轨迹不一致，影响到后续的预测与控制效果。因此，如何建立高精度、高保真的数字孪生模型面临着巨大挑战。

2003年Grieves教授于密歇根大学的产品生命周期管理 （PLM）课程上首次提出数字孪生概念[2]，描述为物理系统的数字化信息结构，数字孪生模型是将数字信息嵌入到物理系统中而形成的独立实体。北航陶飞教授团队[3–6]是国内较早深入开展数字孪生理论与技术研究的团队，提出了数字孪生“五维模型”并建立了数字孪生模型评价指标体系。不同的研究者也从各个方面不断丰富数字孪生的概念[7–8]。Li等[9]构建面向全生命周期动态演化的数字孪生驱动信息架构。Xia等[10]从制造结构层面建立数字孪生驱动的智能工厂，进而在此基础上实现智能感知[11]、预测[12]、过程优化[13]和控制策略[14]等功能。但这些研究都集中在工厂、车间、产线等宏观生产数据的孪生可视化，针对实际物理装备的数字孪生建模技术研究报道较少。

在物理装备孪生建模方面，Kuts等[15]基于AR/VR技术构建了工业机器人虚拟模型和虚拟操控界面，可以通过AR/VR头盔实现对虚拟机械臂的运动控制。王春晓等[16]利用MWorks对三轴数控铣床构建了综合机械、电气等多学科领域的数字孪生模型。Guo等[17]采用Focas数据采集方法，研发了智能制造车间数控机床数据采集与监控系统。这些研究虽然从装备的运动机理和运行状态等进行了数字孪生建模研究，但忽略了数字孪生建模的精度问题。

在面向工业机器人数字孪生建模中，通常采用理论参数来描述机器人内部结构的关联，但由于制造、装配、零部件磨损、重力及负载、温度等因素都会导致机器人本体与理论模型之间存在偏差[18]，严重影响机器人数字孪生建模精度。其中，机器人各连杆的参数误差以及造成的累积误差属于几何误差[19]；齿轮间隙、关节柔性变形等难以量化和衡量的误差属于非几何误差；还存在一些动态误差，会随着机器人使用时间和使用强度的增大而发生变化。

无论是几何误差、非几何误差还是动态误差，它们的存在都是无法避免的，在建模过程中需要对各误差造成的精度损失进行补偿。对于几何参数误差，通常可以采用基于运动学模型的参数标定法[20–22]，通过修正机器人关节参数，来准确反映工业机器人在物理空间的位置，Dennis[23]、Judd[24]、Renders[25]等提出了基于几何误差识别最大似然法的机器人标定方法，用来纠正连杆和关节参数误差。此外，许多学者建立了非运动学参数标定方法[26]，如Zeng等[27–28]建立了基于误差相似度的精度补偿方法，王东署等[29]采用前馈神经网络建立了机器人实际位姿和相应的关节角误差之间的关系，以确定位姿处补偿转角误差值，从而提高机器人的绝对定位精度。

为了减少传动误差等非几何误差，任永杰等[30]利用激光跟踪仪进行末端位姿的在线测量与修正，实现了机器人末端作业全闭环控制。在航空大部件加工中，Moeller等[31]利用激光跟踪仪对工业机器人系统进行实时位姿控制，该方法能够给数字孪生模型提供高精度映射，但需要增加激光跟踪仪等外部辅助设备，不仅提高了作业复杂度，同时只能对机器人末端数据采集，难以获得机器人各个活动关节的精确位置。国外学者将光栅尺安装在机器人各关节轴上，在机器人工作过程中读取光栅尺传感数据获得关节转角[32–35]，计算末端位姿与误差修正量，并在控制器中进行反馈控制[36]。该方法可有效减小关节转角误差，提高机器人的建模精度。

本文从工业机器人数字孪生建模机理出发，对数字孪生建模精度进行了理论定义。通过在机器人上安装光栅传感器，建立数字孪生测量系统，并对测量系统进行标定，从而减少了数字孪生建模中关节传动引起的非几何误差。采用L - M算法[37]对工业机器人建模参数进行辨识修正，补偿了数字孪生建模中的几何误差。试验结果表明了数字孪生建模方法的正确性和建模参数辨识对建模精度补偿的可行性。

1 数字孪生建模

1.1 数字孪生模型初始化

数字孪生模型可以定义为由传感器数据驱动的、具有物理实体运动约束的可视化表征。而模型初始化是在不考虑建模准确度的情况下获得与物理实体结构和功能近似的数字孪生模型。首先，利用OpenGL（开放式图形库）建立数字孪生可视化环境，并对模型显示的材质、场景、光照条件进行全局设置；其次，鉴于工业机器人的运动是各主要运动部件相互协调运动的结果，对工业机器人机械结构进行拆分，确定工业机器人运动所需的各关键运动部件，将其三维模型导入数字孪生可视化环境下；最后，由于OpenGL读取的各个模型位置和姿态是不固定的、任意的，所以需要参考实际机器人各部件模型的位置和装配方式，在可视化环境中进行调整并建立部件模型之间的从属关系。部件模型之间的关系可划分为相对运动与刚体固连关系，机器人各部件模型的从属关系如图1所示。

图1 铣削工业机器人部件从属关系图Fig.1 Subordination diagram of milling industrial robot moving parts

在定义机器人各部件从属关系的基础上，采用MD - H模型[38]表达工业机器人各运动部件之间的关联关系，对各关键运动部件三维模型进行运动约束。为了避免D – H建模中两相邻连杆平行时到达奇异位置难以求解的情况，引入关节平行度参数βi：当关节平行时，定义βi ≠0，di= 0；关节不平行时，定义βi= 0。此时，通过MD - H建立的机器人的第i节连杆相对i– 1节连杆的齐次变换矩阵可以表示为

式中，x、y和z分别为连杆坐标系的坐标轴；ai为连杆长度；αi为连杆扭转角；di为连杆偏距；θi为关节转角；βi为关节平行度。

针对工业机器人转动关节数量的不同，可以推导相对应的齐次变换矩阵。若转动关节数为n，根据式 （1）和理论设计参数，可得出其末端位姿在基坐标系下的齐次变换矩阵[39]，并定义为

式中，F(θS，X)为数字孪生驱动模型，用以后续建模参数标定；θS为关节驱动角度；X为建模结构参数组成的向量；T为旋转矩阵；n为关节数。

最后，利用传感器读取物理实体各关键部件的相对运动量，从而在可视化环境中驱动三维模型运动，保证与物理空间实际机器人具有完全一致的运动行为。

此时完成了数字孪生模型的初始化工作，但由于物理空间的工业机器人并不是理论参数构建的完全复刻体，存在着诸多误差，因此还需要对机器人数字孪生建模精度进行分析。

1.2 数字孪生建模精度

为了更好地表达数字孪生模型和物理实体之间的映射关系，将数字孪生模型与实际物理装备的运动映射误差定义为建模误差，孪生模型对物理空间映射的精准度评价指标定义为建模精度。在建模中应尽量按照工业机器人实际尺寸进行绘制，以保证孪生模型与物理装备运动的一致性。

图2描述了数字孪生的位置建模精度，Oc为物理空间中工业机器人末端测点到达的位置点；G点为孪生模型计算所到达的位置点；Oc点与G点之间误差值为该测点的建模误差AP，表达为

图2 数字孪生位置建模精度定义图Fig.2 Position of digital twin modelling accuracy definition diagram

式中，xj、yj、zj为末端测点到达位置的空间笛卡尔坐标，xi、yi、zi为孪生模型中计算位置点的空间笛卡尔坐标。

在笛卡尔空间中各个方向的建模误差定义为APX、APY、APZ，即

单一的测点难以表征整个数字孪生模型的建模精度，需在物理空间中取N个测点，分别计算出各单一测点的建模误差APN，并求出N个测点的建模误差平均值即

此时则标准偏差SAP值为

2 数字孪生实时测量系统

2.1 构建数字孪生测量系统

通常在无任何辅助测量设备的情况下，数字孪生模型可以通过伺服电机自带的编码器进行传动比换算，间接获得关节的转角值，从而驱动数字孪生模型运动。但由于伺服电机、减速器内部传动的齿轮存在间隙，同时机器人各连杆自重、负载和运动过程中自身的惯性力作用会导致关节处产生柔性变形，使得编码器间接得到的关节转角值与实际关节转角值之间存在角度偏差。因此，在工业机器人各个关节处安装高精度光栅尺传感器作为数字孪生实时测量反馈系统，即通过光栅尺直接测量各关节转角，将获取的光栅数据导入数字孪生模型，进而获得机器人末端的位置和轨迹信息，并通过实时映射来驱动机器人数字孪生模型。工业机器人转动关节部位驱动结构图如图3所示。

图3 工业机器人转动关节部位驱动结构图Fig.3 Diagram of drive structure of rotating joint part of an industrial robot

虽然圆光栅对于测量转角类数据更为合适，但是安装圆光栅需要拆卸机器人的转动轴，再次安装时难以保证光栅和转动轴的同轴度，从而影响传感器的信号采集精度，并且还会破坏机器人自身的传动精度。因此，本文采用将直线光栅粘贴到转动轴外环表面，通过“以直代曲”的方式对机器人关节转角数据进行测量采集，并根据工业机器人各个关节常用转角范围确定光栅尺长度。工业机器人共有6个旋转关节，从基座到末端将各关节定为A1～A6，机器人各个关节安装情况如图4所示。

图4 工业机器人各关节光栅尺安装图Fig.4 Installation diagram for each joint of an industrial robot

2.2 测量系统传感器标定

“以直代曲”的方式需要根据工业机器人各个关节的转角范围将直线位移测量转换成角度测量，为保证关节角映射的准确性，要对光栅尺的读数系统进行重新标定。在光栅尺进行标定时，必须先确定机器人各关节在Home位置、负极限位置、正极限位置的光栅读数，避免在标定过程中出现超限位等安全问题，并将Home位置定为光栅尺新零位。具体的标定过程如下（以A1轴光栅尺标定为例）。

（1）使用示教器控制A1轴关节以3%的均匀速度分别从不同转角返回Home点，重复10次，并依次记录相应的光栅尺读数，求出读数的平均值，定为A1轴的光栅尺新零位。

（2）控制A1轴关节从正极限位置向负极限位置方向以3%的速度匀速转动，以每2°为一个间隔记录光栅尺的读数。到达负极限位置后，以同样的速度和间隔向正极限位置匀速转动A1轴，并记录光栅尺读数。同时在A1关节处安装靶标，使用激光跟踪仪测量靶球当前位置。

（3）通过空间测量分析软件Spatial analyzer（SA）对所采集到的靶标点拟合得到圆A1'和A1''，依次连接各靶标点位与相应A1'、A1''圆心，可得到相邻靶球位置对应的圆心角，即为A1轴两次转动所对应的实际转角。

（4）将所得到的实际关节转角增量与光栅尺读数进行拟合，可得到关节转角与光栅尺读数之间的对应关系。

（5）A2 ～ A6轴光栅尺标定按步骤（1）～（4）分别进行，即可完成机器人各关节轴光栅传感器标定。

3 基于参数辨识的建模方法

3.1 建模误差模型

机器人结构参数误差是导致数字孪生建模过程中物理空间末端实际位姿与孪生模型期望位姿之间误差的来源，数字孪生建模过程中，采用机器人理论参数进行建模可以获得孪生模型的末端理论位姿OT，但由于机器人加工装配中各连杆参数存在误差，使得机器人末端运动过程并未到达理论位置，而当前到达的实际位姿记为OT'，在机器人基座安装等不确定因素影响下，机器人在物理空间的实际基坐标系OB很难确定，但可以利用测量设备对机器人基坐标系进行拟合，并将测量所得的基坐标系定为OB'，后续工作将围绕测量坐标系OB'来完成机器人结构参数标定。各坐标系构建如图5所示。

图5 参数标定坐标系Fig.5 Parametric calibration of the coordinate system

OB相对于OB'存在微动误差B'TB，矩阵B'TB表示为

式中，OB相对于OB'分别在x、y和z轴方向上产生的微小移动量为dxB、dyB和dzB；OB相对于OB'假设中心点的微小转动量为δxB、δyB和δzB。

OT相对OB产生微分运动，表示为

式中，BTT'和BTT分别为OT'和OT相对于OB的齐次变换矩阵；dT为微动齐次变换矩阵。

末端位姿的微小变化是由OT相对于OB产生的微分变换矩阵，表示为

由式（9）和（10）可知，OT'相对于OB'的齐次变换矩阵B'TT'为

式中，I4为4×4的单位矩阵。对式（11）进行变换，并忽略其中高阶微小量，则

式中，B'TT'（BTT）–1–I4为系统的线性标定误差矩阵；B'TB–I4为OB相对于OB'的微分变换。

式（12）表示机器人末端的位姿精度，是由测量基坐标系与实际基坐标系存在的构造误差和机器人运动学参数误差两者引起的。

其中，B'TB–I4可简化为矩阵形式，即

式中，dB和δB分别为OB相对于OB'的微分平移和旋转矢量；I6为6×6的单位矩阵；Δδ为OB相对于OB'的微小变动量所构成的矩阵。

用两个相邻连杆的5个运动学参数误差来表示连杆i的微分变换δAi和位姿误差dAi，多关节机器人的关节数n代表着机器人的自由度，由于机器人运动学具有5个关节参数，故此机器人末端到基坐标系的位姿误差可以用5n个运动学参数误差表示。

设机器人相邻两连杆之间的理论齐次变换矩阵为AiN时，由于机器人结构几何参数误差 （Δθi、Δdi、Δai、Δαi和Δβi）的存在，两相邻连杆之间的实际齐次变换矩阵AiR会与理论齐次变换矩阵AiN存在一个误差dAi。机器人末端实际位姿OT'相对于OB的齐次变换可表示为

将上式展开，忽略高阶无穷小，dT可化简为

式中，BTTi–1和BTTn分别为连杆i–1和n的理论位姿相对于OB的齐次变换矩阵。

由式（10）和（15）可知，机器人末端实际位姿相对于OB的微分变换矩阵Δ'为

机器人各连杆坐标系相对于上一连杆坐标系的微分变换Δ'可以由微分矢量dΔ和δΔ表示为

式中，Jdθ、Jda、Jdβ、Jδθ、Jδα和Jδβ均为雅克比矩阵，是由已知的机器人理论几何参数构建的，令向量ΔP为

将式（10）代入式（18）可以计算出系统线性标定误差向量ΔP、线性标定微分平移矢量d和旋转误差矢量δ；式 （11）可以计算出相对于OB′的微小平移矢量dB和旋转矢量δB，令矩阵J为

式中，I3为3×3的单位矩阵。令参数辨识误差向量ΔX'为

则由式（16）到（18）可得

在机器人参数标定中需要采集大量的点位坐标，用以提高标定的精度，假设存在m个采样点，上式可写为

式中，ΔPm为采样点的位置误差列向量；Jm为雅克比矩阵。

3.2 建模误差参数辨识

机器人建模参数误差辨识最优化过程是将非线性模型中高阶微分量忽略以按照线性分析，参数辨识误差ΔX'越来越小时，输入到孪生模型中获得的末端位置将越接近机器人实际位置，在辨识参数时，将非线性模型线性化的最常用方法是最小二乘法。最小二乘法取近似的中间值会造成计算误差，收敛过程中会使得目标函数离极小值越来越远，最后拟合的效果难以接受，而L-M最小二乘法利用阻尼系数μ限制迭代运算可以避免出现奇异性，提高了最小二乘法拟合效果。采用L-M最小二乘法的参数辨识基本原理如图6所示。

步骤1：初始化数字孪生模型，首次计算建模参数向量Xk，并作为机器人理论结构参数；

步骤2：将一定数量点位作为参数辨识采样点，通过激光跟踪仪测量机器人到达采样点的实际位置P；

步骤3：通过辅助传感器采集该采样点对应的关节转角向量θs；

步骤4：代入建模参数向量Xk下数字孪生模型，计算当前辨识的雅克比矩阵J（Xk），获得末端位置P'；

步骤5：求解P'与P之间存在误差ΔP（Xk）；

步骤6：求解出建模参数误差修正量ΔXk，即

步骤7：更新第k+ 1次迭代时的建模参数向量Xk+1和迭代次数k为

步骤8：更新第k+ 1次迭代时的阻尼系数μk+1为

式 中，||ΔP（Xk+1）||和||ΔP（Xk）||分别为第k+1和k次迭代的位置误差的二范数；δ为迭代参数。

步骤9：当||ΔP（Xk+1）|| – ||ΔP（Xk）|| >ε时，更新第k+1次建模参数向量；

步骤10：迭代运行步骤4到步骤9直 至||ΔP（Xk+1）||–||ΔP（Xk）||≤ε时，两次迭代的位置误差范数之差趋近于ε= 0.0001，证明已经收敛，即得出最优参数误差。

此时，用L-M算法辨识出的机器人结构参数误差ΔX修正机器人结构参数，并代入数字孪生驱动模型中，从而提高数字孪生建模精度。

4 试验验证

4.1 试验平台

为了验证本文提出的数字孪生建模方法的有效性，以及建模参数辨识对建模精度补偿的可行性，本文以工业机器人铣削系统作为试验对象，该系统为搭载铣削末端的KR500 M2830MT型工业机器人，利用FARO vantage E型激光跟踪仪作为机器人末端位姿测量设备，搭建试验平台如图7所示。在末端TCP测量过程中，将安装有激光跟踪仪靶球的特制刀柄固联在机器人末端铣削主轴上，并将靶球测量位置看作机器人末端TCP的位置，对比激光跟踪仪测量数据与数字孪生模型计算数据，完成建模精度试验验证。

图7 试验平台Fig.7 Test platform

4.2 数字孪生交互平台

数字孪生模型并不能独立存在，计算机软件为数字孪生模型的可视化表征提供了载体。综合数字孪生平台的通用性考虑，本文采用Qt编程软件进行数字孪生系统交互平台GUI开发，能够准确地描述工业机器人铣削系统中组件间的复杂相互作用，也可以连接控制器与物理空间建立永久、实时的交互模式。其主要界面如图8所示，包含末端位姿显示区、系统功能区、加工控制区、NC程序操作区。

图8 数字孪生系统交互平台界面Fig.8 Digital twins interactive platform interface

（1）末端位姿显示区。该区域显示由数字孪生模型计算得出的机器人铣削系统末端TCP的位姿，以及执行任务中所要到达的目标点位信息。

（2）系统功能区。该区域主要承担整体系统的初始化、项目工程的管理、加工任务执行进程的应急处置等功能。

（3）加工控制区。该区域是整个交互平台的核心区域，包含了铣削、主轴、电机、基准、法向、加工辅助、刀库、工艺、坐标系、孪生虚拟模型显现窗口等模块，在孪生模块子窗口内，用户可以实时观察到孪生模型对机器人铣削装备任务执行过程的实时映射，而且可以通过平移、旋转、缩放等操作切换不同的视角和细节，能够全面观察机器人的加工运行状态。

（4）NC程序操作区。该区域主要负责NC加工程序的编辑和运行。在加工过程中，查阅编辑框可跟随加工进程自动定位到正在执行的程序，方便用户掌握加工的进度。

4.3 试验验证

为了验证本文提出的数字孪生模型参数修正方法的有效性、正确性，对试验平台KR500机器人的A1～A6轴增加光栅尺，根据上文提到的标定方法对各个关节进行标定，确定光栅尺零点，标定结果如表1所示，后续将测量数据用来驱动数字孪生可视化模型运动。

表1 光栅尺标定结果Table 1 Calibration results for optical scale

在KR500工业机器人的运动可达范围内随机生成100个采样点和150个目标点的位置坐标，将100个采样点使用L-M算法进行参数辨识，另外用150个目标点进行建模误差和建模精度的试验验证。试验过程可分为以下6步。

（1）利用激光跟踪仪建立World坐标系、Base坐标系和末端TCP坐标系，然后在SA软件中建立Base坐标系作为测量基准。

（2）将100个采样点写入NC程序控制机器人运动，使用激光跟踪仪完成建模参数辨识前采样点位置采集，并将各光栅尺读数代入理论建模参数构建的数字孪生模型中以计算末端位置。

（3）将激光跟踪仪测量值与理论参数构建的数字孪生模型计算值对比，计算理论参数建模精度。并利用L-M算法对采样点进行参数辨识，机器人建模参数修正结果如表2所示。

表2 机器人建模参数修正结果Table 2 Robot modelling parameter correction results

（4）将辨识后的参数输入数字孪生交互平台，以修正数字孪生模型建模参数。

（5）将150个目标点写入NC程序控制机器人运动，使用激光跟踪仪测量实际到达位置，以及将各光栅尺读数代入修正的数字孪生模型中计算末端位置。

（6）将激光跟踪仪测量值与修正参数构建的数字孪生模型计算值对比，计算修正参数建模精度。

4.4 试验结果

将激光跟踪仪采集到的点位坐标数据分别与数字孪生模型中理论建模参数和修正建模参数所计算的点位坐标进行对比分析，根据式（3）、（4）和（7）计算综合建模误差、各方向建模误差和建模精度，如表3所示。

表3 修正前后的建模精度Table 3 Modelling accuracy before and after calibration

图9 ～12分别为根据理论建模参数、修正建模参数建立的数字孪生模型反映的点位信息和激光跟踪仪采集到的点位信息的差值。对比发现，采用机器人理论建模参数构建的数字孪生模型与修正建模参数构建的数字孪生模型的差别较大，两者差值最大可达到1.49 mm；最大综合建模误差从修正前的1.4920 mm减少到了0.2962 mm；最小建模误差从0.2590 mm减少到了0.0078 mm；平均建模误差从0.8672 mm减少到了0.1452 mm，映射精度提高了4.97倍；建模精度从± 1.6905 mm提高到± 0.3304 mm，提高了4.12倍，从而更能够准确映射出机器人运动过程中的状态。虽然建立的数字孪生模型与物理模型依旧存在建模误差，但已经在数字孪生模型要求范围内。

图9 X方向上建模误差对比Fig.9 Comparison of coordinate deviations in X-direction

图10 Y方向上建模误差对比Fig.10 Comparison of coordinate deviations in Y-direction

图11 Z方向上建模误差对比Fig.11 Comparison of coordinate deviations in Z-direction

图12 综合建模误差对比Fig.12 Comparison of combined coordinate deviations

5 结论

（1）实现了数字孪生测量系统构建。搭建高精度光栅尺测量系统实时采集铣削机器人在铣削加工时机器人各关节的转角，避免了齿轮间隙、编码器丢码等因素造成的所获关节转角数据不准确的问题，从而提高了数字孪生的建模精度。

（2）提出了数字孪生建模参数辨识方法。本文采用MD-H模型建立连杆坐标系，确定工业机器人铣削系统各连杆之间相互运动关系，基于L-M算法实现多源误差影响下的数字孪生建模参数修正，减小了数字孪生模型与物理空间实体之间差距，提高了数字孪生模型对装备映射的准确性。

（3）开发了数字孪生交互系统平台。基于Qt编程软件开发了与工业机器人铣削系统相匹配的数字孪生人机交互平台，包含了控制、数据采集、可视化等多个模块。

通过激光跟踪仪进行数字孪生建模试验验证，试验结果表明，数字孪生模型建模参数修正后较修正前对工业机器人末端点位的建模精度从± 1.6905 mm提高到± 0.3304 mm，提高了4.12倍，验证了数字孪生模型建立的正确性和建模参数修正的必要性。