面接触应力分布与磨损分析

黄文强 李卫 毛飞宇 黄平

(华南理工大学 机械与汽车工程学院，广东 广州 510640)

在机械工作领域，机器零件接触表面通常从静态转变为动态。此时，接触表面的性质(如应力场、润滑层、磨损等)改变，会影响表面摩擦性能以及机器零件的使用寿命。可以认为，研究表面过渡过程中这些性质随着接触状态的变化有着重要的意义，尤其是接触应力场的演变，因为该因素决定工作表面的损伤趋势[1]。

在接触应力分布的研究中，通常强调接触过程的应力集中。寻找应力集中区域的目的是阐明接触的机械性能，并解释在零件上某些位置产生故障的原因。对于经典的Hertz接触，圆柱与弹性半空间相接触的静态出现最大剪切应力集中在表面下深度约0.78a处，a为一半的接触长度，该区域可能是疲劳磨损的起点[2- 3]。在滑动摩擦状态下，表面下的应力有偏移趋势[4- 5]。通过光弹实验法，对于圆柱与平板接触，可观察到平面应力场从聚焦于接触区域中间偏转到接触处前缘[6- 7]。对于冲头与弹性半空间相接触理论研究，在无摩擦下，压力轴对称分布，集中在冲头的两边[8]，而在有摩擦下压力轴不再对称分布[9]。在二维滑块摩擦模型的研究中，假设接触压力沿滑块长度方向为线性分布，则在水平力作用下接触压力形成三角形分布，而大边缘在运动前方[10]。这些研究主要探讨在单独静态或动态接触条件下的二维模型产生的应力场。

对于三维接触应力分布的研究，人们常用数值分析或实验方法来探讨表面上的接触压力分布以及表面下的应力场。Stepanov等[11]采用数值方法研究了钢球与粘弹性半空间相接触模型，分析了滑动速度、摩擦系数、材料特性等因素对表面压力分布的影响。Bhushan等[5,12- 13]研究了粗糙表面的弹塑性接触模型，发现接触压力分布在不同的接触峰。有些研究采用商业模拟软件(如Ansys、Abaqus等)分析三维接触模型，从仿真结果可直观看到不同场力的分布情况[14- 17]。黄健萌等[14]采用Ansys软件分析了盘式制动器表面上的接触压力分布。Londhe等[15]采用Abaqus软件进行三维接触仿真，扩展了接触力学的Hertz理论，并研究了其对轴承疲劳寿命的影响。Koumi等[16]采用Abaqus软件仿真分析了钢球滚动在粘弹性半空间上的三维压力分布。Yin等[17]采用数值方法分析速度、时间和钢球形状对接触压力分布的影响，发现压力分布明显不是中心线对称的。对于实验方法，可应用Tekscan压力传感器[18- 20]、感压胶片法[21- 22]等直接测量接触表面上的压力分布，但这些实验法存在如下限制：除剪切应力外，用Tekscan压力传感器可测量所有压力；感压胶片只能测量静态最大接触压力分布，不能测量动态接触压力以及剪切应力。这些研究主要涉及静态、动态或在不同时间的接触应力场研究，往往不涉及面接触应力随着状态的变化如何变化。

目前，光弹实验法被广泛应用于平面应力场分析，可以直接观测内部主应力差分布情况。文献[23]对圆柱与平面接触模型采用光弹实验观测了涂层/基体的接触应力分布。对于面接触磨损的研究，销盘试验机作为最普遍的面接触试验机[24]。一般考虑摩擦副、载荷、速度、温度等因素对磨损率的影响，很少关注表面破坏情况。

为研究面接触应力分布以及磨损发生情况，本文对三维滑块摩擦模型(该模型基于机械装配中的楔键连接[25]，与工程实际比较接近)进行了有限元仿真、光弹实验和磨损实验，利用Abaqus软件通过有限元法对干摩擦面接触进行仿真，得到主应力差分布、接触压力分布、摩擦应力分布，同时分析了载荷和摩擦系数对接触应力场的影响；采用光弹实验法观测平面主应力差分布，该实验应力场相当于模拟模型对称平面处的平面应力；通过面接触磨损实验观察接触压力集中所对应的磨损，并观测表面损伤扩展情况。

1 分析模型

在摩擦学中，滑块与平面相对滑动是描述摩擦现象的经典模型。对于该模型，滑块通常受到垂直和水平外力的作用，因此接触表面上出现正压力及切向应力(摩擦应力)，滑块与平面两接触表面应力的值相等，而方向相反。

在接触力学中，接触应力分布问题已被广泛的研究并形成了许多理论，基本上集中在理论方法、实验方法和数值方法3个方面。最近，随着计算机计算速度的提高，有限元法在模拟接触应力分布问题中的应用非常广泛，也取得了不少成果。为了便于研究，本文利用有限元法模拟接触应力分布问题，采用Abaqus软件进行仿真。

