考慮摩擦生熱的扭動微動熱力耦合數值分析

譚升旭，劉 娟，沈火明，陳揚帆

(西南交通大學 力學與工程學院， 成都 610031)

0 引言

扭動微動是指在交變載荷[1]作用下，接觸面之間發生小幅度往復扭轉而引起的相對運動，它發生在人體髖關節、膝關節等部位，也發生在許多工程部件中，如核電設備中的部件[2-4]、高速列車轉向架的中心板[5-7]等。微動損傷是諸多復雜因素共同作用造成，許多問題尚未深入研究以了解損傷的本質和規律。扭動微動已有一些研究，如蔡振兵等[8-10]系統地研究了金屬材料和復合材料的扭動微動磨損，使用潤滑、涂層和氣氛控制等條件減輕微動磨損。半解析方法[11]和ABAQUS數值分析[12]已應用于扭動微動的數值模擬，但關于微動的研究仍集中在切向模式上，其他模式的研究很少，尤其是扭動微動。目前，扭動微動磨損的力學行為和磨損機制仍有待揭示。

研究表明，影響材料微動磨損的參數有50多個[13]。法向載荷、摩擦因數、相對位移幅值等對微動過程的影響已得到廣泛關注[12,14-15]，溫度已被認定為影響微動磨損行為和損傷機制的重要因素[16-17],微動磨損過程中的溫度會受環境溫度和微動接觸面間產生的摩擦耗散能的影響。Attia[18-19]的研究表明，摩擦引起的溫升會顯著影響微動行為。此外，Wen等[20]的研究表明，由摩擦耗散引起的溫升，會對微動過程中氧化物碎片的形成、三體層的發展變化，以及表面的材料性質造成影響。微動磨損速度的變化是高溫促進的氧化物形成而造成的[21]。試驗中很難測量接觸面之間的溫度，所以需要預測技術來促進微動過程中溫度作用機理的研究。許多研究人員(例如Jaeger[22],Archard[23],Tian[24]等)通過引入Péclet數的分析模型，證明了表面溫度與材料的熱性能、接觸面的相對移動速度有關。文獻[25-27]建立了數值模型，分析了溫升和材料變形之間的關系。

綜上所述，溫升在微動機理研究中具有重要意義。因此，本文基于數值方法和摩擦生熱機理，建立預測扭動微動熱力耦合行為數值分析模型，研究表面Von Mises 應力分布、表面溫升和次表面溫升，以及表面塑性變形與溫升的關系，進而探索溫升對扭動微動磨損的影響機制，有助于進一步揭示扭動微動磨損損傷機理。

1 模型和方法

1.1 數值模型

圖1 數值模型結構示意圖

采用主-從接觸算法來滿足非穿透接觸約束，在迭代過程中，定義最大的穿透范圍，且只允許主表面上的節點侵入從表面，因此不會出現導致結果不準確的過度穿透。

法向載荷直接施加會導致局部節點剛體位移，影響模型收斂，因此定義一個參考點，將其和上接觸體的頂面綁定約束。在參考點上施加恒定的法向載荷和水平的正弦扭轉循環角位移，使半球面與平面發生扭轉，實現微動循環的加載[28]，平面的下底部則進行完全約束。載荷加載歷程如圖2所示，設置較短時間的加載和卸載分析步使接觸狀態平穩過渡，增強模型的收斂性。

圖2 加載歷程曲線

此外，基于所有克服摩擦所做的功都以熱源形式消散的假設，熱量將從接觸區域流入接觸體，則可以表示為式(1)的熱功率密度，該項與局部接觸壓力、位置和時間的函數成正比。

Q(t)=Ff|v(t)|

(1)

式中：Ff為摩擦力；v(t)為接觸節點的滑動速度。描述運動的位置函數s(t)可以表示為：

s(t)=rφ(t)=rθsin(ωt)

(2)

式中：ω=2πf，f為微動頻率；θ為角位移幅度；r為距接觸中心的距離。速度如下：

(3)

摩擦力可以表示為：

Ff=Pμ

(4)

故總熱源功率可定義為：

Q(t)=Pμrθω|cos(ωt)|

(5)

1.2 模型驗證

為驗證數值模型的準確性，與Johnson的解析解[28]進行了比較。在接觸載荷的作用下，接觸區壓力分布表示為：

(6)

式中：

(7)

(8)

式中：a0為接觸面積的半徑；p0為最大接觸壓力；E1和E2為兩材料的彈性模量；μ1和μ2為兩材料的泊松比；R1和R2為兩試樣接觸表面的曲率半徑；x為距離接觸中心距離。對于切向接觸問題，基于切向力與壓力沒有相互作用的假設，剪切力可以表示為：

(9)

式中：H( )為Heaviside的階躍函數；μf為摩擦因數，其值為0.6；c為黏著半徑。

摩擦副GCr15/2024合金的材料屬參數如表1所示，所列分別是彈性模量、抗拉強度、屈服強度、泊松比、密度，以及熱膨脹系數和比熱容。

表1 摩擦副的力學和物理性能參數

圖3(a)和(b)分別表示壓力和剪切應力的解析解與數值結果。可以看出，數值結果和解析結果之間有很好的一致性。

