一种基于RCSA的二叉树机器人铣削系统频响预测方法

摘要

本发明涉及机器人铣削加工领域，并公开了一种基于RCSA的二叉树机器人铣削系统频响预测方法，其包括：(a)将铣削系统划分为子结构B、子结构A和柔性结合部；(b)获取子结构B频响函数并对其模态参数标准化，求解子结构B在待求解姿态下的模态参数和耦合函数；(c)建立子结构A有限元模型并优化其刀具材料参数，获取其频响函数和响应矩阵；(d)实验测得待求解姿态下整体结构的频响函数，并通过IRCSA方法反算得到柔性结合部的响应矩阵；(e)根据RCSA方法计算待求解姿态下的铣削系统刀具端频响函数。本发明提出的预测方法，准确获得机器人铣削系统不同姿态下的刀具端频响函数，实现机器人不同姿态下铣削稳定性的精确预测。

说明书

技术领域

本发明涉及铣削装备动态特性测试领域，更具体地，涉及一种基于RCSA的二叉树机器人铣削系统频响预测方法。

背景技术

随着工业自动化与智能化的飞速发展，工业机器人因其具有性价比高、工作空间范围大以及灵活度高等诸多优点，正被越来越多的应用在船用螺旋桨、飞机蒙皮、航空发动机叶片以及火箭壁筒等大型复杂曲面零部件的磨削、铣削、抛光、钻孔、镗削加工。相对于通用CNC装备，机器人加工系统具有刚度低且刚度位姿依赖显著的特点，这严重限制了机器人在机加工领域性能的发挥与提高。

受机器人刚度低的影响，在铣削过程中机器人铣削系统较容易发生颤振现象，这不仅会降低加工效率、影响加工表面质量，而且还将显著加剧刀具磨损，缩短机器人铣削系统使用寿命。刀具端频响是分析铣削系统颤振稳定性极限的重要参数，刀具端频响的准确度直接关系到最终铣削工艺参数确定的合理性，进一步影响机器人铣削质量与效率的提高。

目前预测铣削系统刀具端频响的方法主要针对数控机床，由于机床刀具端频响随位姿变化不大，因此大多基于RCSA耦合法或组件模态综合法分析不同刀具、装夹长度与主轴转速下的刀具端频响变化，而鲜有人对不同姿态、不同位置的机床刀具端频响进行预测。在机器人铣削系统刀具端频响预测方面，主要方法在于建立机器人有限元模型，并通过实验来辨识相关动力学参数，最后来预测不同位姿下的机器人铣削系统刀具端频响，如文献“Mousavi S,Gagnol V,Bouzgarrou B C,etal.Stability optimization in roboticmilling through the control of functional redundancies[J].Robotics andComputer-Integrated Manufacturing,2018(50)：181-192”介绍的预测方法，由于机器人铣削系统需要辨识的参数众多且操作较复杂，因此该方法的预测精度无法很好保证。

针对上述问题，机器人铣削加工领域亟待需要一种操作方便、预测精度高的刀具端频响快速预测方法。

发明内容

针对现有技术的以上缺陷或改进需求，本发明提供了一种基于RCSA的二叉树机器人铣削系统频响预测方法，通过对机器人32种不同姿态下的机器人主轴-刀柄子结构进行锤击实验，利用基于RCSA的二叉树频响预测方法预测出任意给定姿态下的机器人主轴-刀柄子结构响应矩阵，结合刀具子结构有限元分析模型，能够准确的辨识出机器人铣削系统不同姿态下的刀具端频响函数，由此解决刀具端频响预测操作复杂和预测精度低的技术问题。

为实现上述目的，按照本发明，提供了一种基于RCSA的二叉树机器人铣削系统频响预测方法，其特征在于，该方法包括下列步骤：

(a)将待处理六轴铣削系统划分为子结构A、子结构B和连接二者的柔性结合部，其中，子结构A是刀具，子结构B包括机器人、主轴和刀柄，在子结构A上任意设定三个点，分别为点1、点2和点3a，在子结构B上任意设定两点，分别为点4和点3b，其中，点3a和点3b为子结构A与子结构B的交点，在机器人的轴2～轴6的每个关节处分别选择0和θ两个角度，以此获得机器人的32个姿态；

