金锋，肖宏，王连俊，姜子清，杨冀超

(1. 北京交通大学 土木建筑工程学院，北京，100044；2. 北京交通大学 轨道工程北京市重点实验室，北京，100044；3. 中国铁道科学研究院集团有限公司 铁道建筑研究所，北京，100081；4. 中国铁路北京局集团有限公司高铁工务段，北京，100038)

钢轨波浪形磨损是长期以来困扰铁路行业的难题，其发生、发展过程涉及车辆、轨道结构、钢轨材质等多个因素，形成及发展机理复杂，诸多学者从自激振动、材料磨损、耦合动力作用等多个方面开展了系列研究[1-4]，获得了一系列阶段性重要成果，然而，目前尚无法完全避免波磨病害，其依然在部分线路上作为主要的病害类型，严重影响着钢轨使用寿命。尤其在某重载铁路线上，随近年来轴质量与通过总质量的不断增加，较为严重的钢轨波磨病害越来越多，严重影响钢轨使用寿命。现场调研发现，波磨地段主要集中于小半径曲线地段下股钢轨，波长λ集中于160～200 mm。

图1所示为典型的重载铁路小半径曲线地段塑流型波磨病害照片，所在线路为重载铁路重车线，现场铺设75 kg/m 钢轨、Ⅲ型混凝土枕及弹条Ⅱ型扣件系统。由图1可知：钢轨顶部光带边缘两侧存在由于塑性变形所导致的向上的“凸起”，波磨波谷处光带宽度比波峰处光带宽度大。

图1 重载小半径曲线地段典型波磨病害Fig. 1 Typical heavy-haul corrugation in small radius curve

在重载钢轨波磨的发展过程中，钢轨塑性累积作用占主导地位[5-6]，其导致的轨面周期性不平顺将增大车辆通过时的轮轨动态响应，进一步增大波磨区域轮轨动态荷载。更加剧烈的轮轨相互作用将进一步促进轨面的塑性流动，形成波磨发展的正反馈循环，使波磨发展越发迅速。然而，目前对于重载波磨波峰、波谷区域塑性累积特性研究存在不足，尚无法完全掌握这一正反馈循环的作用机理，阻碍了对重载波磨的预防与整治工作的有效开展。

探究轮轨荷载作用下钢轨表面的塑性演变规律是长久以来轨面损伤相关研究中的重要课题之一，主要通过使用弹塑性力学理论结合数值方法开展研究。在轨面伤损模拟过程中，需要建立较为精细的轨面网格，导致求解方程组规模较大，并且由于需要考虑反复加载的情况，使整体计算流程耗费大量的计算资源。受计算效率的限制，早期相关研究主要通过建立二维模型开展研究，如ANDERSSON等[7]通过建立二维移动荷载模型分析了钢轨隐伤(squat)的伤损累积特性。温泽峰等[8]基于简谐特性的赫兹接触压力作为轮轨荷载计算了二维状态下的轮轨塑性累积变形，模型中使用了JIANG等[9-10]提出的弹塑性本构模型，较好地还原了钢轨的塑性累积历程。潘周平[11]使用ANSYS 商用有限元软件建立了二维条件下塑流型波磨仿真计算模型，分析了谐波型动载作用下钢轨塑流型波磨的发展过程。李伟等[12]采用二维轮轨接触斑移动模型，分析了钢轨温度场的特性。

随着计算能力的快速发展，近年来，国内外学者逐步开展了三维钢轨弹塑性分析模型研究，如郭俊等[13]基于法向赫兹接触作为简化的车轮荷载，通过有限单元法计算了30 次循环加载时轨面应力-应变分布特性，并对比了不同钢轨屈服强度对轨面塑性变形的影响。LI 等[14]建立了车轨耦合模型与三维钢轨弹塑性有限元钢轨模型，分析了曲线地段钢轨接头处塑性变形特征及单次荷载作用下接头不平顺大小对塑性变形的影响。REIS等[15]以赫兹法向接触应力作为循环荷载条件，建立了三维钢轨弹塑性分析模型，施加了法向赫兹接触分布应力，研究了多轮次荷载条件下钢轨表面塑性累积变形特性，并结合疲劳分析理论，分析了塑性累积状态下钢轨轨头疲劳寿命。

