重载铁路有砟道床动态行为的离散元模拟与试验研究

张 徐，赵春发，翟婉明，张大伟,3

(1.西南交通大学 牵引动力国家重点实验室，四川 成都 610031；2.广东工业大学 土木与交通工程学院，广东 广州 510006；3.长安大学 汽车学院，陕西 西安 710064)

重载铁路运输具有运能大、效率高、运输成本低等优势，已逐渐成为世界各国煤炭、矿石能源以及大宗货物运输的重要方式。增大列车轴重是提高重载铁路运能的一种有效途径，但同时也会明显增加列车对线路的动力作用，加快轨道结构伤损和道床劣化，最终对行车安全造成不利影响。为保障大轴重条件下的运输安全性，近年来我国开展了轴重25～30 t重载铁路运输技术综合试验，研究了多种措施对线路基础设施进行加固和改善的效果[1-2]。其中，针对桥上和隧道内基础刚度偏大、道床厚度不足的问题，在桥上试验铺设弹性轨枕有砟轨道，在隧道内增设砟下弹性垫层，试验结果表明这些措施能够有效改善相应区段的轨道弹性[1,3]。一些室内模型试验也表明采用弹性轨枕和道砟垫有助于减小道床累积变形，降低道砟的破碎粉化等[4-5]。此外，文献[6]还介绍了采用聚氨酯固化道床技术改善重载铁路轨道的刚度和稳定性，并在山西中南部铁路通道进行了试铺试验。上述有砟道床强化技术的研究主要依赖现场和室内模型试验[3-9]。目前试验还很难获得道砟颗粒之间的接触、移动、磨损等细观力学特征，因此，仅开展试验研究难以深入揭示这些工程措施改善道床服役性能的机理以及强化后道床的性能劣化规律。若在试验研究的同时开展有砟道床力学行为的数值模拟分析，将可以更全面、更细致地研究各种措施对道床力学性能的影响，从而更好地指导重载铁路有砟轨道参数优化与强化技术研究。

事实上，国内外学者很早就尝试采用散体力学、有限单元法等开展道砟力学行为的数值模拟，但基于连续介质力学的理论模型本质上不能反映道砟作为散粒体介质的细观力学特性，其模拟结果与工程实际差别较大。近年来，一些学者将离散单元法引入铁路道砟的力学行为研究中，试图通过将离散元数值模拟与试验研究结合，从宏观与细观层面上更全面地揭示道砟的力学行为。例如，文献[10-12]建立离散元模型，模拟洁净和脏污道砟的直剪行为，分析细小粒径颗粒(煤渣、细沙、粉末等)影响道砟直剪行为的细观机制。文献[13]通过室内试验标定离散元模型中道砟颗粒的细观接触参数，模拟重复轮载作用下的道床累积变形。文献[14]建立具有真实几何外形的道砟颗粒离散元模型，研究道砟静态压碎行为及其破碎机制。这些研究表明，离散单元法能够从细观层面上较好地揭示道砟的力学行为及其变形规律，可以为有砟道床的强化改造技术研发提供理论支撑。但是，目前的研究工作主要考虑了道砟的准静态力学行为，极少关注列车动荷载作用下有砟道床的动态行为，而重载铁路工程实践表明，随着行车速度的提高和轴重的增加，重载列车动力作用引起的道床劣化越来越不能忽视，需要关注。

本文建立重载铁路有砟道床的三维离散元模型，在朔黄铁路现场开展有砟轨道的动力学试验，并利用实测结果验证道床离散元模型的合理性。采用三维道床离散元模型数值模拟重载货车动荷载作用下道床的宏细观力学行为，着重分析不同轴重、不同速度条件下道砟的振动位移、加速度及其空间分布特性，并对轨枕-道砟接触力进行统计分析。这些细观尺度的道床力学行为模拟结果可为重载铁路有砟道床的设计、改造及维护提供理论指导。

1 道床力学行为的数值模拟方法

1.1 道床三维离散元模型

