CN108614905B - 一种基于极大似然法的机车粘着性能参数估计方法 - Google Patents

Abstract

本发明公开了一种估计机车粘着性能参数的方法，可用于机车的粘着控制。通过机车传感器获取的数据，结合极大似然估计方法，构建出机车粘着性能模型的对数似然函数；以似然函数为基础，构建出机车粘着性能参数估计算法的实现流程；在上述流程中加入遗忘因子以及递推实现方法，构建出应用于机车粘着性能参数估计的递推算法。本发明采用极大似然原理和递推估计方法，应用于机车的粘着性能参数估计。

Description

一种基于极大似然法的机车粘着性能参数估计方法

技术领域

本发明属于系统辨识领域，具体地说，涉及对重载机车粘着性能参数的估计。

背景技术

轨面辨识是指通过机车的一些状态量来估计机车轮轨间关系的一种模型辨识方法。正确实时的轮轨面关系模型建立，可以实现在线调整粘着控制器的参数，提高机车粘着利用率并降低轮对出现空转、打滑现象的几率。

由于机车暴露在开放环境中，轮轨间的粘着行为是一个难以建立精确模型的复杂行为。轮轨间的粘着力是依赖于蠕滑运动而产生的，纯滚动产生不了粘着力。除此之外，机车轴重、第三介质、环境温度、湿度等对粘着力的影响也较大。尤其是水膜、冰雪、油膜、撒砂等第三介质，对轮轨间的粘着力有较大程度的影响。例如水介质下的粘着系数相比干态下降低了40％。轮轨间的粘着力还受到一些“连续变化量”如湿度、温度等的影响，粘着力会随着这些连续量的变化而连续变化。

轮轨间的粘着模型除了受到轮对与轨面状态的影响，还受到温度、湿度、轨面清洁程度等不确定性因素的影响，这导致了很难得到精确描述的粘着模型。目前关于蠕滑速度和粘着力的机理模型有线性模型、非线性模型，有基于二维滚动接触、三维滚动接触的模型。有些模型理论值和实际值非常吻合但要求可测的变量较多计算复杂。

目前，对粘着性能参数估计的方法有使用搜索加递推的最小二乘法和模糊逻辑方法。最小二乘法对适应模型动态变化具有一定适应能力，但对噪音的考虑有所欠缺。模糊逻辑方法对粘着性能的考量仅考虑了干燥、潮湿等对粘着性能影响大的因素，对性能参数的估计精准度不够。

极大似然估计方法对于动态系统的参数估计而言是一种具有良好统计性质的估计方法。极大似然算法通过构造一个似然函数，并极大化似然函数获得未知参数的估计值，满足一致性和渐进性等统计性质。

若能实时建立当前轨面条件下的粘着性能模型，得到对应的最佳蠕滑速度和粘着系数，然后控制电机使机车运行在最佳蠕滑点附近，这样就避免了微分运算对干扰的敏感。最优蠕滑的关键在于获取列车当前轮轨状态下的粘着性能参数。在仅考虑列车纵向牵引力前提下，设列车的轮轨间粘着模型为u＝f(λ,θ)，u为粘着系数，λ为蠕滑率，θ为描述模型的参数。如果能得到列车当前轮轨下的参数θ，就能对模型u＝f(λ,θ)求极值从而获得当前轮轨下的粘着峰值点，进而应用最优蠕滑方法对机车进行粘着控制。

发明内容

为了解决现有技术中粘着性能模型难以精确描述的问题，使用固定的参数难以描述轮轨间的粘着特性，例如干燥轨面切换到潮湿轨面时，粘着特性参数将会发生突变。空气湿度变化时，粘着特性参数将会发生缓慢的变化。因此，对粘着性能参数的精确获取需要算法能对轮轨环境的变化进行实时调整。

本发明的目的在于构建一种重载机车粘着性能参数的估计方法，用以解决现有技术中参数估计不准确以及不能快速实时适应参数变化的缺点，通过准确的参数估计获得最佳蠕滑速度和粘着系数，从而建立准确的粘着性能模型。

本发明提供了一种重载机车粘着性能参数的估计方法，利用机车传感器的在线测量数据，使用极大似然方法，达到估计列车粘着性能参数的目的。

本发明的技术解决方案如下：

一种基于极大似然法的机车粘着性能参数估计方法，包括以下步骤：

轮轨模型建立：

采用kiencke粘着模型描述轮轨间的粘着行为；