由以上相关研究可知，既有研究大多采用二维数值模型对轮轨作用行为进行等效分析，建立的部分三维弹塑性模型中，也仅考虑了法向荷载作用对于钢轨表面塑性累积的影响，忽略了横向、纵向接触应力分布的影响。然而，在重载波磨频发的小半径曲线地段，轮轨接触作用伴随较强的横向、纵向蠕滑现象，且实际轨面塑性屈服面的形成与轮轨的切向应力分布密切相关，使得既有的相关研究成果难以指导实际重载地段波磨病害整治。

针对以上不足，本文作者通过建立重载车辆-轨道耦合模型及钢轨弹塑性分析模型，结合半赫兹型轮轨接触算法，对小半径曲线地段钢轨波磨波峰、波谷区域在波磨演化过程中的塑性累积变形特性进行分析，对比波磨发展过程中波峰、波谷区域的塑性累积的差异，并从钢轨整治的角度，探讨同线路条件下重载塑流型波磨的发展趋势特征，以便为重载铁路钢轨塑流型波磨治理提供理论参考。

1 模型建立

本文主要通过建立小半径曲线地段重载列车动力学模型获取轮轨接触几何及法向应力、蠕滑力分布等特征，将其施加于钢轨轨头弹塑性有限元分析模型中，研究车轮循环荷载作用下的钢轨应力、应变特征。有限元仿真结果通过有限单元法进行计算，整体流程控制及数据后处理过程均通过自编模块实现，模型计算过程如图2所示。

图2 计算流程图Fig. 2 Flow chart of simulation

1.1 敞车动力学模型

1) 车辆模型。根据实际调研结果，线路开行敞车种类较多，包括C80、C70A、C64K等，其中C80为线路主要货车种类，本文选取轴质量为25 t的C80型货车建立动力学模型并进行计算。由于3节敞车车辆所组成的列车模型可以反映长编组重载列车中的单个车辆的动力学行为[16]，故本文建立3节敞车动力学模型，主要分析中间车辆的动力学行为及轮轨接触状态，以避免长编组过程中不同车辆所导致的差异，如图3(a)所示。

图3 敞车动力学模型Fig. 3 Dynamics models of wagons

根据实际货车装配的转向架类型及力学属性，构建转K6 型转向架模型，如图3(b)所示。在模型建立过程中，通过建立非线性力元的方式充分考虑了侧架与摇枕、轮对轴箱间的复杂力学关系。车辆车轮踏面LM 型磨耗型踏面，滚动圆半径为420 mm，模型中车辆速度设置为50 km/h。

2) 线路条件。根据实际线路中典型的小半径曲线地段线路参数，动力学计算模型中建立了长为600 m 的线路区段，其中包括长度为200 m、半径为400 m 的圆曲线地段，2 段长度均为100 m 的缓和曲线地段，2 段长度均为100 m 的直线地段。线路实设超高为105 mm，轨底坡度为1/40。

1.2 轮轨接触模型

从计算效率考虑，使用赫兹法向接触结合FASTSIM 切向接触方法是比较高效的手段，但在小半径曲线地段，其接触状态较为复杂，椭圆形接触斑将导致较大的计算误差，在处理小半径曲线地段复杂轮轨接触状态时存在一定不足。而从轮轨接触计算精度角度，使用CONTACT程序获取Kalker轮轨非赫兹型轮轨接触关系精度较高，但计算量大。为了兼顾计算精度及计算效率，本文使用semi-Hertzian(半赫兹)轮轨接触算法进行分析。

