初始粗糙面對切向微動行為影響數(shù)值研究

黃帥宗，劉 娟，沈火明

（西南交通大學力學與工程學院，成都 610031）

微動（（Fretting）是指機械振動、疲勞載荷、電磁振動或熱循環(huán)等交變載荷作用下，接觸界面間發(fā)生的振幅極?。ㄒ话阄灰品禐槲⒚琢考墸┑南鄬\動［1］。切向微動大量存在于實際工程中，其運動幅值微小，且產(chǎn)生的影響不易觀察，但其引起的磨損和疲勞損傷所造成的零部件磨損失效不可忽視。

在機械工程中，表面輪廓形貌在微動磨損機理研究中扮演著非常重要的角色，通過對接觸體表面輪廓形貌分析，可以確定材料的磨損狀態(tài)，并能從磨損表面形貌信息中推測出材料損傷失效形式［2］。近年來，一些研究學者發(fā)現(xiàn)表面形貌在金屬微動磨損機理研究中有著重要的作用［3-5］。

Kubiak K J［6］（2010）等人利用通過銑削和光滑磨料拋光過程所形成的表面，在切向微動中研究了全滑移和部分滑移條件下初始表面粗糙度與加工過程對微動的影響。隨后，Kubiak K J［7］（2011）等人引入磨損活化能的概念，通過試驗研究了初始表面粗糙度對摩擦和磨損過程的影響。觀察到較高的表面初始粗糙度會使得摩擦系數(shù)降低并提高表面材料磨損率。對于越光滑的表面，磨損激活能量越大。在全滑移過程中，較低的表面初始粗糙度可以延緩磨損過程的激活并降低磨損率，但是也會增大摩擦系數(shù)。

由于試驗研究很難發(fā)現(xiàn)微動磨損過程中接觸表面的演變規(guī)律，因此可以使用有限元方法進行預測。Yue Tongyan［8］等人（2014）對微動磨損中部分滑移和全滑移狀態(tài)下應(yīng)力奇異性進行有限元分析。發(fā)現(xiàn)應(yīng)力奇異性與微動機制密切相關(guān)，當接觸狀態(tài)處于部分滑移時，應(yīng)更加注意調(diào)整接觸界面的網(wǎng)格尺寸來減少應(yīng)力奇異性的產(chǎn)生；Pereira K［9］（2017）等人使用多尺度分析的方法，基于有限元模型研究了粗糙度對微動磨損的影響。結(jié)果表明，隨著粗糙度的增加，有效接觸壓力峰值減小而表觀接觸寬度增加，對于較高的法向載荷條件下，有效接觸壓力峰值更接近赫茲解析值；Qin W［10］（2018）等人通過W-M函數(shù)構(gòu)建分形粗糙表面，研究在鋼-鋼微動接觸中表面粗糙度對局部摩擦和溫度分布的影響。研究發(fā)現(xiàn)接觸表面溫升幅值隨著粗糙度的增大而增大，材料導熱性能越好，接觸表面溫升幅值越小。

考慮到表面形貌對切向微動磨損有重要的影響，但是以上有限元分析都沒有綜合考慮表面初始粗糙度和磨損對切向微動的影響。因此本研究可以為研究工程機械中的微動磨損基礎(chǔ)力學特性和磨損行為提供參考價值。

本文考慮表面初始粗糙度對切向微動的影響，利用二維Weierstrass-Mandelbrot函數(shù)（簡稱W-M分形函數(shù)）建立分形粗糙表面，由UMESHMOTIOM用戶子程序引入Archard磨損法則建立二維切向微動數(shù)值模型，對微動磨損中的磨損程度進行分析研究。

1 粗糙表面生成

1983年Mandekbrot B提出分形理論［11］。經(jīng)實驗表明，分形理論可以表征許多零件粗糙表面微凸體的形狀與分布，并具有連續(xù)性、不可微性、尺度不變性和自相關(guān)性特征。二維表面輪廓由W-M函數(shù)給出［12］。

式（1）中，z（x）為隨機表面輪廓高度；x表示輪廓的位置坐標；G是反映表面輪廓振幅的分形尺度系數(shù)；D是決定高頻和低頻分量在表面輪廓中的相對比重的分形維數(shù)（1＜D＜2）；γ控制表面輪廓中頻率的密度，對于大多數(shù)表面輪廓，γ取1.5較為合適［12］；n是分形尺度指數(shù) ，n1和n2代表最低和最高截止頻率指數(shù)。令ω為輪廓空間頻率且ω=γn，那么ωL為起始頻率，由樣本長度L確定且ωL=1/L，ωU是頻率上限，由輪廓分辨率δ確定且ωU=δ/2。若一條被模擬的粗糙輪廓曲線的均方根為Rq，則特征尺度系數(shù)G可由方程（2）［2］確定，故可以根據(jù)D和G求出粗糙面的輪廓均方根Rq。

