8Cr4Mo4V轴承钢微观裂纹萌生与扩展机制的分子动力学模拟

Molecular Dynamics Simulation of Microscopic Crack Initiation and Extension Mechanism in 8Cr4Mo4V Bearing Steels

摘要

为研究8Cr4Mo4V轴承钢中渗碳体对基体力学性能及微裂纹萌生与扩展的影响，采用分子动力学方法系统分析了渗碳体的几何参数（如形状、尺寸、位置）对裂纹萌生和扩展的影响机制，并结合内聚力理论研究了界面裂纹扩展特性。研究结果表明：渗碳体显著提高了bcc-Fe基体的力学性能，渗碳体尺寸越小，它对基体力学性能的增强效果越显著；渗碳体的形状和位置对力学性能影响较小，但尖锐的夹杂加速了裂纹扩展，且夹杂的位置决定了裂纹的扩展路径；在bcc-Fe基体与渗碳体间的界面和错向角较大的孪晶界面，裂纹更难萌生。

Abstract

To investigate the influences of cementite on the mechanics properties of the matrix and the initiation and propagation of microcracks in 8Cr4Mo4V bearing steels， molecular dynamics models were used to systematically analyze the effects of cementite's geometric parameters （such as shape， size， and position） on crack initiation and extension mechanism. And combined with cohesive force theory， the characteristics of interface crack propagation were studied. The results indicate that cementite significantly enhances the mechanics properties of the bcc-Fe matrix， with smaller cementite particles providing a more pronounced strengthening effectiveness. While the shape and position of cementite exert a relatively minor impact on overall mechanics performance， sharper inclusions accelerate crack propagation， and the position of inclusions determines the crack propagation path. Furthermore， interfaces between the bcc-Fe matrix and cementite， as well as twin boundaries with larger misorientation angles， exhibit increased resistance to crack initiation and propagation.

Graphical abstract

关键词

8Cr4Mo4V轴承钢 / 渗碳体 / 分子动力学 / 裂纹扩展 / 内聚力参数

Key words

8Cr4Mo4V bearing steel / cementite / molecular dynamics / crack extension / cohesive parameter

0 引言

航空发动机是现代航空工业的核心，其性能和可靠性对飞行的安全和效率起着至关重要的作用。航空发动机的可靠性与其零部件所使用的材料性能息息相关。8Cr4Mo4V高温轴承钢（ASTM牌号为M50）因优异的力学性能和热稳定性而被用作航空发动机的主要轴承材料^［1］。滚动接触疲劳^［2-3］是航空发动机轴承失效的主要模式之一，这种失效通常源于材料内部的微观缺陷，如微裂纹或碳化物夹杂。为了避免或延缓微观缺陷的产生，国内外学者对提高轴承钢材料的性能进行了大量研究。YUE等^［4］利用超声喷丸技术实现8Cr4Mo4V轴承钢表面纳米晶化从而提高材料的机械强度，同时加深了对轴承钢纳米级机械强化的理解。孙玉凤等^［5］发现激光冲击会使8Cr4Mo4V轴承钢表面的碳化物碎化和次表面二次碳化物析出，抑制疲劳裂纹的萌生，显著提高8Cr4Mo4V轴承钢的旋转弯曲疲劳性能。YANG等^［6］研究了稀土对8Cr4Mo4V轴承钢滚动接触疲劳性能的影响，结果表明稀土元素显著提高了其滚动接触疲劳寿命。还有学者通过建立轴承宏观分析模型，与赫兹接触理论、断裂力学和连续介质力学等理论相结合，探索轴承钢的失效机制。GUO等^［7］采用Voronoi有限元方法对轴承钢损伤-寿命本构关系进行了数值分析。然而，上述对航空发动机轴承材料变形和破坏的研究尚未涉及微观机理，尤其是微裂纹的起源与扩展机制并未有深入研究。

在微观裂纹扩展研究方面，分子动力学模拟已成为近些年的常用手段^［8-11］。目前针对8Cr4Mo4V轴承钢基体bcc-Fe的裂纹模拟已经有了大量研究。MA等^［11］研究了bcc-Fe单晶在循环载荷下疲劳裂纹的扩展以及位错的发射机理。曹莉霞等^［12］进一步指出，裂纹面的晶体取向对微裂纹的扩展有强相关作用。对于8Cr4Mo4V中的碳化物，GHAFFARIAN等^［13］对含渗碳体多晶模型在温度和晶粒尺寸影响下的变形行为进行了分子动力学模拟研究，发现在拉伸载荷作用下通过晶界会形成孔洞进而萌生损伤，并最终发生破裂。NAKAMURA等^［14］研究了Cr₂₃C₆在bcc-Fe晶体内以及晶界处析出对临界剪切应力的影响。然而，在渗碳体对bcc-Fe基体裂纹的萌生和扩展影响机制方面，包括从两种晶体界面裂纹的扩展这一角度的研究仍然较少。