在半赫兹轮轨接触算法计算过程中，沿纵向滚动方向上依然采用赫兹接触理论，而横向上将可能的接触区域划分为N个条带，并分别计算每个条带的相对曲率计算法向应力分布，可真实反映轮轨接触斑的非椭圆、非对称特性。后续切向应力分布则通过KALKER[17]所提出的FASTSIM 算法计算，具体算法推导过程参照文献[18-19]。在本文计算过程中，将条带间隔Δy设为0.5 mm。根据半赫兹轮轨接触模型计算中所获得的横向蠕滑力分布与纵向蠕滑力的结果，绘制接触班内切向应力向量分布情况，如图4所示，图中箭头长度代表轮轨切向应力的大小，箭头方向表示切向应力方向。

图4 基于半赫兹算法的轮轨接触切向应力计算结果Fig. 4 Calculation results of tangential stress of wheel-rail contact by semi-Hertzian method

从图4可以看出：在小半径曲线地段轮轨切向力应力分布中横向作用较为明显，最大切应力方向指向曲线中心，与TORSTENSSON 等[20]中对小半径地段钢轨进行金相分析的结果相符。此外，小半径曲线地段切向作用明显，且轮轨应力分布呈较为明显的非对称特征，这说明已有研究中使用赫兹法向力作为轨面荷载的方法，对小半径曲线地段钢轨真实塑性累积模拟存在一些不足。

1.3 弹塑性轨头子模型建立

1.3.1 重载波磨几何学模拟

在钢轨有限元模型中，本文通过调节轨头模型单元坐标的方式，对钢轨波磨几何特征进行模拟。考虑到钢轨波磨几何的三维特性，使用三维波磨几何描述方法[21]，编制了自动化调整有限单元模型单元节点程序，对钢轨波磨几何特征进行描述，本文使用的波磨三维几何描述方程如下：

式中：x0为波磨施加的纵向起始位置；x为当前节点纵向坐标；l为波磨波长；d为波磨波深；w为波磨宽度的一半。根据现场测试结果，波磨波长选为160 mm，波磨宽度选为60 mm，波深选为0.5 mm。

弹塑性累积过程的模拟需要往复多次进行，将耗费大量计算成本，然而，在移动荷载过程中，若模型长度较短，则存在较为显著的边界效应，将对求解精度、收敛性均造成一定的影响。考虑到关注区段钢轨波磨波长为160 mm，结合模型计算效率与边界效应的影响[22]，建立长度为480 mm的钢轨子模型。本文通过C3D8R 型三维8 节点减缩积分单元对钢轨求解区域进行离散，并对钢轨前后及底部单元节点自由度进行约束。为了确保弹塑性变形区域计算的准确性，对轨面轮轨接触区域进行加密，加密区网格长度约为1 mm。

图5所示为三维钢轨波磨不平顺施加后几何与网格划分示意图，为了更好地展示波磨波峰、波谷区域几何特点，将波磨幅值设置为实际值的5倍(图中幅值为2.5 mm)。

图5 钢轨波磨几何与网格划分示意图Fig. 5 Schematic diagram of rail corrugation geometry and mesh division

1.3.2 钢轨弹塑性本构模型

重载小半径曲线地段一般采用高强钢轨，本文弹塑性模型的钢轨材料基本参数[15]如表1所示。

表1 弹塑性模型参数Table 1 Parameters of elastoplastic model

从图4可以看出，法向荷载外切向轮轨作用也较为显著，轮轨接触区域处于较为复杂的三维荷载作用中。因此，分析钢轨表面塑性累计过程中除考虑三向荷载作用外，应选择合适的三维弹塑性累积模型。本文弹塑性分析过程中，弹塑性准则使用屈服函数表达为

式中：F为屈服函数；σeq为等效应力；σy为屈服应力。等效应力通过Mises等效应力描述为

式中：S为偏应力张量；αdev为背应力的偏应力分量。在对应的塑性变形过程中，塑性流动准则可以表示为

式中：为塑性流动速率；σ为二阶应力张量；为等效塑性流动速率，

使用随动强化模型考虑多次碾压过程中钢轨表面材料特性的强化过程，选取用Ziegler 强化法则对其进行描述，其强化过程可描述为

式中：C为动态强化模量；α为背应力；σ0为弹性极限。

