随机轮轨力作用下基于2.5维离散支撑模型的轨道垂向振动分析

朱志辉，王 凡，罗思慧，李 奇，蒋丽忠

(1.中南大学 土木工程学院，湖南 长沙 410075；2.中南大学 高速铁路建造技术国家工程实验室，湖南 长沙 410075；3.同济大学 土木工程学院， 上海 200092)

在以往的钢轨振动研究中，多用格林函数法结合Euler-Bernoulli 梁或Timoshenko 梁理论表示钢轨的垂向动力学行为，通常假设钢轨受等效连续支撑。Remington[1]采用位移导纳及传递率方法求解轨道作为Euler梁在弹性地基支撑上的解。徐志胜等[2]使用模态叠加法建立了车辆-轨道垂向耦合振动频域分析模型，分析了轮轨高频振动的基本特征。孙文静等[3-4]用格林函数法分别建立了基于两层连续支撑Timoshenko 梁模型和单层离散支撑Euler梁模型的车辆-轨道耦合频域分析模型，孙宇等[5]基于格林函数法建立了车辆-轨道耦合时域分析模型，但因未考虑扭力，梁模型不能很好地预测钢轨的横向响应。格林函数法所建梁模型只能用中性轴代替整个钢轨截面的位移，且无法准确模拟高频时的横截面变形。传统的3维有限元方法可以考虑剪切变形和横截面变形的影响，但3维模型的规模较大，对于需要高频分析的问题，计算效率低。而短实体模型研究的钢轨长度有限，不能全面分析钢轨在各频率的振动传递特性和振型。刘林芽等[6]、吴天行等[7]建立了钢轨的实体有限元模型，发现纵向单元长度为0.05 m时才能满足计算精度。

近年来，一些学者提出了2.5维有限元方法，用于表示具有二维几何形状且在第三方向不变的结构。采用二维有限元网格对截面进行建模和常规形状函数假设，在垂直截面方向采用波数解。这种方法已被用于研究轨道、薄壁梁[8]、肋板加筋板[9]和汽车轮胎[10]等结构。早期相关工作由Knothe等[11]和Gavric[12]完成，已获得自由轨道上传播的自由波。Ryue等[13]利用2.5维有限元和边界元耦合模型确定了在80 kHz内连续支撑轨道上传播的波；Nilsson等[14]计算了在一个点力激励下连续支撑轨道的振动和声辐射。但在这些研究中，边界条件只包括单层轨道垫支撑。Zhang等[15]建立了多层支撑轨道模型，轨枕采用柔性梁模拟，一定程度上降低了计算效率。

本文基于弹性理论，建立无限长2.5维有限元钢轨模型，采用荷载等效法模拟离散支撑并求得钢轨的振动响应。模型考虑了钢轨几何断面和载荷的非对称性，且只对钢轨结构进行二维离散，计算效率高。通过与格林函数结果进行对比，验证了本文模型的正确性。结合分步求解法和柔度矩阵法求解轨道高低不平顺引起的垂向动态轮轨力，运用虚拟激励法及叠加法快速求得轮轨力功率谱；将虚拟轮轨力施加到离散支撑的2.5维轨道模型上，计算得到轮轨多点接触工况下的轨道系统随机振动响应。利用该模型对钢轨横截面不同位置、沿钢轨纵向不同距离、不同垫板刚度下的钢轨振动响应进行了分析。

1 随机轮轨力功率谱

1.1 轮轨力计算公式

单节移动高速车辆模型考虑为具有二系悬挂的10自由度多刚体模型，如图1所示，沿轨道以恒定速度v行驶。其中：mc和Jc分别为车体质量和转动惯量；mt和Jt分别为转向架质量和转动惯量；mw为轮对质量；zc和θc分别为车体沉浮和点头自由度；zti(i=1,2,…,4)和θti(i=1,2)分别为前后转向架沉浮和点头自由度；zwi为4个轮对沉浮自由度；kpz和cpz分别为一系悬挂刚度和阻尼；ksz和csz分别为二系悬挂刚度和阻尼；lc为车辆定距之半；lt为转向架轮对定距之半。