本文采用分子动力学方法构建了含渗碳体的8Cr4Mo4V轴承钢基体bcc-Fe模型，研究渗碳体的几何参数如形状、大小、位置对材料力学性能和裂纹扩展特性的影响。同时研究不同界面下的双晶模型对晶界裂纹扩展的影响，通过计算界面模型的内聚力-位移关系获得内聚力参数，为材料的微观结构提供了更深入的理解，并可以将其作为之后多尺度计算模型中的关键参数。

1 模拟方法

本文使用由美国Sandia国家实验室开发的大规模原子/分子并行模拟器（LAMMPS）^［15］进行模拟，研究8Cr4Mo4V轴承钢材料中渗碳体在拉伸过程中对bcc-Fe基体的影响，并采用OVITO软件实现模拟过程中原子构型的可视化，其中晶体结构通过共近邻（CNA）法进行分析^［16］。使用LAMMPS 进行分子动力学模拟的步骤大致可以分为建立模型、设置边界条件、施加力场、弛豫、施加载荷、输出结果。建立图1所示的拉伸模型。

图1a为基体及含渗碳体拉伸模型示意图。该模型由基体bcc-Fe和渗碳体组成，其中bcc-Fe为体心立方结构，空间群为Im3̅m，晶格参数a=b=c=0.286 64 nm，渗碳体为斜方晶格结构，空间群为Pnma，晶格参数a=0.449 nm、b=0.503 nm、c=0.674 nm。初始模型尺寸为20.3 nm×20.3 nm×1.43 nm，厚度设置为bcc-Fe晶格常数的四倍以上，当晶格常数足够大时，可以消除沿厚度方向的额外尺寸效应。为了研究边裂纹的扩展情况，在模型前端预制长3 nm、宽0.4 nm的裂纹，裂纹平面的晶面指数为（-110），即bcc-Fe的取向为x-［011］、y-［0-11］、z-［100］。模型中包含块状或球状的渗碳体夹杂，渗碳体取向为x-［100］、y-［010］、z-［001］，在2.2节中，块状渗碳体的xy平面尺寸为3 nm×3 nm，球状渗碳体的半径r=1.5nm，z方向厚度均与整体模型厚度相同。为了将渗碳体插入基体中，从具有预制裂纹的基体中切出球形或块状空隙，并放置在y方向与裂纹面相同的位置。之后将渗碳体从完美的渗碳体块中切出放置在空隙内。如图1a所示，模型在y方向上下边界设置边界层。拉伸过程中，为避免控压导致模型尺寸变化，采用正则系综（NVT）调控体系状态，使体系温度始终控制在300 K。模型在x方向设置为固定边界条件，y和z方向均设置为周期性边界条件。在y方向上施加应变条件用以实现模型的拉伸，应变率为3×10⁹ s{}^{-\mathrm{1}}。模型拉伸之前，将模型在等温等压系综（NPT）下以300 K的温度弛豫30 ps使模型达到最稳定状态。

图1b为界面拉伸模型示意图。该模型沿x方向可看作两个晶体并排得到，中心便形成了沿y方向的晶界。研究三种不同晶界情况下的中心裂纹扩展情况：①不同晶向排列导致的晶界，其中，左侧bcc-Fe取向为x-［100］、y-［010］、z-［001］，右侧bcc-Fe取向为x-［011］、y-［0-11］、z-［100］；②孪晶晶界，其中一侧晶体与晶界的错向角为α，两晶体间夹角即为2α；③渗碳体与bcc-Fe两相晶界，为了满足渗碳体与bcc-Fe之间的Bagaryatskii晶体取向，重新构造bcc-Fe单个晶胞，通过坐标变换将按x-［100］、y-［010］、z-［001］排列的bcc-Fe晶胞转换成按x-［11-2］、y-［1-10］、z-［111］排列。三种晶界模型的局部图也在图1b中展示，其中渗碳体与bcc-Fe晶界模型为yz平面的渗碳体及重构后bcc-Fe晶胞的平面图。三种晶界模型的初始尺寸（单位为nm）分别为20.2×20.1×1.43、21.1×21.1×1.43、22.5×20.1×3.37，中心均预制一个长（x方向）0.4 nm、宽（y方向）2 nm的裂纹。晶界模型的模拟过程与拉伸模型相同，只是拉伸沿x方向。如图1b所示，为了获得裂纹扩展过程中分离位移与界面应力（内聚力）之间的关系，将晶界左右两侧各1 nm宽的区域定义为内聚力区，通过计算该区域原子的平均位移与平均原子应力，将不同时刻得到的分离位移、内聚力值绘制于同一坐标系中得到内聚力-位移（T-S）曲线。在经典分子动力学模拟中，结果的准确性在很大程度上取决于原子间势函数的选择，原子间势函数用于描述不同原子之间的原子相互作用。对仅含bcc-Fe基体的模型采用文献［17］拟合的嵌入原子势（EAM），对包含渗碳体或渗碳体界面模型采用文献［18］拟合的Tersoff/zbl势函数。为了验证所选势函数的合理性，计算bcc-Fe及渗碳体的弹性常数，并与文献［19-20］中实验和密度泛函理论（DFT）计算所获得的参数进行对比，结果如表1和表2所示，相对误差基本在±10%以内。其中渗碳体弹性常数C₄₄较小，计算误差比较大。由于采用的是拉伸模型，C₄₄弹性常数对整体的影响可以忽略，因此认为所选势函数适用于后续的分子动力学计算。

