Application of Good Lattice Point with Power Generator Method in Stochastic Dynamic Analysis of Vehicle-Bridge System
-
摘要:
针对难以精确选取具有代表性的超高维随机相位角问题，采用方幂好格子点法生成代表性轨道不平顺样本，并将该样本作用于车-桥系统，得到轮轨力的均值与标准差；然后，通过对比虚拟激励法、确定性时程法和蒙特卡洛法的计算结果来探究该方法的计算精度与计算效率；最后，采用线性与非线性轮轨接触关系研究考虑列车日运营量的脱轨系数阈值. 以和谐号通过桥梁为例，计算结果表明：与蒙特卡洛法相比，采用方幂好格子点法生成不同方向的轨道不平顺样本之间的均匀性较好；方幂好格子点法求得的随机动力响应的概率特征参数与其他方法相比，具有较高的计算精度，其计算效率较蒙特卡洛法提高了近5倍；分别采用线性与非线性轮轨接触关系时所得的脱轨系数阈值相差达4.68%，方幂好格子点法具有较广泛适用性.
Abstract:Since it is difficult to accurately select representative super-high-dimensional random phase angles, the good lattice point with power generator method (GLPPGM) was utilized to generate samples of representative track irregularities, which were applied to the vehicle-bridge system to obtain the mean and standard deviation of random dynamic responses. Then, the calculation accuracy and efficiency of this method were explored by comparing the results of the pseudo-excitation method, deterministic time history method, and Monte Carlo method (MCM). Finally, the threshold value of the derailment factor considering the daily operation volume of trains was studied by using linear and nonlinear wheel-rail contact relationships. The Harmony train passing over a bridge was studied, and the results show that compared with that by the MCM, the uniformity between the samples of track irregularities generated by the GLPPGM in different directions is better. The probability characteristic parameters of the random dynamic response obtained by the GLPPGM have higher calculation accuracy than different methods, and its calculation efficiency is nearly five times higher than that of the MCM. When linear and nonlinear wheel-rail contact relationships are considered, the threshold value of the derailment factor differs by 4.68%, and the GLPPGM has a wider applicability.
-
高速铁路的快速发展使车-桥耦合振动下列车运行安全受到人们普遍重视，尤其是考虑随机荷载[1-2]. 轨道不平顺作为主要随机外荷载，直接影响列车的运行安全. 因此，探究基于随机轨道不平顺的列车运行安全具有较大的工程意义.
随机振动研究方法众多，传统方法有蒙特卡洛法、拉丁法、增强法与子集模拟法[3]等. 随着随机理论的发展，学者们又提出了随机摄动法、虚拟激励法[4]、概率密度演化法. 蒙特卡洛法理论比较简单，Xu等[5-6]采用此方法研究了车-轨-桥系统与车-轨-路面系统的动力响应. 但是，蒙特卡洛法计算量大，计算成本高，因此，上述研究基于给定或少量的轨道不平顺样本，研究结果具有随机性. 林家浩等[4]提出虚拟激励法；朱艳等[7-9]采用该方法研究车-桥系统随机振动，并求得桥梁结构安全限值；Zhang等[10]采用此理论研究行车安全指标的概率分布，为行车安全研究提供新思路. 虚拟激励法计算效率高，但只能应用于线性系统，而采用线性轮轨关系研究车桥系统随机振动存在一定的局限性. 对于非线性随机问题，李杰等[11]提出概率密度演化理论. 基于该理论，Yu等[12]考虑桥梁参数的随机性与随机轨道不平顺，研究桥梁结构的安全；Xu等[13]采用概率密度演化法建立轨道不平顺的概率模型，并采用该理论进行了车体、轨道结构等的可靠性研究. 概率密度演化法依赖于数论选点法，因此，适用范围受选点法优劣的影响.
行车安全指标由轮轨力得到，而轮轨力是高频响应，其值与轨道不平顺样本密切相关，对样本选取要求十分严格. 对于轨道不平顺随机相位角的降维方法主要采用加强正交基条件来实现. 该方法对于仅考虑一个方向轨道不平顺而言，有较好的降维效果[14]. 但是，对于空间轨道不平顺，正交基之间相互独立. 因此，需要一种可直接选取超高维随机变量的方法. 方幂好格子点法是近些年提出的均匀设计方法，可以高效地生成大量随机数，其随机数相互独立、具有较高精度且均匀性较好. 该方法可以处理超高维随机变量间的均匀性问题. Liu等[15]采用此方法模拟了多元随机平稳过程.
本文采用方幂好格子点法选取超高维随机相位角，生成具有代表性的轨道不平顺样本； 然后，采用竖向密贴与横向简化Kalker线性蠕滑轮轨接触关系求得随机车-桥系统动力响应的均值与标准差，并通过比较方幂好格子点法与虚拟激励法、蒙特卡洛法和确定性时程法计算结果的拟合程度，验证该方法的精确性与高效性； 最后，采用竖向密贴与沈式修正横向Kalker线性蠕滑轮轨接触关系研究考虑列车日运营量的脱轨系数阈值，验证该方法的适用性.
1. 车-桥系统振动模型
1.1 车-桥系统模型
每节列车由1个车体、2个转向架与4个轮对构成，均被视为刚体，分别由线性弹簧和黏滞阻尼连接. 车体（转向架）分别有5个自由度：竖向zc（zt）、横向yc（yt）、扭转θc（θt）、点头φc（φt）和摇头ψc（ψt）. 当采用线性轮轨关系时，轮对仅考虑3个自由度：竖向zw、横向yw、扭转θ w. 当采用非线性轮轨关系时，轮对考虑4个自由度：竖向、横向、扭转和摇头（ψw）. 车辆模型如图1所示，其详细描述可参考文献[13].
本文采用直接刚度法建立桥梁动力学方程. 车辆和桥梁分别视为独立子系统，二者通过轮轨力实现耦合. 车-桥耦合振动方程如下：
\left\{ {\begin{array}{*{20}{c}} {{{\boldsymbol{M}}_{\text{v}}}{{\ddot {\boldsymbol{X}}}_{\text{v}}} + {{\boldsymbol{C}}_{\text{v}}}{{\dot {\boldsymbol{X}}}_{\text{v}}} + {{\boldsymbol{K}}_{\text{v}}}{{\boldsymbol{X}}_{\text{v}}} = {{\boldsymbol{F}}_{\text{v}}}}, \\ {{{\boldsymbol{M}}_{\text{b}}}{{\ddot {\boldsymbol{X}}}_{\text{b}}} + {{\boldsymbol{C}}_{\text{b}}}{{\dot {\boldsymbol{X}}}_{\text{b}}} + {{\boldsymbol{K}}_{\text{b}}}{{\boldsymbol{X}}_{\text{b}}} = {{\boldsymbol{F}}_{\text{b}}}}, \end{array}} \right. (1) 式中：M、K和C分别为系统的质量、刚度和阻尼矩阵，{\boldsymbol{X}}、\dot {\boldsymbol{X}}、\ddot {\boldsymbol{X}}和F分别为系统的位移、速度、加速度和外荷载向量，下标v和b分别表示车辆和桥梁子系统.
车-桥系统的动力响应采用全过程迭代法求解[16]. 对车-桥系统做如下假设：1） 忽略轨道模型，将轨道结构与桥梁视为整体，二者间无相对位移；2） 轮轨力荷载只分配到相邻的桥梁节点，按照静力等效法分配；3） 列车匀速直线运行.
1.2 轮轨接触关系
线性轮轨关系采用竖向密贴与横向简化Kalker线性蠕滑理论；非线性轮轨关系采用竖向密贴与沈式修正横向Kalker线性蠕滑理论.
1.3 车-桥系统动力响应
由于方幂好格子点法选取超高维代表性点的优势和计算车桥系统随机响应的计算精度与计算效率，本文采用此方法研究桥梁中加速度动力响应和轮轨力动力响应. 最后，讨论线性与非线性轮轨关系对脱轨系数阈值的影响来体现此方法的优势[17].
2. 车-桥系统动力响应的概率特征参数及脱轨系数阈值的求解
车-桥系统随机因素包含系统参数和轨道不平顺等. 通常，不平顺对列车影响较为明显. 因此，本文仅考虑随机轨道不平顺.
2.1 虚拟激励法
车-桥系统是时变系统，因此，随机动力响应具有时变特性. 根据虚拟激励法，随机轮轨力{F_r}为
{F_r}(w,t) = {\boldsymbol{\varGamma}} (t){\boldsymbol{V}}(t){\boldsymbol{G}}(t)\sqrt {{{\boldsymbol{S}}_r}(w)} {{\mathrm{e}}^{{\text{i}}wt}}, (2) 式中：w、t分别为频率与时间，r\in \{ {\textit{z}}、y, 0 \} ， {\boldsymbol{\varGamma}} (t) 、{\boldsymbol{G}}(t)、{\boldsymbol{V}}(t)分别为轮轨力指示向量、调制函数与时滞矩阵，S为不平顺功率谱.
将式（2）代入式（1）即可得到系统的虚拟响应{\boldsymbol{u}}，进而求得随机动力响应的功率谱{{\boldsymbol{S}}_{\text{u}}}.
则方差计算式为
{\boldsymbol{\sigma}} _{}^2(t) = \sum\limits_{i = 1}^{\Delta m} {{{\boldsymbol{S}}_{\text{u}}}({w_i},t)} {\text{d}}w\text{,} (3) 式中：σ为响应的标准差，wi为第i分量的频率大小，{\text{d}}w为频率增量，\Delta m 为频率增量数量.
2.2 方幂好格子点法[18]
假设a为生成元，npt为随机变量代表点总数，则有
{a^{m + 1}} \text{mod} {n_{{\text{pt}}}}= 1, (4) 式中：mod(·)为取余数运算符号， m为a对npt的次数.
假设超立方体为s维，当a与npt互质且a对npt的次数不小于s−1时，方幂生成向量ha可以写为式（5），且ha中元素均互不相等.
{{\boldsymbol{h}}_{\text{a}}} = ({a^0},{a^1}, \cdots ,{a^{s - 1}})(\text{mod} {\text{ }}{n_{{\text{pt}}}}). (5) 由于本文随机变量的维数较高，即s较大，因此建议a与npt分别取质数来减小随机数之间的偏差. 随机点xj,k为
{\begin{split} & {{x_{j,k}} = \dfrac{{{h_{{\text{a}}j}}k}}{{{n_{{\text{pt}}}}}} - {{\mathrm{int}}} \left(\dfrac{{{h_{{\text{a}}j}}k}}{{{n_{{\text{pt}}}}}}\right)}, \\ &\quad {k = 1,2, \cdots ,{n_{{\text{pt}}}},{\text{ }}j = 1,2, \cdots ,s}, \end{split}} (6) 式中：haj为ha中的第j个元素.
构造的随机点都是单位超立方体{[0,1]^{\mathrm{s}}}内均匀分布的散布随机点，即
{x_{j,k}} \in [0,1]. (7) 随机长度的超立方体[−l, l]s内均匀分布的散布随机点 {\theta _{j,k}} 为
{\theta _{j,k}} = 2({x_{j,k}} - 0.5) {\text{π}}. (8) 通常而言，{n_{{\text{pt}}}}会很大，不利于计算. 因此，可以采用超球体法[12]对随机点进行筛选，如式（9）.
{\sum\limits_{j = 1}^s {{{({x_{j,q}} - 0.5)}^2} \leqslant {{\left(\frac{{{r_0}}}{2}\right)}^2}} },\quad {q = 1,2, \cdots ,{N_{{\text{pt}}}},{\text{ }}} (9) 式中：r0为有界半径系数， 对于高维随机点而言，分布主要集中在超球体范围内，因此，通常假设r0=1；Npt为经过筛选后的随机变量代表点总数.
将式（9）代入式（8）可得
\begin{split} & {\theta _{j,q}} = 2({x_{j,q}} - 0.5) \lambda . \end{split} (10) 当常数 \lambda 等于π时，可以求得随机相位角代表点. 本文采用三角级数法，通过轨道不平顺功率谱生成轨道不平顺时域样本{r_{r}}{\boldsymbol{(\varTheta}} ,x) ，如式（11）.
{r_r}({\boldsymbol{\varTheta}} ,x) = \sum\limits_{i = 1}^m {\sqrt {2{\text{π }} {{\boldsymbol{S}}_r}({\varOmega_i})} } \cos ({\varOmega _i} x + \phi _{_{ir}}({\boldsymbol{\varTheta}} )), (11) 式中：{\boldsymbol{\varTheta}} 为轨道不平顺的随机变量，{\boldsymbol{\varTheta}} = ({\xi _1},{\xi _2}, \cdots , {\xi {{_{q}}}},\cdots,{\xi _{{N_{{\text{pt}}}}}})；x为轨道不平顺的长度，即x向的距离；\xi_{q} 为随机相位角的离散点集；Sr(·)为轨道不平顺功率谱；{\varOmega _{i}}为第i分量的空间频率；\phi _{_{ir}}({\boldsymbol{\varTheta}} )为第i分量频率的随机相位角.
假设{\xi _q} = \{ {\phi _{q,1}},{\phi _{q,2}}, \cdots ,{\phi _{q,s}}\} \in {\boldsymbol{\varTheta}} ，表示随机相位角的代表性点集，每个点集的概率为Pq，则车-桥系统的概率{ P} _\varTheta (·)可由式（12）求得
\sum\limits_{q = 1}^{{N_{{\text{pt}}}}} {\int_{{\xi _q}} {{{{P}}_\varTheta }} } ({\xi _q}){\text{d}}\xi = \sum\limits_{q = 1}^{{N_{{\text{pt}}}}} {{P_{q}}} = 1. (12) 式中：{ P} _\varTheta (·)为车-桥系统的概率.
根据概率守恒定律[12]，每个点集获得的随机动力响应概率也为Pq. 车-桥耦合随机动力方程为
\left\{ {\begin{array}{*{20}{l}} {\boldsymbol{M}}{_{\text {v}}} \ddot X_{\text {v}}({\boldsymbol{\varTheta}} ,t) + C{_{\text {v}}}\dot X_{\text {v}}({\boldsymbol{\varTheta}} ,t) + K_{\text {v}}X_{\text {v}}({\boldsymbol{\varTheta}} ,t) = {\boldsymbol{F}}_{\text {v}}({\boldsymbol{\varTheta}} ,t) , \\ {\boldsymbol{M}}{_{\text {b}}} \ddot {\boldsymbol{V}}_{\text {b}}({\boldsymbol{\varTheta}} ,t) + C{_{\text {b}}}\dot {\boldsymbol{V}}_{\text {b}}({\boldsymbol{\varTheta}} ,t) + K_{\text {b}}X_{\text {b}}({\boldsymbol{\varTheta}} ,t) = {\boldsymbol{F}}_{\text {b}}({\boldsymbol{\varTheta}} ,t) . \end{array}} \right. (13) 假设 {\boldsymbol{U}}(\varTheta ,t) 为车-桥系统随机状态向量，如式（14）.
\begin{split} & {\boldsymbol{U}}({\boldsymbol{\varTheta}} ,t) = G[\ddot {\boldsymbol{X}}({\boldsymbol{\varTheta}} ,t),\dot {\boldsymbol{X}}({\boldsymbol{\varTheta}} ,t),{\boldsymbol{X}}({\boldsymbol{\varTheta}} ,t)]=\\ &\quad [U_1(\xi_1,t),U_2(\xi_2,t),\cdots,U_q(\xi_q,t)), \end{split} (14) 式中：G(·)为传递函数.
随机动力响应的均值与标准差分别如式（15）和式（16）.
u(t) = {P_{q}} \sum\limits_{q = 1}^{{N_{{\text{pt}}}}} {{{E}}(U_q({\xi _q},t))} , (15) \begin{split} & \sigma^2(t)=P_q\sum\limits_{q=1}^{N_{\text{pt}}}(U_q(\xi_q,t)-E(U_q(\xi_q,t)))^2= \\ &\quad P_q\sum\limits_{q=1}^{N_{\text{pt}}}(U_q(\xi_q,t)-u(t))^2, \end{split} (16) 式中：E(·)为随机序列期望.
随机响应的协方差与相关系数分别为
\begin{split} & {{\mathrm{cov}}} (t) = \sum\limits_{q = 1}^{{N_{{\text{pt}}}}} {(E({U_{\text{1}q}}({\xi _q},t) - {u_1}(t){\text{)}} {\text{(}}{U_{2q}}({\xi _q},t) - {u_2}(t){\text{)}}} {\text{) }} {\text{ = }} \\&\quad {P_{q}} \sum\limits_{q = 1}^{{N_{{\text{pt}}}}} (({U_{1q}}({\xi _q},t) - {u_1}(t){\text{)}} {\text{(}}{U_{2q}}({\xi _q},t) - {u_2}(t)) ) , \end{split} (17a) {\rho _{12}}(t) = \frac{{{{\mathrm{cov}}} (t)}}{{{\sigma _1}(t) {\sigma _2}(t)}} , (17b) 式中：{\sigma _1} (·)、{\sigma _2} (·)为2个具有相关性的随机动态响应，U_{1q} 、U_{2q} 分别为2个不同状态向量中的元素.
忽略轨道不平顺求得车-桥系统随机动力响应的均值（确定性时程法），计算仅需一次，结果具有确定性. 蒙特卡洛法依据大数定理，当样本足够多时，计算结果会趋于精确解. 本文采用虚拟激励法、确定性时程法和蒙特卡洛法作为对照组.
2.3 脱轨系数阈值的求解
文献[13]通过试验与理论计算发现轮轨力近似服从正态分布. 为了研究需要，本文假设每时刻的轮轨力服从高斯分布，竖向力恒大于0. 根据文献[19]，当横、竖向轮轨力均服从高斯分布，且竖向轮轨力恒大于0时，脱轨系数服从：
\begin{split} & {{{F}}_X}{\text{(}}x{{) = P}}\left(\frac{{{F_{{y}}}}}{{{F_{\textit{z}}}}} < x\right)\approx \\&\quad \varPhi \left(\frac{{{u_{{F_{\textit{z}}}}}x - {u_{{F_{{y}}}}}}}{{\sqrt {{\sigma _{{F_{{y}}}}} + {\sigma _{{F_{\textit{z}}}}}{x^2} - 2x{\rho _{\text{F}}}{\sigma _{{F_{{y}}}}}{\sigma _{{F_{\textit{z}}}}}} }}\right) \end{split}, (18) 式中：FX(·)为概率分布函数； \varPhi (·)为标准正态分布函数；Fy、Fz分别为横、竖向轮轨力； {u_{{F_{{y}}}}} （ {u_{{F_{\textit{z}}}}} ）、 {\sigma _{{F_{y}}}} （ {\sigma _{{F_{\textit{z}}}}} ）、 {\rho _{\text{F}}} 分别为横（竖）向轮轨力的均值、标准差与相关系数.
已知随机变量底分布为F(x)，则底分布的极大值概率分布函数为
{{{F}}_{{\max }}}(y) = F^n(x) , (19) 式中：{{{F}}_{{{\max }}}}为极大值概率分布函数，n为极大值样本数量.
将式（18）代入式（19），可得脱轨系数最大值概率分布函数为
{{{F}}_{{\text{Dmax}}}}(x) = {\varPhi^n \left(\frac{{{u_{{F_{\textit{z}}}}}x - {u_{{F_{{y}}}}}}}{{\sqrt {{\sigma _{{F_{{y}}}}} + {\sigma _{{F_{\textit{z}}}}}{x^2} - 2x{\rho _{\text{F}}}{\sigma _{{F_{\textit{z}}}}}{\sigma _{{F_{\textit{z}}}}}} }}\right)} . (20) 根据脱轨系数的极值概率分布函数和设定的置信度可获得脱轨系数的阈值.
3. 算例验证
本文采用和谐号通过5跨简支梁为例进行方法验证. 车辆参数见表1；桥梁参数见文献[20]. 采用德国低干扰谱，车速v假设为270 km/h. 列车上桥前先运行100 m，出桥后再行驶50 m. 高低、水平与方向不平顺的随机相位角数量均为32个，采用方幂好格子点法得到代表性相位角. 当考虑线性轮轨接触时，生成199条轨道不平顺样本，采用虚拟激励法、确定性时程法和蒙特卡洛法生成
1 000 条轨道不平顺样本作为对照组. 采用非线性轮轨接触时，采用方幂好格子点法生成4 003 条轨道不平顺样本，确定性时程法和蒙特卡洛法生成5 000 条轨道不平顺样本作为对照组. 对比不同轮轨接触关系时的脱轨系数阈值的差异. 工况示意如图2所示.表 1 车辆参数Table 1. Train parameters名称 拖车 动车 名称 拖车 动车 转向架轴距/m 2.5 2.5 x 向转向架转动惯量/（t•m2） 3.2 1.6 车辆定距/m 18 18 y 向转向架转动惯量/（t•m2） 7.2 1.7 一系悬挂横向跨距/m 0.33 0.34 z 向转向架转动惯量/（t•m2） 6.8 1.7 转向架导一系悬挂/m 2.05 2.05 一系纵向阻尼（单侧）/（kN•s•m−1） 960 960 二系悬挂横向跨距/m 2.05 2.05 一系横向阻尼（单侧）/（kN•s•m−1） 960 960 车体中心到二系悬挂/m 0.36 0.83 一系竖向阻尼（单侧）/（kN•s•m−1） 1040 700 二系悬挂到转向架/m 0.24 0.15 二系纵向阻尼（单侧）/（kN•s•m−1） 240 210 轮对滚动园半径/m 0.4575 0.4575 二系横向阻尼（单侧）/（kN•s•m−1） 240 210 轮对质量/t 1.9 2.2 二系竖向阻尼（单侧）/（kN•s•m−1） 400 350 轮对转动惯量/（t•m2） 1.067 1.63 一系纵向阻尼（单侧）/（kN•s•m−1） 0 0 转向架质量/t 3.4 1.7 一系横向阻尼（单侧）/（kN•s•m−1） 0 0 x 向车体转动惯量/（t•m2） 101.5 74 一系竖向阻尼（单侧）/（kN•s•m−1） 30 38 y 向车体转动惯量/（t•m2） 1064.4 1370 二系纵向阻尼（单侧）/（kN•s•m−1） 120 150 z 向车体转动惯量/（t•m2） 867.2 1370 二系横向阻尼（单侧）/（kN•s•m−1） 30 15 车体质量/t 42.4 44 二系竖向阻尼（单侧）/（kN•s•m−1） 33 40 3.1 方幂好格子点法生成轨道不平顺样本的偏差
标准差是一种用来衡量数据偏离平均值程度的统计量，时程样本间均匀性越好，其标准差时程曲线越平稳. 通常假设轨道不平顺是零均值高斯过程，因此，本文采用其标准差来验证代表性样本的均匀性.
采用方幂好格子点法生成高低、水平、方向轨道不平顺样本各199条； 然后，随机对各个方向随机样本加减，求其每计算步的标准差； 最后，通过观察标准差的平稳性来判定方幂好格子点法的优劣性. 3个方向轨道不平顺样本的标准差、一阶导数标准差、二阶导数的标准差为 \sigma_1 ({{a}} {r_{{{\textit{z}}}}} + {{b}} {r_{{y}}} + {{c}} {r_{{u}}}) 、 \sigma_2 ({{a}} {\dot r_{{{\textit{z}}}}} + {{b}} {\dot r_{{y}}} + {{c}} {\dot r_{{u}}}) 、 \sigma_3 ({{a}} {\ddot r_{{{\textit{z}}}}} + {{b}} {\ddot r_{{y}}} + {{c}} {\ddot r_{{u}}}) . 其中：a、b、c均为常数，本文随机取值；rz、ry、ru分别表示高低、水平与方向不平顺样本. 图3是随机轨道不平顺样本随机加减后标准差的计算结果.
从图3可以看出，方幂好格子点法生成199条轨道不平顺样本标准差计算结果是一条光滑直线，满足轨道不平顺平稳性假定[21]. 即使采用蒙特卡洛法生成1 000条轨道不平顺样本计算标准差，其时程曲线仍具有明显波动，其二阶导数的波动最显著. 该方法得到轨道不平顺、一阶导数、二阶导数的标准差近似等于蒙特卡洛法计算结果的均值. 因此，方幂好格子点法可以较好地选取代表性随机相位角，并获得均匀性良好的样本.
3.2 方幂好格子点法的计算精度与计算效率
采用蒙特卡洛法求得随机动力响应的均值与标准差；确定性时程法与虚拟激励法分别求得随机动力响应的均值与标准差. 对比方幂好格子点法的计算结果，通过三者的拟合程度来判断方幂好格子点法的计算精度与计算效率. 由于虚拟激励法仅适用于线性系统，因此，本节采用线性轮轨接触关系.
3.2.1 标准差的拟合
标准差的计算结果如图4所示.
图 4 车-桥系统动力响应的标准差Figure 4. Standard deviation of dynamic response of vehicle-bridge system由图4可知：采用蒙特卡洛法、虚拟激励法与方幂好格子点法求得的桥梁跨中加速度动力响应拟合效果较好；对于轮轨力而言，方幂好格子点法和虚拟激励法计算结果拟合效果相当，蒙特卡洛法得到的结果收敛于二者方法，虚拟激励法求得的结果是精确的. 因此，方幂好格子点法求解车-桥系统动力响应的标准差具有很高精度.
3.2.2 均值的拟合
均值拟合曲线如图5所示.
从图5可知：采用蒙特卡洛法、时程法与方幂好格子点法求得桥梁跨中加速度动力响应的拟合精度高；时程法和方幂好格子点法计算结果几乎吻合，与蒙特卡洛法存在较小差异；对于轮轨力而言，方幂好格子点法和时程法计算结果的拟合精度仍然较高，蒙特卡洛法得到的结果收敛于二者方法，时程法求得结果是精确的，因此，方幂好格子点法求解车-桥系统动力响应的均值具有很高的精度.
3.2.3 计算时间比较
计算采用的电脑配置为：英特尔®酷睿i7-10700 K处理器、32 G镁光DDR4内存和2 TB西数硬盘. 蒙特卡洛法计算所需时间最长，为152 min，虚拟激励法和方幂好格子点法所需时间相近，分别为34 min和31 min. 方幂好格子点法具有较高的计算效率.
3.3 方幂好格子点法的适用性
脱轨系数是单侧轮对横向力与竖向力之比. 采用Kalker蠕滑理论，横向力会受竖向力的影响，因此二者存在相关性且横向力存在非线性. 由于传统随机方法很难考虑轮轨力的相关性和非线性轮轨力标准差的求解，所以常采用线性轮轨接触关系研究其概率问题[10]. 脱轨系数的研究可以检验方幂好格子点法的适用性. 当采用虚拟激励法时，采用线性轮轨接触关系；当采用蒙特卡洛法和方幂好格子点法时，采用非线性轮轨接触关系.
由图3可知，方幂好格子点法得到轨道不平顺样本的均匀性很好. 因此，采用线性轮轨接触关系时，可以采用少量样本求得较高精度的结果. 当采用非线性系统时，计算得到的随机动力响应样本间的均匀性会发生改变. 对于低频响应（如车体、桥梁动力响应），少量样本仍然可以获得较高精度解；对于高频响应（如轮轨力），少量样本则难以获得较高精度解. 因此，计算轮轨力随机动力响应时，需要增加样本数量. 基于本文算例，采用
4003 条代表性轨道不平顺样本研究横向轮轨力的标准差与横、竖向轮轨力的相关系数.3.3.1 脱轨系数的概率特征参数
根据式（20）可知，脱轨系数的概率特征参数包括均值、标准差和相关系数. 由于不同轮轨接触关系均采用竖向密贴接触假定，因此，本节仅讨论横向轮轨力， 图6是横向力的均值与标准差，图7是横、竖向力的相关系数.
图 7 横、竖向轮轨力的相关系数Figure 7. Correlation coefficient of lateral and vertical wheel-rail forces从图6可知，方幂好格子点法计算横向力的均值与标准差是一条光滑曲线，且蒙特卡洛法的计算结果收敛于该方法. 因此，考虑复杂轮轨接触时，方幂好格子点法能够准确的求得横向力的均值与标准差. 从图7可知，方幂好格子点法求得的相关系数存在明显波动. 相关系数的计算需要先求得横向力的绝对值，此计算过程极大地削弱了响应样本间均匀性，因此，计算结果存在明显波动. 但是就结果而言，精度依旧高于蒙特卡洛法.
3.3.2 考虑列车日运营量的脱轨系数阈值
假设桥梁的使用寿命是100年，每天运营40辆列车，根据式（20）可求得极大值概率分布函数. 本文假设阈值的置信度为99.87%. 车-桥系统是时变系统，每时刻的极大值概率分布函数不一样，因此，本节选取桥梁跨中处来研究脱轨系数的阈值. 当采用蒙特卡洛法时，相关系数波动较大，因此，采用整个时刻的平均值. 图8是脱轨系数极大值分布函数.
从图8可知，采用蒙特卡洛法计算的结果为0.340，而采用虚拟激励法与本文方法计算结果分别为0.358、0.342；采用蒙特卡洛法求得脱轨系数阈值与方幂好格子点法几乎一样，而采用虚拟激励法计算结果明显大于其他工况，差异达4.68%. 根据图6与式（20）可知，二者差异主要由相关系数引起. 采用沈氏修正的Kalker线性蠕滑轮轨接触关系更接近实际情况. 因此，行车安全研究需要采用非线性轮轨接触关系.采用非线性轮轨接触关系时，方幂好格子点法仍然能准确计算出轮轨力的概率特征参数. 该方法适用性更广泛.
4. 结 论
本文通过方幂好格子点法生成代表性的轨道不平顺样本，进而得到车-桥系统随机响应的均值与标准差. 通过对比不同方法的计算结果，验证本文方法具有较高的计算精度与计算效率. 采用非线性轮轨接触关系研究考虑列车行日运营量下的脱轨系数阈值，验证该方法的适用性. 通过案例计算得到如下结论：
1） 方幂好格子点法生成的不同方向轨道不平顺样本具有良好的均匀性. 该方法生成轨道不平顺样本标准差的时程曲线与蒙特卡洛法相比，平稳性较好.
2） 方幂好格子点法具有较高的计算精度与计算效率. 通过不同方法求得的随机动力响应均值与标准差验证了其计算精度. 计算效率与蒙特卡洛法相比，提高了近5倍.
3） 方幂好格子点法可以应用于非线性轮轨接触关系. 当考虑列车日运营量时，采用线性轮轨接触关系与非线性轮轨接触关系计算的脱轨系数阈值存在一定差异. 就本文工况而言，差异达4.68%. 因此，行车安全研究需要采用轮轨非线性接触关系.
-
图 4 车-桥系统动力响应的标准差
Figure 4. Standard deviation of dynamic response of vehicle-bridge system
图 7 横、竖向轮轨力的相关系数
Figure 7. Correlation coefficient of lateral and vertical wheel-rail forces
表 1 车辆参数
Table 1. Train parameters
名称 拖车 动车 名称 拖车 动车 转向架轴距/m 2.5 2.5 x 向转向架转动惯量/（t•m2） 3.2 1.6 车辆定距/m 18 18 y 向转向架转动惯量/（t•m2） 7.2 1.7 一系悬挂横向跨距/m 0.33 0.34 z 向转向架转动惯量/（t•m2） 6.8 1.7 转向架导一系悬挂/m 2.05 2.05 一系纵向阻尼（单侧）/（kN•s•m−1） 960 960 二系悬挂横向跨距/m 2.05 2.05 一系横向阻尼（单侧）/（kN•s•m−1） 960 960 车体中心到二系悬挂/m 0.36 0.83 一系竖向阻尼（单侧）/（kN•s•m−1） 1040 700 二系悬挂到转向架/m 0.24 0.15 二系纵向阻尼（单侧）/（kN•s•m−1） 240 210 轮对滚动园半径/m 0.4575 0.4575 二系横向阻尼（单侧）/（kN•s•m−1） 240 210 轮对质量/t 1.9 2.2 二系竖向阻尼（单侧）/（kN•s•m−1） 400 350 轮对转动惯量/（t•m2） 1.067 1.63 一系纵向阻尼（单侧）/（kN•s•m−1） 0 0 转向架质量/t 3.4 1.7 一系横向阻尼（单侧）/（kN•s•m−1） 0 0 x 向车体转动惯量/（t•m2） 101.5 74 一系竖向阻尼（单侧）/（kN•s•m−1） 30 38 y 向车体转动惯量/（t•m2） 1064.4 1370 二系纵向阻尼（单侧）/（kN•s•m−1） 120 150 z 向车体转动惯量/（t•m2） 867.2 1370 二系横向阻尼（单侧）/（kN•s•m−1） 30 15 车体质量/t 42.4 44 二系竖向阻尼（单侧）/（kN•s•m−1） 33 40 
下载: 导出CSV
-
[1] 雷虎军,黄江泽. 考虑桩土相互作用的车-轨-桥系统地震响应分析[J]. 西南交通大学学报,2021,56(2): 229-237.LEI Hujun, HUANG Jiangze. Seismic responses analysis of train-track-bridge system considering pile-soil interaction[J]. Journal of Southwest Jiaotong University, 2021, 56(2): 229-237. [2] 何佳骏,向活跃,龙俊廷,等. 桥隧过渡段高速列车行车抗风安全分析[J]. 西南交通大学学报,2021,56(5): 1056-1064.HE Jiajun, XIANG Huoyue, LONG Junting, et al. Wind-resistant safety analysis of high-speed trains passing through bridge-tunnel transition[J]. Journal of Southwest Jiaotong University, 2021, 56(5): 1056-1064. [3] 王伟,张亚辉,欧阳华江. 考虑轨道不平顺随机性的车轨耦合系统动力可靠度分析[J]. 大连理工大学学报,2019,59(2): 111-117.WANG Wei, ZHANG Yahui, OUYANG Huajiang. Dynamic reliability analysis of vehicle-track coupled systems considering randomness of track irregularities[J]. Journal of Dalian University of Technology, 2019, 59(2): 111-117. [4] 林家浩,张亚辉,赵岩. 虚拟激励法在国内外工程界的应用回顾与展望[J]. 应用数学和力学,2017,38(1): 1-31. doi: 10.1007/s10483-016-2152-6LIN Jiahao, ZHANG Yahui, ZHAO Yan. The pseudo-excitation method and its industrial applications in China and abroad[J]. Applied Mathematics and Mechanics, 2017, 38(1): 1-31. doi: 10.1007/s10483-016-2152-6 [5] XU L, ZHAI W M. Probabilistic assessment of railway vehicle-curved track systems considering track random irregularities[J]. Vehicle System Dynamics, 2018, 56(10): 1552-1576. doi: 10.1080/00423114.2018.1424916 [6] XIA C Y, XIA H, DE ROECK G. Dynamic response of a train-bridge system under collision loads and running safety evaluation of high-speed trains[J]. Computers & Structures, 2014, 140: 23-38. [7] 朱艳,李小珍,强士中. 高速铁路简支梁桥车桥系统随机响应[J]. 西南交通大学学报,2011,46(4): 535-540. doi: 10.3969/j.issn.0258-2724.2011.04.002ZHU Yan, LI Xiaozhen, QIANG Shizhong. Stochastic responses of train-bridge coupling system for simply-supported high-speed railway bridge[J]. Journal of Southwest Jiaotong University, 2011, 46(4): 535-540. doi: 10.3969/j.issn.0258-2724.2011.04.002 [8] HE X H, SHI K, WU T. An efficient analysis framework for high-speed train-bridge coupled vibration under non-stationary winds[J]. Structure and Infrastructure Engineering, 2020, 16(9): 1326-1346. doi: 10.1080/15732479.2019.1704800 [9] 赵岩,张亚辉,林家浩. 列车-结构耦合作用的随机振动分析方法[M]. 北京:科学出版社,2018:156-190. [10] ZHANG N, ZHOU Z J, WU Z Z. Safety evaluation of a vehicle-bridge interaction system using the pseudo-excitation method[J]. Railway Engineering Science, 2022, 30(1): 41-56. doi: 10.1007/s40534-021-00259-6 [11] 李杰,陈建兵. 随机动力系统中的概率密度演化方程及其研究进展[J]. 力学进展,2010,40(2): 170-188. doi: 10.6052/1000-0992-2010-2-J2009-105LI Jie, CHEN Jianbing. Advances in the research on probability density evolution equations of stochastic dynamical systems[J]. Advances in Mechanics, 2010, 40(2): 170-188. doi: 10.6052/1000-0992-2010-2-J2009-105 [12] YU Z W, MAO J F. A stochastic dynamic model of train-track-bridge coupled system based on probability density evolution method[J]. Applied Mathematical Modelling, 2018, 59: 205-232. doi: 10.1016/j.apm.2018.01.038 [13] XU L, ZHAI W M, GAO J M. A probabilistic model for track random irregularities in vehicle/track coupled dynamics[J]. Applied Mathematical Modelling, 2017, 51: 145-158. doi: 10.1016/j.apm.2017.06.027 [14] 刘章军,张传勇,齐东春. 高速列车乘坐舒适度评价的概率密度演化方法[J]. 振动与冲击,2018,37(14): 149-155.LIU Zhangjun, ZHANG Chuanyong, QI Dongchun. Probability density evolution method for the ride comfort evaluation of high speed trains[J]. Journal of Vibration and Shock, 2018, 37(14): 149-155. [15] LIU Z J, LIU Z H, PENG Y B. Simulation of multivariate stationary stochastic processes using dimension-reduction representation methods[J]. Journal of Sound and Vibration, 2018, 418: 144-162. doi: 10.1016/j.jsv.2017.12.029 [16] ZHANG N, XIA H. Dynamic analysis of coupled vehicle–bridge system based on inter-system iteration method[J]. Computers & Structures, 2013, 114/115: 26-34. [17] 张楠,郭薇薇,夏禾. 高速铁路车桥耦合动力学[M]. 北京:北京交通大学出版社,2018:171-177. [18] 方开泰,刘民千,覃红,等. 均匀试验设计的理论和应用[M]. 北京:北京科学出版社,2018. [19] MARSAGLIA G. Ratios of normal variables[J]. Journal of Statistical Software, 2006, 16(4): 1-10. [20] ZHU Y, LI X Z, JIN Z B. Three-dimensional random vibrations of a high-speed-train–bridge time-varying system with track irregularities[J]. Proceedings of the Institution of Mechanical Engineers Part F Journal of Rail and Rapid Transit, 2015, 230(8): 1851-1876. [21] 翟婉明. 车辆-轨道耦合动力学(上册)[M]. 4版. 北京:科学出版社,2015:120-128. -
下载:
下载:
点击查看大图
计量
- 文章访问数: 75
- HTML全文浏览量: 36
- PDF下载量: 8
- 被引次数: 0