(b)在子结构B上建立主轴坐标系，对于所述32个姿态，在主轴坐标系中分别对每个姿态下的子结构B进行锤击实验，获得每个姿态下子结构B两点各自对应的频响函数，同时将该频响函数按照满二叉树进行存储，利用模态理论获得所述频响函数对应的模态参数，建立将θ角度向90度转化的模态参数标准化关系式(一)实现对每个关节的模态参数正交化，以此获得正交姿态下的模态参数，建立正交姿态和待求解姿态模态参数的关关系式(二)，以此获得待求解姿态模态参数，并由此获得该待求解姿态下子结构B在3b处的响应矩阵；

(c)在主轴坐标系中建立子结构A的有限元模型，同时在子结构A上建立刀具坐标系，并在该刀具坐标系下进行锤击实验，以此获得刀具坐标系下子结构A上每个端点的频响函数，建立该频响函数与子结构A的材料参数的关系式(三)，以此获得子结构A的材料参数，将其带入所述有限元模型中，以此获得主轴坐标系中子结构A的频响函数和响应矩阵；

(d)将子结构A装夹在子结构B上形成铣削系统整体结构，并在所述待求解姿态下对其进行锤击实验获得整体结构的频响函数，利用步骤(b)获得的子结构B点3b处的频响函数求解获得子结构B点3b处的响应矩阵，采用RCSA方法建立所述整体结构频响函数、步骤(c)获得的子结构A频响函数、子结构B点3b处的响应矩阵和待求解姿态下柔性结合部在点3a处的响应矩阵之间的关系式(四)，并由此计算获得在所述待求解姿态下柔性结合部在点3a处的响应矩阵部；

(e)根据RCSA理论将步骤(b)获得的在待求解姿态下子结构B在3b处的响应矩阵、步骤(c)获得的子结构A的响应矩阵、以及步骤(d)在待求解姿态下结合部在点3a处的响应矩阵部进行耦合计算，以此获得待求解姿态下整体结构点1处的响应矩阵，由此完成机器人铣削系统频响预测。

进一步优选地，在步骤(b)中，所述表达式(一)优选按照下列进行，

其中，a＝cos²θ,b＝sin²θ，与分别表示在α姿态下的节点1处频响第r阶模态质量、模态阻尼与模态刚度，α＝0°,θ或90°。

进一步优选地，在步骤(b)中，所述系式(二)优选按照下列进行，

其中，分别表示在待求解姿态θ^*下的节点1处频响第r阶模态质量、模态阻尼与模态刚度。

进一步优选地，在步骤(c)中，所述关系式(三)优选按照下列进行，

其中，U(X)是目标函数，是刀具坐标系下子结构A在节点i或j的y方向施加力而在节点i或j的y方向处拾取的频响函数，是刀具坐标系下由子结构A刀具参数仿真获得的在节点i或j的y方向施加力而在节点i或j的y方向处拾取的频响函数，i＝1,2或3a，j＝1,2或3a。

进一步优选地，在步骤(d)中，所述关系式(四)按照下列进行，

其中，R_3b3b是子结构B在点3b处激励，在点3b处拾取的响应矩阵，R_3a3a是子结构A在点3a处激励，在点3a处拾取的响应矩阵,R_J是子结构A与子结构B的柔性结合部响应矩阵，R_3a,J表示R_3a3a和R_J之和的中间变量，没有实际意义，h₁₁、h₁₂和h₂₂H_3a1是整体结构在端点1,2的频响函数，H_3a1、N_3a1、H_3a2、N_3a2、H₁₁、H₁₂和H₂₂是子结构A的频响函数，D是响应矩阵之和，y₁，y₂和y₃是响应矩阵之和中的元素。

进一步优选地，在步骤(e)中，所述待求解姿态下整体结构点1处的响应矩阵G₁₁优选按照下列表达式计算，

G₁₁＝R₁₁-R_13aD^-1R_3a1

其中，R₁₁是子结构A在点1处激励，在点1处拾取的响应矩阵，R_13a是子结构A在点3a处激励，在点1处拾取的响应矩阵，R_3a1是子结构A在点1处激励，在点3a处拾取的响应矩阵。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，能够取得下列有益效果：

1、本发明通过采用二叉树存储各个姿态的频响数据，建立了多轴机器人位姿描述模型，实现了n个关节的耦合；

2、本发明通过采用RCSA的方法，经过合理简化，提出了一种根据正交姿态实验频响结果计算任意姿态频响的数值方法；

3、本发明与背景技术中提到的文献相比，不需要辨识各个机械臂的动力学参数，对实验误差敏感度比较小，因此实现起来更为简单，且精度更高。

附图说明

