风力涡轮机表面感应硬化轴承塑性变形和屈服准则的实证公式 郭晓宇、 赵燕、 娄文静、 王晓波、 郭峰、 刘海超

《Lubricants》：Empirical Formulas for Plastic Deformation and Yield Criterion of Surface Induction-Hardened Bearings for Wind Turbines Xiaoyu Guo, Yan Zhao, Wenjing Lou, Xiaobo Wang, Feng Guo and Haichao Liu

　　摘要：本研究开发了一种弹塑性有限元（FE）接触模型，用于评估风力涡轮机中使用的表面感应淬火锥形滚子轴承的塑性变形。该模型考虑了深度依赖的材料特性和热处理引起的残余应力分布。通过将计算出的塑性变形与静态压缩实验的测量结果进行比较，验证了该模型的有效性。结果表明，感应淬火产生的残余应力对轴承的弹塑性行为有显著影响。基于该模型，进行了参数分析，以研究表面硬化深度（SHD）、接触压力和残余应力对表面塑性变形的影响。开发了经验公式来预测表面硬化锥形滚子轴承的表面塑性变形并评估材料的屈服性能，从而防止服役过程中的过度变形。这使得能够快速估计不同硬化深度下的最大塑性变形，并提供了一种有效的评估屈服风险的方法。

1. 引言
主轴承是确保风力涡轮机安全稳定运行的关键部件。锥形滚子轴承（TRBs）因其能够有效承受径向和轴向载荷而得到广泛应用[1,2]。这些轴承在实际应用中经常面临严峻的工作条件，包括低转速、重载荷和循环冲击载荷[3]。为了提高风力涡轮机轴承的服役可靠性，采用了多种策略，包括表面涂层、先进的润滑技术和高性能润滑剂添加剂。随着风力涡轮机额定功率的不断增加，主轴承的尺寸也在相应增大。对于这些大型轴承，由于材料和工艺的限制，进行硬化处理往往不切实际[4,5]。因此，表面感应淬火在工程实践中得到了广泛应用，因为它提高了轴承滚道的抗塑性变形能力和耐磨性[6,7]。表面感应淬火后，近表面材料的机械性能（如硬度、屈服强度和本构响应）随深度显著变化[8]。根据硬度分布，材料通常可以分为三个区域：硬化层、过渡层和核心材料。先前的研究[9,10,11,12]已经证实，表面硬化部件的屈服起始强烈依赖于表面硬化深度（SHD）。当SHD较小时，塑性变形更可能发生在过渡层与核心材料的界面附近。随着SHD的增加，屈服起始位置向硬化层移动，通常位于最大冯·米塞斯应力区域附近。因此，设计满足服役要求的适当SHD具有重要意义。

许多研究致力于评估表面感应淬火轴承的塑性变形响应。Glode?等人[13]采用有限元（FE）分析研究了具有表面感应淬火的大型回转球轴承的接触行为。他们的发现表明，在相似条件下，表面淬火球轴承的抗塑性变形能力通常低于整体淬火的轴承。他们还发现，在轴向力较小而倾覆力矩较大的载荷条件下，轴承间隙会显著降低抗变形能力。G?ncz等人[14,15]和He等人[16]进一步研究了滚子轮廓和滚子倾角对表面感应淬火的大型三排圆柱滚子回转轴承接触行为和变形的影响。他们的结果显示，随着SHD的增加，抗塑性变形能力增强。当SHD超过某个临界值时，抗变形能力接近整体淬火轴承[15]。此外，优化滚子轮廓被发现可以进一步提高接触性能，而滚子倾角被认为是影响轴承承载能力的重要参数[15,16]。基于FE模拟，Li等人[17]通过线性回归分析建立了表面感应淬火圆柱滚子回转轴承表面塑性变形的预测公式。考虑到材料的屈服标准，Lai等人[18,19]开发了分析模型，用于预测在线接触、点接触和椭圆接触条件下的表面感应淬火轴承的塑性变形。然而，大多数现有的关于大型表面淬火轴承塑性变形的分析没有考虑残余应力，这可能导致理论预测与实际运行条件之间的不一致。