图1 车轨耦合系统模型

车辆系统运动微分方程可表示为

(1)

式中：MV、CV、KV分别为车辆系统的质量、阻尼和刚度矩阵；UV(t)为系统的位移向量；FV(t)为由轨道不平顺随机激励产生的车辆动荷载向量。动荷载、动位移可分别表示为

(2)

将式(2)代入式(1)得

(KV+iωCV-ω2MV)UV(ω)=FV(ω)

(3)

频域内车辆系统车轮位移点的导纳可表示为

(4)

式中：AV(xij)为荷载作用在j车轮，所得到i车轮的传递导纳，AV为一节车厢在轮对处的导纳矩阵。

采用线性赫兹接触模型来模拟轮轨关系[16]，轮轨接触导纳矩阵AC可表示为

(5)

式中：Uδ为轮轨接触变形位移；kH为接触弹簧刚度，由轮轨预加载荷和曲率半径通过赫兹接触理论得到。

车辆-轨道耦合系统的激励包含4个轮轨接触点的激励，采用叠加法[17]得到轮轨接触点处车轮和轨道的振动响应UV和Ur，可分别表示为

UV=-AVFV

Ur=-ArFC

(6)

式中：Ar为轨道导纳矩阵。

单位简谐荷载下，离散支撑轨道位置x处的导纳Ar(xij)为

(7)

式中：Ar(xij)为荷载作用在j车轮，所得到i车轮处钢轨的传递导纳；Ar为在一节车厢四个轮对处钢轨的导纳矩阵。

轮轨接触变形位移可表示为

Uδ=UV-Ur-r

(8)

根据柔度矩阵法，以轨道不平顺r作为输入，根据式(4)、式(5)和式(7)结果，得到轮轨力计算表达式为

(9)

1.2 基于虚拟激励法的轮轨力功率谱

采用高效的虚拟激励法求解轮轨力功率谱。同一轨道不平顺谱下，不同轮轴所受激励的相位不同，激励可表示为

(10)

式中：Δr(t)为轨道不平顺平稳随机激励的时间历程函数；n为假设系统所受不同相位激励点数，n=4m,m为车厢数量。假设Δr(t)的时域自功率谱密度为SX(ω)，则由虚拟激励法可得相应虚拟激励表达式为

(11)

简化可得

(12)

空间域功率谱密度函数SX(ω)和时间域功率谱密度函数SX(k)的转换关系为

SX(ω)=SX(k)/c

(13)

式中：c为波速，m/s；k为激励空间波数，rad/m。

(14)

(15)

2 离散支撑的2.5维轨道动力学模型

2.1 自由边界2.5维钢轨运动方程

假设钢轨在y-z平面上的横截面恒定，沿x方向上材料特性、截面形状不变，在圆频率ω上的简谐运动用隐式时间依赖eiωt假设。将钢轨视为弹性体，由Hamilton原理得

δU(ω)-δT(ω)-δW(ω)=0

(16)

式中：δ表示一阶变分；U和T分别为系统的势能和动能；δW为结构上的虚功。

将保守力产生体系势能、体系总动能和非保守力所做功的一阶变分代入式(16)得

(17)

式中：ε为材料中的应变振幅；u为节点位移；ρm为重力密度；f为钢轨所受外载荷；V对应一个2.5维有限元单元；D为材料刚度矩阵。

应变可用位移表示为

(18)

式中：B0为包含y、z方向的偏导；B1为与x方向偏导有关的参数。位移u和虚位移δu可近似表示为

(19)

式中：N(y,z)为八节点等参四边形单元的二维有限元形状函数。将式(18)和式(19)代入式(17)，结构在圆频率ω的简谐激励下的运动方程可表示为[14]