其中，μ₀是粘着特性曲线初始斜率，λ为蠕滑率，μ_λ为蠕滑率等于λ时的粘着系数，p₁,p₂为描述参数。对于各种轨面，其粘着特性曲线的初始斜率都近似相同，待估参数为p₁,p₂，记θ＝[p₁,p₂]。可求出该模型的最大粘着系数μ_m(λ_m)及其对应的蠕滑率λ_m为：

由(2)式可知，p₁,p₂共同决定了最大粘着系数的大小，最佳蠕滑率则由p₂决定。

若通过在线辨识的方式得到p₁,p₂的值，就可得到当前轨面的粘着性能模型。

建立似然函数

使用极大似然法求取待估计参数θ(k)＝[p₁(k),p₂(k)]。极大似然法是构造一个以测量数据和未知参数有关的似然函数，并通过极大化这个函数来获得模型的参数辨识值。

将(1)式变形如下：

记输出Z(k)＝μ₀·λ-μ(λ)，输入

待估计参数为θ(k)＝[p₁(k),p₂(k)]，考虑高斯测量噪音v(k):N(0，σ)，则(2)式可记为：

Z(k)＝p₁(k)·U₁(k)+p₂(k)·U₂(k)+v(k) (3)

由(3)式结合极大似然算法原理可得对数似然函数如下：

要使上式的对数似然函数取得极小值，需使

取得极小值。由(3)式可得v(k)＝Z(λ)-p₁(k)·U₁(k)-p₂(k)·U₂(k)，代入

可得：

令：

由于z²(k)为观测值，求极值时为一确定常数，则可将(6)式写为：

由(6)式可知，求解J(k)的极小值转变为一个求解二次规划问题，求使得(6)取极小值的x(k)即为欲求解的参数值。

上述(6)可用解二次规划的迭代方法求解，例如信赖域法，乘子法，拟牛顿法，单纯形法等解此二次规划问题。

引入时变遗忘因子

为了及时跟踪到粘着特性的实时变化，可引入遗忘因子如下式：

上式中，η为遗忘因子。

遗忘因子的大小影响了估计方法的敏感度；当遗忘因子较小时，历史数据对估计的影响下降，算法对环境的变化更为敏感，同时也会一定程度上影响算法的精度，此种情况适合于轮轨间环境发生变化下的情形；当遗忘因子较大时，历史数据保留较为完整，算法对环境的变化敏感度下降，但算法的精度上升，此种情况下算法适合于轮轨间环境稳定的情形。

考虑到机车行驶中速度变化和环境变化的不可测，这里使用P控制器来调节遗忘因子。遗忘因子的取值定为以下时变量：

上式中，η_k为k时刻的遗忘因子；ξ_k为算法的实时误差，定义为粘着系数u(k)的|预测值-观测值|；K_p为可调节的比例参数；η_min和η_max为η_k的取值范围；式中，当ξ_k越大，说明轮轨环境发生变化，此时遗忘因子减小，提高对环境的适应能力；当ξ_k接近零时，遗忘因子增大，增强了算法的精度。

算法流程

算法整体如下：

D1.给定初始点以及误差范围，考虑初始点设定为干燥轨面参数；

D2.在k时刻，由

计算H(k),f(k)；

D3.设校正矩阵

并计算在x⁽¹⁾(k)处的梯度g₁；

D4.计算搜索步长l_i，l_i为使得

成立的值；

D5.从x⁽¹⁾(k)出发，沿d⁽ⁱ⁾＝-M_ig_i搜索。令x⁽ⁱ⁺¹⁾＝x⁽ⁱ⁾+l_id⁽ⁱ⁾；

D6.计算梯度

令p⁽ⁱ⁾＝x⁽ⁱ⁺¹⁾-x⁽ⁱ⁾，q⁽ⁱ⁾＝g_i+1-g_i，计算矫正矩阵

D7.重复D3到D4直到满足误差要求，输出x(k)；

D8.在k+1时刻，重复D2到D7。

附图说明

图1是机车粘着性能参数估计方法的结构图

图2为模型对比图

图3为P1值辨识结果图

图4为P2值辨识结果图

图5为最大粘着系数值辨识结果图

图6为轨面切换发生时P1值辨识结果图

图7为轨面切换发生时P2值辨识结果图

图8为轨面切换发生时最大粘着系数值辨识结果图

图9为时变遗忘因子下最大粘着系数值辨识结果图

具体实施方式