2 波磨区段弹塑性演化特性分析

2.1 钢轨波峰、波谷弹塑性行为分析

图6所示为小半径曲线货车内轨非椭圆轮轨接触移动荷载循环作用下钢轨表面应力-应变关系曲线，根据滞回曲线收敛情况，计算了50 次循环加载的情况。

图6 轨面应力-应变关系曲线Fig. 6 Stress-strain relationship curves of rail surface

从图6可以看出：随车轮荷载增大，钢轨表面形成塑性变形，钢轨表面应力-应变曲线开始移动，并且每次加载过程中的应力-应变曲线均不闭合，说明钢轨表面材料尚未进入安定状态；随车轮荷载作用次数增加，钢轨表面应力、应变将继续向新的平衡状态发展，在经过一定次数的车轮荷载作用后，由于材料的强化作用，轮载作用下轨头变形呈现以弹性变形为主的状态，此时，应力-应变曲线呈稳定的闭合曲线。因此，在长编组列车多次施加荷载的条件下，塑性变形导致的钢轨波磨演化的现象将趋于稳定，与实际塑流型钢轨波磨病害演化过程中病害发展后期波磨波深幅值变化趋缓的特性相符。

在多次荷载作用下，钢轨表面垂向塑性累积较为明显，在轮轨作用初期塑性变形变化较大，而后呈现减小的情况。对比波磨波峰及波谷位置垂向累积应变可知，波谷处垂向累积塑性应变明显比波峰处的垂向累积塑性应变大。这说明在波磨产生后，即使在稳态轮轨接触荷载作用下，波磨波峰、波谷区段轨面塑性应变的发展差异依然存在。因此，钢轨波磨病害形成后，即便使用调整轮轨系统共振频率的方式对钢轨波磨进行整治，也难以对已经形成的波磨病害有较好的整治效果。由图6(b)可知：塑流型钢轨波磨演化过程中波磨波谷塑性变形比波磨波峰的塑性变形大，与图1所示的波磨波谷塑流纹路偏向钢轨外侧的实际情况相符。由图6(c)可知：在重载波磨发展过程中，波磨波谷与波峰区域滞回环中心区域基本重合，产生的纵向相位差异较小，说明实际线路中波磨病害的纵向位置并不会随波磨发展的进程产生较大变化。

图7所示为循环移动荷载作用下轨面等效塑性累积应变。

图7 循环移动荷载作用下轨面等效塑性累积应变Fig. 7 Equivalent plastic accumulation strain of rail surface with cyclic moving load

从图7可以看出：等效塑性应变随加载次数的增加呈初期塑性累积速度较快而后减缓的特征。对比波磨波谷、波峰等效塑性累积应变特征可知，波谷区域处于安定状态时，等效塑性累积应变比波峰区域的大，且进入安定状态所需的荷载施加次数也较多。

图8所示为循环加载完成后，波磨轨波峰、波谷钢轨横断面残余等效塑性应变云图。从图8可以看出：无论是波峰还是波谷，最大的等效塑性应变均出现在钢轨表面以下2 mm之内。波磨波谷处塑性应变最大值较大，而且在邻近钢轨表面区域也存在较为明显的塑性累积现象。

图8 波磨轨波峰、波谷累积塑性应变云图Fig. 8 Contours of accumulated plastic accumulation strain at peaks and troughs of corrugation

2.2 钢轨波峰、波谷塑性分布特性分析

图9所示为在现场实际线路条件(半径为400 m，车速为50 km/h，超高为105 mm)下，重载波磨波峰、波谷区域等效塑性累积应变沿钢轨深度的变化。从图9可以看出：波磨波峰区域与波谷区域等效塑性累积应变均随着深度呈现先增加后减小的特点，最大等效塑性应变均出现在轨面以下。此外，在波磨波谷处除表面塑性累积较波峰处显著外，钢轨内部塑性累积也较大。

图9 等效塑性累积应变沿钢轨深度的变化Fig. 9 Variations of equivalent plastic accumulation strain along rail depth