(20)

式中：

采用高斯求积法求解二维截面上的积分。

为了得到波数域内的运动方程，对x坐标进行傅里叶变换

(21)

式中：k为x方向上的波数；～表示波数域中的傅里叶变换性质。将式(20)转换到波数域可表示为

(22)

钢轨在单频率荷载下的位移可通过对k进行积分等到

(23)

当在xe=0处施加幅值为1的单位集中力时，自由钢轨任意位置x处的位移导纳为

(24)

2.2 考虑离散支撑的2.5维轨道运动方程

图2 离散支撑轨道动力学模型

假设轨道受圆频率为ω的简谐激励作用，支撑点i的位移可以表示为[15]

(25)

通过阻抗Z将每个支撑点的力与相应的位移联系起来,即

(26)

式中：

(27)

其中：sp为轨道垫线刚度；sb为道砟线刚度。

各点竖向位移可表示为

(28)

在所有支撑点处xj(j=1,…,N)按式(28)进行计算，则可以得到支撑位置处的一组线性方程，求解方程可得支撑点位移，再次利用式(28)求解得到任意位置x处的位移。由此得到钢轨在外荷载和支撑荷载共同作用下任意位置的位移响应。

联立式(15)和式(7)，2.5维离散支撑轨道系统响应功率谱可表示为

(29)

2.3 计算流程图

根据上述理论，基于Matlab平台编制了高速列车引起的2.5维有限元轨道高频振动特性求解程序，主要计算流程如图3所示。

图3 计算流程图

2.4 2.5维离散支撑轨道动力学模型验证

研究表明[18]，当轮轨接触点与轨道边界的距离不小于10倍轨枕距时，轨道长度对频响结果影响较小。本文选取131个轨枕离散支撑的轨道模型，以第65个轨枕上方钢轨及第64个与第65个轨枕跨距中间钢轨的轨顶中心为单位简谐荷载作用点，利用2.5维有限元法求解该轨道模型，分别计算得到轨枕上方钢轨及轨枕跨距中间的原点位移导纳。本文选用CHN60钢轨进行分析，所取参数与文献[7]一致，钢轨及轨道支撑参数如表1所示。

表1 CHN60有砟轨道参数

其中钢轨横截面包括20个平面8节点单元，共有99个节点和297个自由度，钢轨横截面网格划分情况如图4所示。

图4 横截面网格划分图

图5为2.5维有限元法与格林函数法计算结果对比。

图5 2.5维有限元法与格林函数法结果对比

由图5可知，频率80 Hz处出现第一个峰值，此频率是轨枕和道砟组成的系统的共振频率；频率530 Hz处出现第二个峰值，此频率是轨枕和扣件组成的系统的共振频率，此时钢轨和轨枕受垫板刚度影响异相振动；频率250 Hz处出现1个反共振峰值，这是由于钢轨受双层弹簧阻尼结构支撑，轨枕充当了钢轨的动力吸振器，此处轨枕振幅较大。频率1 090 Hz和4 210 Hz两处出现1阶和2阶“pinned-pinned”共振现象，力作用在轨枕正上方时响应出现最小值；力作用在两轨枕中间时响应出现峰值，此时钢轨中的垂向弯曲振动波波长λ分别等于2倍和1倍轨枕间距。

3 离散支撑2.5维轨道随机振动分析

3.1 计算参数

采用Sato高低短波不平顺及我国高速铁路轨道不平顺联合谱作为输入[19]。不平顺波长小于2 m时采用Sato谱，大于2 m时采用我国高速铁路轨道谱。轨道结构计算参数如表1所示，车辆采用CRH3型高速动车，计算参数如表2所示。

表2 CRH3车辆参数

3.2 轮轨系统响应

