Identification of Wheel-Rail Vertical Forces of Rail Vehicles Based on Square Root Cubature Kalman Filter Algorithm
-
摘要:
轮轨作用力是评估轨道车辆运行品质的关键指标,为实现对轨道车辆轮轨垂向力的在线监测,提出一种基于平方根容积卡尔曼滤波(SRCKF)算法的识别方法. 以考虑悬挂元件非线性的车辆-轨道垂向耦合动力学模型为例,建立包含轮轨垂向力和车辆部件状态变量的非线性过程函数,将车体、构架和轮对垂向加速度作为观测量,基于SRCKF算法递推识别轮轨垂向力;在此基础上,建立整车动力学模型及其对应的17自由度轮轨垂向力估计模型,对车辆在实际不平顺激扰下的左右侧轮轨垂向力进行识别. 仿真结果表明:所提方法识别垂向车辆模型在随机不平顺、钢轨波磨不平顺和钢轨焊缝不平顺作用下的轮轨垂向力时,轮轨垂向力识别值在时域和频域同仿真值均有较高的吻合度,相关系数分别为0.988、0.999和0.969;在识别整车模型的轮轨垂向力时,左、右侧轮轨垂向力的相关系数最低分别为0.747和0.720,左、右侧轮轨垂向力之和的相关系数为0.999.
Abstract:The wheel-rail interaction force is the key index to evaluate the operation quality of rail vehicles. An identification method based on the square root cubature Kalman filter (SRCKF) algorithm was proposed for online monitoring of the wheel-rail vertical force of rail vehicles. By taking the vehicle-track vertical coupled dynamics model considering the nonlinearity of suspension elements as an example, a nonlinear process function of wheel-rail vertical forces and state variables of vehicle components was established. The vertical accelerations of the vehicle body, frame, and wheelset were adopted as observations, and the SRCKF algorithm was employed for the recursive estimation of wheel-rail vertical forces. Furthermore, a vehicle dynamics model and its corresponding wheel-rail vertical force estimation model of 17 degrees of freedom were established to identify the vehicle’s left and right wheel-rail vertical forces under actual irregularity excitation. Simulation results indicate that the proposed method can precisely identify wheel-rail vertical forces of vehicles excited by the random irregularity, rail corrugation irregularity, and rail weld irregularity, with the identified value of wheel-rail vertical forces is in good agreement with the simulation value in both time domain and frequency domain. The correlation coefficients are 0.988, 0.999, and 0.969, respectively. When the proposed method is used to identify the vehicle’s wheel-rail vertical force, the lowest correlation coefficients of the left and right wheel-rail vertical forces are 0.747 and 0.720, respectively, and the correlation coefficient of the sum of the left and right wheel-rail vertical forces is 0.999.
-
轮轨相互作用力是轨道车辆运行过程中车辆与轨道相互联系、耦合的核心纽带. 列车运行过程中的轮轨力难以直接测量,却又直接关乎车辆运营中的列车运行安全性、牵引制动、轮轨磨损等问题[1-3]. 因此,实现轮轨力的在线识别对提高列车运维品质、保障列车安全运行有着不可或缺的作用.
研究人员从载荷识别的角度开展了大量轮轨力识别研究:Xia等[4]建立基于车辆逆动力学模型的灰盒模型以识别轮轨力,该方法只适用于测量低频轮轨力;Wei等[5]提出基于轮对准静态方程和轴箱加速度的轮轨力反演方法;Zhu等[6]提出一种基于杜哈梅积分和轴箱加速度的时域轮轨力识别方法;Sun等[7]提出一种基于动力学逆模型的轮轨力识别方法,但应用该方法需要已知线路的不平顺激励;王明猛等[8]提出一种以车辆结构响应为输入,通过逆结构滤波识别轮轨力的方法. 以上研究尚未考虑轨道车辆中止挡、减振器等结构的非线性特性.
Ward等[9]提出将轮轨力视作系统的内在状态,引入状态识别方法并采用卡尔曼滤波(Kalman filter,KF)对半车模型中轮轨力进行识别. 对于非线性过程状态识别问题,基于卡尔曼滤波已经发展出扩展卡尔曼滤波(extended Kalman filter,EKF)、无迹卡尔曼滤波(unscented Kalman filter,UKF)和容积卡尔曼滤波(cubature Kalman filters,CKF)等改进算法,这些算法在目标跟踪[10-12]、电池荷电状态估计[13-14]和导航[15-16]等领域的非线性状态识别问题中已有较为成熟的应用. 其中,容积卡尔曼滤波[17]的改进思路是基于三阶球面径向容积准则逼近贝叶斯滤波理论框架中的后验均值和方差. 与扩展卡尔曼滤波相比,容积卡尔曼滤波对非线性函数均值和方差的近似避免了计算非线性系统的雅可比矩阵,有效降低了估计误差和计算量. 此外,容积卡尔曼滤波较无迹卡尔曼滤波在高维度和强非线性问题上有更好的估计效果.
平方根容积卡尔曼滤波(square root cubature Kalman filter,SRCKF)[18]是一种将协方差矩阵的平方根代替误差协方差矩阵,以增强滤波稳定性的容积卡尔曼滤波方法. 本文建立考虑悬挂元件非线性的车辆-轨道垂向耦合模型,导出包含轮轨垂向力的车辆状态非线性递推函数,采用SRCKF算法结合车体、构架和轮对振动加速度识别轮轨垂向力. 为进一步验证所提出算法在实际工程中的可行性,建立整车动力学模型及其对应的17自由度轮轨垂向力估计模型,分析车辆模型在实际不平顺激扰下的左右侧轮轨垂向力识别结果.
1. 动力学建模与识别算法
1.1 车辆-轨道垂向耦合动力学模型
车辆-轨道垂向耦合动力学模型是铁路动力学系统研究中的基础模型[19]. 如图1所示,模型中车辆部分由车体、构架和轮对组成,考虑部件的浮沉和点头运动. 图中:Fsi为构架i的二系悬挂力,i=1,2;Fpj为轮对j的一系悬挂力,j=1,2,3,4;FNj为轮对j的轮轨垂向力;zwj为轮对j的浮沉自由度;zti为构架i的浮沉自由度;zc为车体的浮沉自由度;βc为车体的点头自由度;βti为构架i的点头自由度.
车辆模型中各部件之间的悬挂元件采用弹簧-阻尼单元模拟. 实际上,车辆系统中部分悬挂元件具有如图2所示的非线性特性,因此,模型中采用此类非线性阻尼特性模拟一系悬挂元件的阻尼特性. 图中:δ为间断点对应速度值,阻尼力${F_{\text{d}}}(\dot x) $如式(1)所示.
Fd(˙x)={C1˙x,|˙x|⩽δ,C1δ+C2(|˙x|−δ)sgn(˙x),|˙x|>δ, (1) 式中:Ci为阻尼系数,$ \dot x $为相对速度.
模型中轮轨垂向力${F_{\text{N}}} $采用Hertz接触理论[20]计算,如式(2)所示.
FN = {(δzwr/G)3/2,δzwr>0,0,δzwr⩽0, (2) 式中:G为轮轨接触常数;$ {\delta _{{{\textit{z}}_{{\text{wr}}}}}} $为轮轨接触点处的法向弹性压缩量,当$ {\delta _{{{\textit{z}}_{{\text{wr}}}}}} $< 0时轮轨分离.
轨道模型由钢轨、轨枕、道床和路基组成,除路基固定不动外其余部件只考虑其垂向运动,各部件之间采用弹簧-阻尼单元模拟,其中,钢轨采用连续弹性离散点支承的简支Euler梁模拟.
1.2 SRCKF算法识别原理
对于上述车辆-轨道垂向耦合动力学模型中的10自由度车辆系统,状态变量如式(3).
{\boldsymbol{x}} = ({{\textit{z}}_{\mathrm{c}}},{\beta _{\mathrm{c}}},{{\textit{z}}_{{\mathrm{t}}i}},{\beta _{{\text{t}}i}},{{\textit{z}}_{{\text{w}}j}},{\dot {\textit{z}}_{\mathrm{c}}},{\dot \beta _{\mathrm{c}}},{\dot {\textit{z}}_{{\mathrm{t}}i}},{\dot \beta _{{\text{t}}i}},{\dot {\textit{z}}_{{\text{w}}j}}). (3) 将车辆振动系统中车轮受到的轮轨垂向力看作系统的内在状态,从而可将待识别的状态变量x扩维为
{\boldsymbol{x}} = ({{\textit{z}}_{\mathrm{c}}},{\beta _{\mathrm{c}}},{{\textit{z}}_{{\mathrm{t}}i}},{\beta _{{\text{t}}i}},{{\textit{z}}_{{\text{w}}j}},{\dot {\textit{z}}_{\mathrm{c}}},{\dot \beta _{\mathrm{c}}},{\dot {\textit{z}}_{{\mathrm{t}}i}},{\dot \beta _{{\text{t}}i}},{\dot {\textit{z}}_{{\text{w}}j}},{F_{{\text{N}}j}}). (4) 此时,车辆系统的状态空间方程为
\left\{ {\begin{array}{*{20}{l}} {\dot {\boldsymbol{x}} = f({\boldsymbol{x}}) + {\boldsymbol{w}}}, \\ {{\boldsymbol{y}} = h({\boldsymbol{x}}) + {\boldsymbol{v}}} , \end{array}} \right. (5) 式中:y为系统观测量,f(x)、h(x)分别为非线性状态转移函数和量测函数,w为协方差矩阵Q的过程噪声,v为协方差矩阵R的量测噪声.
扩维状态变量的一阶微分为
\dot {\boldsymbol{x}} = ({\dot {\textit{z}}_{\mathrm{c}}},{\dot \beta _{\mathrm{c}}},{\dot {\textit{z}}_{{\mathrm{t}}i}},{\dot \beta _{{\text{t}}i}},{\dot {\textit{z}}_{{\text{w}}j}},{\ddot {\textit{z}}_{\mathrm{c}}},{\ddot \beta _{\mathrm{c}}},{\ddot {\textit{z}}_{{\mathrm{t}}i}},{\ddot \beta _{{\text{t}}i}},{\ddot {\textit{z}}_{{\text{w}}j}},{\dot F_{{\text{N}}j}}), (6) 其中,加速度项和轮轨垂向力的微分项可通过式(7)得到.
\left\{\begin{array}{*{20}{l}}m_{\mathrm{c}}\ddot{\text{z}}_{\mathrm{c}}=F_{\text{s1}}+F_{\text{s2}}, \\ I_{\mathrm{cy}}\ddot{\beta}_{\mathrm{c}}=(F_{\text{s1}}-F_{\text{s2}})l_{\mathrm{c}}, \\ m_{\mathrm{t}}\ddot{\text{z}}_{\mathrm{t}i}=F_{\text{p}(2i-1)}+F_{\text{p}(2i)}-F_{\text{s}i}, \\ I_{\mathrm{ty}}\ddot{\beta}_{\text{t}i}=(F_{\text{p}(2i-1)}-F_{\text{p}(2i)})l_{\mathrm{t}}, \\ m_{\text{w}}\ddot{\text{z}}_{\text{w}j}=-F_{\text{p}j}-P_{\text{st}}+F_{\text{N}j}, \\ \dot{F}_{\text{N}j}=0,\end{array}\right. (7) 式中:mc、mt和mw分别为车体、构架和轮对的质量,Icy、Ity分别为车体和构架的点头惯量,lc、lt分别为车辆定距、轴距的一半,Pst为车辆静轴重.
为实现轮轨垂向力的识别更新,将车辆系统状态方程离散为式(8),并使用四阶龙格库塔法实现递推式的迭代更新.
\boldsymbol{x}_{k+1}=F(\boldsymbol{x}_k), (8) 式中:xk为迭代过程第k步的状态变量,F(·)为时间离散下状态变量x的递推函数.
车辆系统中各个部件的加速度相比位移和速度更容易测得,故而选取车体垂向加速度、构架垂向加速度和轮对垂向加速度作为观测量,如式(9).
\boldsymbol{y}=(\ddot{{\textit{z}}}_{\mathrm{c}},\ddot{{\textit{z}}}_{\mathrm{t}i},\ddot{{\textit{z}}}_{\text{w}j}). (9) 结合式(5)~(9)和SRCKF算法,即可设计轮轨垂向力的状态观测器. 应用SRCKF算法时应给出初始化的状态量$ {\hat {\boldsymbol{x}}_{0|0}} $和误差协方差矩阵$ {{{\boldsymbol{P}}}_{0|0}} $,采用第k步识别的状态变量$ {\hat {\boldsymbol{x}}_{k|k}} $递推识别第k + 1步的状态变量$ {\hat {\boldsymbol{x}}_{k + 1|k + 1}} $(下标$ k|k $表示第k步的后验估计,$ k + 1|k $表示第k + 1步的先验估计). 首先,进行时间更新,需要采用$ {\hat {\boldsymbol{x}}_{k|k}} $和$ {{{\boldsymbol{P}}}_{k|k}} $计算容积点$ \boldsymbol{X}_{k|k}^{\left(a\right)} $(式(10)),将容积点代入非线性过程函数F(·)计算传播后的容积点$ \boldsymbol{X}_{k+1|k}^{\left(a\right)} $(式(11));再计算第k步的状态量预测值$ {\hat {\boldsymbol{x}}_{k + 1|k}} $(式(12))和误差协方差矩阵的平方根$ {{{\boldsymbol{S}}}_{k + 1|k}} $(式(13)).
\boldsymbol{X}_{k|k}^{\left(a\right)}=\boldsymbol{S}_{k|k}\boldsymbol{\xi}_a+\hat{\boldsymbol{x}}_{k|k}, (10) \boldsymbol{X}_{k+1|k}^{\left(a\right)}=F(\boldsymbol{X}_{k|k}^{\left(a\right)}), (11) \hat{\boldsymbol{x}}_{k+1|k}=\frac{1}{m}\sum\limits_{a=1}^m\boldsymbol{X}_{k+1|k}^{\left(a\right)}, (12) {{{\boldsymbol{S}}}_{k + 1|k}} = {\text{Tria}}([{{\boldsymbol{\chi}} }_{k + 1|k}^* \;{{{\boldsymbol{S}}}_{Q_k}}]), (13) 式中:$ {{{\boldsymbol{S}}}_{k|k}} $为$ {{{\boldsymbol{P}}}_{k|k}} $乔列斯基分解得到的上三角矩阵;$ {{{\boldsymbol{\xi}} }_a} $为基本的容积点集,a=1, 2, …, m,m为容积点的总数且为状态变量维数的2倍;Tria(·)表示对矩阵进行三角化并获得矩阵的方阵[16];${\boldsymbol{S}}_{Q_k} $为第k步协方差矩阵Qk的平方根;$ {{\boldsymbol{\chi}} }_{k + 1|k}^* $为中心加权矩阵.
然后,进行量测更新,使用预测值$ {\hat {\boldsymbol{x}}_{k + 1|k}} $更新容积点(式(14)),将更新后的容积点$ \boldsymbol{X}_{k+1|k}^{\left(a\right)} $代入量测方程计算传播后的容积点$ \boldsymbol{Y}_{k+1|k}^{\left(a\right)} $(式(15)),再计算第k+1步的量测预测值$ {\hat {\boldsymbol{y}}_{k + 1|k}} $(式(16))、误差协方差矩阵平方根$ {{{\boldsymbol{S}}}_{{\mathrm{yy}},k + 1|k}} $(式(17))和协方差矩阵$ {{{\boldsymbol{P}}}_{{\mathrm{xy}},k + 1|k}} $(式(18)).
\boldsymbol{X}_{k+1|k}^{\left(a\right)}=\boldsymbol{S}_{k+1|k}\boldsymbol{\xi}_a+\hat{\boldsymbol{x}}_{k+1|k}, (14) \boldsymbol{Y}_{k+1\mid k}^{\left(a\right)}=h\left(\boldsymbol{X}_{k+1\mid k}^{\left(a\right)}\right), (15) \hat{\boldsymbol{y}}_{k+1|k}=\frac{1}{m}\sum\limits_{a=1}^m\boldsymbol{Y}_{k+1|k}^{\left(a\right)}, (16) \boldsymbol{S}_{\mathrm{yy},k+1|k}=\text{Tria}([\boldsymbol{\zeta}_{k+1|k}\; \boldsymbol{\ S}_{\boldsymbol{R}_k}]), (17) {{{\boldsymbol{P}}}_{{\mathrm{xy}},k + 1|k}} = {{{\boldsymbol{\chi}} }_{k + 1|k}}{{\boldsymbol{\zeta}} }_{k + 1|k}^{\mathrm{T}}, (18) 式中:${{{\boldsymbol{S}}}_{{\boldsymbol{R}}_k}} $为第k步协方差矩阵Rk的平方根,$ {{{\boldsymbol{\zeta}} }_{k + 1|k}} $和$ {{{\boldsymbol{\chi }}}_{k + 1|k}} $均为中心加权矩阵.
最后,进行状态更新,分别计算第k+1步的卡尔曼增益Kk+1、状态变量识别量$ \hat{{\boldsymbol{x}}}_{k + 1| k + 1} $和误差协方差平方根识别值$ {{{\boldsymbol{S}}}_{k + 1|k + 1}} $,如式(19)~(21)所示.
\boldsymbol{K}_{k+1}=(\boldsymbol{P}_{\mathrm{xy},k+1|k}/\boldsymbol{S}_{\mathrm{yy},k+1|k}^{\mathrm{T}})/\boldsymbol{S}_{\mathrm{yy},k+1|k}, (19) {\hat {\boldsymbol{x}}_{k + 1|k + 1}} = {\hat {\boldsymbol{x}}_{k + 1|k}} + {{{\boldsymbol{K}}}_{k + 1}}({{\boldsymbol{y}}_{k + 1|k}} - {\hat {\boldsymbol{y}}_{k + 1|k}}), (20) \boldsymbol{S}_{k+1\mid k+1}=\mathrm{Tria}\left(\left[\boldsymbol{\chi}_{k+1\mid k}-\boldsymbol{K}_{k+1}\boldsymbol{Y}_{k+1\mid k}^{\left(a\right)}\ \ \boldsymbol{K}_{k+1}\boldsymbol{S}_{\boldsymbol{R}_k}\right]\right). (21) 2. 算法验证
2.1 车辆-轨道垂向耦合动力学模型验证
为使得车辆-轨道垂向耦合动力学模型能够真实反映其动力学行为,动力学仿真时域积分步长设置为
0.0001 s;同时,动力学仿真的采样步长和观测器状态迭代更新的时间步长均设置为0.005 s,以体现实际中在线观测器计算能力的限制. 为减小观测器时间更新过程中的误差,需要将误差协方差矩阵与轮轨垂向力相关的参数设置为较大的值. 观测器中状态变量初值、误差协方差矩阵、系统噪声、量测噪声等参数的详细设置如式(22)所示. 在此条件下分析车辆系统以80 km/h的运行速度在随机不平顺、钢轨波磨不平顺和钢轨焊缝不平顺作用下的轮轨垂向力识别结果.\qquad\qquad\qquad\left\{ \begin{gathered} {{\hat x}_{0|0}} = (\underbrace {0, \cdot \cdot \cdot ,0,}_{20}{P_{{\text{st}}}},{P_{{\text{st}}}},{P_{{\text{st}}}},{P_{{\text{st}}}}), \\ {{{\boldsymbol{P}}}_{0|0}} = {\text{diag}}(\underbrace {1 \times{10^{ - 2}}, \cdot \cdot \cdot ,1\times{10^{ - 2}}}_{24}), \\ {{\boldsymbol{Q}}} = {\text{diag}}(\underbrace {1\times{10^{ - 3}}, \cdot \cdot \cdot ,1\times{10^{ - 3}},}_6\underbrace {1\times{10^{ - 6}}, \cdot \cdot \cdot ,1\times{10^{ - 6}},}_4 \\ \quad \underbrace {1\times{10^{ - 3}}, \cdot \cdot \cdot ,1\times{10^{ - 3}},}_6\underbrace {1\times{10^{ - 6}}, \cdot \cdot \cdot ,1\times{10^{ - 6}},}_41\times{10^9},1\times{10^9},1\times{10^9},1\times{10^9}), \\ {{\boldsymbol{R}}} = {\text{diag}}(\underbrace {1\times{10^{ - 8}}, \cdot \cdot \cdot, 1\times{10^{ - 8}}}_7). \\ \end{gathered} \right. (22) 图3为车辆在随机不平顺激扰下运行时一位轮对轮轨垂向力的时域响应、功率谱密度和均方根误差. 由图3可以看出:轮轨垂向力的识别曲线和仿真采样的轮轨垂向力曲线十分接近,两者的相关系数为0.988;功率谱密度显示轮轨垂向力在高于3 Hz时有较好的识别效果;轮轨垂向力识别值的均方根误差在识别开始后快速减小,稳态均方根误差的均值和方差分别为0.175 kN和0.027 kN2.
钢轨波磨不平顺[21]和钢轨焊缝不平顺[22]作为常见的轮轨确定性激扰,会显著加剧车辆运行时的轮轨相互作用. 其中,钢轨波磨不平顺一般为谐波叠加形成,单一谐波的激励输入可以描述为
{{\textit{z}}_0}(t) = \frac{1}{2}a_{\mathrm{g}}(1 - \cos \omega t),\begin{array}{*{20}{c}} {}&{0 \leqslant t \leqslant \dfrac{\lambda }{v}} , \end{array} (23) \omega = \frac{{2{\text{π}} v}}{\lambda }, (24) 式中:z0(·)为车辆以一定速度v通过一谐波激励在时刻t的激励输入函数,ag为钢轨波磨波深,ω为钢轨波磨频率,λ为钢轨波磨波长.
仿真中钢轨波磨不平顺设置为单个频率的谐波,波深和波长分别为0.1 mm和0.5 m. 在运行速度为80 km/h的条件下,该轮轨激励对应的轮轨垂向力频率f0为44.44 Hz. 图4为车辆在钢轨波磨不平顺下运行时一位轮对轮轨垂向力的时域响应、功率谱密度和均方根误差. 由图4可以看出:轮轨力识别曲线和仿真采样的轮轨力曲线基本一致,两者相关系数达到0.999;轮轨垂向力识别值的功率谱密度可以反映钢轨波磨对应的轮轨垂向力频率和峰值;轮轨垂向力识别值的均方根误差同样在识别初期快速减小,稳态均方根误差值的均值和方差分别为0.142 kN和0.028 kN2.
钢轨焊缝不平顺一般为2个单一谐波叠加形成,其激励输入可以描述为
\begin{gathered} {{\textit{z}}_{{\text{w0}}}}(s) = \\ \left\{ {\begin{array}{*{20}{l}} {\dfrac{1}{2}{a_1}\left(1 - {\mathrm{cos}}\;\dfrac{{2{\text{π}} {{s}}}}{{{\lambda _1}}}\right),\;\; 0 \leqslant s < \dfrac{{{\lambda _1} - {\lambda _2}}}{2}或 \dfrac{{{\lambda _1} + {\lambda _2}}}{2} \leqslant s < {\lambda _1}}, \\ {\dfrac{1}{2}{a_1}\left(1 - {\mathrm{cos}}\;\dfrac{{2{\text{π}} s}}{{{\lambda _1}}}\right) + \dfrac{1}{2}{a_2}\left[1 - {\mathrm{cos}}\;\dfrac{{\text{π}} }{{{\lambda _2}}}(2s - {\lambda _1} + {\lambda _2})\right],\;} \\ \quad {\dfrac{{{\lambda _1} - {\lambda _2}}}{2} \leqslant s < \dfrac{{{\lambda _1} + {\lambda _2}}}{2}} , \end{array}} \right. \end{gathered} (25) 式中:s为轨道方向的里程坐标,λ1、λ2分别为钢轨焊缝不平顺长波和短波的波长,a1、a2分别为钢轨焊缝不平顺长波和短波的波深.
仿真中,钢轨焊缝不平顺长波的波深和波长分别为0.3 mm和1.0 m,短波的波深和波长分别为0.3 mm和0.1 m. 图5为车辆在钢轨焊缝不平顺下运行时一位轮对轮轨垂向力的时域响应、功率谱密度和均方根误差. 由图5可知:轮轨力识别曲线和仿真采样的轮轨力曲线的十分接近,两者的相关系数为0.969;轮轨垂向力识别值在频域内钢轨焊缝不平顺引起的宽频响应,轮轨垂向力在高于6 Hz时有较好的识别效果;车辆通过钢轨焊缝处不平顺时轮轨垂向力识别值均方根误差不大,结果表明:均方根误差值的均值和方差分别为0.101 kN和0.004 kN2.
2.2 整车三维动力学模型验证
实际工程中车辆各个部件均有6个运动自由度且运行时会受到轨道的多向激扰. 为进一步验证本文所提出算法,采用SIMPACK软件建立如图6所示的某型地铁车辆动力学模型,基于本文提出算法对该模型仿真的轮轨垂向力进行估计.
为区分车辆运行时的左右侧轮轨垂向力,在10自由度垂向车辆模型的基础上,考虑车体、构架和轮对的侧滚自由度建立17自由度车辆动力学模型,车辆各部件侧滚振动方程如式(26)所示. 基于此模型建立与1.2节中类似的非线性过程函数,在车体、构架和轮对上选取非质心位置处的垂向加速度作为观测量,如式(26)所示.
\left\{ \begin{array}{l} I_{\mathrm{c}x} \ddot{\phi}_{\mathrm{c}}=\left(F_{\mathrm{szL} 1}+F_{\mathrm{szL} 2}-F_{\mathrm{szR} 1}-F_{\mathrm{szR} 2}\right) d_{\mathrm{s}}, \\ I_{\mathrm{t}x} \ddot{\phi}_{\mathrm{t}i}=\left(F_{\mathrm{pzL} 2 i-1}+F_{\mathrm{pzL} 2 i}-F_{\mathrm{pzR} 2 i-1}-F_{\mathrm{pzR} 2 i}\right) d_{\mathrm{w}}, \\ I_{\mathrm{w}x} \ddot{\phi}_{\mathrm{w} j}=\left(F_{\mathrm{pzL} j}-F_{\mathrm{p} \mathrm{zR} j}\right) d_{\mathrm{w}}+\left(F_{\mathrm{N} \mathrm{l} j}-F_{\mathrm{N} \mathrm{r} j}\right) a_0, \end{array}\right. (26) 式中:$ \ddot \phi $为侧滚自由度,下标c、t和w分别为车体、构架和轮对;Icx、Itx和Iwx分别为车体、构架和轮对的侧滚惯量;FszLi、FszRi分别为架构i左、右侧的二系悬挂力;FpzLj、FpzRj分别为轮对j左、右侧的一系悬挂力;dw、ds分别为一系和二系悬挂横向跨距的一半;a0为滚动圆跨距的一半;FNLj、FNRj分别为左、右侧的轮轨垂向力.
图7展示了整车模型在随机不平顺激扰下运行时的轮轨垂向力估计结果. 图7(a)在时域中对比了一位轮对左侧轮轨垂向力的识别和仿真结果,左侧轮轨垂向力识别值的均方根误差均值和方差分别为0.394 kN和0.067 kN2. 图7(b)在频域中对比一位轮对左侧轮轨垂向力功率谱密度的识别和仿真结果,由图可知轮轨垂向力在高于6 Hz时有较好的识别效果. 从图7(c)中可以看到整车模型中左、右侧轮轨垂向力相关系数最低值分别为0.747和0.720,属于显著相关;左右侧轮轨垂向力之和相关系数最低能达到0.999,属于高度相关,证实了所提出算法的有效性. 分开估计左、右侧轮轨垂向力时,估计结果的相关系数相对低于合并估计轮轨垂向力时的相关系数,其原因是整车模型采用的过程函数相对于整车模型做了一定简化,并未考虑轮轨接触点位置和轮轴横向力波动的影响.
3. 结 论
1) 基于SRCKF算法识别10自由度车辆在随机不平顺激励的轮轨垂向力,识别结果与仿真结果时域曲线的相关系数可达0.988.
2) 对于钢轨波磨不平顺引起的轮轨力大幅波动信号和钢轨焊缝不平顺引起的轮轨力宽频信号,采用本文方法识别轮轨垂向力在时域和频域均有不错识别结果,2种激励下识别结果的相关系数分别为0.999和0.969.
3) 基于17自由度车辆模型可有效识别整车模型左、右侧轮轨垂向力,仿真结果表明左、右侧轮轨垂向力相关系数最低分别能达到0.747和0.720,左、右侧轮轨垂向力之和相关系数能达到0.999.
-
-
[1] 刘永锋,赵旭峰,栗杨等. 轮轨界面低黏着对大功率电力机车车轮滚动接触疲劳的影响[J]. 铁道学报,2022,44(12): 22-29.LIU Yongfeng, ZHAO Xufeng, LI Yang, et al. Effect of adhesion between wheel and rail on wheel rolling contact fatigue of high-power AC locomotives[J]. Journal of the China Railway Society, 2022, 44(12): 22-29. [2] 杨飞,吴细水,孙宪夫,等. 钢轨轧制不平顺激扰下的动车组动力响应特性[J]. 西南交通大学学报,2022,57(2): 267-276,294.YANG Fei, WU Xishui, SUN Xianfu, et al. Dynamic response characteristics of EMU under excitation of rail straightening irregularity[J]. Journal of Southwest Jiaotong University, 2022, 57(2): 267-276,294. [3] 施以旋,戴焕云,毛庆洲等. 基于车轨耦合的地铁车轮多边形形成机理研究[J]. 西南交通大学学报,2024,59(6): 1357-1367,1388.SHI Yixuan, DAI Huanyun, MAO Qingzhou, et al. Study on formation mechanism of metro wheel polygonal based on vehicle-track coupling[J]. Journal of Southwest Jiaotong University, 2024,59(6): 1357-1367,1388. [4] XIA F J, COLE C, WOLFS P. Grey box-based inverse wagon model to predict wheel-rail contact forces from measured wagon body responses[J]. Vehicle System Dynamics, 2008, 46(S1): 469-479. [5] WEI L, ZENG J, WU P B, et al. Indirect method for wheel–rail force measurement and derailment evaluation[J]. Vehicle System Dynamics, 2014, 52(12): 1622-1641. doi: 10.1080/00423114.2014.953180 [6] ZHU T, WANG X R, FAN Y W, et al. A time domain method for wheel-rail force identification of rail vehicles[J]. Vehicle System Dynamics, 2022, 60(3): 790-809. doi: 10.1080/00423114.2020.1838562 [7] SUN Y Q, COLE C, SPIRYAGIN M. Monitoring vertical wheel-rail contact forces based on freight wagon inverse modelling[J]. Advances in Mechanical Engineering, 2015, 7(5): 1-11. [8] 王明猛,朱涛,王小瑞,等. 一种逆结构滤波法的轨道车辆轮轨力识别[J]. 振动工程学报,2019,32(4): 602-608.WANG Mingmeng, ZHU Tao, WANG Xiaorui, et al. An inverse structural filter method for wheel-rail contact forces identification of railway vehicles[J]. Journal of Vibration Engineering, 2019, 32(4): 602-608. [9] WARD C P, GOODALL R M, DIXON R, et al. Adhesion estimation at the wheel–rail interface using advanced model-based filtering[J]. Vehicle System Dynamics, 2012, 50(12): 1797-1816. doi: 10.1080/00423114.2012.707782 [10] 苏义鑫,赵俊,张华军. 带有UKBF的船舶动力定位预测控制器设计[J]. 西南交通大学学报,2018,53(3): 589-595.SU Yixin, ZHAO Jun, ZHANG Huajun. Predictive controller with UKBF for marine dynamic positioning system[J]. Journal of Southwest Jiaotong University, 2018, 53(3): 589-595. [11] 宋义彤,舒红宇,陈仙宝等. 分布式驱动电动汽车双无迹卡尔曼滤波状态参数联合观测[J]. 机械工程学报,2019,55(22): 93-102. doi: 10.3901/JME.2019.22.093SONG Yitong, SHU Hongyu, CHEN Xianbao, et al. State and parameters estimation for distributed drive electric vehicle based on unscented Kalman filter[J]. Jouranl of Mechanical Engineering, 2019, 55(22): 93-102. doi: 10.3901/JME.2019.22.093 [12] 肖业方,郭蕴华,高海波,等. 基于多传感器鲁棒容积卡尔曼融合的半潜式起重平台立柱快速压载舱液位估计[J]. 船舶力学,2023,27(1): 70-80.XIAO Yefang, GUO Yun Hua, GAO Haibo, et al. Level estimation of column side ballast tank of semi-submersible crane vessel based on multi-sensor robust cubature Kalman fusion[J]. Journal of Ship Mechanics, 2023, 27(1): 70-80. [13] 杜常清,吴中意,武冬梅等. 基于KF-EKF算法的动力电池SOC估计[J]. 武汉理工大学学报,2022,44(4): 84-92.DU Changqing, WU Zhongyi, WU Dongmei, et al. State of charge estimation of power battery based on KF-EKF algorithm[J]. Jouranl of Wuhan University of Technology, 2022, 44(4): 84-92. [14] 肖仁鑫,李斌,黄志强,等. 基于SRCKF的电动汽车锂离子电池荷电状态估计[J]. 电源技术,2021,45(11): 1443-1447.XIAO Renxin, LI Bin, HUANG Zhiqiang, et al. State of charge estimation of electric vehicle lithium ion battery based on SRCKF[J]. Chinese Journal of Power Sources, 2021, 45(11): 1443-1447. [15] 乔美英,高翼飞,李宛妮,等. 基于模糊鲁棒自适应CKF算法的MEMS-IMU姿态估计[J]. 中国惯性技术学报,2022,30(3): 296-303.QIAO Meiying, GAO Yifei, LI Wanni, et al. Attitude estimation of MEMS-IMU based on fuzzy robust adaptive CKF algorithm[J]. Journal of Chinese Inertial Technology, 2022, 30(3): 296-303. [16] 章司怡,陈熙源. 运动约束辅助的基于SVD-CKF的组合导航方法[J]. 电子测量与仪器学报,2022,36(4): 82-89.ZHANG Siyi, CHEN Xiyuan. Motion constraint aided integrated navigation method based on SVD-CKF[J]. Journal of Electronic Measurement and Instrumentation, 2022, 36(4): 82-89. [17] ARASARATNAM I, HAYKIN S. Cubature Kalman filters[J]. IEEE Transactions on Automatic Control, 2009, 54(6): 1254-1269. doi: 10.1109/TAC.2009.2019800 [18] 穆静,蔡远利. 平方根容积卡尔曼滤波算法及其应用[J]. 兵工自动化,2011,30(6): 11-13,26.MU Jing, CAI Yuanli. Square root cubature Kalman filter algorithm and application[J]. Ordnance Industry Automation, 2011, 30(6): 11-13,26. [19] 翟婉明. 车辆-轨道耦合动力学理论的发展与工程实践[J]. 科学通报,2022,67(32): 3793-3807.ZHAI Wanming. Development of vehicle-track coupled dynamics theory and engineering practice[J]. Chinese Science Bulletin, 2022, 67(32): 3793-3807. [20] 翟婉明. 车辆-轨道耦合动力学下册[M]. 4版. 北京: 科学出版社,2015: 41-42. [21] 陈光雄. 钢轨波磨预测模型验证工况的研究[J]. 西南交通大学学报,2022,57(5): 1017-1023,1054. doi: 10.3969/j.issn.0258-2724.20200842CHEN Guangxiong. Study on validation conditions of rail corrugation prediction models[J]. Journal of Southwest Jiaotong University, 2022, 57(5): 1017-1023,1054. doi: 10.3969/j.issn.0258-2724.20200842 [22] 陈美,翟婉明,閤鑫,等. 高速铁路多边形车轮通过钢轨焊接区的轮轨动力特性分析[J]. 科学通报,2019,64(25): 2573-2582. doi: 10.1360/N972018-01288CHEN Mei, ZHAI Wanming, GE Xin, et al. Analysis of wheel-rail dynamic characteristics due to polygonal wheel passing through rail weld zone in high-speed railways[J]. Chinese Science Bulletin, 2019, 64(25): 2573-2582. doi: 10.1360/N972018-01288 -