本文的研究模型如图1所示，主要由夹具、滑块和块体组成。夹具用来给滑块加载，沿oz轴垂直移动，而限制ox和oy轴的水平方向移动。夹具上面受到压力p的作用，因此在一般情况下夹具向下走。当滑块承受夹具的垂直力时，滑块底面被压而与块体接触。当块体从左向右移动时，滑块同时也将承受夹具的悬出部分施加的水平力。块体代表无限弹性半空间(足够大)，其底部设置限制上下运动边界条件，静态下块体不运动，而在运动下其从左向右运动。

根据滑块的对称性，在仿真分析中只对模型的一半进行建模，将对称边界条件设在对称截面，如图1(b)所示。模型各部分的尺寸分别为：夹具、滑块和块体的尺寸分别为12 mm×4 mm×4 mm、4 mm×4 mm×2 mm、80 mm×50 mm×40 mm，夹具的悬出部分长2 mm。模型的材料性质如下：夹具是钢制的，其密度ρ1=7.85×103kg/m3，弹性模量E1=2.10×105MPa，泊松比ν1=0.29；滑块由聚碳酸酯制成，ρ2=1.20×103kg/m3，E2=1.78×103MPa，ν2=0.36；环氧树脂块体的密度ρ3=1.19×103kg/m3，E3=2.10×103MPa，ν3=0.40。

为了保证精度，同时优化计算时间，在特定区域(滑块与块体接触区)进行更精细的网格划分，最小网格尺寸为0.05 mm。网格划分的计算结果如下：夹具、滑块、块体的单元数分别为633、173 800、606 253。相互作用区域符合库仑摩擦定律，假设每个接触点的摩擦系数不变，滑块与夹具的摩擦系数为0.25，滑块与块体的摩擦系数为0.35。这些摩擦系数从销盘式摩擦磨损实验机测得[7]。施加在夹具上面的压力p为1.29、1.42和1.56 MPa，分别对应垂直载荷为92.8、102.6和112.4 N。设块体的两侧表面的运动边界条件速度为0.5 mm/s，沿ox轴方向。为防止加载过程中产生应力波，压力p和速度v的幅值从0稳步增加到目标值1，并稳定在1，如图2所示。

根据以上描述，将仿真分为100个时间步长，总仿真时间为0.3 s。采用Intel(R)Core(TM)i3- 7100 CPU、RAM 8.00 GB的计算机进行计算仿真，此过程必须连续运行17 h以上才能完成一次仿真。

图1 模型示意图Fig.1 Schematic diagram of the model

图2 压力和速度的平滑阶跃振幅Fig.2 Smooth step amplitudes of pressure and velocity

2 仿真结果与分析

2.1 主应力差分布

主应力差、Tresca应力或最大剪应力的分布常用于判断材料在载荷作用下的强度、破坏和应力集中趋势。为了直观地描述在块体接触区主应力差的分布随着接触状态的改变而变化，在此采用三维等值面图显示，如图3所示。从图中可以看出，在静态接触状态下，主应力差(S)集中在滑块边缘周围，由于滑块具有对称性，所以应力场呈对称分布。当改变状态时，滑块与块体接触表面的剪切应力增大，然后达到最大值(此时发生相对滑动)，在这过程中主应力差逐渐偏移，集中在滑块的相对运动方向前缘，而尾缘周围的应力最小。此外，主应力差的最大值也随着增加，可观察到图中等值面层的数量在增多，并且外等值面层(3.6 MPa)更大。当载荷从92.8 N增加到102.6 N再达到112.4 N时，主应力差的最大值也随之增大，等值面的范围越来越大，颜色越来越深。

在载荷为102.6 N、滑块与块体摩擦系数f分别为0.20、0.35和0.50时主应力差分布如图4所示。从图中可以看出，在滑动状态下，随着摩擦系数的增大，主应力差越趋集中在运动前方，最大值分别为7.2、14.4 和24.3 MPa。

2.2 接触压力分布

在摩擦学中，表面接触压力分布对于磨损计算有着重要的意义，承受高压力的接触区损伤更严重[1]。因此，从观察接触压力分布可确定摩擦过程中哪些接触区发生损伤大。仿真结果采用带状等高线图显示块体表面的接触压力分布，如图5所示。从图中可知：当静态接触时，接触压力(S22)集中在滑块下的边缘周围，最大值(绝对值)在边角处，而最小值在中间区；当从静态向滑动转变时，接触压力偏移集中在滑块相对运动方向前缘，后缘的接触压力更小，可达到0(见图6(c))；当完全滑动时，前缘的接触压力最大，而在侧缘向后缘逐渐变小形成U形分布区(模型中得到半U形)。由于磨损量与接触压力成正比，因此可判断在滑动摩擦过程中，滑块前缘可能获得最大损伤。