图1是按照本发明的优选实施例所构建的频响预测方法流程图；

图2是按照本发明的优选实施例所构建的频响预测方法流程详尽图；

图3是按照本发明的优选实施例所构建的机器人铣削系统刀具端频响预测模型子结构划分示意图；

图4是按照本发明的优选实施例所构建的子结构A、B和柔性结合部的示意图；

图5是按照本发明的优选实施例所构建的多轴加工系统频响二叉树排列结构图；

图6是按照本发明的优选实施例所构建的单转动关节机械臂模态参数标准化示意图；

图7是按照本发明的优选实施例所构建的多轴加工系统正交姿态频响二叉树排列结构图；

图8是按照本发明的优选实施例所构建的空间梁单元受力示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

图1是按照本发明的优选实施例所构建的频响预测方法流程图，图2是按照本发明的优选实施例所构建的频响预测方法流程详尽图，如图1和2所示，一种基于RCSA的二叉树机器人铣削系统频响快速预测方法，其特征在于，具体包括如下步骤：

(a)图3是按照本发明的优选实施例所构建的机器人铣削系统刀具端频响预测模型子结构划分示意图，如图3所示，整个铣削系统划分为2个子结构，即：机器人-主轴-刀柄子结构(子结构B)与刀具子结构(子结构A)，子结构A与子结构B通过柔性结合部连接，如图4所示，在所述刀具子结构A上标记节点1、2与3a，在所述机器人主轴-刀柄子结构B上标记节点3b与4，其中，节点1与2距离为L，节点3b与4距离为l，且节点3a与3b在同一条直线上；

(b)在机器人Axis2～Axis6每个关节处选择2个角度，对确定的32个姿态进行锤击实验获得子结构B在点3b、4处对应的直接与跨点频响函数，具体地，所述机器人Axis2～Axis6每个关节处选取2个角度(k＝2,3,...,6)中共可以确认2⁵＝32个不同的姿态，对每个姿态进行锤击实验获得所述子结构B对应的直接与跨点频响函数H_3b3b,xx,H_3b4,xx,H_3b3b,yy,H_3b4,yy,H_3b3b,zz，其中，H_3b3b,xx表示在坐标系x方向锤击节点3b在节点3b上x方向拾取振动信号得到的频响函数，H_3b4,xx表示在坐标系x方向锤击节点4在节点3b上x方向拾取振动信号得到的频响函数而，H_3b3b,yy表示在坐标系y方向锤击节点3b在节点3b上y方向拾取振动信号得到的频响函数，H_3b4,yy表示在坐标系y方向锤击节点4在节点3b上y方向拾取振动信号得到的频响函数，H_3b3b,zz表示在坐标系z方向锤击节点3b而在节点3b上z方向拾取振动信号得到的频响函数；

(c)将32组频响数据按照满二叉树的排列方式进行存储，根据模态理论，拟合出每组频响数据的模态参数，并按照下列步骤计算机器人正交姿态下的模态参数：

(1)按照满二叉树的结构将每个姿态的频响进行排列，如图5所示，共有n层，每一层节点数为2^k个，叶子节点数为2ⁿ，从最上层开始，依次在每一层中取一个，就可以得到一个具有n个自由度的姿态，即每个叶子节点对应一个姿态；

(2)根据所述模态理论，结构频响函数是多个独立模态的叠加，可以用下式表示：

式中，N表示模态阶数，由锤击实验结果确定，λ＝ω/ω_r，ω_r为第r阶模态自然频率，ξ_r表示第r阶模态阻尼比，K_er表示第r阶模态刚度。频响函数H(ω)实部的最大最小值对应的频率为和虚部的最小值对应频率为自然频率ω＝ω_r，对应谷值为A^r，频响函数确定后，可获得上式中的各个参数值，由此可得到每组频响函数模态参数拟合方法为

(3)对于只含有一个转动关节的机械臂，如图6所示，所述模态参数标准化计算如下

其中，a＝cos²θ,b＝sin²θ，与分别表示机械臂在α姿态下的节点1处频响第r阶模态质量、模态阻尼与模态刚度，且有

与分别为机械臂在α姿态下的节点1处频响第r阶模态刚度、模态阻尼比与自然频率，于是，90°姿态下的节点1处频响第r阶模态参数可以根据与由MP广义逆计算得到；

(4)对于含有n个转动关节的机械臂，所述正交姿态模态参数标准化计算方法为：