铁路道砟是一种典型的粗粒径颗粒材料。为建立有砟道床的离散元模型，首先需要构建单个道砟颗粒的离散元模型。考虑到道砟颗粒具有不规则的几何形状和尖锐的棱角，文献[15-16]给出了一种采用三维激光扫描技术获取道砟颗粒真实几何形态并据此构建道砟簇颗粒离散元模型的方法。图1为采用该方法获取的部分道砟颗粒样本三维扫描图像及构建的簇颗粒离散元模型。从图1可以看出，所构建的簇颗粒模型能够较逼真地模拟道砟样本的真实不规则几何形态。采用这样的簇颗粒模型，可在道床力学行为模拟中避免道砟发生过度转动，从而在细观尺度上较好地模拟道砟颗粒体系的接触互锁效应。

图1 道砟颗粒的真实几何形状及其簇颗粒离散元模型

对一批碎石道砟样本进行簇颗粒建模后，构建道砟簇颗粒模型样本库[16]。按照一级道砟的级配标准，从样本库中随机选择不同粒径的道砟簇颗粒模型。在离散元软件PFC中对这些簇颗粒进行堆积、压实，建立如图2所示的重载铁路有砟道床离散元模型。由于三维离散元计算比较耗时，同时考虑到轨道的对称性，图2中的离散元模型用21 371个簇颗粒模拟道砟，用长度为1 300 mm的近似长方体刚性簇颗粒模拟轨枕，构建了轨道中心线一侧的半轨枕-道床结构。道床顶面宽度为1.8 m，边坡比为1∶1.75，砟肩堆高为15 cm。在轨道纵向平面内，设置道床离散元模型长度为600 mm，模拟枕盒道砟以及轨枕下方道砟接触力在道床中的扩散现象。考虑到后文现场试验中道床的实际厚度约为500～600 mm，离散元模型的道床厚度设为550 mm。

此外，为了便于向道床离散元模型施加重载列车引起的枕上动压力，在轨枕顶面的承轨台位置设置一个刚性加载颗粒。道床离散元模型的边界设置如图3所示，即顶面和边坡为自由边界，而底面和侧面为固定平面边界。

图2 重载铁路有砟道床三维离散元模型(单位：mm)

图3 道床离散元模型的边界

采用线性弹簧-阻尼模型计算道砟颗粒间的接触力，其中，切向接触力满足Coulomb摩擦定律。道砟与轨枕、边界的接触力也采用线性接触模型计算。道床离散元模型的细观接触力学参数取值见表1。参考文献[10,12-13]以及文献[15]开展的有砟轨道室内实尺模型试验研究结果，选取道砟颗粒法向接触刚度为15 MN/m，切向接触刚度为10 MN/m。考虑到现场试验中道砟颗粒表面比较粗糙、呈现突出的表面纹理，取道砟颗粒之间的摩擦系数为0.7。另外，为了模拟路基对道床的弹性支承，离散元模型中道砟与底面边界的法向接触刚度取值较小，设为0.2 MN/m。

表1 道床离散元模型参数

1.2 离散元模型输入荷载

为了获得重载列车运行时有砟道床承受的动荷载，采用文献[16]建立的车辆-轨道耦合动力学模型计算重载货车通过时的钢轨支点动压力。参照现场试验地段的轨道结构，按75 kg/m钢轨、Ⅲ型混凝土轨枕等建立轨道动力学模型，其中，扣件垂向刚度取为120 MN/m，扣件间距0.6 m。参照现场试验中KM98型漏斗车(轴重30 t)和C80型敞车(轴重25 t)分别建立货车动力学模型，以朔黄重载铁路实测轨道随机不平顺为激励，模拟计算车辆分别以60、70、80、90和100 km/h速度运行时的钢轨支点动压力。

图4为C80型和KM98型货车以80 km/h速度运行时钢轨支点动压力的时程曲线。从图4可以看出，在C80和KM98货车作用下，最大钢轨支点压力分别为51.6 kN和64.8 kN，后者比前者大25.6%。文献[3]指出，在朔黄铁路现场试验中，当行车速度达到75 km/h时，在C80和KM96货车作用下路基地段钢轨支点压力幅值分别为52.7 kN和67 kN。可见，在相似工况下，本文车辆-轨道耦合动力学数值模拟得到的钢轨支点压力幅值与现场试验结果较接近。