车辆轨道定点激励，轨道不平顺以列车运行速度在列车和钢轨之间移动，形成相对位移激励，采用虚拟激励法，在车速为300 km/h时对离散支撑轨道模型进行随机振动分析。当轮轨接触点位于轨枕上方及轨枕跨距中间时，两者除 “pinned-pinned”共振频率外结果基本相同。以下分析中均假设轮轨接触点位于轨枕上方，求解车辆-轨道系统垂向随机振动响应。

图6为轮轨耦合系统导纳幅频曲线，10～35 Hz低频段内车轮导纳幅值大于钢轨和轮轨接触导纳幅值，说明在该范围内车轮导纳在轮轨耦合系统的振动响应中起主要影响；35～1 820 Hz中频段内主要由钢轨导纳决定轮轨耦合系统导纳；1 820～5 000 Hz高频段内轮轨耦合系统导纳主要由钢轨与接触导纳决定。

图6 轮轨系统导纳幅频曲线

图7为车辆的车体中心、转向架、车轮和轨道轨顶中心的振动加速度功率谱。车轮和钢轨的振动加速度功率谱在低于50 Hz时基本一致，在高频振动时差别逐渐明显。由于一系、二系悬挂的减振作用，轮轨力及车轮振动在传递至转向架及车体的过程中存在一定程度上的衰减，并且频率越高减振效果越好，衰减速率越快。车体在中高频处均有多个峰值和谷值，其中峰值点与车体的激励频率有关，谷值点与车辆的几何滤波现象有关。

图7 车-轨系统加速度功率谱

3.3 钢轨横截面不同点响应

为分析高频随机振动对钢轨横截面变形的影响，选取钢轨横截面5个采样点进行分析，采样点选取示于图4。单位谐荷载作用下钢轨位移导纳结果如图8(a)所示，共振频点包括250、530、1 120、3 000、4 210 、4 860 Hz。低于250 Hz时钢轨截面各部位同相位整体垂向振动，未产生振动变形；250～3 000 Hz内轨底翼缘振幅小于其他部位；高于3 000 Hz时因根部翘曲运动导致轨底弹性变形逐渐显著，轨底翼缘振幅最大。高于1 200 Hz时钢轨截面出现全面的变形，不同频段各采样点振幅的大小关系为：1 200～3 000 Hz内轨顶>轨底>轨腰>轨头翼缘>轨底翼缘；3 000～4 210 Hz内轨底翼缘>轨顶>轨底>轨头翼缘>轨腰；4 210～5 000 Hz内轨底翼缘>轨顶>轨头翼缘>轨底>轨腰。

图8(b)为一节车厢运行下钢轨横截面不同点的振动加速度功率谱。低于1 000 Hz时5个采样点加速度功率谱密度曲线变化一致，高于2 800 Hz时轨脚与其他四点出现明显差别。4 800 Hz处存在明显的谷值，各采样点达到谷值的顺序为轨脚、轨顶中心、轨头脚部、轨腰、轨底中心。45 Hz处加速度功率谱密度最大，此时钢轨振动最为剧烈。

图8 x=0处钢轨横截面不同位置处响应

3.4 沿钢轨方向不同距离处响应

车辆-轨道耦合系统的激励包含4个轮轨接触点的激励，因而需计算沿轨道纵向不同距离的钢轨垂向位移幅值特性，分析钢轨垂向振动沿钢轨纵向传递过程中的变化。

图9(a)为单位激励点一侧沿钢轨纵向不同距离(分别距激励点0、2.4、17.4、19.8 m)轨顶中心的垂向位移导纳，100 Hz以下各截面振动的幅值差别明显；70～160 Hz和300～1 090 Hz间钢轨传递导纳与原点导纳的差值逐渐变小，160～300 Hz间差值又显著增大。“pinned-pinned”频率之后钢轨上出现波的振动，幅值衰减且衰减幅度增大。振动传递到距激励点17.4、19.8 m处时幅值曲线基本吻合，证明已经衰减稳定，但是在100～600 Hz因为轨道垫、轨枕和道砟的系统的影响，导致幅值出现较为明显的振动。低频振动在传递到距激励点2.4 m幅值已经衰减了20倍左右，传递到距激励点17.4 m衰减了107倍左右，在荷载作用点附近几跨衰减更明显。