本文中，取G=2.85E-4μm，D=1.5，δ=0.5μm，L=［50，170］（μm），則n=［-11，0］，生成的粗糙表面輪廓曲線如圖1所示，輪廓均方根Rq為0.1302μm。

圖1 粗糙表面輪廓曲線

2 有限元模型

2.1 有限元模型建立

二維柱-平面切向微動有限元模型如圖2所示。考慮到分形函數(shù)生成的粗糙表面有很多微凸體，施加載荷后，可能會造成接觸失敗與網(wǎng)格嚴重畸變，故利用python中的scipy模塊對輪廓曲線進行插值優(yōu)化。為了保證計算精度，對接觸區(qū)網(wǎng)格進行細分，細分后的網(wǎng)格尺寸為0.25μm，單元類型選擇CPE4。半圓柱半徑為1200μm，平板尺寸為2400μm×1200μm。為了減少計算量，降低磨損計算時間，僅在平板接觸表面中心區(qū)域（長120μm）生成粗糙面。平板和半圓柱材料分別是AISI1034和AISI52100，材料參數(shù)如表1所示，其中硬度（H）為屈服強度的2.8倍［13］。

圖2 有限元模型

表1 材料參數(shù)

為使半圓柱體均勻下壓，故將半圓柱上表面與參考點建立剛體綁定約束，對參考點施加法向載荷F和切向周期位移載荷u。加載過程如圖3所示。其中0-1表示壓緊階段，在參考點上施加Y軸負向法向載荷F；1-2表示切向微動階段，參考點上施加法向載荷F和X方向周期性位移載荷u；2-3表示卸載階段。

圖3 加載曲線

為了接觸過程中同時考慮主面和從面的形狀變化，不會因為過度滲透而導致結(jié)果不準確，故半圓柱和平板接觸部分設(shè)置為主面-從面接觸，其中粗糙平板為從面，半圓柱光滑面為主面。離散方法采用面對面離散。法向接觸使用硬接觸，約束方法使用增廣Lagrange法，并允許接觸后可以分離；切向接觸使用各項同性的庫倫摩擦公式，摩擦系數(shù)為0.5。

本文磨損模型如式（3）［14］所示：

其中 τ為時間，Δh（x，τ）為磨損深度增量，p（x，τ）為接觸應(yīng)力，δ（x，τ）為相對滑動，ki為磨損系數(shù)，ΔN為循環(huán)跳躍系數(shù)。本文中為不類似的兩個金屬接觸，因此磨損系數(shù)ki為6.7e-5［15］。UMESHMOTION將式（3）所計算出的節(jié)點磨損量作為節(jié)點偏移量進行網(wǎng)格調(diào)整。通過重復此過程便可實現(xiàn)對磨損程度的有限元預測，具體流程如圖4所示。

2.2 有限元模型驗證

在Y-K分形粗糙表面接觸模型［13］中，存在一個臨界接觸面積a′c，如式（4）所示：

當最大微凸體截斷面積a′L大于a′c時：（a）Ds=2.5（D=1.5），彈性部分作用力Fe與塑性部分作用力Fp分別如式（5）與式（7）所示；（b）Ds≠2.5（D≠1.5），彈性部分作用力Fe與塑性部分作用力Fp如式（6）與式（7）所示。當最大微凸體截斷面積a′L小于等于a′c時，彈性部分作用力Fe為0，塑性部分作用力Fp如式（7）所示。

對總作用力F和接觸面積A進行無量綱化處理：復合彈性模量E′=E/（1-v2），無量綱接觸面積A′=A/A0，無量綱應(yīng)力P′=F/（A′×E′）。模型驗證中，試樣長度L=120μm，G=0.12μm，Ds=2.5（D=1.5），γ=1.5，變形體材料參數(shù)為E=210 GPa，ν=0.3，H=1834 MPa。

如圖5所示，將無量綱化載荷與接觸面積有限元結(jié)果和解析解對比，兩者具有較好一致性，故驗證了本模型的合理性。

圖4 微動磨損有限元分析流程圖

圖5 有限元結(jié)果與解析解對比

3 結(jié)果分析

3.1 粗糙面對磨損的影響

圖6可以看出，微動循環(huán)后的粗糙表面與光滑表面磨損區(qū)域和磨損形貌基本一致。兩者磨損區(qū)域外的表面由于受到相鄰微凸體的擠壓均發(fā)生塑性變形而突起。