2 结果与讨论

2.1 基体裂纹扩展机理

首先对bcc-Fe基体在边裂纹情况下进行拉伸模拟，模拟后应力数值通过统计拉伸方向上的原子应力获得，应变数值通过计算拉伸方向上的模型尺寸变化率获得，每100步输出一组应力-应变数值，最终绘制曲线如图2所示，图3为边裂纹bcc-Fe基体模型随应变增加的构型演变图及应力云图。其中，图2中标注的序号①~④对应于图3中的①~④和①’~④’。可以看出，模型的拉伸断裂过程可以大致分为以下几个阶段：线弹性阶段、初始裂纹扩展、脆性断裂、裂纹尖端钝化直至扩展断裂。由图2可以看出，在拉伸过程中，模型起初处于弹性变形阶段，晶格发生弹性变形，如果此时卸掉外力则模型会恢复原貌。当应变ε=0.0465时，裂纹前端开始出现少量的应力集中，并且可以看到在裂纹前端出现原子构型转变，即出现Fcc构型，此时应力-应变曲线仍处于线弹性阶段。当ε=0.0741时，应力达到峰值（11.95 GPa），此时可以在构型演变图中看到裂纹前端出现沿中心向上下两边对称扩展的相变区域，扩展方向与水平方向夹角大约为78°，并且在应力云图中也明显看出该区域的应力较大。这是由于在拉伸过程中，bcc-Fe晶格发生晶体滑移，裂纹尖端上下晶格中的原子沿（110）面滑移。需要指出的是，本文所建模型近似于二维模型，没有考虑原子在z方向的滑移，更准确的结果应通过三维模型及透射电镜试验辅助获得。当ε>0.0741时，应力开始逐渐下降，此时模型的裂纹前端较为尖锐，表现为脆性断裂，如图3中的③所示，此时ε=0.1。之后裂纹进一步扩展，当应变ε=0.12时，可以看到裂纹尖端发生钝化，裂尖出现部分无规则构型，裂纹尖端形貌从V形转变为U形，这表明此时裂纹扩展更具韧性而不是纯脆性断裂，随后裂纹扩展贯通。这与之前部分学者研究的bcc-Fe裂纹扩展形态与机理相似^{［12，21］}，说明该拉伸模型适用于后续研究。

2.2 渗碳体对基体裂纹萌生与扩展的影响分析

为了研究拉伸过程中渗碳体对基体裂纹萌生与扩展的影响，建立了4种含渗碳体夹杂的bcc-Fe基体模型，分别为边裂纹球状夹杂模型、无裂纹球状夹杂模型、边裂纹块状夹杂模型和无裂纹块状夹杂模型，图4、图5分别为这4种模型拉伸过程中的应力-应变曲线、构型演变图和应力云图。4种模型的应力-应变曲线与边裂纹bcc-Fe基体模型的（图2）相似，均为开始经历一段线弹性阶段，到达应力峰值点，之后应力开始下降至平稳。从图5构型演变图中可以看出，在拉伸的初始阶段，含边裂纹的模型（模型一和模型三）与不含边裂纹的模型（模型二和模型四）呈现不同的规律。模型一和模型三在此阶段，裂纹尖端原子沿裂纹扩展面（-110）向y方向两边对称滑移，发生从Bcc到Fcc的相变，这与前文分析的bcc-Fe基体模型拉伸过程相同。当模型含有夹杂时，在裂纹扩展初期应力达到峰值前，夹杂区域随应力的增大逐渐产生微孔洞。从图中可以明显看到当模型一ε=0.0699以及模型三ε=0.0663时，渗碳体夹杂下方出现孔洞，此时应力达到峰值。在对应的应力云图中也可以看到此时裂纹尖端周围的原子滑移区以及渗碳体夹杂区域边界应力较大，而产生孔洞的区域应力得到释放。模型二ε=0.0684和模型四ε=0.0711时应力集中于渗碳体夹杂边界区域，这也导致了在渗碳体两边的基体沿（110）晶面发生了从Bcc相到Fcc相和Hcp相的转变，同时在夹杂区域的下方生成初始裂纹。

