CN115470691A - 一种滚珠丝杠元动作单元进给***参数辨识方法 - Google Patents

Abstract

本发明属于数控机床精细化控制技术领域，公开了一种滚珠丝杠元动作单元进给***参数辨识方法，首先，将数控机床按照结构化分解为机械元动作单元并选取滚珠丝杠元动作单元为研究对象，再建立滚珠丝杠进给***动力学模型，构造***的最小二乘估计，然后将滚珠丝杠动力学方程写为递推最小二乘的形式，确定该形式下所需要的初始参数以及观测值长度，计算此时刻的增益矩阵和协方差矩阵，然后使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子，之后计算出当前时刻的参数估计值，通过不断迭代输出最终的辨识结果。本发明通过引入动态调节遗忘因子的递推最小二乘法进行参数的估计，提高了参数的估计精度和收敛性能。

Description

一种滚珠丝杠元动作单元进给***参数辨识方法

技术领域

本发明属于数控机床精细化控制技术领域，特别是涉及一种滚珠丝杠元动作单元进给***参数辨识方法。

背景技术

数控机床作为装备制造业的“工作母机”，其技术水平和产品质量直接体现国家在机械加工领域的水平。目前，由于滚珠丝杠副具有高速、高精度、高刚度等优点被广泛用作数控机床伺服进给***的关键部件。滚珠丝杠进给***一般由伺服电机、联轴器、丝杠、螺母、轴承等部件通过不同类型的结合部联结而成，而结合部子结构间在机床运行过程中由于复杂的动态特性产生多方向自由度的微幅振动。因此，建立滚珠丝杠动力学模型有助于分析滚珠丝杠的动态特性，对提高数控机床的可靠性有重要作用。

动力学建模是进行精准控制的基础，建模精度和参数辨识精度直接决定了动态特性分析的精度，从而影响实时控制的效果。数控机床是一个非常复杂的非线性***，建立的动力学模型中包含有大量的参数，其中，刚度、阻尼系数、摩擦系数等受机床***的装配状态、运行工况的影响，导致参数随着时间的变化而变化，如何精准地辨识出这些动态参数对于模型的精确性有重要意义。

近年来，在机床参数辨识领域，基于传统智能优化算法的辨识算法存在收敛速度缓慢或容易陷入局部最优等缺点，在实际应用中存在一定的局限性。并且大多数方法仅适用于单个组合部的动态参数识别，缺乏对模型参数的整体分析。

发明内容

为解决现有技术存在的不足，本发明的目的在于提供一种机床元动作单元进给***参数辨识方法。

为实现以上目的，本发明采用如下技术方案：

一种机床元动作单元进给***参数辨识方法主要包括以下步骤：

步骤1.建立数控机床中典型的元动作单元-滚珠丝杠进给***的动力学模型；

步骤2：构造***最小二乘估计；

步骤3：建立滚珠丝杠的递推最小二乘形式，定义观测值长度L，初始化参数P(0)以及θ(0)；

步骤4：计算当前时刻的增益矩阵K(k)和协方差矩阵P(k)；

步骤5：使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子λ；

步骤6：更新当前时刻的参数估计值

步骤7：迭代次数N后输出最终的辨识结果。

进一步地，步骤1，将进给***的柔性件等效为***的刚度和阻尼，所述柔性原件包括联轴器和螺母；将工作台和滚珠丝杠转化为***的等效转动惯量，建立滚珠丝杠动力学模型如下：

其中，T_M为伺服电机的输出转矩；J_M为电机转子惯量；J_L为工作台的等效转动惯量；K、C分别为柔性元件的等效刚度和等效阻尼；θ₁、θ₂分别为电机轴和工作台的输出转角；

进一步地，步骤2，所述构造最小二乘法估计包括：

步骤2.1、建立***的回归方程为：

y(i)＝θ₁u₁(i)+θ₂u₂(i)+···+θ_nu_n(i)i＝1，2，3···L (2)

式中，u₁(i),u₂(i),···,u_n(i)为***在t_i时刻的输入，y(i)为***在t_i时刻的输出值，i＝1,2,3,···,L为观测次数；θ＝{θ₁,θ₂,···,θ_n}为***的待估计参数；n为***阶数；

步骤2.2、将(2)式写成矩阵的形式，具体如下：

Y＝φΘ (3)

式中，

由于在实际运行过程中各类误差的存在，输出值Y往往不等于ФΘ，因此假设误差矩阵为：

E＝[e(1),e(2),...e(L)]^T