图4 行车速度为80 km/h时钢轨支点压力时程曲线

将计算得到的钢轨支点动压力施加到图2所示道床离散元模型的加载颗粒质心上，即可模拟重载货车动荷载作用下有砟道床的力学行为。

2 有砟道床现场试验

2.1 试验工况

分别于2014年12月和2015年8月在朔黄铁路开展两次重载铁路轨道动力学现场试验。试验列车由神八机车、KM98型和C80型货车等混编而成。试验中行车速度为60～100 km/h，各种速度工况下的行车试验均重复进行了3次。

2.2 现场测点布置

现场测试断面位于朔黄铁路上行线K118+500和K70+940附近，均为土质路基地段有砟轨道，铺设75 kg/m重型钢轨、弹条Ⅱ型扣件和Ⅲ型混凝土轨枕，扣件间距0.6 m。道床采用一级碎石道砟，厚度约为500～600 mm，道床边坡1∶1.75，砟肩堆高约为15 cm。

图5为现场试验中轨枕垂向位移和道砟振动加速度测点的位置示意。其中，在混凝土轨枕两端的上表面设置2个垂向位移计，在两股钢轨正下方、轨枕底面以下15 cm处埋设2个振动加速度传感器。

图5 轨枕位移和道砟加速度测点布置示意(单位：cm)

3 道床动态行为分析

3.1 轨枕位移响应

图6给出了行车速度在60～100 km/h范围内轨枕位移幅值的离散元模拟和现场实测结果，每个速度工况下均有3个实测值。从图6可以看出，相同速度工况下重复测试结果的一致性较好，离散元模拟得到的轨枕位移幅值与实测值接近；模拟和试验结果均显示，轨枕位移幅值随行车速度的提高略有增加，5种速度工况下其最大值约为1 mm。

图6 轨枕位移幅值离散元模拟值与实测值对比

图6还表明，在相同的行车速度下，轴重30 t货车作用下的轨枕位移幅值略大于轴重25 t货车。例如，当速度为80 km/h时，两种轴重下的轨枕位移幅值分别为0.92 mm和0.86 mm，可见，增大轴重对轨枕位移幅值影响较小。

3.2 道砟振动响应

为了对比分析道砟振动响应的离散元模拟结果与试验结果，在道床离散元模型中选择4个大小和形状接近的簇颗粒(图7)作为道砟振动响应观测点。测点道砟颗粒位于钢轨支点位置正下方、轨枕底面以下15 cm，与现场试验中道砟振动加速度传感器的埋深一致。

图7 离散元模型中测点颗粒位置

取4个测点道砟颗粒位移和速度幅值的平均值，表2列出不同轴重和行车速度下道砟振动位移和速度的模拟计算结果。由表2可知，在轴重25 t和30 t货车作用下，测点道砟颗粒的位移和速度幅值均随行车速度的提高呈增大趋势，但变化幅度较小。在相同行车速度下，轴重30 t货车比25 t货车作用下的道砟颗粒位移和速度幅值略大。

表2 不同轴重和行车速度下道砟颗粒的位移和速度幅值

图8给出了行车速度在60～100 km/h范围内道砟振动加速度幅值的离散元模拟和现场实测结果。其中，每种行车速度工况均给出了3个实测值和4个测点道砟颗粒加速度幅值的离散元模拟值。从图8可以看出，道砟振动加速度的离散元模拟值与实测值较吻合。轴重30 t货车以60 km/h速度运行时，离散元模拟的道砟振动加速度幅值平均值为0.54g；当行车速度提高至80 km/h和100 km/h时，道砟加速度幅值平均值分别为1.24g和1.65g，较速度60 km/h时增大了约1.3倍和2倍。这说明重载列车的行车速度对道砟振动加速度的影响较明显。至于轴重对道砟振动加速度的影响，图8显示当行车速度为60 km/h和70 km/h时，两种轴重下的道砟加速度幅值较接近；速度为80～100 km/h时，轴重30 t货车作用下道砟振动加速度幅值略大于轴重25 t货车。