表面感应淬火不仅改变了近表面材料的机械性能，还改变了部件内的残余应力状态。通常，在硬化层产生压缩残余应力，而在过渡层和核心材料产生拉伸残余应力[20,21,22,23,24]。大量研究表明，残余应力显著影响轴承的变形行为和裂纹演化。压缩残余应力可以抑制裂纹的起始并延缓裂纹扩展，从而提高抗塑性变形能力并延缓损伤累积[25,26,27,28,29,30,31]。相比之下，拉伸残余应力通常会促进变形并加速损伤演化[32,33,34]。Torkamani等人[35]发现，与传统淬火轴承相比，表面感应淬火轴承产生的残余应力提高了轴承的抗微点蚀和表面耐磨性。Areitioaurtena等人的研究[36]表明，表面感应淬火引起的残余应力可以将轴承的滚动接触疲劳寿命提高1.56倍。Londhe等人[37]报告称，表面感应淬火可以显著降低轴承的疲劳裂纹增长率，从而提高其断裂韧性。Dechant等人[38]强调，在高载荷滚动-滑动接触条件下，部件的接触响应和损伤演化由次表面微观结构和残余应力的共同作用决定。总体而言，表面感应淬火带来的性能提升源于表面机械性能的改善以及诱导残余应力的综合效应。

以往关于残余应力的研究主要关注其对小尺度、具有浅层感应硬化层轴承的变形行为的影响。为了了解残余应力对风力涡轮机大型表面感应淬火TRBs的弹塑性接触行为的影响，本研究进行了数值模拟和实验研究。首先，开发了一个用于表面感应淬火大型TRBs接触问题的弹塑性FE模型，并在其中模拟了残余应力场。其次，分析了残余应力、SHD以及硬化层和核心材料的本构行为对弹塑性接触问题的影响。最后，基于数值模拟结果提出了预测表面感应淬火轴承表面塑性变形和材料屈服的经验公式。

2. 不均匀材料的弹塑性接触问题的FE模型
本研究调查的大型TRB的几何形状如图1a所示。内外圈的直径分别为3.25米和2.65米。它可以作为18 MW风力涡轮机的前主轴承或20 MW风力涡轮机的后主轴承使用。内外圈的滚道经过表面感应淬火处理，而滚子则进行了整体淬火。滚子母线和内外圈滚道之间的接触简化为线接触问题，如图1b所示，可以等效地建模为滚子-平面接触配置。从另一个角度来看，滚子-滚道接触的椭圆率通常超过35，因此接触可以合理简化为平面应变问题。因此，在平面应变条件下建立了弹塑性FE模型，如图1b的放大视图所示。在此模型中，表示滚子与轴承滚道之间的等效接触半径。滚子-内圈接触和滚子-外圈接触的等效接触半径分别使用公式（1）和（2）确定[39]：
(1)
(2)
其中为滚子半径，为内圈滚道半径，为外圈滚道半径。

图1. FE模型的示意图：(a) 锥形滚子轴承的示意图；(b) 平面应变问题的示意图；(c) FE模型中使用的网格。

图2显示了实验室中表面感应淬火后轴承材料的深度依赖硬度分布。在本研究中，SHD定义为从轴承滚道顶面到维氏硬度达到表1中指定阈值（HV1）的位置的距离。过渡层的厚度随SHD变化，根据参考文献[18,19]，大约是SHD值的0.1倍。

表1. FE模型中使用的材料参数。表面感应淬火引起的硬度梯度导致材料具有深度依赖的本构行为，包括应力-应变响应和屈服强度的变化。在本研究中，硬化层材料（硬度为，屈服强度为）和核心材料（硬度为，屈服强度为）的应力-应变曲线是通过单轴拉伸试验获得的，如图3所示。这些拉伸数据用于校准本构模型以及FE分析中采用的冯·米塞斯等效应力屈服标准。根据图2中的硬度分布，然后通过参考文献[19]中的公式（3）确定的插值方法确定过渡层的应力-应变关系。由于滚子的硬度和屈服强度高于滚道，因此可以认为滚子在接触过程中的塑性变形可以忽略不计[13,14,15,16,17,18,19]。此外，作者还研究了滚子的弹塑性建模的影响，发现其对滚道最终塑性变形的贡献可以忽略不计。因此，在分析中将滚子建模为线性弹性材料，以确保计算效率的同时保持足够的准确性。

图3. 硬化层、过渡层和核心材料的应力-应变图。在本研究中，采用Voce各向同性硬化模型[40]来表征轴承材料的塑性变形行为。其相应的本构关系表示如下：
(4)
其中表示从拉伸试验中获得的初始屈服应力，表示饱和流应力，表示饱和指数。这两个参数是通过拟合图3中的应力-应变曲线确定的。在本框架中，拉伸屈服应力作为各向同性硬化律的输入参数，并在多轴接触应力状态下与冯·米塞斯标准一起进行评估。FE模型中采用的材料参数列在表1中。

表面感应淬火过程后轴承滚道的残余应力分布是通过数值模拟确定的。参考文献[41]中呈现的感应淬火试样的XRD测量结果观察到相似的残余应力分布。如图4所示，在硬化层产生压缩残余应力，而在过渡层和核心材料产生拉伸残余应力。在本研究中，残余应力分布作为初始应力场纳入FE模型，分别在x和z方向上。