圖6 100次循環(huán)后光滑表面與粗糙表面磨損形貌對比：F=60 mN，u=4μm，t=3 s

圖7 100次循環(huán)后光滑表面與粗糙表面磨損量對比：F=60 mN，u=4μm，t=3 s

在圖7中，光滑表面磨損量明顯比粗糙表面磨損量小，兩者均呈現(xiàn)出自中心區(qū)域逐漸減小的趨勢。初始表面為粗糙面，則粗糙峰會引起接觸應(yīng)力的不均勻分布，并產(chǎn)生強烈的應(yīng)力集中，導致磨損量的分布具有離散性。初始光滑面與粗糙面下的接觸應(yīng)力分布具體如圖8所示。粗糙表面微凸體的分布對磨損區(qū)域和磨損形貌也有一定影響。在圖6與圖7中，x=138μm附近左側(cè)的微凸體明顯較高，因此在法向載荷作用下，左側(cè)表面發(fā)生磨損，磨損后的表面低于磨損前的表面；右側(cè)表面沒有發(fā)生磨損，在相鄰微凸體的擠壓下，微動循環(huán)后的表面高于磨損前的表面。

圖8 光滑表面與粗糙表面接觸應(yīng)力分布對比：F=60 mN，u=4μm，t=1 s

3.2 微動工況對磨損的影響

圖9和圖10分別表示不同法向載荷下的磨損量和等效塑性應(yīng)變對比。結(jié)果表明，隨著法向載荷的增加，接觸中心區(qū)域磨損量和等效塑性應(yīng)變差別較小，接觸區(qū)域邊緣磨損量和等效塑性應(yīng)變差別較大。這是因為，粗糙峰在微動過程中會產(chǎn)生強烈的應(yīng)力集中，在一定量的微動循環(huán)后，在應(yīng)力集中明顯的接觸中心區(qū)域，其磨損程度較高，微凸體高度趨于一致，在一定范圍內(nèi)法向載荷的增加對接觸中心區(qū)域的應(yīng)力集中影響較小。而接觸邊緣趨于磨損程度較低，微凸體之間差異較大，在一定范圍內(nèi)法向載荷的增加對接觸中心區(qū)域的應(yīng)力集中影響較大。

圖11和圖12分別表示不同切向位移幅值下的磨損量和等效塑性應(yīng)變對比。由圖可見，位移幅值在一定范圍內(nèi)變化時，等效塑性應(yīng)變變化不大，而磨損量顯著增加。這是因為同一法向載荷下，接觸區(qū)域的應(yīng)力水平基本不變，而切向位移幅值的增加導致半圓柱和微凸體相對位移增加，使得磨損量顯著增加。

圖9 100次循環(huán)后不同法向載荷下磨損量對比：u=4μm，t=3 s

圖10 100次循環(huán)后不同法向載荷下等效塑性應(yīng)變對比：u=4μm，t=3 s

圖11 100次循環(huán)后不同位移幅值下磨損量對比：F=60 mN，t=3 s

圖12 100次循環(huán)后不同位移幅值下等效塑性應(yīng)變對比：F=60 mN，t=3 s

4 結(jié)論

建立分形粗糙表面數(shù)值模型，對二維切向微動下磨損程度進行分析，主要結(jié)論如下：

（1）光滑表面磨損量明顯比粗糙表面磨損量小，兩者均呈現(xiàn)出自中心區(qū)域逐漸減小的趨勢，且粗糙表面磨損量分布有很大的的離散性。這是因為，初始表面為粗糙面時，粗糙峰會引起接觸應(yīng)力和摩擦切應(yīng)力的不均勻分布，并產(chǎn)生強烈的應(yīng)力集中，導致磨損量顯著變大并具有離散性。

（2）初始表面為粗糙面時，由于微凸體的高度與分布具有不均勻性與隨機性，接觸時產(chǎn)生強烈的應(yīng)力集中，較大的塑性變形使得粗糙表面和光滑平面磨損形貌差異不大。

（3）在一定載荷范圍內(nèi)，法向載荷的增加對接觸中心區(qū)域磨損程度影響較小，而對接觸邊緣區(qū)域磨損程度影響較大。這是因為隨著微動磨損的進行，磨損嚴重的接觸中心區(qū)域微凸體高度趨于一致，一定范圍的載荷增加對其應(yīng)力水平影響較小。

（4）切向位移載荷的增加使得磨損程度明顯變大。這是因為切向位移幅值的增加導致半圓柱和粗糙表面微凸體相對滑動變大，導致磨損量增加。