圖3 解析解與數值結果曲線

2 結果和討論

2.1 Von Mises 應力

普遍認為，第一屈服準則對材料磨損機制具有重要指導作用。Johnson[28]已經證明了剪切應力、Von Mises 應力分布以及最大應力的確切位置。在扭動微動的半解析法(SAM)研究中，Liu等[12]已經證明，在扭動微動磨損機制中，表面的首次屈服主要表現為黏著磨損，接觸面材料損傷加劇表現為明顯的塑性變形。

為揭示表面塑性變形和磨損，在法向載荷P=100 N作用下，分別取角位移幅值θ為0.15°，1°，8°，呈現出扭動微動的3種不同滑移狀態，如圖4所示。Von Mises 應力分布的數值結果如圖5所示：在部分滑動狀態下整個接觸面的應力較小，沒有發生屈服[圖5(a)]，而對于混合滑移狀態和完全滑移狀態，高應力幾乎遍布整個表面，且發生了明顯的屈服[圖5(b)和5(c)]；在混合滑動狀態，高應力主要表現在發生磨損的環形滑動區域[30]；在完全滑移狀態，整個接觸區域的顯著屈服意味著嚴重的表面磨損，這與Cai等[31]的實驗結果一致。在混合滑移狀態的微動循環期間，Von Mises 應力值并不總是大于屈服強度，這表明表面磨損并不是主要的損壞機制；事實上，磨損機制主要是與疲勞相關的微裂紋，因為屈服主要發生在次表面[12]。

圖4 3種不同滑移狀態微動曲線

圖5 微動周期內 Von Mises 應力的分布曲線

2.2 溫升分布分析

角位移幅值是影響摩擦生熱的主要因素之一，圖6表示混合滑移狀態和完全滑移狀態下表面和次表面溫升的預測結果，表面溫度的距離接觸中心x=0.2 mm的位置，次表面是距離表面監測點h=0.2 mm深度的位置。顯然，溫升受位移幅值的影響較大，表面和次表面溫升隨著扭轉角位移幅值的增加而顯著增加，且完全滑移狀態的溫升大于混合滑移狀態的溫升；次表面溫升較接觸表面溫升顯著降低；溫升在很短的微動時間內顯著增加。圖7表示高頻率(f=200 Hz)和低頻率(θ=100 Hz)下10 000 s的溫升隨時間的變化。結果表明，在這2種情況下，前幾個周期的溫升都急劇增加，而上升速度隨時間的增加迅速下降[如圖7(a)]。圖7(b)、(c)和(d)分別顯示了3個不同階段的表面溫升，可以看出3個階段的溫升增加速度差異。在高頻率f=200 Hz的扭動作用下，最大增長速度出現在A階段，B階段的增速大于C階段。經過數千秒后，溫升幾乎達到一個穩定值。

圖6 表面(x=0.2 mm)和次表面(h=0.2 mm)溫升隨時間變化的預測曲線(P=100 N, f=100 Hz)

圖7 較長時間的溫升曲線 (x=0.2 mm,P=100 N,θ=5°)

2.3 塑性變形

2.1節的結果表明，表面屈服引起了大范圍的表面塑性變形，表面磨損是完全滑移狀態的主要損傷機制。且在2.2節的分析中表明，微動頻率高時表面溫升大，特別是在完全滑移狀態下，表面溫升明顯。相關研究認為，接觸中的溫度會影響磨屑形成穩定層和從接觸區排出的方式[20]。

為了研究溫升對扭動微動表面磨損的影響，本文預測了表面塑性變形與溫升之間的關系，在P=100 N,θ=5°的工況下，探究溫度變化對材料塑性變形的影響。如圖8所示，溫升影響了表面塑性變形行為，在循環一周時，塑性變形變化尚微小，而隨著微動周期的增加，扭動微動中的塑性變形隨著溫度的升高顯著增大。此外，圖9顯示了在距離接觸中心0.2 mm處監測到的扭動微動中的塑性變形隨隨溫度升高的變化。隨著溫度的升高，表面塑性變形顯著增加，但增加的速度緩慢下降。

圖8 不同微動循環下塑性應變分布曲線(P=100 N，θ=5°)

圖9 表面塑性變形隨溫升的變化曲線(x=0.2 mm,P=100 N,θ =5°)

3 結論

基于摩擦耗散的摩擦生熱機理，建立了預測扭動微動的熱力耦合數值分析模型。在模型驗證基礎上，研究了表面Von Mises 應力分布、表面溫升和次表面溫升，以及表面塑性變形與溫升間的關系。結果表明：Von Mises 應力分析結果顯示表面磨損是完全滑移狀態的主要損傷機制，而與疲勞相關的微裂紋是混合滑移狀態的主要磨損機制，該結論與以往研究中的實驗現象一致；表面和次表面溫升均對角位移幅值和微動頻率敏感，但表面溫升比次表面溫升更顯著，且溫升變化速度先增大后急劇減小，在數千秒后便達到穩定狀態；表面溫升顯著影響表面塑性變形，進而影響表面磨損，因此摩擦生熱引起的表面溫升明顯影響扭動微動中的表面磨損行為。