根据满二叉树排列顺序，从1～n(自顶到底)依次对每个转动关节进行单个转动关节模态参数标准化计算，当对第k个关节进行模态参数标准化时，选择两组除第k个关节不同外其余关节均相同的姿态，即姿态与(m₁,m₂,…,m_n∈enum{1,2})，其中根据单个转动关节模态参数标准化方法可以计算姿态对应的频响模态参数。对n个关节完成标准化后即可得到n个关节正交的模态参数，其满二叉树排列如图7所示。

(d)建立刀具子结构有限元模型，根据锤击实验实际频响测试结果，并优化得到刀具材料参数，以此进一步计算刀具子结构A的频响函数，其中，所述刀具子结构A材料参数优化方法按照下列步骤：

(1)根据Timoshenko梁理论建立考虑剪切变形的刀具子结构A有限元分析模型。将m个Timoshenko梁单元质量矩阵M^e、刚度矩阵K^e按照有限元方法进行组集，得到总体刚度矩阵K与总体质量矩阵M

根据动力学方程的构造方法，可得具有m个Timoshenko梁单元结构的频域动力学方程，其中H为6×(n+1)维方阵，

式中，C为Rayleigh阻尼，即C＝αM+βK(α与β分别为质量与刚度系数)。

(2)将刀具子结构利用细绳悬挂起来通过锤击法测得刀具坐标系下的端点频响函数与如图8所示。把材料参数E、μ、ρ、α、β设为未知量X，则应找到一组X使得仿真得到的端点频响与实验测量结果相差最小，即求解刀具材料参数问题可变换成一个最优化问题

利用遗传算法或粒子群算法即可以进行刀具材料参数优化。

(3)将优化的刀具材料参数代入刀具子结构A的有限元模型中，计算主轴坐标系下频响函数H_3a1,xx、N_3a1,yx、H_3a2,xx、N_3a2,yx、H_11,xx、H_12,xx、H_22,xx、H_11,zz、H_3a1,zz，其中H_3a1,xx、N_3a1,yx表示刀具子结构A在节点1的主轴坐标系x方向施加力而在节点3a的x方向处拾取的平动、y方向拾取转动位移频响函数，其他频响函数含义可类推。；

(e)实验测得机器人-主轴-刀柄-刀具端频响，并通过IRCSA方法反算得到刀柄与刀具结合部的响应矩阵，该响应矩阵计算过程包括：

(1)将刀具子结构A与机器人主轴-刀柄子结构B装夹到一起形成铣削系统整体结构，选择32个频响测试姿态中的一个姿态分别在x方向测量频响函数h_11,xx、h_12,xx、h_22,xx，在y方向测量频响函数h_11,yy、h_12,yy、h_22,yy，在z方向测量频响函数h_11,zz，其中，h_12,xx表示在机器人铣削系统整体结构刀具端节点2的x方向施加力而在节点1的x方向处拾取的平动位移频响函数，其他频响函数含义可类推。

(2)根据RCSA理论，所述机器人铣削系统刀具端响应矩阵

G₁₁＝R₁₁-R_13aD-¹R_3a1

其中，G₁₁表示机器人铣削系统刀具端在1点激励，在1点拾取的响应矩阵，R₁₁表示子结构A在节点1处激励，在节点1处拾取的响应矩阵别，而式中

D＝R_3a3a+R_3b3b+R_J

其中，R_J表示柔性结合部响应矩阵。R_3b3b表示子结构B在节点3b处激励，在节点3b处拾取的响应矩阵，R_3a3a是子结构A在点3a处激励，在点3a处拾取的响应矩阵,。

Ⅰ.在主轴坐标系的x或y方向，有

对于R_3b3b

则由差分公式可得

因此，根据测得的H_3b3b,xx,H_3b4,xx,H_3b3b,yy,H_3b4,yy，可以分别计算机器人主轴-刀柄子结构在x与y向的响应矩阵R_3b3b。

Ⅱ.在主轴坐标系的z方向，G_ef与R_uv简化为

G_ef＝[h_ef,zz](e,f∈enum{1,2})R_uv＝[H_uv,zz](u,v∈enum{1,2,3a,3b,4})

对于R_3b3b

R_3b3b＝[H_3b3b,zz]＝[X_3b3b,zz/F_3b3b,zz]

(3)所述IRCSA方法求解结合部响应矩阵步骤为：

由于机器人变换姿态后，只有R_3b3b发生变化，将B改写为