图6为4种模型渗碳体区域在萌生孔洞或微裂纹时的原子应力云图，图中虚线圈或框表示渗碳体球或块，实线圆圈表示孔洞或微裂纹区域。可以观察到，随着应变的增大，渗碳体夹杂区域周围由于原子滑移产生了应力集中，即图中红色原子区域。而孔洞或微裂纹则生成在两处应力集中区域的中心位置，可以认为是由于原子滑移位错堆积挤压从而导致夹杂物边缘形成孔洞或微裂纹。并且它们倾向于在夹杂物垂直于拉伸方向的边界中间位置形成（图6a、c、d））。

本文仿真的拉伸过程中Bcc晶格转变示意图见图7，可以看到。在渗碳体夹杂两边发生的Fcc相转变和Hcp相转变是交替出现的，这是由两种晶格的密排面结构相同、密排面之间的堆叠形式不同导致。当应力到达峰值点后，裂纹开始逐渐扩展。图5中模型一在ε=0.075和模型三在ε=0.08时边裂纹发生明显扩展，渗碳体夹杂下方的孔洞也发展为新的裂纹。对应应力云图中也能观察到两裂纹之间以及新裂纹后方的相变区应力显著提高，而已经断裂区域的应力明显得到释放。模型二在ε=0.08和模型四在ε=0.0966时裂纹垂直于拉伸方向沿（-110）面向两边扩展，并观察到在由基体原子滑移产生的Fcc相和Hcp相交替排布的相变区，应力显著高于其他区域。文献［22］研究了纯铁形变诱导马氏体相变晶内和晶界成核的主要因素，在多晶铁的拉伸模型中也发现了交替排布的相变区域，之后含边裂纹模型中的相变区逐渐向模型右侧移动，当到达边界后裂纹迅速扩展。图5中模型一在ε=0.09和模型三在ε=0.1时，裂尖周围原子滑移方向与夹杂处新产生孔洞时周围的滑移方向不同，如画圈处所示，此时裂尖原子的滑移沿（001）面，并且在对应的应力云图中可以观察到该区域的应力偏低。同样，该规律也适用于不含边裂纹模型，图5中模型二在ε=0.1和模型四在ε=0.11的图中画圈处也展现了这种规律，说明渗碳体的夹杂会导致裂纹扩展时裂尖原子周围滑移区的滑移方向发生改变。

图8所示为上述4个模型与边裂纹bcc-Fe基体模型共5个模型的应力-应变曲线对比，可以看出，边裂纹bcc-Fe基体模型的刚度与最大应力值最低，含渗碳体夹杂模型的刚度与最大应力值均有明显提高。这说明渗碳体的夹杂对bcc-Fe基体的力学性能有强化作用，这种强化作用是通过其对位错运动的阻碍作用而表现出来的，当位错遇到夹杂在基体中的渗碳体时，位错运动阻力上升，使材料强度提高。在渗碳体夹杂模型中，无裂纹模型的刚度与最大应力值相较于含裂纹模型有小幅度的上升，这也是因为裂纹这一缺陷的存在导致材料力学性能的降低。对比块状和球状渗碳体夹杂，尽管在无裂纹模型中块状夹杂的最大应力值略高于球状夹杂，但是整体而言两者不论是刚度还是最大应力值都相差不大，可以认为在微观尺度下，这两种夹杂物形状对材料力学性能的影响较小。结合图5可以观察到，块状夹杂下裂纹扩展速率更快。图8中块状（含裂纹）曲线在下降段比球状（含裂纹）曲线更早到达低点也说明了这个观点。这可能是由于相比于球状夹杂，块状夹杂边缘更尖锐，导致应力更加集中。

2.3 渗碳体尺寸大小对裂纹扩展的影响分析

为了研究渗碳体尺寸大小对8Cr4Mo4V基体裂纹扩展的影响，建立5个含不同尺寸大小球状渗碳体的基体裂纹模型。渗碳体的半径r分别为0.5、1.0、1.5、2.0、2.5 nm，所含原子个数分别为106、493、1054、1913、2937。图9为5个模型在拉伸过程的应力-应变曲线，可以看出，渗碳体夹杂的半径越小，最大应力值越大，并且最大应力值对应的应变也越大。这说明在微观尺度下，半径越小的渗碳体夹杂对材料基体中裂纹扩展的阻碍能力越强。