下面结合具体实施方式对本发明作进一步的说明。其中，附图仅用于示例性说明，表示的仅是示意图，而非实物图，不能理解为对本专利的限制；为了更好地说明本发明的实施例，附图某些部件会有省略、放大或缩小，并不代表实际产品的尺寸；对本领域技术人员来说，附图中某些公知结构及其说明可能省略是可以理解的。

实施例1

学者kiencke给出了以下的模型描述轮轨间粘着特性：

其中，μ₀是粘着特性曲线初始斜率，λ为蠕滑率，μ(λ)为蠕滑率等于λ时的粘着系数，p₁,p₂为描述参数。对于各种轨面，其粘着特性曲线的初始斜率都近似相同，待估参数为p₁,p₂，记θ＝[p₁,p₂]。可求出该模型的最大粘着系数μ_m(λ_m)及其对应的蠕滑率λ_m为：

由上式可知，p₁,p₂共同决定了最大粘着系数的大小，最佳蠕滑率则由p₂决定。

Kiencke模型具有精度高、参数少的优点，非常适应于在线辨识的实时性要求。

如图2所示，Kiencke模型和著名的Burckhardt模型的对比，可以看出两个模型很接近。对于最受关注的粘着系数峰值而言，两个模型的差距小于2％，对应的蠕滑率差距约为0.1。

由上述分析可知，若通过在线辨识的方式得到p₁,p₂的值，就可得到当前轨面的粘着性能模型。如图1所示，机车粘着性能参数估计方法的结构图；

通过在线辨识方式获得p₁,p₂的值包括以下方法：

使用极大似然法求取待估计参数θ(k)＝[p₁(k),p₂(k)]。极大似然法是构造一个以测量数据和未知参数有关的似然函数，并通过极大化这个函数来获得模型的参数辨识值。

将(1)式变形如下：

记输出Z(k)＝μ₀·λ-μ(λ)，输入

待估计参数为θ(k)＝[p₁(k),p₂(k)]，考虑高斯测量噪音v(k):N(0，σ)，则(2)式可记为：

Z(k)＝p₁(k)·U₁(k)+p₂(k)·U₂(k)+v(k) (3)

由(3)式结合极大似然算法原理可得对数似然函数如下：

要使上式的对数似然函数取得极小值，需使

取得极小值。由(3)式可得v(k)＝Z(λ)-p₁(k)·U₁(k)-p₂(k)·U₂(k)，代入

可得：

令：

由于z²(k)为观测值，求极值时为一确定常数，则可将(5)式写为：

由(7)式可知，求解J(k)的极小值转变为一个求解二次规划问题，求使得(6)取极小值的x(k)即为欲求解的参数值。

该式可用解二次规划的迭代方法求解，例如拟牛顿法解此二次规划问题。

实施例2

如图3所示，模拟机车在单一轨面行驶的模型辨识结果，进行单一轨面情形下的参数辨识，并对最重要的当前轨面粘着峰值进行估计。输入信号考虑到列车实际的行驶状况，大部分时间都在蠕滑区行驶，只有少部分数据点会出现空转，因此设计的仿真实验中数据只包含了蠕滑区数据，取λ_umax＝0.1。噪音选择v(k):N(0，0.01)，估计算法的初值选择真值的20％，以验证估计算法在实际应用中的适应性，

由图3、图4可以看出，估计算法在蠕滑率从0到0.03时，便已估计出正确值。实验说明算法对数据不完全的适应度较高，只收集列车行驶在正常蠕滑区间内的数据，便能估计出列车的性能参数。例如通过上述的两个参数的估计，可以获知列车当前所在的轮轨面间的最大粘着系数，

由图5可以看出，当列车从蠕滑率为0加速到蠕滑率为0.03时，估计算法所预报的粘着系数峰值已经收敛到真值。

实施例3

如图6所示，主要考察估计算法对轨面环境突变的适应能力。仿真实验设计了列车从干燥轨面切换到潮湿轨面的情形。取λ_umax＝0.1，即输入输出数据仅选择了0到0.1的蠕滑区数据，噪音选择v(k):N(0，0.01)。考虑在列车的运行环境中，轮轨间的粘着性能是一个慢时变的过程，因此遗忘因子的值选为0.99，

图6和图7中的实验结果表明，在轨面环境发生突变的情形下，估计算法依然及时跟踪上了轨面的变化，正确辨识出了新轨面环境下的粘着性能参数。

图8的结果展示了估计算法在轨面切换发生时估计算法对最大粘着系数的预报。实验结果说明，估计算法能适应轮轨间环境的变化，及时地跟踪上新的轮轨环境下的性能参数。