D＝R_3b3b+R_3a,J(R_3a,J＝R_3a3a+R_J)

由于x、y方向与z方向的响应矩阵维度不同，结合部响应矩阵求解也存在一定差别。

Ⅲ.在主轴坐标系的x或y方向，总体响应矩阵与子结构响应矩阵满足下列关系

设

根据互易定理及向量化算子，上式可以变换为

在x方向，将实验得到的机器人铣削系统刀具端频响h_11,xx、h_12,xx、h_22,xx与有限元计算的刀具子结构频响H_3a1,xx、N_3a1,yx、H_3a2,xx、N_3a2,yx、H_11,xx、H_12,xx、H_22,xx代入上式即可求得D^-1，进一步代入所述机器人姿态θ^*下的R_3b3b,xx即可计算得到x方向的响应矩阵R_3a,J,xx；同样的，在y方向，将实验得到的机器人铣削系统刀具端频响h_11,yy、h_12,yy、h_22,yy与有限元计算的刀具子结构频响H_3a1,xx、N_3a1,yx、H_3a2,xx、N_3a2,yx、H_11,xx、H_12,xx、H_22,xx代入上式即可求得D^-1，进一步代入所述机器人姿态θ^*下的R_3b3b,yy计算得到y方向的响应矩阵R_3a,J,yy。

Ⅳ.在主轴坐标系的z方向，总体响应矩阵与子结构响应矩阵满足下列关系

G₁₁＝R₁₁-R_13aD^-1R_3a1

将实验得到的机器人铣削系统刀具端频响h_11,zz与有限元计算的刀具子结构频响H_11,zz、H_3a1,zz代入上式即可求得D^-1，进一步代入所述机器人姿态θ^*下的R_3b3b,zz即可计算得到z方向的响应矩阵R_3a,J,zz。

R_3a,J包括三个方向的分量R_3a,J,xx，R_3a,J,yy和R_3a,J,zz，R_3a,J表示柔性结合部的响应矩阵，与机器人的姿态无关。

(f)根据基于RCSA的二叉树机器人铣削系统刀具端频响预测方法，利用步骤(c)获得的正交姿态下的模态参数，计算任意待求解姿态θ_obj下的子结构B在点3b、4处的频响函数，并进一步计算子结构B在3b处的响应矩阵，该响应矩阵计算过程如下：

(1)对于单转动关节的机械臂，任意角度θ处的模态参数可以由所述标准化后的0°与90°姿态的模态参数根据MP广义逆得到

(2)对于含有n个转动关节的机械臂，所述频响计算方法为：

给定机械臂待求解的姿态θ_obj＝[θ₁,θ₂,L,θ_n]，根据满二叉树排列顺序，从n-1(自底到顶)依次对每个转动关节进行单个转动关节模态参数计算，当对第k个关节进行模态参数标准化时，选择两组除第k个关节不同外其余关节均相同的姿态，即姿态与(m₁,m₂,L,m_k-1∈enum{1,2})，这里第k个以后的关节角与给定关节角一致，第k个关节角为0°与90°，而第k个之前的关节角则为正交姿态，因此具有2^i-1对，反复利用单个转动关节模态参数计算方法对n个关节完成计算后可以得到给定姿态对应的频响模态参数，并进一步计算得到相应的频响函数H_3b3b,xx,H_3b4,xx,H_3b3b,yy,H_3b4,yy。

(3)根据所述差分公式，可以计算所述待求解姿态θ_obj＝[θ₁,θ₂,L,θ_n]下的机器人主轴-刀柄子结构B在3b处的响应矩阵R_3b3b,xx、R_3b3b,yy与R_3b3b,zz。

(g)将子结构B在3b处的响应矩阵、子结构A的响应矩阵(频响函数)以及刀柄与刀具结合部的响应矩阵根据RCSA方法进行耦合，计算得到待求解姿态θ_obj＝[θ₁,θ₂,L,θ_n]下的机器人铣削系统刀具端频响函数，具体包括下列步骤：

(1)将所述由IRCSA计算得到的柔性结合部响应矩阵R_3a,J与待求解姿态θ_obj＝[θ₁,θ₂,L,θ_n]机器人主轴-刀柄子结构B在3b处的响应矩阵R_3b3b代入下式

D＝R_3b3b+R_3a,J

(2)根据RCSA理论，机器人铣削系统刀具端节点1的响应矩阵与各子结构的响应矩阵关系如下

G₁₁＝R₁₁-R_13aD^-1R_3a1

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。