在滑动时，摩擦系数对接触压力分布的影响非常显著，如图6所示，其中载荷为102.6 N。从图中可知：随着摩擦系数的增大，接触压力更明显集中在运动前缘，而尾缘更小，如摩擦系数过大，滑块后缘可不再与块体接触(见图6(c)，当f=0.50时，尾缘处的接触压力为0)；接触压力最大值分别达到15.8、22.9和27.2 MPa，表明接触压力最大值随着摩擦系数的增大而增大。

图3 由静态到滑动的演变过程中主应力差分布的三维等值面图Fig.3 Three-dimensional isosurface map of principal stress difference distribution in transition from static to sliding

图4 载荷为102.6 N时摩擦系数对主应力差分布的影响Fig.4 Influence of friction coefficient on distribution of principal stress difference under load of 102.6 N

2.3 摩擦应力分布

在仿真中，块体接触表面上的剪切应力是摩擦应力，该场力在xoy面，摩擦应力分布的等高线图如图7所示。从图中可知，静态下摩擦应力(S12)集中在前缘和后缘，而中间区的摩擦应力趋于0(叫中性区)。该结果与摩擦学的粘滑现象以及微动磨损的研究一致[3]。

在静载荷作用下，由于弹性膨胀的作用，滑块的右侧面与夹具的悬出部分接触，从而将滑块向左移动导致接触表面产生一定的残余应力，因此在静态块体接触表面上的应力不完全对称。这也导致了主应力差分布不对称(见2.1节)，同时在xoy面上的接触压力分布也会偏移(见2.2节)。

图5 由静态到滑动的演变过程中接触压力分布的带状等高线图Fig.5 Strip contour map of contact pressure distribution in transition from static to sliding

图6 载荷为102.6 N时摩擦系数对接触压力分布的影响Fig.6 Influence of friction coefficient on distribution of contact pressure under load of 102.6 N

图7 由静态到滑动的演变过程中摩擦应力分布的等高线图Fig.7 Contour map of frictional stress distribution in transition from static to sliding

当从静态转变到动态时，中性区逐渐向尾部移动，然后在滑动状态下完全消失。此外，摩擦应力分布集中在运动的前缘，与接触压力分布相同。当载荷为92.8、102.6 和112.4 N时，最大摩擦应力增加分别达到6.5、7.1和7.8 MPa。

当摩擦系数增大时，摩擦应力越趋集中在前缘，其最大值也随之增大，如图8所示。当摩擦系数过大时，滑块后缘可不与块体接触，从而导致摩擦应力为0，这结果与前面的接触压力分布分析相似。

图8 载荷为102.6 N时摩擦系数对摩擦应力分布的影响Fig.8 Influence of friction coefficient on distribution of frictional stress under load of 102.6 N

3 平面应力场实验

由于所研究模型的对称性，将对称截面的应力场视为二维平面应力(见图1(b)中的xoz面)，可利用光弹实验法观察平面应力分布。光弹实验法的本质是通过光系统的圆偏振光场来反映材料内部的主应力差分布(可参考文献[23]的光弹实验方法与理论)。光弹实验原理图见图9，主应力差计算式为[23]

(1)

式中：σ1与σ2分别为两个主应力，MPa；fσ为材料应力条纹系数，N/(mm·条)；n为等差线条纹次数；d为试样厚度，mm。

实验中材料的性质与仿真所设的材料一致(见第1部分)，在此平板环氧树脂作为光弹材料，其光弹应力条纹系数为14.6 N/(mm·条)[7]。实验各部件的尺寸如下：平板和滑块的尺寸分别为80 mm×50 mm×4 mm、4 mm×4 mm×10 mm，夹具厚10 mm，

图9 光弹实验原理Fig.9 Principle of photoelasticity experiment

夹具悬出部分的长度为2 mm。在室温、干摩擦条件下进行实验。所有实验结果由CCD摄像机记录，记录速度设为120 f/s。光弹实验结果如图10所示，图中显示了从静态接触转变到滑动接触状态的干涉条纹，可观察到主应力差的分布由聚焦于两点变为集中在滑块的相对运动方向前缘。根据式(1)可算出每条条纹的主应力差的值(见图10)。在仿真结果中，xoz平面上主应力差的分布(Tresca应力)如图11所示。图中的等高线间隙为3.6 MPa，与实验图相同。图10和图11中的曲线较相似，这验证了仿真结果的正确性。

图10 光弹实验结果Fig.10 Photoelasticity experimental results

图11 对称截面主应力差分布的仿真结果Fig.11 Simulation results of principal stress difference distribution of symmetrical section

4 表面磨损实验