实施例4

如图9所示，主要考察时变遗忘因子对估计算法的影响。(交代参数)仿真实验主要考察轨面切换下最大粘着系数的估计结果相比固定遗忘因子下的提升。遗忘因子的取值范围设为0.971～1，P控制器的参数设置为0.1，其余参数和实验2一致。

比较图8和图9，在轨面环境发生切换前，估计算法对最大粘着系数的估计在收敛之后的波动消失了。在轨面环境发生切换后，算法相比固定遗忘因子提前了500个数据点收敛到真值附近。实验结果说明，可变遗忘因子的引入，增强了算法对环境的适应能力，能够更快更精确地估计粘着性能参数。

显然，上述实施例仅仅是为清楚地说明本发明的技术方案所作的举例，而并非是对本发明的实施方式的限定。对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动。凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明权利要求的保护之内。

Claims (7)

1.一种基于极大似然法的机车粘着性能参数估计方法，其特征在于，包括以下步骤：

S1.通过机车传感器获取的实时数据；

所述数据包括机车粘着系数μ(λ)和蠕滑率λ；

S2.建立轮轨间模型；

所述模型采用kiencke粘着模型：

其中，μ₀是粘着特性曲线初始斜率；

λ为蠕滑率；

μ(λ)为蠕滑率等于λ时的机车粘着系数；

p₁,p₂为待估计的描述参数；

S3.建立似然函数；

对步骤S2中kiencke粘着模型变形；构造参数的辨识模型；所述参数辨识模型构造过程如下：

其中，记输入为U，输出为Z，

令输入

令输出Z(k)＝μ₀·λ(k)-μ(λ)，考虑高斯测量噪音v(k)～N(0,σ)，构造出参数辨

识模型如下：Z(k)＝p₁(k)·U₁(k)+p₂(k)·U₂(k)+v(k) (3)

S4.结合极大似然算法原理构造对数似然函数；

对S3中参数辨识模型(3)构造对数似然函数并对其求极小值；

S5.以步骤4中的机车粘着模型的对数似然函数为基础，构造目标函数，转化为二次规划问题；

由(3)式可得，v(k)＝Z(K)-P₁(k)·U₁(k)-P₂(k)·U₂(k)，代入

可得：

其中，令

上式变形为：

S6.采用拟牛顿法求解步骤S5的二次规划问题，得到待估参数p₁,p₂。

2.根据权利要求1所述的参数估计方法，其特征在于，还包括引入时变遗忘因子，如下式：

f(k)＝ηf(k-1)+[-2·Z(k)·U₁(k) -2·Z(k)·U₂(k)]。

3.根据权利要求2所述的参数估计方法，其特征在于，所述时变遗忘因子采用比例控制器来调节。

4.根据权利要求2所述的参数估计方法，其特征在于，所述遗忘因子取值包括以下时变量：

上式中，η_k为k时刻的遗忘因子；ξ_k为算法的实时误差，定义为粘着系数的|预测值-观测值|；K_p为可调节的比例参数；η_min和η_max为η_k的取值范围。

5.根据权利要求2所述的参数估计方法，其特征在于，所述拟牛顿法构造出递推方法得到参数θ(k)＝[p₁(k),p₂(k)]的算法流程为：

D1.给定初始点以及误差范围，考虑初始点设定为干燥轨面参数；

D2.在k时刻，由

计算H(k),f(k)；

D3.设校正矩阵

并计算在x⁽¹⁾(k)处的梯度g₁；

D4.计算搜索步长l_i，l_i为使得

成立的值；

D5.从x⁽¹⁾(k)出发，沿d⁽ⁱ⁾＝-M_ig_i搜索， 令x⁽ⁱ⁺¹⁾＝x⁽ⁱ⁾+l_id⁽ⁱ⁾；

D6.计算梯度

令p⁽ⁱ⁾＝x⁽ⁱ⁺¹⁾-x⁽ⁱ⁾，q⁽ⁱ⁾＝g_i+1-g_i，计算矫正矩阵

D7.重复D3到D4直到满足误差要求，输出x(k)；

D8.在k+1时刻，重复D2到D7。

6.根据权利要求1所述的参数估计方法，其特征在于，还包括步骤S7：仿真结果分析。

7.根据权利要求1所述的参数估计方法，其特征在于，所述二次规划问题求解方法还包括信赖域法、乘子法和单纯形法。

Images