图9给出行车速度为80 km/h时，轴重30 t货车作用下道砟振动加速度模拟与实测结果的功率谱。从图9可以看出，离散元模拟的道砟振动加速度功率谱幅值和频带分布特性均与试验结果基本吻合。模拟和试验结果均表明，道砟振动加速度的优势频率在60～150 Hz范围内，与文献[17-18]的理论和试验研究结论一致。

图8 道砟加速度幅值离散元模拟值与实测值对比

图9 道砟振动加速度功率谱

3.3 道砟振动的空间分布形态

为了分析道砟振动响应在轨道纵向和横向平面内的分布规律，图10和图11分别给出了KM98型货车以80 km/h速度通过(第二位轮对作用于轨枕上方)时道砟的位移和速度矢量分布，箭头方向表示道砟颗粒的运动方向，箭头颜色表示矢量大小。从图10、图11可以看出，在轨道纵向平面内，道砟主要呈现垂直向下的振动位移和速度，道砟颗粒的位移和速度均在临近轨枕底面区域最大，并沿道床深度方向逐渐衰减。

在轨道横向平面内，枕下道砟以垂直向下运动为主，较大的道砟位移和速度主要出现在钢轨支点以下位置，且自轨枕底面向下，随着道床深度的增加，道砟的位移和速度幅值逐渐衰减，位移和速度等值线基本呈扇形。这体现了动荷载由钢轨支点集中传递至轨枕，再由散体道床逐渐分散的荷载分布与传递特征。

图10 道砟颗粒位移分布

图11 道砟颗粒速度分布

从图10和图11的结果还可以看出，尽管离散元模拟时仅施加了垂向钢轨支点压力荷载，但位于砟肩范围的道砟不仅存在垂直方向振动，还出现了沿水平方向的振动位移和速度，这是由轨枕下方颗粒的横向挤压导致的。

上述结果表明，重载货车通过时，位于钢轨支点正下方的道砟振动响应最明显。该区域的道砟颗粒在重载列车长期循环作用下将会被逐渐压实，从而导致有砟道床发生累积变形，形成轨道沉降。因此，为改善轨道的几何平顺性，有砟轨道的捣固作业应主要针对钢轨下方及附近区域的道砟进行。作为道砟振动最明显的区域，位于钢轨支点下方的道床也是有砟道床强化改造的重点关注区域。此外，尽管在图10和图11中砟肩仅出现了较小的垂向和水平振动响应，但是可以推测，由于受到重载列车引起的动荷载或冲击荷载的长期作用，枕下道砟的振动和挤压效应可能使砟肩发生横向流动甚至边坡松塌，这可能是导致有砟道床流变的主要原因之一。

3.4 轨枕-道砟接触压力

由于目前的测试技术还很难获得轨枕与道砟之间的接触信息，当前对轨枕与道砟接触行为的认识仍较缺乏。本文基于离散元模拟结果分析轨枕底面与道砟颗粒的接触压力特征。表3列出行车速度为80 km/h、车辆第一位轮对运行至轨枕上方时，轨枕底面与道砟颗粒的接触压力统计值。其中，平均接触力f0通过接触力总和除以接触点数计算。由表3可知，轴重25 t和30 t货车作用下，轨枕底面(半根轨枕)与道砟颗粒接触点数分别为158和164，最大接触力均超过平均接触力f0的4倍。若定义小于f0的接触力为弱力链，介于f0和3f0之间的接触力为强力链，超过3f0的接触力为极强力链，那么，在所有接触力中，超过55%为弱力链，强力链数约占30%～40%，而极强力链数不足10%，说明轨枕-道砟接触力呈现出明显的不均匀性。此外，轴重为25 t时，弱力链、强力链和极强力链承担的荷载分别为10.0、28.3和12.0 kN，三者占总荷载的比例分别为20%、56%和24%；而轴重为30 t时，三者承担的荷载比例分别为20%、63%和17%。可见，尽管弱力链数量超过55%，但承担的荷载比仅为20%，而80%的荷载由强力链和极强力链(65～70个接触点)承担。