(5)

式中，e(i)为误差，i＝1,2,3,···,L为观测次数。

因此，(3)式可以写成：

Y＝φΘ+E (6)

步骤2.3、设定性能指标为：

式中，J为代价函数，也称目标函数；E为***的误差矩阵；Y为***的输出矩阵；Φ为***输入矩阵。

步骤2.4、为使步骤2.3中的目标函数J取最小值，需要J对θ求取一阶导数，并令其为0，得到J取最小值时的参数

从而得到的

为最小二乘算法的最优辨识参数：

上述也即最小二乘估计的过程。

进一步地，步骤3，建立滚珠丝杠的递推最小二乘形式，具体包括以下步骤：

步骤3.1、对(1)式建立的动力学模型进行拉普拉斯变换得到传递函数；

式中，G(s)代表***的传递函数，是复变数s的有理代数分式；

代表输出转角的一阶导数；T_M为伺服电机的输出转矩；J_M为电机转子惯量；J_L为工作台的等效转动惯量；K、C分别为柔性元件的等效刚度和等效阻尼。

步骤3.2、采用零极点匹配法将连续的传递函数进行离散化处理，

离散后的传递函数为：

式中，b₁、b₂、b₃、a₁、a₂、a₃为需要辨识的六个参数，z为时移算子；

步骤3.3、将离散后的传递函数转化为差分方程的格式；

其中，

式中，y(k)为***第k时刻的输出值，u(k)为***第k时刻的输入值，e(k)表示为***的噪声；

步骤3.4、将步骤3.3中的方程转化为最小二乘的格式；

式中，y(k)为***第k时刻的输出值，ψ(k)为***的输出量的观测矩阵，e(k)表示为***的噪声，θ为***需要辨识的参数向量；

其中，

至此，滚珠丝杠元动作单元最小二乘形式构造完成。

进一步地，步骤5，具体包括如下步骤：

定义适应度函数为：

式中，y_t(k)为记录的输出实测值，

为估计出来的输出值。

选择λ作为优化变量，以J最小化为目标进行迭代寻优，迭代终止

时的遗忘因子即为当前时刻合适的遗忘因子值。

进一步地，使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子。

进一步地，最优遗忘因子λ传送给递推最小二乘法的算法中，完成相应参数的辨识。

进一步地，采用带遗忘因子的递推最小二乘法算法进行参数的辨识，递推公式为:

式中，

表示为k时刻辨识参数的估计向量，

表示第k+1时刻***输出的观测值，y(k+1)表示为***的实际输出，ψ(k+1)表示***第k+1时刻的观测矩阵，K(k+1)表示为***的估计增益矩阵，P(k)表示为***的误差协方差矩阵，λ为遗忘因子。

进一步地，步骤6具体为，将计算得到的P(k)、K(k)以及遗忘因子λ代入到式(16)中，得到当前时刻的参数估计值

进一步地，步骤7具体为，如若迭代次数达到N次，则输出参数估计值；若次数不足，则k＝k+1返回步骤5继续迭代，直到满足要求。

与现有技术相比，本发明具有以下有益效果：

(1)本发明将联轴器、螺母等柔性件等效为***的刚度和阻尼，将工作台和滚珠丝杠转化为***的转动惯量，建立动力学模型后，为解决传统最小二乘算法在辨识过程中会随数据量增大而出现“参数饱和”从而造成参数跟踪精度差的问题，通过在递推最小二乘算法中引入遗忘因子来提高该算法的收敛速度和跟踪精度，从而精准地辨识出动力学模型中的动态参数。

(2)在以往的带遗忘因子最小二乘算法中，遗忘因子一般都是固定值，无法使得算法的效果达到最优，因此，本发明中通过使用遗传算法改进的粒子群算法优化遗忘因子的选择，从而提高滚珠丝杠元动作单元进给***参数辨识的精准性。

附图说明

通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：

图1为本发明中的参数辨识方法的流程图；

图2为本发明建立的滚珠丝杠物理结构示意图；

图3为本发明建立的滚珠丝杠双惯量动力学模型；

图2中标记：1-电机，2-联轴器，3-滚珠丝杠，4-左轴承，5-右轴承，6-螺母，7-工作台，8-导轨，9-底座。

具体实施方式

下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进。这些都属于本发明的保护范围。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

如图1所示，一种滚珠丝杠元动作单元进给***参数辨识方法，包括以下步骤：