3 线路条件对波磨塑性累积的影响

3.1 曲线半径对波磨塑性累积的影响

根据调研结果，波磨主要产生于小半径曲线地段，为了分析曲线半径R对钢轨塑流型波磨病害发展的影响，本节选取平衡超高状态下，300～1 000 m曲线半径作为典型计算工况。

不同半径条件下钢轨表面等效塑性累积应变如图10所示。

图10 不同半径下钢轨表面等效塑性累积应变Fig. 10 Equivalent plastic accumulation strain of rail surface with different radii

从图10 可以看出：在不同曲线半径下，钢轨表面塑性累积特征存在显著差异；曲线半径越小，安定状态下等效塑性累积应变越大，到达安定状态时所需要的荷载次数也越多；此外，在曲线半径为500 m及以上区段，钢轨表面均未产生塑性应变。

表2所示为安定条件下钢轨表面最终等效塑性累积应变随曲线半径变化的对比结果。从表2可以看出：钢轨表面的塑性累积特性对曲线半径极其敏感。

表2 钢轨表面等效塑性累积应变随曲线半径的变化Table 2 Variation of equivalent plastic strain of rail surface with curve radii

为了进一步分析塑性累积应变在深度方向上的分布情况，对比了不同曲线半径条件下钢轨表面波谷区域等效塑性累积应变沿深度方向的分布情况，如图11所示。

图11 不同曲线半径下轨面等效塑性累积应变分布特性Fig. 11 Equivalent plastic accumulation strain distributions of rail surface with different curve radii

从图11 可以看出：在深度方向上，小半径曲线地段塑性累积效应明显比大半径地段塑性累积效应大；当曲线半径较小时，等效塑性应变最大值出现在钢轨表面，而后随深度增加，等效塑性应变逐渐减小；当曲线半径较大时，钢轨等效塑性应变最大值出现在距离轨面2～3 mm处；随深度增加，等效塑性应变呈先增大后减小的情况，并在曲线半径大于或等于500 m时，钢轨表面未产生塑性累积，主要塑性应变产生在钢轨内部。这说明随曲线半径减小，钢轨塑性累积及相关伤损将更易在钢轨表面产生，与实际小半径曲线地段塑流型波磨、钢轨剥离掉块现象较为严重的工务实际情况较为契合。

为进一步分析不同曲线半径条件导致钢轨塑性区分布差异的原因，绘制了钢轨表面最大轮轨横向蠕滑力所在纵断面图，如图12所示。

图12 不同曲线半径下钢轨表面横向蠕滑力分布Fig. 12 Lateral creep force distributions of rail surface with different curve radii

从图12 可以看出：随曲线半径减小，最大横向蠕滑力呈明显增大的趋势；当半径为300 m 时，最大横向蠕滑力为640.12 MPa，较半径为500 m时最大横向蠕滑力372.73 MPa增加71.7%。

图13 所示为在不同曲线半径条件下，当钢轨表面处于安定状态时，轨面垂向、横向、纵向的应力-应变曲线。

图13 不同半径下轨面应力-应变曲线Fig. 13 Stress-strain curves of rail surface with different radii

从图13(a)可以看出：由于处于安定状态，各工况应力-应变曲线均闭合。曲线半径越小，应力-应变曲线包围面积越大，说明单次加载所耗散的能量越大；当曲线半径在500 m 及其以上情况时，应力-应变曲线均呈弹性特征。从图13(b)可以看出：由于小半径曲线地段横向冲切作用显著，随曲线半径减小，安定状态下钢轨表面横向塑性应变显著增加，且单次耗散的能量也显著增加。从图13(c)可以看出：纵向应力-应变曲线所包围面积随曲线半径的减小而增大。

综上可知，小半径地段由于轮轨横向蠕滑作用加剧，在钢轨表面复杂的三维荷载作用下，塑性累积变形在钢轨波磨病害发展中起重要作用，从理论角度解释了所调研山区小半径曲线地段重载波磨严重的原因。因此，在重载铁路选线过程中，应尽量避免使用半径小于500 m的曲线线型。