图10为5种模型在渗碳体夹杂区域出现孔洞时的构型演变图及此刻的应变值。从图10中也能看到，随着渗碳体夹杂的半径增大，出现孔洞时的应变逐渐减小。在渗碳体夹杂半径较小的模型中（r=0.5、1.0、1.5 nm），在边裂纹与渗碳体夹杂之间以及夹杂的后方出现Fcc相和Hcp相交替排布的相变区，根据2.2节中的分析可知，这也带来了高的应力区。在渗碳体夹杂半径较大的模型中（r=2.0、2.5 nm），渗碳体边界的孔洞出现时并未伴随这些，这也说明当渗碳体夹杂的半径过大时，在拉伸过程更易产生孔洞。值得注意的是，半径为0.5 nm的渗碳体夹杂模型在拉伸过程中，孔洞出现在渗碳体夹杂的中间，而其他大小夹杂模型的孔洞出现在渗碳体夹杂的上下边界。这是因为当夹杂物尺寸较小且可以变形时，位错会切过夹杂物使其同基体一起变形，此时会出现新的表面，导致界面能升高，并且产生弹性应力场阻碍位错运动。这也说明了较小尺寸的渗碳体夹杂对基体的强化作用更大。

2.4 渗碳体与初始裂纹相对位置对裂纹扩展的影响分析

进一步研究渗碳体与初始裂纹相对位置对裂纹扩展的影响。图11为渗碳体分别位于模型中心、右侧中心、右侧下方和中心下方4个位置的模型在拉伸过程的应力-应变曲线。从图11中可以观察到，4个模型在线性阶段几乎相同，仅有夹杂位于中心下方模型的最大应力值偏低。这说明在微观尺度下，单一渗碳体与初始裂纹的相对位置对基体的力学性能影响不大，只对裂纹的后续扩展有影响。图12显示了这4个模型在裂纹扩展阶段的构型演变，可以看到，裂纹沿与拉伸方向垂直的（110）平面扩展，在渗碳体两侧边缘均形成了2.3节分析的Fcc相和Hcp相交替排布区域。区别是，当渗碳体位置处于模型的中心及右侧中心时，也就是在y方向与预制的边裂纹位置相同时，渗碳体边缘产生的孔洞也近似在这个高度。这就使得裂纹更易与孔洞连接，最终贯通模型。并且在扩展过程中裂纹尖端周围原子的滑移方向也与上文分析相同，即沿（001）面滑移。而当渗碳体位于模型的右侧下方和中心下方时，由于渗碳体边缘产生的孔洞在y方向和预制的边裂纹存在高度差，裂纹和新生成的孔洞都会先沿垂直于拉伸方向的平面进行扩展，直至二者在x方向扩展到相同位置，如图12中渗碳体位于中心下方的模型构型演变所示。此时，在边裂纹与孔洞扩展形成的裂纹之间形成一个平行四边形的相变区，仍是由基体原子沿（001）面滑移得到的Fcc相和Hcp相交替排布形成。之后，随着模型在y方向的应变增大，相变区逐渐缩小，直至裂纹贯通。

2.5 双晶界面裂纹扩展分析

建立3种双晶模型来研究双晶界面裂纹的扩展。图13所示为3种双晶模型在裂纹扩展时3个应变状态的构型演变。

图13a所示为两种bcc-Fe晶向排列的双晶界面，左侧bcc-Fe取向为x-［100］、y-［010］、z-［001］，右侧bcc-Fe取向为x-［011］、y-［0-11］、z-［100］。当ε=0.054时，晶界上裂纹周围作为应力集中位置，开始激活局部的滑移系统，裂纹两端开始发生Fcc相变，值得注意的是相变区域均处于左侧晶体内，这说明左侧的晶向在相同的拉伸应力下更容易发生原子滑移。当应变增大到0.064时，随着应力的增加可以明显观察到裂纹两端周围的原子在左侧晶体中沿｛110｝晶面族滑移，形成Fcc相变区。此时裂纹开始逐渐沿晶界扩展，模型发生塑性变形。当ε=0.075时，左侧晶体内Fcc相变区域进一步增大，已经滑移至模型上下边界。而右侧晶体中仅在裂纹两端出现少量沿（110）面滑移形成的Fcc相变区。这表明右侧bcc-Fe的取向在拉伸过程中保持了更好的结构稳定性，表现出更高的抗裂纹扩展能力。同时，裂纹扩展的过程是通过晶界附近滑移系统的激活和原子重排进行，表明在受拉作用下，材料的塑性变形主要通过晶界处的滑移机制来调节。

图13b所示为bcc-Fe孪晶界面在拉伸过程的裂纹扩展构型演变情况，其中两侧晶体与孪晶面的夹角α=35°。当ε=0.060时，晶界与裂纹两端应力相对集中，原子开始滑移，发生少量的Fcc相和Hcp相转变。随着应力的进一步增大，应变到达0.080时，在裂纹下方可以明显观察到沿（110）面滑移而产生的Fcc相变区，而在裂纹上方的原子沿晶界同样发生滑移，形成Fcc相与Hcp相交替排布的相变区，这与上文渗碳体夹杂模型中所见情况类似，说明此处应力更大，并且裂纹还未扩展。当ε=0.090时，可以看到此时裂纹沿-y方向扩展，相变区域进一步扩大，这也进一步证明了上面分析结论的正确性。

