1.本发明属于铣削加工相关技术领域，更具体地，涉及一种辨识铣削机器人低频频响函数的增量自激励方法。
背景技术：
2.铣削机器人工作空间大、结构灵活且可重构性强，在航空领域大型构件的加工中应用广泛。但机器人的串联结构导致了更低的动刚度，由此引发的铣削低频颤振极大限制了加工效率。机器人频响函数是预测铣削低频颤振的重要基础，因此，机器人姿态依赖的频响函数辨识对提升机器人大范围加工效率十分关键。
3.机器人频响函数辨识方法一般可分为建模法和实验法。在建模法方面，huynh等人2020年建立了基于拉格朗日方程的动力学模型，其对大部分模态的频率、阻尼比及振型有较好的预测效果。chen等人2020年考虑刀柄与刀具之间的耦合效应，采用位移导纳耦合的子结构分析方法计算铣削机器人任意姿态下的刀尖频响。wu等人2020年基于有限元分析提出了一种预测机器人低阶动力学的半解析方法。在实验法方面，sun等人2020年用关节角作为输入变量，使用偏最小二乘方法得到了六自由度机器人固有频率的回归模型。li等人2022年提出了一种利用拉丁超立方体采样技术，基于盲克里金模型的工业机器人固有频率预测方法。wang等人2022年结合机器学习方法和模态识别理论，提出了一种基于随机森林的多模态预测方法。
4.然而，对于复杂的机电系统，一般很难建立全面有效的动力学模型，因为偏差无处不在且难以预计，例如机器人关节的刚度、阻尼及摩擦特性往往随关节角呈现复杂的非线性变化趋势。因此，建模法虽然对任意机器人姿态都是通用的，但模型及其参数的偏差限制了其辨识精度。实验法虽然在实验空间附近具有较高的精度，但由于实际动力学数据的获取往往依赖人工试验，要想在大空间范围取得精确的预测结果意味着需要相当的时间和人力成本。
技术实现要素：
5.针对现有技术的以上缺陷或改进需求，本发明提供了一种辨识铣削机器人低频频响函数的增量自激励方法，解决机器人任意姿态的频响函数的辨识问题。
6.为实现上述目的，按照本发明，提供了一种辨识铣削机器人低频频响函数的增量自激励方法，该方法包括下列步骤：
7.s1在机器人预设标准姿态下利用质量块进行增量自激励实验，计算在该预设标准姿态下机器人末端的质量块的冲量和冲量矩，然后在所述机器人的工具中心点进行模态测试获得机器人在所述预设标准姿态下的频响函数，利用该频响函数和在该预设标准姿态下机器人末端的质量块的冲量和冲量矩求解机器人在所述预设标准姿态下的惯性自激励；
8.s2再次利用所述质量块在机器人的目标姿态下进行增量自激励实验，计算在该目标姿态下机器人末端的质量块的冲量和冲量矩；
9.s3利用步骤s2获得的在目标姿态下的所述质量块的冲量、冲量矩和步骤s1中在预设标准姿态下机器人的惯性自激励，求解在所述目标姿态下的机器人频响函数。
10.进一步优选地，在步骤s1和s2中，所述增量自激励实验按照下列步骤进行：
11.(a)机器人末端关节在预设速度下恒速运行，然后突然停止，测量获取突然停止时机器人末端的工具中心点处的加速度时域信号；
12.(b)在机器人末端附加一个质心位于所述工具中心点处的质量块，重复步骤(a)，以此获得带质量块后所述工具中心点处的加速度时域信号。
13.进一步优选地，在步骤s1中，所述机器人在所述预设的标准姿态下的惯性自激励按照下列步骤获得：
14.s11根据步骤s1中机器人在所述预设标准姿态下的频响函数，确定模态总阶数e，第e阶模态角频率ωe；
15.s12将增量自激励实验中步骤(a)和(b)获得的加速度时域信号，分别转化为加速度的频域信号，将两个加速度的频域信号相减以此获得相减后的加速度频域信号在第e阶模态角频率ωe下的幅值，即f
mea,e
；
16.s13计算在在所述预设标准姿态下的频响函数的频谱密度函数和所述f
mea,e
之间的残差平方和；
17.s14利用所述残差平方和计算获得标准差的辨识结果；
18.s15将所述标准差的辨识结果和步骤s1中获得在预设标准姿态下的质量块的冲量代入惯性自激励函数中，以此获得机器人在所述预设标准姿态下的惯性自激励。
19.进一步优选地，在步骤s13中，所述残差平方和按照下列关系式进行：
[0020][0021]
其中，r(σ)是残差平方和，σ为f(t)的标准差，f(t)为增量自激励，e表示模态角频率阶数，ωe表示第e阶模态角频率，e为用来辨识σ的模态数量，ii为关节i运行时质量块产生的冲量，f
mea,e
为根据所测响应和频响计算的自激励在ωe处的幅值。
[0022]
进一步优选地，在步骤s14中，所述标准差的辨识结果为：
[0023][0024]
其中，r(σ)是残差平方和。
[0025]
进一步优选地，在步骤s15中，所述惯性自激励函数f(t)按照下列关系式进行：
[0026][0027]
其中，t是时间，e表示指数底数，δii为关节i运行时质量块产生的冲量，|
·
|表示求向量模。
[0028]
进一步优选地，在步骤s1和s2中，所述质量块的冲量和冲量矩均按照下列关系式进行：
[0029][0030]
[0031]
其中，δii是机器人末端质量块的冲量，δi
mi
是机器人末端质量块的冲量矩，δji和δj
mi
分别是机器人末端质量块的动量和动量矩，m、i分别是质量块的质量和惯性张量，vi是质量块随关节i运动时质心的速度，ωi为预设的关节i的角速度，i是机器人关节的编号。
[0032]
进一步优选地，在步骤s3中，所述在目标姿态下的机器人的频响函数按照下列步骤求解获得：
[0033]
s31将步骤s1中获得的机器人在所述预设标准姿态下的惯性自激励转化为频域信号，获得δf
qi
(ω)，根据δf
qi
(ω)和δf
il
的关系，获得δf
il
，以此获得该δf
il
的方向
[0034]
s32采用传感器测量为关节i自激励时，加装质量块前后tcp处的平动(即x、y、z方向)响应的差值的频域表达，根据获得δγ
xi
、δγ
xi
和δγ
zi
，对δf
il
与δγ
xi
[0035]
、δγ
xi
、δγ
zi
做模态分析，以此获得三个不同方向上的第e阶模态的频响函数h
fi,x,e
(ω)、h
fi,y,e
(ω)和h
fi,z,e
(ω)；
[0036]
s33根据δf
il
和末端频响之间的关系，构建第e阶模态下机器人工具中心的模态振动方向的关系式；
[0037]
s34利用所述和构建在目标姿态下的机器人的频响函数的关系式，以此获得所需的在任意姿态下的机器人的频响函数。
[0038]
进一步优选地，在步骤s33中，所述第e阶模态下机器人工具中心的模态振动方向的关系式按照下列进行：
[0039][0040]
其中，分别为以δf
il
为输入，以δγ
xi
、δγ
yi
、γ
zi
为输出，通过模态分析拟合的第e阶模态的频响函数。
[0041]
进一步优选地，所述在目标姿态下的机器人的频响函数按照下列关系式进行：
[0042][0043]
其中，he(ω)是第e阶模态的频响函数he(ω)，是第e阶模态下机器人工具中心的模态振动方向，是质量块激励的方向，是质量块激励的方向，分别是三个方向上第e阶模态的频响函数。
[0044]
总体而言，通过本发明所构思的以上技术方案与现有技术相比，具备下列有益效果：
[0045]
1.本发明通过在机器人工具中心处附加质量块获得完全可知和可控的激励增量，基于附加质量块带来的响应增量与激励冲量增量之间的线性关系，绕开复杂的机器人惯性参数，使用激励增量和响应增量辨识机器人末端频响，辨识精度高；
[0046]
2.本发明根据动量定理建立了增量自激励冲量(矩)的模型，针对频响计算过程对噪声十分敏感的问题，将增量自激励假设为高斯脉冲，给出了其频谱的简化解析式，以此为基础提供了一种使用实测频响函数辨识增量自激励的方法；
[0047]
3.本发明建立了脉冲响应增量与增量自激励的关系模型，提供了一种基于铣削机器人末端模态方向性的频响函数计算方法，结果包含所有的直接项和交叉项，实例验证结果表明，本方法可以有效辨识铣削机器人tcp处的低频频响函数。
附图说明
[0048]
图1是按照本发明的优选实施例所构建的一种辨识铣削机器人低频频响函数的增量自激励方法的流程图；
[0049]
图2是按照本发明的优选实施例所构建的机器人增量自激励方法中的装置示意图，其中，(a)是机器人自激励的结构示意图，(b)是有质量增量和无质量增量的机器人末端示意图；
[0050]
图3是按照本发明的优选实施例所构建的圆柱-球形工具和刀具的几何信息；
[0051]
图4是按照本发明的优选实施例所构建的使用姿态1～3下实验结果计算的误差平方和曲线；
[0052]
图5是按照本发明的优选实施例所构建的姿态4的tcp频响辨识结果，其中，(a)是x向激励x向响应的频响，(b)是y向激励x向响应的频响，(c)是z向激励x向响应的频响，(d)是x向激励y向响应的频响，(e)是y向激励y向响应的频响，(f)是z向激励y向响应的频响，(g)是x向激励z向响应的频响，(h)是y向激励z向响应的频响，(i)是z向激励z向响应的频响；
具体实施方式
[0053]
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
[0054]
如图1所示，机器人增量自激励方法的步骤如下：
[0055]
(1)在机器人末端关节按照预设速度恒速运行突然停止，测量突然停止时带来的机器人末端工具中心点的加速度a1；
[0056]
(2)在机器人末端附加了一个质心位于工具中心点tcp处的质量块；机器人末端携带质量块m按照预设的速度v恒速运行突然停止，测量突然停止时带来的机器人末端工具中心点的加速度a2；
[0057]
对于频响函数的获取方法，具体步骤如下：
[0058]
s1、增量自激励冲量(矩)计算
[0059]
关节运行是一种宏观运动，在末端关节预设的速度恒速运行突然停止时，由于铣削机器人柔性集中在关节，宏观运动上述的运动转换为关节内部的振动动态变形的微观运动。这一过程发生在极短的时间内，因此可以将停止运行时关节内部变形的受激力视作脉冲，根据动量定理，其冲量大小正是关节停止运行前一瞬机器人运动部件的动量。
[0060]
首先进行机器人增量自激励实验，以关节i为例(i＝1,2,
…
,n，n为机器人驱动关节数量)展开说明。按denavit-hartenberg方法建立关节坐标系，记关节坐标系i为csi(o
i-x
i-y
i-zi)，机器人基坐标系为cs0，刀具坐标系为cs
t
，其原点o
t
位于tcp处。当只有关节i运行时，以和分别表示机器人运动部件的动量和动量矩。在机器人末端附加了一个质心位
于工具中心点tcp处的质量块，以δji和δj
mi
分别表示其的动量和动量矩。记总的动量和动量矩分别为ji和j
mi
，则有：
[0061][0062]
其中的变量都以cs
t
为基准(后文中没有特殊说明也是如此)。如前所述，和与复杂的机器人惯性参数相关，难以获取，但δji和δj
mi
可由下式计算：
[0063][0064]
其中m、i分别为质量块的质量和惯性张量，vi为质量块随关节i运动时质心的速度，ωi为预设的关节i的角速度，满足：
[0065]
vi＝ωi×oit
ꢀꢀ
(3)
[0066]
其中，o
it
表示关节坐标系的原点点oi到点o
t
的矢量。
[0067]
根据动量(矩)定理，质量块的动量(矩)会转换成冲量(矩)，于是有：
[0068][0069][0070][0071]
其中δii、ii、δi
mi
、i
mi
分别表示相应的冲量和冲量矩。
[0072]
s2、增量自激励辨识方法
[0073]
增量自激励的频域表达式是通过响应增量计算频响的必要条件。显然，机器人的(增量)自激励是无法直接测量的。这里提出一种结合所测频响函数辨识增量自激励的方法。
[0074]
结合响应增量和沿着增量自激励方向锤击所得的频响，可以获得实际的增量自激励的频谱。具体操作是响应增量的频谱与所测频响在每个频率处相除，但由于频响在大多数频率处幅值非常小，这样所得的结果对噪声非常敏感，是无效的。因此，本发明将增量自激励假设为形式固定的脉冲信号，再通过响应增量和频响在模态频率处所提供的少量离散数据来辨识脉冲的参数。
[0075]
在机器人铣削加工中，一般只考虑力而不考虑力矩，因此下面只关注自激励中力的辨识。将惯性自激励记为f(t)，假设其为高斯脉冲。标准高斯脉冲g(t)表示为：
[0076][0077]
其中，t为时间、σ为标准差。
[0078]
f(t)只存在于t≥0。由于g(t)主体部分在[-3σ,3σ]范围内，这里以g(-3σ)为f(t)的起点，即：
[0079]
[0080]
其中δii为关节i运行时质量块产生的冲量，可根据式(4)获得。
[0081]
f(t)的频谱密度函数f(ω)为：
[0082][0083]
上式可化为：
[0084][0085]
当|jσω|≤3，即ω≤3/σ时，有：
[0086][0087]
然后有：
[0088][0089]
而当ω＞3/σ时，f(ω)衰减到f(0)的1％左右，此时增量自激励难以造成有效响应，因此在使用响应增量辨识增量自激励时，式(10)是可用的。实际实施时，可以先根据响应增量在最低频模态处的幅值初步辨识σ，再根据所有小于3/σ的模态角频率处的幅值来提升辨识精度。模态频率可以使用最小二乘复指数法从响应增量中识别。
[0090]
下面介绍f(t)的辨识方法。由式(10)可知，f可以看作是关于ω和σ的函数。基于最小二乘原理可以辨识σ。设用来辨识σ的模态频率共有e个，则残差平方和(sse)表示为：
[0091][0092]
其中，ωe表示第e阶模态角频率，f
mea,e
表示(a2-a1)的频域信号在ωe处的幅值和频响函数在ωe处的幅值(由频响函数获得)的比值。式中
“×
2”的原因为模态频率一般只考虑正频率，而|f(ω)|关于纵轴对称，只考虑正频率时幅值需要乘2。
[0093]
这里只关注增量自激励的辨识。具体计算方法为：在tcp处测量沿着增量自激励方向锤击的三向频响；由于机器人末端模态响应具备方向性，以计算模态频率处三向频响幅值的二范数作为相应模态的合成动柔度；按照同样的方法计算响应增量在各阶模态频率处的合成振幅；然后合成振幅除以合成动柔度就获得了各阶模态频率处增量自激励的幅值f
mea,e
。
[0094]
辨识结果为：
[0095][0096]
式(12)是一个非线性最小二乘问题，可以使用高斯-牛顿法、列文伯格-马夸尔特法等方法来求解。由于不同关节的减速性能可能不一样，其增量自激励需要分别辨识。
[0097]
s3、机器人频响函数辨识方法
[0098]
将机器人tcp处的q方向激励p方向响应的单位脉冲响应函数记为h
pq
(t)(p,q＝x,y,z,ω
x
,ωy,ωz，分别表示cs
t
的坐标轴方向和绕坐标轴旋转的方向)，为了表示方便，以δf
qi
(t)第i个关节运动停止时质量块在六个方向上引起的激励力，即δfi＝[δf
xi δf
yi δf
zi
]
t
，使用上标
“ⅰ”
和
“ⅱ”
分别表示没有装和装有质量块，则使用关
节i自激励时tcp处p方向的响应和为：
[0099][0100][0101]
将式(13)与式(14)相减并变换到频域：
[0102][0103]
其中类似，表示傅里叶变换，ω表示角频率。式(15)可以变换为如下形式：
[0104][0105]
其中后者表示加装质量块之后tcp处频响的变化。
[0106]
铣削机器人本体质量往往达1吨以上，质量块的质量相对非常小，其对tcp频响的影响也往往可以忽略。而另一方面，常见的工业机器人的结构组成形式往往是靠近基座的连杆占据绝大多数质量，而末端连杆质量占比很小。这就使得使用末端关节采用运动-停止的方式进行自激励时，质量块带来的自激励是必须考虑的，其引起的响应增量是可以被有效捕捉的。也就是说，加装质量块后，tcp处频响的变化率远小于自激励的变化率。以和分别表示加装质量块之后tcp处频响和自激励的变化比，根据上面的论述，有：
[0107][0108]
式(16)就可以简化为：
[0109][0110]
上式即为脉冲响应增量与增量自激励的关系模型。
[0111]
下面计算机器人末端频响函数。
[0112]
式(18)可以表示成如下矩阵形式：
[0113]
γi＝hⅱδf
i (19)
[0114]
其中δγi、δfi是分别以δγ
pi
、δf
qi
为元素的6维列向量，hⅱ是以为元素的6
×
6维矩阵。将式(19)中的向量和矩阵分割成表示线位移/力和角位移/力矩的子向量和子矩阵：
[0115][0116]
在机器人铣削研究中，一般关注的是末端线位移，末端角位移的影响可以忽略。因此接下来只考虑而不考虑其满足：
[0117][0118]
在铣削机器人动态特性研究中，末端频响一般指的是hⅱll
，即线位移与力之间的传
递函数。在机器人自激励实验中，可以通过三向加速度传感器测得。hⅱla
是未知的，但它和hⅱll
在同一量级(这是由于机器人的柔性部位在关节，而作用在机器人末端的力和力矩投射形成的致使关节变形的扭矩显然是同一量级的)。由于质量块自激励中的冲量矩的大小应远小于冲量大小：
[0119]
δf
ia
＜＜δf
il
ꢀꢀ
(22)
[0120]
式(21)可以简化为：
[0121][0122]
下面提供一种基于铣削机器人末端模态方向性的频响计算方法。
[0123]
铣削机器人末端交叉频响本质上是末端模态振动方向性在不同激励和测量方向的体现。记δf
il
的方向和大小分别为单位矢量和δf
il
，满足显然与vi方向相同。对δf
il
与δγ
xi
、δγ
yi
、δγ
zi
做模态分析，记拟合的第e阶模态的频响函数分别为模态分析方法可使用最小二乘复指数法。基于式(23)，第e阶模态下机器人tcp的模态振动方向可由下式计算：
[0124][0125]
其中ωe为第e阶模态角频率。
[0126]
则方向上第e阶模态的频响函数he(ω)可由下式计算：
[0127][0128]
且有：
[0129][0130][0131]
其中h
..e
为he的元素，
[0132]
下面将结合具体的实施例进一步说明本发明。
[0133]
以abb公司的irb6660型重载铣削机器人为对象验证本方法。
[0134]
实验姿态如表1所示。通过式(12)，使用姿态1～3来共同辨识增量自激励，以减小辨识误差。使用姿态4来验证所提的频响辨识方法。
[0135]
实验现场如图2所示。图2中(a)展示了机器人自激励的示意图，即机器人从某姿态开始，匀速转动某一末端关节(这里是6关节)，然后在历经一定角度后突然停止。要注意的是，所关注的只是关节停止后的振动信号，之所以要遍历一定的角度，是为了使关节启动所引起的振动在关节停止前就衰减完毕，从而不影响关节停止所引起的振动信号。使用圆柱-球形工具作为质量增量，使用刀具作为对照，它们可以在机器人末端放大图图2中(b)看到。圆柱-球形工具由不锈钢做成，其密度与刀具相近，其圆柱部分的直径与刀具相同。质量增量为圆柱-球形工具质量减去刀具质量，其质心近似在球心。刀具为sandvik生产的600-025a25-10h面铣刀。圆柱-球形工具和刀具的几何信息如图3所示，两者的质量分别为
4.50kg和0.48kg，因此质量增量为4.02kg。使用dytran 3263a2三向加速度计来测量响应，其放置在接近球心的地方。结合dytran 5802a脉冲式大力锤测试频响。
[0136]
刀具坐标系cs
t
的x
t
轴与机器人6轴同向，z
t
轴与刀轴同向，如图2中(b)所示。使用6关节实施自激励，关节运行速度为16.7
°
/s。为了充分衰减6关节启动带来的振动，自激励的遍历角度设置为-150
°
。锤击方向为 y
t
，与自激励方向相同。
[0137]
增量自激励辨识的sse曲线如图4所示，辨识结果为σ＝0.0198s。图5中分别展示了姿态4下根据式(26)的频响辨识结果。实验结果表明本方法对机器人末端固有频率的辨识是准确的，且对于拥有更高幅值的频响来说，其幅值辨识的精度更高。然而，对于超过20hz的模态，如20.5hz和43hz模态，在本次实验中并没能被辨识出来。这是由于增量自激励的有效频带比锤击更低，无法激发起这些模态。但这并不影响本方法的实用性，因为总体上铣削机器人的动柔度会在频率升高时下降，发生本体颤振的频率一般低于20hz，在绝大部分情况下增量自激励的有效频带已经足够了。总的来说，实验结果表明，本方法可以有效辨识铣削机器人tcp处的直接和交叉低频频响函数。
[0138]
表1机器人姿态
[0139]
序号关节角1[0,0,30,0,-30,90]2[0,-41,19,0,-19,90]3[-38,-20,27,-60,-45,140]4[0,47,16,0,-16,90]
[0140]
本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。