3. 静态压缩试验
为了验证建立的弹塑性FE模型，对轴承环样品进行了静态压缩试验，本节详细介绍了实验设置和程序。仅为本研究定制的大型TRB的内外圈被切割成多个部分。在500吨液压试验机（Tianduan 500T，GENECTEC TIANDUAN，天津）上进行了静态压缩试验，以测量在不同载荷下的表面塑性变形。图5a显示了实验装置的局部放大视图，包括试验机、夹具和加载装置。图5. 静态压缩试验和塑性变形的轮廓测量：(a) Tianduan 500T液压试验机的照片；(b) 探针型表面轮廓仪和测量位置的照片；(c) 测量的表面轮廓示例。在每个内外圈部分的左侧、中间和右侧三个位置进行了静态加载试验。在每个试验位置，载荷均匀施加10秒，然后以恒定速率卸载，即一个完整的加载-卸载循环。测试结束后，使用了一种笔式表面轮廓仪来测量样品上部、中部和下部的滚道轮廓。滚子与内圈的接触位置以及轮廓测量位置在图5b中以示意图形式展示。测量得到的轮廓显示在图5c中。结果表明，样品的中部位置（4-2、5-2和6-2）承受了最大的接触压力，并且发生了最大的表面塑性变形。因此，这些位置代表了有限元（FE）模型进行塑性变形分析时最关键的区域。用于后续模型验证的外圈和内圈的几何参数分别列在表2和表3中。表2. 外圈的几何参数。表3. 内圈的几何参数。

4. 有限元模型的边界条件与验证
4.1. 左右边界条件对塑性变形的影响
在进行表面感应硬化处理后，内外圈的残余应力在圆周方向上大致是均匀的。然而，当圈被分割时，左右边界处的圆周约束被消除，导致这些自由边缘附近的残余应力部分释放并发生局部重新分布。为了研究这种重新分布对弹塑性接触的影响，考虑了两种类型的边界条件。在边界条件1（图6a）下，底部边缘被完全约束（狄利克雷边界条件，位移约束），而左右边界则受到与预定残余应力一致的应力边界条件（诺伊曼边界条件，应力约束）。相应的残余应力场显示在图7a中。这种配置可能代表了未分割轴承的实际服役条件，其中保持了圆周残余应力的连续性。在边界条件2（图6b）下，底部边缘也被完全固定，而左右边界不受约束，允许残余应力在这些边缘自由释放。相应的残余应力分布显示在图7b中。这种情况反映了压缩测试中使用的分割样品的实验条件。

4.2. 有限元模型的验证
用于FE模型的内外圈滚道的材料属性已在表1中列出。几何参数分别列在表2和表3中，模拟中使用的计算参数列在表6和表7中。表6. 外圈情况的计算参数。表7. 内圈情况的计算参数。图8和图9分别将FE模拟结果与内外圈的实验结果进行了比较。当将残余应力场纳入模型时，数值模拟结果与实验测量结果非常吻合，最大误差小于15%。相反，当忽略残余应力时，预测值与实验结果的偏差显著超过50%。这些结果表明残余应力边界条件会影响预测的表面塑性变形。由于边界条件2更接近实验中使用的实际段配置，因此与测量数据的吻合度更好。

5. 考虑残余应力的塑性变形与经验公式
5.1. 表面硬化深度和接触压力对弹塑性变形的影响
图10显示了表面塑性变形随表面硬化深度（SHD）的变化情况。图中，表示等效滚子直径，表示用于标准化的参考接触压力，设定为4 GPa [17]。结果表明，无论是否考虑残余应力，表面塑性变形都会随着比的增加而减小。然而，两种情况下的变化趋势有显著差异。当忽略残余应力时，表面塑性变形随的变化表现出明显的两个阶段特征。在第一阶段，随着的增加，变形迅速减小。超过某个特定临界值后，进入第二阶段，进一步增加几乎不会引起表面塑性变形的变化。此外，这个临界值随着的增加而增大。这种行为与先前的研究结果一致[13,14,15,16,17,18,19]，表明当SHD足够大时，轴承的塑性变形接近全硬化轴承的塑性变形。相比之下，当考虑残余应力时，表面塑性变形随着的增加而单调减小，没有明显的平台期。值得注意的是，考虑和不考虑残余应力得到的表面塑性变形曲线在较小的值处相交。然而，随着的增加，考虑残余应力时的预测变形小于或不考虑残余应力时的变形。