图13c所示为渗碳体与bcc-Fe双晶界面拉伸过程的构型演变，左侧为渗碳体结构，由于它并不属于常见的3种排列结构，因此在Ovito软件CNA分析中呈现为灰色原子。当ε=0.041时，晶界处的bcc-Fe晶体出现少量Fcc相变，随着应力的增加，当ε=0.056时，晶界周围原子的滑移区域变大。与先前两种模型不同的是，此时相变区并不是集中在裂纹两端或沿裂纹两端延伸，而是较为平均地分布在晶界上，并且呈现沿裂纹所在x方向对称排布于晶界的构型。这也导致随后裂纹的扩展并不是由裂纹尖端发展，而是在裂纹两端的前方分别先产生新的孔洞，之后再沿晶界贯通，如图13c中的③所示。

图14分别给出了3种双晶模型在裂纹扩展过程中3个应变时刻的平均原子拉伸应力沿y轴的分布曲线。其中，模型的下方定义为原子位置为0 nm。从图14中可以明显看出预置裂纹区域，即原子位置处于模型中间（10 nm左右）的区域不论在哪种应变状态下均为低应力区，且随着应变的增加，尤其是当裂纹开始扩展后，裂纹区域的应力进一步降低。在两种晶向排列的双晶界面模型拉伸过程中，当应变到达0.064之后，裂纹开始扩展，如图14a所示，此时裂纹两侧到模型边缘的平均原子应力值基本不会随应变增加到0.075而发生变化，仅是低应力区域分别向+y和-y方向呈现出对称性的扩大，这也说明此时裂纹沿晶向向两边扩展的速率相同。该规律同样出现于图14c渗碳体与bcc-Fe双晶界面模型中，区别是渗碳体与bcc-Fe双晶界面中裂纹扩展得更快，且平均原子应力值也更高。而图14b孪晶（α=35°）界面模型随着应变的增大，-y方向（下方）的应力值更低，说明裂纹向这个方向扩展的速率更快。同时可以观察到在原子位置0~5 nm之间又出现一个应力峰值，这可能也是由于晶界周围原子滑移而导致的。

2.6 双晶界面内聚力区参数获取与分析

内聚力区的概念来自于内聚力模型，内聚力模型是基于弹塑性断裂力学，提出在裂纹尖端存在一个微小的内聚力区，在内聚力区内，裂纹两侧物质原子或分子之间存在相互作用的内聚力。内聚力区的尺寸是一个相对很小的常数，且与加载方式无关。内聚力模型中裂纹面的应力为张开位移的函数，称为开裂界面上的内聚力-位移关系，也就是内聚力-位移曲线。临界断裂能G_IC、初始损伤刚度K₀、损伤应力T_max、损伤位移S_in、及开裂位移S_fail等是内聚力模型计算时的参数，而这些参数均可以通过内聚力-位移曲线获得，如图15所示。

图16分别展现了2.5节所描述的3种双晶界面模型的内聚力-位移曲线，其中，内聚力通过对模型中内聚力区原子的维里应力进行平均获得，分离位移通过界面两侧不同晶体的平均原子位移的差获得。将所得结果用双线型内聚力模型拟合，不同双晶界面模型的内聚力区参数如表3所示，模型1代表两种晶向排列的双晶界面，模型2代表孪晶（α=35°）界面，模型3代表渗碳体与bcc- Fe双晶界面。从表3中可以看出，渗碳体与bcc-Fe双晶界面的损伤应力值最大，为12.5011 GPa，说明相较于bcc-Fe基体间的界面，存在渗碳体的界面裂纹扩展需要更大的应力，这与文献［23］的研究结果相似。孪晶（α=35°）界面具有最大的临界断裂能，为9.914N/m，这说明这种界面在裂纹扩展过程中需要吸收更多的能量。而两种晶向排列的双晶界面的损伤应力和临界断裂能都最小，表明最易发生裂纹的扩展。

为了研究孪晶中晶界错向角的大小对内聚力区的影响，分别建立了晶界错向角从10°到40°的7个模型，并绘制其内聚力-位移曲线，如图17所示。表4给出了各个模型的损伤应力与临界断裂能。从图17和表4中可以看出，晶界错向角较大时（α=35°和α=40°），损伤应力与临界断裂能都较大，分别为11.0262 GPa、9.914 N/m和11.5953 GPa、10.851 N/m。晶界错向角为10°、15°、25°和30°时的损伤应力值非常接近，临界断裂能除了10°以外也都比较接近，这是因为晶界错向角为10°的模型在拉伸过程中裂纹后期向一侧晶体扩展，导致内聚力区的原子平均应力急速下降。而当晶界错向角为20°时，具有最小的损伤应力，并且临界断裂能也显著小于其他模型，这可能是因为在该晶界错向角下模型的晶界附近原子及位错的不同滑移形式导致的。