将数控机床按照FMA(“功能-运动-动作”)结构化分解为最小的单元—机械元动作单元，并选取滚珠丝杠元动作单元为研究对象。

步骤1：建立滚珠丝杠动力学模型；

首先构建滚珠丝杠进给***物理结构模型，如图2所示，一般数控机床进给***由电机1、联轴器2、滚珠丝杠3以及轴承(包括左轴承4和右轴承5)、螺母6、工作台7、导轨8、底座9组成。借助滚动体在丝杠滚道与螺母6滚道间循环往复的螺旋运动，实现将丝杠的旋转运动转化为螺母的直线运动，螺母带动工作台7实现直线进给运动。

如图3所示，将联轴器、螺母等柔性件等效为***的刚度和阻尼，将工作台和滚珠丝杠转化为***的转动惯量。

建立滚珠丝杠动力学模型如下：

式中，T_M为伺服电机的输出转矩；J_M为电机转子惯量(对应图3中驱动惯量)；J_L为工作台的等效转动惯量(对应图3中执行惯量)；K、C分别为柔性元件的等效刚度和等效阻尼；θ₁、θ₂分别为电机轴和工作台的输出转角。

其中，由于数控机床复杂的非线性和时变性等特点，在机床运行情况下由于摩擦、工作台重量发生变化从而导致建立的数学模型的精准性不高，滚珠丝杠进给***的精确性降低。因此，需要对数学模型中某些时变参数进行辨识。

步骤2：构造***最小二乘估计，具体为；

步骤2.1、建立***的回归方程为：

y(i)＝θ₁u₁(i)+θ₂u₂(i)+···+θ_nu_n(i)i＝1，2，3···L (2)

式中，u₁(i),u₂(i),···,u_n(i)为***在t_i时刻的输入，y(i)为***在t_i时刻的输出值，i＝1,2,3,···,L为观测次数；θ＝{θ₁,θ₂,···,θ_n}为***的待估计参数；n为***阶数；

步骤2.2、将(2)式写成矩阵的形式，具体如下：

Y＝φΘ (3)

式中，

由于在实际运行过程中各类误差的存在，输出值Y往往不等于ФΘ，因此假设误差矩阵为：

E＝[e(1),e(2),...e(L)]^T

(5)

式中，e(i)为误差，i＝1,2,3,···,L为观测次数。

因此，(3)式可以写成：

Y＝φΘ+E (6)

步骤2.3、设定性能指标为：

式中，J为代价函数，也称目标函数；E为***的误差矩阵；Y为***的输出矩阵；Φ为***输入矩阵。

步骤2.4、为使步骤2.3中的目标函数J取最小值，需要J对θ求取一阶导数，并令其为0，得到J取最小值时的参数

从而得到的

为最小二乘算法的最优辨识参数：

上述也即***最小估计的过程。

步骤3、建立滚珠丝杠的递推最小二乘形式，定义观测值长度L，初始化参数P(0)以及θ(0)；

步骤3.1、将建立的动力学模型进行拉普拉斯变换得到如下的传递函数：

式中，G(s)代表***的传递函数，是复变数s的有理代数分式；

代表输出转角的一阶导数；T_M为伺服电机的输出转矩；J_M为电机转子惯量；J_L为工作台的等效转动惯量；K、C分别为柔性元件的等效刚度和等效阻尼。

如上所述，现对求得的传递函数进行辨识。

由于本发明中利用光栅尺采集工作台的输出位移信号，但采集到的信号中包含有不可测量的有色噪声，因此，若是采用最小二乘法或者是递推最小二乘法来辨识***的参数，将无法识别到***的无偏估计值，所以，此发明采用带有遗忘因子的最小二乘法进行改进，继而对上述求得的传递函数进行估计。

步骤3.2、采用时域辨识的方法，需要对连续的传递函数进行离散化，因此，本发明中采用零极点匹配法将连续的传递函数进行离散化处理。

处理后的离散传递函数如下：

式中，b₁、b₂、b₃、a₁、a₂、a₃为需要辨识的六个参数，z为时移算子。

步骤3.3、将离散后的传递函数转化为差分方程的格式；

其中，

式中，y(k)为***第k时刻的输出值，u(k)为***第k时刻的输入值，e(k)表示为***的噪声。

步骤3.4、将上式转化为最小二乘的格式；

式中，y(k)为***第k时刻的输出值，ψ(k)为***的输出量的观测矩阵，e(k)表示为***的噪声，θ为***需要辨识的参数向量。

步骤3.5、确定待辨识θ以及观测矩阵