5.2. 大型表面感应硬化滚子轴承表面塑性变形的经验公式
目前用于计算表面感应硬化滚子轴承表面塑性变形的经验公式存在一定的局限性。Li等人[17]提出的公式没有包括屈服准则，因此即使在不应发生屈服的低载荷下也预测出非零的塑性变形。Lai等人[18]通过线性叠加硬化层和核心材料的塑性贡献来预测表面塑性变形，并根据屈服条件分配相应的贡献。然而，他们的屈服准则假设了一个恒定的最大冯·米塞斯应力因子（对于线接触为0.56），这在塑性开始和应力场重新分布后变得可疑。此外，这些研究中使用的FE模型没有考虑残余应力。正如本研究所示，残余应力对表面塑性变形有不可忽视的影响。因此，基于一系列包含残余应力的弹塑性接触模拟，本研究提出了一个特别适用于大型表面感应硬化滚子轴承的表面塑性变形公式。

滚子母线中点的接触压力（对应于表2中的OR-4-2–5-2和表3中的IR-5-3–5-2）是所有测量位置中最高的，因此代表了潜在失效的最关键位置。因此，后续数值模拟中使用的几何参数选为与OR-4-2–5-2和IR-5-3–5-2相匹配的参数。此外，为了更好地模拟完整轴承的实际服役条件，使用了边界条件1（图6a），该条件保持圆周连续性。图12和图13分别显示了表面塑性变形与和的依赖性。结果表明，表面塑性变形随的增加而非线性增加，而随着的增加则非线性减小。根据这些趋势，表面塑性变形表示为以自然常数为底的指数函数，和为主要控制变量：(5) 其中，，，和是通过FE模拟得到的表面塑性变形数据确定的拟合参数。所提出的公式进一步包含了屈服准则，其中??表示麦克劳林括号。该公式中的项表示在核心材料中最大残余拉应力位置计算的无量纲冯·米塞斯应力。其数学表达式如下：(6) 其中，，和是通过FE模拟得到的冯·米塞斯应力值确定的常数。

对于方程(5)，为了确保足够的拟合精度，在和间隔进行FE采样进行了回归分析。内外圈的决定系数都达到了0.9994，表明拟合表达式与数值数据之间有极好的一致性。图14比较了FE结果与公式拟合结果：(a) 外圈；(b) 内圈。对于方程(5)中表面塑性变形表达式的回归分析，FE采样在和间隔进行。内外圈的决定系数都达到了0.9995，表明拟合精度很高。如图15所示，FE结果与拟合表面非常吻合，进一步验证了所提公式的准确性。表面塑性变形的FE结果与公式拟合结果比较：(a) 外环；(b) 内环。图16展示了所提公式预测的表面塑性变形与静态压缩试验所得结果的对比。图中，渐变色表面代表计算值，红色点对应实验测量值。橙色区域表示纯弹性行为，其余渐变色区域表示弹塑性区域。在该区域内，材料开始屈服，表面塑性变形随接触压力和硬度差异（SHD）非线性变化。比较结果显示，公式预测值与实验测量值的偏差在大多数情况下小于10%，只有少数数据点的偏差接近15%。考虑到第4.1节讨论的边界条件对表面塑性变形的影响，实际工程应用中的差异预计将进一步减小。总体而言，这些结果证明了所提公式的有效性和预测能力。

结论与展望
本研究探讨了大型表面感应淬火锥形滚子主轴承在风力涡轮机中的弹塑性接触行为，特别关注残余应力的作用。通过将FE模拟与实验测量结果进行比较，评估了残余应力对弹塑性接触响应的关键影响。基于参数化FE分析，提出了明确考虑残余应力的表面塑性变形和材料屈服的预测公式。主要结论如下：
(1) 考虑残余应力的弹塑性接触模拟与静态压缩实验结果的一致性显著提高，这证实了在大型表面感应淬火滚子轴承的弹塑性接触分析中应考虑残余应力，以提高模拟精度。
(2) 残余应力改变了表面感应淬火轴承内部的初始应力状态，导致表面塑性变形随接触压力和硬度差异的变化模式不同，这取决于是否考虑了残余应力。因此，轴承的预测应力场和塑性变形行为受到残余应力是否包含的影响。
(3) 提出了表面感应淬火轴承表面塑性变形和材料屈服的实用公式。前者能够快速估计最大表面塑性变形，后者为评估材料屈服的开始提供实用标准。
本研究中的实验、数值模拟和推导出的公式均基于单个滚子-滚道接触对。在未来的研究中，所提出的实用公式可以结合轴承系统的全局载荷分布分析，考虑多个滚子与滚道之间的相互作用，从而扩展其在评估整个轴承每个接触位置的塑性变形中的应用。此外，本研究中的实验和数值研究都是在单次加载-卸载循环下进行的。后续研究将进行循环加载-卸载试验和数值模拟，以进一步明确材料强化对塑性变形演变的影响。还将研究不同的表面感应淬火工艺参数（例如加热温度、冷却速率）对残余应力分布和塑性变形行为的影响，为表面感应淬火工艺优化提供更全面的理论基础。