值得说明的是，通过分子动力学计算得出的内聚力-位移曲线能有效揭示裂纹尖端的复杂变形机制，如位错发射和滑移等，因此可将其作为构建多尺度计算模型的桥梁。例如，将内聚力模型与有限元方法结合，可有效分析微观缺陷如何影响宏观层面的应力强度因子和断裂韧性等断裂特性，进而为更加复杂的断裂机制的探索提供新的视角。

3 结论

1）基体裂纹行为的进展： 在含边裂纹的情况下，8Cr4Mo4V轴承钢基体的裂纹扩展过程呈现明确的阶段性，包括线弹性阶段、初始裂纹扩展、脆性断裂以及裂纹尖端钝化等。这一发现对理解材料在实际应用中的裂纹扩展机理具有重要意义，尤其是在工程设计中评估裂纹扩展风险时具有实际应用价值。

2）渗碳体对力学性能的影响：渗碳体在8Cr4Mo4V轴承钢基体中的强化作用，尤其在没有初始裂纹时，显著提高了材料的刚度和最大应力值。研究发现，在拉伸应力达到峰值时渗碳体夹杂物周围的原子沿（110）面滑移形成Fcc相与Hcp相交替排布的区域，在垂直于拉伸方向的表面会生成新的孔洞或裂纹；在裂纹扩展阶段，裂纹周围原子沿（001）面滑移，形成的新缺陷与边裂纹贯通或向模型两端扩展。

3）渗碳体的形状在微观尺度对8Cr4Mo4V基体的力学性能没有明显的影响，但块状相比于球状增大了裂纹扩展速率；渗碳体的尺寸大小对8Cr4Mo4V基体的力学性能和断裂特性具有一定的影响，碳化物尺寸越大，模型的最大应力值越低，越容易在夹杂区域出现新的孔洞，当尺寸较小（r=0.5 nm）时，孔洞出现在夹杂物中，裂纹会贯穿渗碳体夹杂。在实际应用中控制渗碳体的尺寸和位置对提高材料的抗裂纹能力至关重要。

4）不同类型的双晶界面对裂纹扩展具有不同的影响，特别是孪晶（α=35°）界面和渗碳体与bcc-Fe双晶界面模型的差异，从而提供了对裂纹扩展方向和速度的深刻理解。

5）通过对双晶界面内聚力区的计算发现，不同类型的双晶界面对裂纹萌生的抵抗力有所不同，研究表明，孪晶界面对裂纹扩展的抑制作用更强，尤其在晶界错向角较大的情况下。

6）本研究还提供了8Cr4Mo4V轴承钢基体与渗碳体界面的内聚力参数，进一步揭示了内聚力参数在评估微观缺陷对宏观力学行为影响中的重要作用，为未来的多尺度材料设计提供了新的理论支持。结合有限元分析方法，未来可以更深入地探讨微观缺陷对宏观应力强度因子、断裂韧性等力学性能的影响。这不仅对提高轴承钢等材料的可靠性和耐久性具有重要意义，也为相关工业应用中材料的优化设计提供了理论支持。

参考文献

[1] BHADESHIA H K D H. Steels for Bearings［J］. Progress in Materials Science， 2012， 57 （2）：268-435.

[2] PANDKAR A S， ARAKERE N， SUBHASH G. Microstructure-sensitive Accumulation of Plastic Strain due to Ratcheting in Bearing Steels Subject to Rolling Contact Fatigue［J］. International Journal of Fatigue， 2014， 63：191-202.

[3] XU Funing， DING Ning， LI Nan， et al. A Review of Bearing Failure Modes， Mechanisms and Causes ［J］. Engineering Failure Analysis， 2023， 152：107518.

[4] YUE Xiong， HU Shan， WANG Xiaokang， et al. Understanding the Nanostructure Evolution and the Mechanical Strengthening of the M50 Bearing Steel during Ultrasonic Shot Peening［J］. Materials Science and Engineering A—Structural Materials Properties Microstructure and Processing， 2022， 836：142721.

[5] 孙玉凤， 刘伟军， 张宏伟， 等. 8Cr4Mo4V钢激光冲击残余应力的演化仿真及其对疲劳性能的影响 ［J］. 材料研究学报， 2023， 37 （12）：933-942.