表3 轨枕底面与道砟的接触压力统计

为分析轨枕-道砟接触力的空间分布特征，图12给出了轴重25 t货车轮对作用在轨枕上方时轨枕底面与道砟的接触力分布，图中箭头的长度与矢量大小成正比。从图12可以看出，接触力空间分布呈现明显的区域集中现象，大部分强力链和极强力链(平均接触力f0=318.2 N)分布在钢轨支点下方及附近区域。受道砟颗粒不规则形状和排列位置的影响，在靠近轨道中心和轨枕边缘处也出现了个别极强力链。

图12 25 t轴重下轨枕-道砟接触力分布

根据轨枕-道砟接触力模拟结果进一步分析接触面的应力分布。文献[19]将压力盒安装在轨枕底面，测得轴重25 t货车作用下钢轨正下方枕下压应力约为250 kPa。为了与该试验结果进行对比，考虑到压力盒直径约为12 cm，将离散元模型中轨枕和道砟接触面划分为13 cm×10 cm的网格，计算各网格面积内的平均压应力，如图13所示。从图13可以看出，钢轨支点下方区域的压应力较大，最大压应力达到247 kPa，与文献[19]的测试结果接近；枕中位置附近的压应力较小，最小压应力仅为12 kPa，约为最大压应力的1/20，说明轨枕-道砟接触面的压应力分布极不均匀，特别是钢轨支点下方与轨道中心附近区域的压应力差异较大。

图13 25 t轴重下轨枕-道砟接触面的平均应力分布

当货车轴重为25 t时，模拟计算得到的枕下压力为48.3 kN，将其除以枕底有效支承面积(0.386 m2)，得到道床顶面平均应力为125 kPa，由此可知应力不均匀系数约为2.0。当轴重为30 t时，离散元模拟得到的13 cm×10 cm区域最大压应力为290 kPa，道床顶面平均应力为158 kPa，应力不均匀系数约为1.8。可见，两种计算工况下的道床顶面应力不均匀系数均略大于工程设计常用值1.6。因此，我国在发展大轴重铁路货运技术时，可考虑适当提高道床顶面应力不均匀系数的取值。

4 结论

建立重载铁路有砟道床三维离散元模型，模拟分析动轮载作用下道床的宏细观动态行为，并与重载铁路现场试验结果进行对比分析，得到以下结论：

(1)在轴重25 t和30 t重载货车作用下，离散元模拟得到的轨枕位移幅值、道砟振动加速度幅值及功率谱均与实测结果吻合较好，表明本文建立的有砟道床三维离散元模型是合理的，离散元模拟结果基本可靠。

(2)在相同行车速度下，车辆轴重由25 t增大至30 t以后，轨枕和道砟振动响应幅值均略有增大。行车速度在60～100 km/h范围内，轨枕和道砟颗粒位移幅值均随速度的增加而略有增大，而道砟振动加速度幅值快速增大，这会加剧道砟颗粒之间的接触挤压、磨损和破碎，加速道床服役状态和性能的劣化。

(3)在车辆动荷载作用下，枕下道砟主要发生垂直运动，其振动位移和速度在钢轨支点下方最大，随着道床深度的增加，道砟位移和速度幅值逐渐衰减，位移和速度等值线基本呈扇形。砟肩由于受到枕下道砟的挤压作用同时出现垂向和水平振动。道床不同区域道砟振动响应的差异从一定程度上揭示了有砟道床发生流变、边坡松塌的细观机制，同时，道砟振动响应的空间分布特征也可为弹性轨枕枕下胶垫的参数设计和道砟胶固化方案设置提供理论参考。

(4)在动荷载作用下，轨枕底面与道砟的接触力呈现明显的不均匀性，模拟结果显示80%的枕下压力由约70个接触点(半根轨枕)的强力链承担，且主要分布在钢轨支点下方区域。模拟结果还表明，枕下最大压应力约为道床顶面平均应力的1.8～2倍，应力不均匀系数略大于其工程设计系数1.6，这一点还有待进一步的理论研究和试验验证。