和

其中，a₁、a₂、a₃、b₁、b₂、b₃为滚珠丝杠进给***利用最小二乘法直接辨识出的六个机械参数。

至此，滚珠丝杠元动作单元最小二乘形式构造完成，进一步地，确定该算法中观测值长度L，初始化参数P(0)以及θ(0)；

步骤4、计算当前时刻的增益矩阵K(k)和协方差矩阵P(k)。

步骤5、使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子λ；

由于在传统的最小二乘法(RLS)辨识中，遗忘因子λ的取值往往为固定值，对于该算法的收敛速度和跟踪精度以及抵抗噪声的能力都是有很大的影响。

本发明引入动态遗忘因子，在遗忘因子最小二乘法中引入遗传算法改进的粒子群算法来优化遗忘因子的选择，从而提高参数辨识的精准性。

由于粒子群算法在没有搜索到群体最优解的情况下容易收敛到个体最优解，因此，本发明提出依靠遗传算法改进粒子群算法。具体解释为引入基因的交叉和变异过程来增加群体中个体的多样性，另一方面，加入淘汰机制，随机引入新的粒子，增强粒子的全局搜索能力。

步骤5.1、定义适应度函数为：

式中，y_t(k)为记录的输出实测值，

为估计出来的输出值。

选择λ作为优化变量，λ变化影响增益矩阵K(k)与协方差矩阵P(k)，以J最小化为目标进行迭代寻优。

步骤5.2、遗传算法改进粒子群算法过程如下：

初始化粒子群：

设置种群大小i，迭代次数it，学习因子C₁，C₂，惯性因子ω₁，交叉概率P_c，变异概率P_m；

种群位置初始化：λ₁，λ₂......λ_i；

种群速度初始化：V(λ₁)，V(λ₂)......V(λ_i)；

计算种群各粒子初始化适应度函数值，比较得出个体最优位置Pi和当前群体最优位置Pgi；

更新惯性权重值，计算当前迭代次数，再更新粒子的速度和位置；

计算适应值并排序，比较更新个体最优Pi；

将每个粒子的适应值Pi与经历过的最好位置Pgi进行比较，若较好，则更新群体最优Pgi；

用遗传算法对粒子个体当前位置进行交叉和变异；

若是达到迭代次数it，则停止算法的计算，输出群体最优位置Pgi，故而得到当前的全局最优解λ；

步骤5.3、将优化得到的遗忘因子λ传送给递推最小二乘法的算法中，完成相应参数的辨识。

本发明采用带遗忘因子的递推最小二乘法算法(FFRLS)进行参数的辨识，递推公式为:

式中，

表示为k时刻辨识参数的估计向量，

表示第k+1时刻***输出的观测值，y(k+1)表示为***的实际输出，ψ(k+1)表示***第k+1时刻的观测矩阵，K(k+1)表示为***的估计增益矩阵，P(k)表示为***的误差协方差矩阵，λ为遗忘因子。

步骤6、更新当前时刻的参数估计值

将步骤5-6中计算得到的P(k)、K(k)以及遗忘因子λ代入到式(16)中，得到当前时刻的参数估计值

步骤7、迭代次数N后输出最终的辨识结果；

具体为，如若迭代次数达到N次，则输出参数估计值，若次数不足，则k＝k+1返回步骤5继续迭代，直到满足要求。

最终可由输出的

计算得到滚珠丝杠进给***的可变参数。

以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变形或修改，这并不影响本发明的实质内容。

Claims (10)

1.一种机床元动作单元进给***参数辨识方法，其特征在于，包括以下步骤：

步骤1.建立数控机床中典型的元动作单元-滚珠丝杠进给***的动力学模型；

步骤2：构造***最小二乘估计；

步骤3：建立滚珠丝杠的递推最小二乘形式，定义观测值长度L，初始化参数P(0)以及θ(0)；

步骤4：计算当前时刻的增益矩阵K(k)和协方差矩阵P(k)；

步骤5：使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子λ；

步骤6：更新当前时刻的参数估计值

；

步骤7：迭代次数N后输出最终的辨识结果。

2.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤1，将进给***的柔性件等效为***的刚度和阻尼，所述柔性原件包括联轴器和螺母；将工作台和滚珠丝杠转化为***的等效转动惯量，建立滚珠丝杠动力学模型如下：

式中，T_M为伺服电机的输出转矩；J_M为电机转子惯量；J_L为工作台的等效转动惯量；K、C分别为柔性元件的等效刚度和等效阻尼；θ₁、θ₂分别为电机轴和工作台的输出转角。