根据磨损理论，磨损微观机理十分复杂，磨损的影响因素有许多，如接触压力、材料性能、环境温度等。其中接触压力为主要因素，磨损量与接触压力成正比[1]。因此，在相同接触条件下，如整个表面上的任何区域磨损较大，则表明接触应力集中在该位置。为了验证和分析接触应力集中以及应力分布的规律，本文建立磨损实验并观察表面滑块的损伤情况。

磨损实验装备如图12所示。上下面非常平整光滑的透明聚碳酸酯样品作为滑块(见图13(a))，其尺寸为6 mm×6 mm×6 mm。采用砂纸2000Cw贴合在钢平面上用来研磨。通过砝码加载(总载荷为P)将滑块与砂纸接触，然后通过以角速度ω旋转螺纹杆，将水平运动v传递给滑动部分。当滑块与砂纸相对滑动时，滑块表面发生损伤。

磨损实验前后滑块表面的磨损情况和三维形貌(由TALYSURF CLI 1000轮廓仪扫描滑块表面形貌测得)分别如图13、14所示。在磨损前，滑块的表面上下都是透明光滑的，见图13(a)、图14(a)。经过磨损实验后，滑块下方出现损伤形成划痕，见图13(b)。从图14(b)可明显看出，磨损区形成斜面。从实验中可观察到，磨损主要集中在黄色边缘(见图12和图13(b)、13(c))。该边缘是运动方向的前缘。这表明，在运动过程中，接触压力主要集中在运动前缘，进一步证实了上述仿真结果的正确性。

图13(c)给出了载荷为33.9 N时表面滑块损伤随着相对滑动距离的变化情况。从图中可以看出，随着相对距离的增加，表面滑块损伤从运动方向前缘向后缘扩展，磨损面积不断增加，最后覆盖整个接触面。这可解释如下：开始磨损时，接触压力集中在运动方向前缘，故受到最大损伤，磨平后，该区域的接触压力减小，为承受载荷，接触区向后缘扩展，导致损伤也向后缘展开。这表明，在接触磨损过程中，随着磨损距离的增加，损伤面不断向后增大以达到最佳接触(均匀接触压力)，从而避免应力集中。

图12 磨损实验装置Fig.12 Wear experimental setup

在不同载荷作用下，磨损面积占名义接触面积的比例(简称磨损面积比例)随磨损距离的变化如图15所示。从图中可知：随磨损距离的增加，磨损面积比例呈上升趋势，最后达到100%；当磨损距离相同时，磨损面积比例随着载荷的增大而增加。因此，若接触面所施加的载荷越大，则损伤程度越高，与磨损理论一致。

图13 滑块的磨损情况Fig.13 Worn of slider

图15 磨损面积比例随磨损距离的变化Fig.15 Variation of wear area ratio with wear distance

在静载荷作用下，接触应力集中在边缘周围是因为端部效应，这已在文献[3，26]中提出。本文采用的滑块是方形，较圆形端部的应力集中更显著[3]。机械零件端部形状一般采用圆形[25]。

由静态向滑动转变时，接触应力逐渐偏移到运动方向前缘。这与多位学者的研究结果一致。例如，对于二维滑块摩擦模型的理论研究，Huang等[10]提出了接触压力分布形成三角形，而大边缘在运动方向前缘；对于球体滑动在平面上，Hamilton等[4]绘制了接触应力图像，发现应力集中于运动前方；对于盘式制动器，在摩擦力作用下其界面接触压力分布明显集中在入口处，也是运动前缘[14]。此外，在动态线接触光弹实验研究中，也观察到接触应力偏移到接触区域前缘的现象[6- 7]。另外，本文的模型虽然与文献[4，10]的几何模型不同，但摩擦系数的影响规律一致。

本文的面接触磨损实验结果表明，运动前缘的应力集中会导致表面接触受到破坏，呈不均匀。开始滑动磨损时首先前缘被破坏，然后损伤区随磨损距离的增大逐渐向后扩散。这一规律对研究及控制摩擦磨损具有重要的意义。

5 结论

本文建立三维滑块摩擦模型，利用有限元法分析了面接触应力分布，并进行了平面应力场的光弹实验以及磨损实验，得出如下结论：

(1)仿真结果表明，随着接触状态由静态向滑动转变，接触主应力差、表面压力和摩擦应力从集中在滑块边缘周围逐渐偏移集中在运动方向前缘，载荷及摩擦系数对接触应力分布的影响比较显著。

(2)平面应力场的光弹实验结果表明，接触应力集中在运动前缘，在滑动状态下，运动前缘的干涉条纹级数最大，越向后越少，这表明接触几乎只发生在运动前方。

(3)表面接触磨损实验结果表明，开始接触磨损时，滑块前缘出现最大损伤，对应的接触应力最大；随着磨损距离的增加，表面损伤向后扩展，最后覆盖整个接触面；接触面所受的载荷越大，磨损程度越高，这与磨损理论一致。