图9 轨顶中心沿钢轨纵向不同距离处响应

图9(b)为一节车厢运行下四个轮对处的振动加速度功率谱，低于20 Hz时轮对二和轮对四加速度功率谱略大于轮对一和轮对三，20～50 Hz时轮对一、二、三加速度功率谱略大于轮对四。中频段四个轮对处的加速度功率谱密度基本相同，1 200～4 100 Hz时轮对二、三振动较轮对一、四剧烈。

3.5 不同垫板刚度的影响

为分析不同垫板刚度对钢轨垂向振动的影响，选取三种不同刚度的垫板(柔性板180 MN/m、半刚性板350 MN/m、刚性板600 MN/m)。

单位谐荷载作用下轨顶中心位移导纳结果如图10(a)所示，190 Hz以下低频段垫板刚度的影响较小，该频段特性主要由道砟刚度决定。190～1 090 Hz时随频率增加，位移幅值先由一阶谷值达到二阶峰值，然后逐渐衰减到“pinned-pinned”频率处，该频段为垫板刚度主要影响频段。垫板刚度影响一阶谷值和二阶峰值大小及频率，刚度越大，一阶谷值和二阶峰值越小，两处频率越高。二阶峰值随着垫板刚度增加分别为340、480、610 Hz，此时钢轨和轨枕在垫板刚度作用下异相振动。在“pinned-pinned”频率处垫板刚度越大，钢轨的振动幅值越小，但是对该频率大小没有影响。1 090 Hz以上高频范围保持较好的一致性，在1 110、3 000、4 210、4 860 Hz共振频率处，垫板刚度越大，对应频率略大，幅值越大，但是随频率增加垫板刚度的影响逐渐减小。垫板刚度对导纳在全频域内产生影响。

图10(b)为一节车厢运行下轨顶中心振动加速度功率谱，170 Hz以上垫板刚度对加速度产生影响。170～300 Hz时刚度越大，幅值越大；300～1 090 Hz时刚度越小幅值越大，但是差值逐渐减小；1 090～3 000 Hz时刚度越大，钢轨振动越剧烈；3 000 Hz后垫板刚度影响变小。

图10 x=0处不同垫板刚度下轨顶中心响应

4 结论

(1)本文将2.5维有限元法所建立的离散支撑轨道应用于车轨耦合系统随机振动过程，将钢轨视为实体单元，建立离散支撑弹性轨道模型，采用2.5维有限元法及荷载等效法得到全频域上离散支撑钢轨任意位置的频响特性，结合高速CRH3型车辆动力学模型，以真实轨道谱为激励，使用虚拟激励法和荷载叠加法快速精确地对轨道结构进行随机振动分析。

(2)单位谐响应分析中，1 200～5 000 Hz主要产生三种不同形式横截面变形；振动传递到距激励点17.4、19.8 m处时，在100～600 Hz因为轨道垫、轨枕和道砟的系统的影响，导致幅值出现较为明显的振动；垫板刚度对导纳在全频域内产生影响。

(3)随机振动分析中，高于2 80 0Hz轨道结构的横截面开始明显变形；低于50 Hz时四个轮对处的振动加速度功率谱密度差别显著；在170～3 000 Hz垫板刚度对轨道结构振动加速度影响较大。

(4)2.5维有限元法解决了格林函数法不能考虑横截面变形的问题，包含2.5维离散支撑弹性轨道的模型更加符合实际情况；与3维有限元法相比具有更高的计算效率。随着列车速度提高，轮轨激振频率增加，更加符合实际情况且为今后轨道结构高频振动研究及耦合系统减振降噪提供依据。