[6] SUN Yufeng， LIU Weijun， ZHANG Hongwei， et al. Simulation of Residual Stress Evolution of 8Cr4Mo4V Steel Induced by Laser Shock and Its Influence on Faticue Performance ［J］. Chinese Journal of Materials Research， 2023， 37 （12）：933-942.

[7] YANG Liqi， XUE Weihai， GAO Siyang， et al. Rolling Contact Fatigue Behavior of M50 Bearing Steel with Rare Earth Addition ［J］. International Journal of Fatigue， 2023， 177：107940.

[8] GUO Wei， MA Tianyu， CAO Hongrui， et al. Numerical Analysis of Rolling Contact Fatigue Crack Initiation Considering Material Microstructure ［J］. Engineering Failure Analysis， 2022， 138：106394.

[9] LIU Yazhou， LUO Yun， SU Shuo， et al. Molecular Dynamics Simulation of Phase Transition and Crack Propagation in Metastable High Entropy Alloy ［J］. Materials Today Communications， 2022， 33：104642.

[10] MOLAEI F. Molecular Dynamics Simulation of Edge Crack Propagation in Single Crystalline Alpha Quartz［J］. Journal of Molecular Graphics & Modelling， 2022， 111：108085.

[11] XING Zheyuan， FAN Haidong， KANG Guozheng. Molecular Dynamics Simulations on the Intergranular Crack Propagation of Magnesium Bicrystals［J］. Computational Materials Science， 2022， 210：111058.

[12] MA Lei， XIAO Shifang， DENG Huiqiu， et al. Molecular Dynamics Simulation of Fatigue Crack Propagation in bcc Iron under Cyclic Loading［J］. International Journal of Fatigue，2014， 68：253-259.

[13] 曹莉霞， 王崇愚. α-Fe裂纹的分子动力学研究［J］. 物理学报， 2007（1）：413-422.

[14] CAO Lixia， WANG Chongyu. Molecular Dynamics Simulation of Fracture in α-iron［J］. Acta Physica Sinica， 2007 （1）：413-422.

[15] GHAFFARIAN H， TAHERI A K， KANG K， et al. Molecular Dynamics Simulation Study of the Effect of Temperature and Grain Size on the Deformation Behavior of Polycrystalline Cementite［J］. Scripta Materialia， 2015， 95：23-26.

[16] NAKAMURA K， KUMAGAI T， OHNUMA T. Atomistic Simulation of Shear Deformation at bcc-Fe Grain Boundary and Precipitation Strengthening by Cr₂₃C₆ ［J］. Materials Today Communications， 2022， 33：104711.

[17] PLIMPTON S. Fast Parallel Algorithms for Short-range Molecular-dynamics［J］. Journal of Computational Physics，1995， 117 （1）：1-19.

[18] STUKOWSKI A. Visualization and Analysis of Atomistic Simulation Data with Ovito—the open Visualization Tool［J］. Modelling and Simulation in Materials Science and Engineering， 2009， 18（1）：015012.

[19] MENDELEV M I， HAN S， SROLOVITZ D J， et al. Development of New Interatomic Potentials Appropriate for Crystalline and Liquid Iron［J］. Philosophical Magazine，2003， 83 （35）：3977-3994.

[20] HENRIKSSON K O E， BJÖRKAS C， NORDLUND K. Atomistic Simulations of Stainless Steels：a Many-body potential for the Fe-Cr-C System［J］. Journal of Physics-condensed Matter， 2013， 25（44）：445401.

[21] LEE B J， BASKES M. Second Nearest-neighbor Modified Embedded-atom-method Potential［J］. Physical Review B，2000， 62 （13）：8564-8567.

[22] JIANG C， SRINIVASAN S G， CARO A， et al. Structural， Elastic and Electronic Properties of Fe₃C from First Principles［J］. Journal of Applied Physics， 2008， 103（4）：2127.

[23] ZHU Jiaqi， HE Xiaoqiao， YANG Dong， et al. A Peridynamic Model for Fracture Analysis of Polycrystalline bcc-Fe Associated with Molecular Dynamics Simulation［J］. Theoretical and Applied Fracture Mechanics， 2021， 114：102999.

[24] SHIMOKAWA T， FUJII K， NIIYAMA T. Atomic Simulation Study of the Factors Affecting Nucleation in Deformation-induced Martensitic Transformation in Grains and at Grain Boundaries in Pure Iron［J］. Acta Materialia， 2024：265：119629.

[25] WANG Huiling， WANG Feng， QIAN Dongsheng， et al. Investigation of Damage Mechanisms Related to Microstructural Features of Ferrite-cementite Steels via Experiments and Multiscale Simulations ［J］. International Journal of Plasticity， 2023， 170：103745.

基金资助

国家科技重大专项(J2019-IV-0004-0071)

材料服役安全研究评价设施开放项目(MSAF-2023-004)