3.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤2，所述构造***最小二乘估计包括：

步骤2.1、建立***的回归方程为：

y(i)＝θ₁u₁(i)+θ₂u₂(i)+···+θ_nu_n(i)i＝1，2，3···L (2)

式中，u₁(i),u₂(i),···,u_n(i)为***在t_i时刻的输入，y(i)为***在t_i时刻的输出值，i＝1,2,3,···,L为观测次数；θ＝{θ₁,θ₂,···,θ_n}为***的待估计参数；n为***阶数；

步骤2.2、将(2)式写成矩阵的形式，具体如下：

Y＝φΘ (3)

式中，

由于在实际运行过程中各类误差的存在，输出值Y往往不等于ФΘ，因此假设误差矩阵为：

E＝[e(1),e(2),...e(L)]^T (5)

式中，e(i)为误差，i＝1,2,3,···,L为观测次数；

因此，(3)式可以写成：

Y＝φΘ+E (6)

步骤2.3、设定性能指标为：

式中，J为代价函数，也称目标函数；E为***的误差矩阵；Y为***的输出矩阵；Φ为***输入矩阵；

步骤2.4、为使步骤2.3中的目标函数J取最小值，需要J对θ求取一阶导数，并令其为0，得到J取最小值时的参数

从而得到的

为最小二乘算法的最优辨识参数：

上述也即最小二乘估计的过程。

4.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤3，建立滚珠丝杠的递推最小二乘形式，具体包括以下步骤：

步骤3.1、对(1)式建立的动力学模型进行拉普拉斯变换得到传递函数如下：

式中，G(s)代表***的传递函数，是复变数s的有理代数分式；

代表输出转角的一阶导数；T_M为伺服电机的输出转矩；J_M为电机转子惯量；J_L为工作台的等效转动惯量；K、C分别为柔性元件的等效刚度和等效阻尼；

步骤3.2、采用零极点匹配法将连续的传递函数进行离散化处理，

离散后的传递函数为：

式中，b₁、b₂、b₃、a₁、a₂、a₃为需要辨识的六个参数，z为时移算子；

步骤3.3、将离散后的传递函数转化为差分方程的格式；

其中，

B(z^-1)＝b₁z^-1+b₂z^-2+b₃z^-3

A(z^-1)＝1+a₁z^-1+a₂z^-2+a₃z^-3 (12)

式中，y(k)为***第k时刻的输出值，u(k)为***第k时刻的输入值，e(k)表示为***的噪声；

步骤3.4、将步骤3.3中的方程转化为最小二乘的格式；

式中，y(k)为***第k时刻的输出值，ψ(k)为***的输出量的观测矩阵，e(k)表示为***的噪声，θ为***需要辨识的参数向量；

其中，

至此，滚珠丝杠元动作单元最小二乘形式构造完成。

5.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤5，具体包括如下步骤：

定义适应度函数为：

式中，y_t(k)为记录的输出实测值，

为估计出来的输出值；

选择λ作为优化变量，以J最小化为目标进行迭代寻优，迭代终止时的遗忘因子即为当前时刻合适的遗忘因子值。

6.根据权利要求5所述机床元动作单元进给***参数辨识方法，其特征在于，使用遗传算法改进的粒子群算法迭代出当前时刻的最优遗忘因子。

7.根据权利要求6所述机床元动作单元进给***参数辨识方法，其特征在于，最优遗忘因子λ传送给递推最小二乘法的算法中，完成相应参数的辨识。

8.根据权利要求7所述机床元动作单元进给***参数辨识方法，其特征在于，采用带遗忘因子的递推最小二乘法算法进行参数的辨识，递推公式为:

式中，

表示为k时刻辨识参数的估计向量，

表示第k+1时刻***输出的观测值，y(k+1)表示为***的实际输出，ψ(k+1)表示***第k+1时刻的观测矩阵，K(k+1)表示为***的估计增益矩阵，P(k)表示为***的误差协方差矩阵，λ为遗忘因子。

9.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤6具体为，将计算得到的P(k)、K(k)以及遗忘因子λ代入到式(16)中，得到当前时刻的参数估计值

10.根据权利要求1所述机床元动作单元进给***参数辨识方法，其特征在于，步骤7具体为，如若迭代次数达到N次，则输出参数估计值；若次数不足，则k＝k+1返回步骤5继续迭代，直到满足要求。

Images