3.2 曲线超高对波磨塑性累积的影响

线路超高是重要的线型参数，对于曲线内外轨的受力均衡、横向冲切程度均有较大影响。本节将讨论在不同超高状态下，波磨轨表面塑性变形特征及欠超高、过超高对于塑流型波磨病害发展演化的作用。

大轴重万吨长编组列车在实际运营过程中，由于受小半径曲线较多、线路纵坡较大等客观条件的限制，在部分困难地段无法达到预期运营速度。尤其在小半径曲线地段，普遍存在线路实设超高比平衡超高大的现象，以线路常见的400 m小半径曲线为例，在平均通过速度为50 km/h时，其平衡超高为75 mm，然而，现场实设超高为105 mm，对实际运营重载货运列车的过超高量约为30 mm。

图14 所示为在不同超高H条件下，曲线内轨波磨波谷处等效累积塑性应变随车轮荷载的变化。

从图14 可以看出：随曲线超高增加，钢轨波磨波谷区域塑性累积应变明显增大。表3所示为等效塑性累积应变随超高变化百分比，由表3 可见：对比不同超高状态下钢轨塑性等效变形可知，过超高状态下轨面塑性应变累积趋势增加明显，说明超高是重载铁路小半径曲线地段内轨波磨发展的重要原因之一。

表3 等效塑性累积应变随超高的变化Table 3 Variation of equivalent plastic accumulation strain with superelevation

图15所示为不同超高状态下轨面波谷应力-应变关系曲线。

图15 不同超高条件下轨面波谷应力-应变关系曲线Fig. 15 Stress-strain relationship curves of rail surface at toughs of corrugation with different superelevations

从图15(a)可以看出：随超高增加，塑流型波磨波谷垂向应变将明显增加，整体垂向应力-应变曲线呈现向平移的特征，且滞回环面积明显增加。从图15(b)可以看出：随超高增加，横向塑性应变的累积逐渐增加，这说明在过超高情况下，曲线内轨额外的轮轨接触力将明显增加钢轨表面的横向塑性累积变形。从图15(c)可以看出：对于钢轨波磨的纵向累积应变整体较小，线路超高状态对于钢轨纵向塑性累积影响较小。

以上计算结果表明，对于小半径曲线波磨地段，在过超高状态下，由于曲线内外轨受力不均衡，导致钢轨在车轮荷载作用下垂向、横向塑性累积变形明显增加，曲线内轨更易产生钢轨波磨病害，为了减缓钢轨波磨的发展趋势，应尽量避免线路过超高现象。

4 结论

1) 在小半径曲线地段重载列车多轴移动轮轨循环荷载作用下，曲线内轨波峰、波谷区域等效塑性应变呈初期累积速度快而后逐步进入安定状态的趋势。在进入安定状态后，波磨轨波谷处的塑性累积应变速率比波峰处的塑性累积应变速率大；波磨波峰、波谷的塑性累积特性差异将导致钢轨波磨的塑性演化，使波磨波深增加；在波磨演化过程中，波峰与波谷区域产生的纵向相位差异较小，在波磨发展过程中，线路产生的波磨病害位置不会产生较大变化。

2) 在不同曲线半径下，当曲线半径较小时，等效塑性应变最大值出现在钢轨表面，并随深度增加，等效塑性应变逐渐减小；当曲线半径较大时，随深度增加，等效塑性应变呈现先增加后减小的情况，钢轨等效塑性应变最大值出现在轨面以下2～3 mm 区间内，并在曲线半径大于等于500 m时，钢轨表面未产生塑性累积，因此，建议在重载选线过程中尽量选择500 m 及以上曲线半径。

3) 波磨轨表面塑性累积状态对曲线超高较为敏感，当存在30 mm 过超现象时，轨面等效塑性累积应变相较于平衡超高情况增加明显，其中垂向塑性累积占主导因素。这说明在过超高情况下，额外的塑性累积将使波磨现象更易发展，因此，应保证列车运营速度，避免在曲线地段出现过超高现象。