脂介質(zhì)條件下輪軌摩擦因數(shù)數(shù)值模擬分析

劉思磊 許玉德 劉一鳴

1.同濟大學(xué)道路與交通工程教育部重點實驗室，上海201804；2.深圳市綜合交通設(shè)計研究院有限公司，廣東深圳518003

近年來，重載鐵路向著大軸重、長編組、高密度方向發(fā)展，輪軌狀態(tài)惡化加劇。運用輪軌摩擦控制劑來控制輪軌間摩擦是延長輪軌使用壽命、降低輪軌間磨耗與損傷的有效手段之一。

早在2004年，就有日本、加拿大、西班牙等國家開始對線路進行潤滑試驗，并獲得了較大的經(jīng)濟效益［1-3］。對于輪軌摩擦因數(shù)的變化規(guī)律，日本的伴巧、松本耕輔、深貝晉也等［4-6］通過在二圓筒試驗機上的試驗，發(fā)現(xiàn)摩擦因數(shù)的改變量與轉(zhuǎn)動周數(shù)呈指數(shù)關(guān)系。近年來，國內(nèi)學(xué)者通過試驗進行了相關(guān)研究。劉騰飛［7］利用滾動摩擦磨損試驗機研究了不同介質(zhì)對輪軌黏著的影響。李偉［8］發(fā)現(xiàn)，進行曲線摩擦控制可改善軌距面的磨耗，增加車輛的曲線通過能力，大大降低輪軌橫向力。張信等［9］設(shè)計了試驗，分析丙烯酸樹脂乳液的加入量對涂覆層的影響。宋靖東等［10］分析了固體潤滑劑在輪軌間的摩擦性能和有效工作時間情況。王珊珊［11］在點接觸試驗機上通過球盤接觸試驗，研究了滑滾比對潤滑脂纖維團在接觸區(qū)分布的影響。

既有文獻中關(guān)于摩擦控制條件下摩擦因數(shù)的影響機理研究較少，且缺少對摩擦因數(shù)改善效果的定量研究。本文在調(diào)研國內(nèi)外鐵路摩擦控制及摩擦學(xué)相關(guān)研究的基礎(chǔ)上，進行脂介質(zhì)條件下的輪軌摩擦因數(shù)分析與計算。

1 計算模型

輪軌摩擦控制屬于彈性流體動壓潤滑問題，即在流體動壓潤滑的基礎(chǔ)上考慮了彈性變形，而彈性變形因素的比重與膜厚比有關(guān)。在求解時，一般將Hertz接觸壓力分布作為初始壓力分布。本文引用彈流理論基本方程，引入脂介質(zhì)的流變本構(gòu)模型并考慮接觸中熱效應(yīng)等問題。

1.1 脂介質(zhì)下彈性流體動壓潤滑理論基本假設(shè)

對流體的摩擦狀態(tài)進行研究時，所用的基本方程是Reynolds方程和膜厚方程［12］。其中Reynolds方程是從黏性流體力學(xué)的基本方程出發(fā)，作了一定的假設(shè)而導(dǎo)出的微分方程。這些基本假設(shè)是：①液膜厚度很小，僅有數(shù)千納米，因此沿液膜厚度方向，液膜壓力為常數(shù)；②輪軌摩擦表面的曲率半徑遠(yuǎn)大于液膜厚度，因此接觸面可展為平面；③摩擦控制劑視為脂類，剪應(yīng)力與剪應(yīng)變率呈非線性關(guān)系；④液膜中流體作層流運動；⑤流體的慣性力忽略不計。

1.2 基本方程

1.2.1 Reynolds方程

Reynolds方程是基于動量方程和連續(xù)性方程推導(dǎo)出來的流體動力潤滑理論，于1886年由Reynolds提出，是彈流問題中最基本的方程［11］。加入Ostwald模型后，可推導(dǎo)出脂介質(zhì)的一維Reynolds方程，即

其中，Ostwald模型為

式中：n為流體指數(shù)；ρ為流體密度；h為液膜厚度；Φ為塑性黏度；p為流體壓力；us為上下接觸表面的平均速度；τ為剪應(yīng)力；γ˙為剪應(yīng)變率。

入口區(qū)的邊界條件為p(x0)=0；出口區(qū)的邊界條件為其中，x0表示入口坐標(biāo)，xe表示出口坐標(biāo)。

1.2.2 膜厚方程

膜厚是Reynolds方程的一個重要子項［12］。兩個半徑分別為R1和R2的圓柱體作線接觸時，若兩圓柱未變形時的間距為h0，考慮到表面的彈性變形，實際的線接觸膜厚公式為

式中：R可由得到；v(x)為彈性變形。

如圖1所示，v(x)由彈性力學(xué)中表面上的法向分布荷載在半無限大平面上引起的表面變形公式得到，即

式中：p(s)為荷載分布函數(shù)；s是x軸上的附加坐標(biāo)，表示任意線荷載p(s)ds與坐標(biāo)原點的距離；實際計算中，c為常數(shù)，可以同式（3）中的h0合并。

1.2.3 荷載平衡方程

求解彈流問題是在給定荷載條件下進行的，所以求出的壓力必須滿足荷載平衡條件。當(dāng)輪軌接觸處于部分膜彈流狀態(tài)時，荷載由固體粗糙峰荷載We和液膜荷載Wa共同承擔(dān)，荷載方程［13］為

式中：pA為固體A上粗糙峰的壓力。

根據(jù)Greenwood-Tripp接觸理論，pA的計算方程為

當(dāng)h∕σ在0.000 1～0.010 0之間變化時，kA的取值為0.000 3～0.003 0。由于F5∕2(h∕σ)中積分的計算較復(fù)雜，因此采用擬合公式進行計算，即

1.2.4 摩擦因數(shù)計算方程

本文彈流問題的荷載由摩擦控制劑的動壓力We和液膜荷載Wa共同承擔(dān)，當(dāng)求出We和微凸體峰接觸所承擔(dān)的荷載Wa后，總的摩擦因數(shù)f的表達(dá)式為［14］

式中：μe為脂類摩擦控制劑全膜彈流狀態(tài)下的摩擦因數(shù)，軌頂使用的脂類摩擦控制劑μe約為0.30；μa為微凸體峰接觸摩擦因數(shù)，一般由試驗得到，本文取μa=0.10。

2 數(shù)值計算

為方便數(shù)值求解，改善方程的適用性，須采用多重網(wǎng)格法求解彈性變形，然后用壓力迭代修正。

2.1 多重網(wǎng)格法求解彈性變形

彈流問題的關(guān)鍵是將彈性變形方程疊加到膜厚方程上，彈性變形方程的求解要解偏微分方程，求解過程中選擇合適的網(wǎng)格常常是困難的，使用稀疏網(wǎng)格得到的解精度太低，而且對于非線性問題常得不到收斂的解。20世紀(jì)60年代初提出了一種求解代數(shù)方程組的高效率迭代方法——多重網(wǎng)格法［15］，可以有效地解決上述困難。

利用多重網(wǎng)格法計算彈性變形主要包含以下步驟：

1）壓力下傳?？紤]壓力變化，將細(xì)網(wǎng)格上的壓力插值后傳遞到粗網(wǎng)格上。

2）積分系數(shù)下傳。采用相應(yīng)計算公式將細(xì)網(wǎng)格的積分系數(shù)下傳給粗網(wǎng)格。

3）粗網(wǎng)格積分。

4）積分值插值回代。部分未計算積分值的細(xì)網(wǎng)格節(jié)點上，如果與粗網(wǎng)格節(jié)點重合則作映射運算，如果不重合則作插值運算。

5）細(xì)網(wǎng)格修正。

2.2 壓力迭代修正

在初始壓力分布下（本文為Hertz接觸壓力）應(yīng)用多重網(wǎng)格計算方法得到彈性變形后，可以計算膜厚和黏度值，并帶入Reynolds方程求解得到新的壓力分布，對前一次的壓力分布進行迭代修正，然后帶入能量方程求解溫度。利用新的溫度修正黏度值，再迭代求解壓力。重復(fù)此過程，直至達(dá)到精度要求，迭代結(jié)束。最終得到所求的壓力分布和含彈性變形的膜厚分布。

壓力迭代修正同時使用Gauss-Seidel迭代和Jacobi雙極子迭代方法［13］。在低壓力區(qū)，選擇Gauss-Seidel迭代法；在高壓力區(qū)，為避免發(fā)散，選擇Jacobi雙極子迭代法。

1）Gauss-Seidel迭代法

修正時須在本點加上修正量。表達(dá)式為

2）Jacobi雙極子迭代法

修正時須在本點加上修正量和在前點減去修正量，表達(dá)式為

式中：c2為松弛因子。

3 計算結(jié)果分析

按照上述計算流程求解基本方程，依據(jù)計算結(jié)果討論脂類摩擦控制劑的黏度、列車軸重及列車速度對摩擦因數(shù)的影響。

3.1 脂介質(zhì)黏度對摩擦因數(shù)的影響

基于經(jīng)典彈流理論和脂介質(zhì)本構(gòu)方程，計算黏度為10～500倍水黏度時的摩擦因數(shù)，模擬中等黏度的脂類摩擦控制劑。

表面粗糙度為6 μm時，選取脂介質(zhì)黏度分別為0.01、0.10、0.30、0.50 Pa·s（依次對應(yīng)10、100、300、500倍水黏度）計算得到輪軌接觸區(qū)壓力分布情況及膜厚變化情況見圖2。其中，坐標(biāo)進行了無量綱化處理；橫坐標(biāo)-1.0～1.0為接觸區(qū)。

由圖2可知：

圖2 不同脂介質(zhì)黏度下接觸壓力和液膜厚度變化曲線

1）從接觸區(qū)寬度看，隨著脂介質(zhì)黏度增加，液膜和粗糙峰的壓力接觸區(qū)入口逐漸向前移動，導(dǎo)致接觸區(qū)寬度增加。當(dāng)脂介質(zhì)黏度為10倍水黏度（0.01 Pa·s）時，液膜壓力和粗糙峰壓力的接觸區(qū)寬度與Hertz接觸區(qū)寬度基本一致；當(dāng)脂介質(zhì)黏度增加至500倍水黏度（0.50 Pa·s）時，液膜壓力的接觸區(qū)寬度約為1.5倍Hertz接觸區(qū)寬度，粗糙峰壓力的接觸區(qū)寬度約為1.2倍Hertz接觸區(qū)寬度。

2）從液膜接觸壓力看，隨著黏度增加，在出口區(qū)附近，液膜接觸壓力會出現(xiàn)逐漸明顯的二次尖峰。

3）從液膜厚度看，在膜厚出口處均有頸縮現(xiàn)象，隨著黏度增加液膜厚度逐漸增加，使得液膜承擔(dān)的荷載比重上升，粗糙峰承擔(dān)的荷載比重下降。

不同輪軌表面粗糙度情況下，黏度變化對摩擦因數(shù)的影響見圖3。

圖3 不同表面粗糙度下黏度變化對摩擦因數(shù)的影響

由圖3可知，隨著脂介質(zhì)黏度增加，摩擦因數(shù)呈增加趨勢，且在輪軌表面比較光滑時，摩擦因數(shù)更早地逼近臨界值0.3。產(chǎn)生這一現(xiàn)象的原因是，當(dāng)黏度增加時，脂介質(zhì)所形成的液膜將承擔(dān)更多的荷載，根據(jù)式（7），接觸區(qū)的摩擦因數(shù)由液膜承擔(dān)的荷載比重決定，而脂介質(zhì)的摩擦因數(shù)大于粗糙峰的摩擦因數(shù)，因此引起摩擦因數(shù)增加。而摩擦因數(shù)的上限是全膜彈流狀態(tài)下的0.3，因此摩擦因數(shù)隨黏度的增加趨勢將趨于平緩。

3.2 軸重對摩擦因數(shù)的影響

表面粗糙度為6 μm時，選取軸重分別為10、14、18、22 t，計算得到輪軌接觸區(qū)壓力分布情況及膜厚變化情況見圖4。其中，坐標(biāo)進行了無量綱化處理；橫坐標(biāo)-1.0～1.0為接觸區(qū)。

圖4 不同軸重下接觸壓力和液膜厚度變化曲線

由圖4可知：

1）從接觸區(qū)寬度看，隨著軸重增加，液膜和粗糙峰的壓力接觸區(qū)入口逐漸向后移動，導(dǎo)致接觸區(qū)寬度縮小。當(dāng)軸重為10 t時，液膜壓力的接觸區(qū)寬度約為1.9倍Hertz接觸區(qū)寬度，粗糙峰壓力的接觸區(qū)寬度約為1.4倍Hertz接觸區(qū)寬度；當(dāng)軸重為22 t時，液膜壓力的接觸區(qū)寬度與Hertz接觸區(qū)寬度基本一致，粗糙峰壓力的接觸區(qū)寬度約為1.25倍Hertz接觸區(qū)寬度。

2）從液膜接觸壓力看，隨著軸重增加，液膜接觸壓力的二次尖峰逐漸向出口區(qū)移動，當(dāng)軸重大于22 t以后，尖峰不再明顯，總應(yīng)力分布已經(jīng)接近Hertz應(yīng)力分布，說明此時摩擦效果接近失效。

3）從液膜厚度看，在膜厚出口處均有頸縮現(xiàn)象，隨著軸重增加，液膜厚度逐漸減小，使得液膜承擔(dān)的荷載比重下降，粗糙峰承擔(dān)的荷載比重上升。

不同輪軌表面粗糙度情況下，軸重變化對摩擦因數(shù)的影響見圖5。可知，隨著列車軸重增加，摩擦因數(shù)呈下降趨勢。這是因為隨著軸重增加，輪軌間脂介質(zhì)形成的液膜厚度將被壓縮，使得粗糙峰承擔(dān)的荷載比重增加，而粗糙峰的摩擦因數(shù)小于脂介質(zhì)的摩擦因數(shù)，因此引起摩擦因數(shù)下降。

圖5 不同表面粗糙度下軸重變化對摩擦因數(shù)的影響

3.3 速度對摩擦因數(shù)的影響

表面粗糙度為6 μm時，選取列車速度分別為60、70、80、90 km∕h，計算得到輪軌接觸區(qū)壓力分布情況及膜厚變化情況見圖6。其中，坐標(biāo)進行了無量綱化處理；橫坐標(biāo)-1.0～1.0為接觸區(qū)。

圖6 不同速度下接觸壓力和液膜厚度變化曲線

由圖6可知：

1）從接觸區(qū)寬度看，列車速度由60 km∕h增加至90 km∕h，液膜壓力的接觸區(qū)寬度由1倍Hertz接觸區(qū)寬度增加至約1.08倍Hertz接觸區(qū)寬度，粗糙峰壓力的接觸區(qū)寬度基本保持在1.25倍Hertz接觸區(qū)寬度，不隨速度變化而改變。

2）從液膜接觸壓力看，隨著列車速度增加，液膜接觸壓力的二次尖峰位置和峰值大小變化不大。

3）從液膜厚度看，在膜厚出口處均有頸縮現(xiàn)象，隨著速度增加，液膜厚度逐漸增加，使得液膜承擔(dān)的荷載比重上升，粗糙峰承擔(dān)的荷載比重下降。

不同輪軌表面粗糙度情況下，列車速度變化對摩擦因數(shù)的影響見圖7?？芍S著列車速度增加，摩擦因數(shù)呈增加趨勢，但增加的幅度很小，當(dāng)輪軌表面粗糙度處于一般水平10 μm時，速度由50 km∕h增加至100 km∕h，摩擦因數(shù)僅由0.260增加至0.263。

圖7 不同表面粗糙度下速度變化對摩擦因數(shù)的影響

4 試驗驗證

為驗證模型的可靠性以及數(shù)值計算結(jié)果的正確性，通過脂介質(zhì)條件下的滾動接觸試驗，驗證軸重因素對摩擦因數(shù)的影響，總結(jié)不同軸重下摩擦因數(shù)的變化規(guī)律。為了能夠穩(wěn)定地將脂介質(zhì)涂覆在接觸面上，試件的轉(zhuǎn)速設(shè)定為10 r∕min。按照上述模擬準(zhǔn)則，試驗所采用的試驗力及所對應(yīng)的模擬軸重見表1，其他主要參數(shù)見表2。試驗現(xiàn)場如圖8所示。

表1 模擬軸重與試驗力

表2 試驗其他參數(shù)

圖8 試驗現(xiàn)場

對軸重分別為10.0、12.5、15.0、17.5、20.0、22.5 t時試驗得出的摩擦因數(shù)進行分析。其中軸重15 t時的摩擦因數(shù)變化曲線見圖9?？芍孩偻扛仓橘|(zhì)后，輪軌摩擦因數(shù)迅速爬升，之后在0.2附近波動。涂覆階段摩擦因數(shù)有時會出現(xiàn)較大波動的原因是，滾動接觸的過程中須不斷手工涂刷脂介質(zhì)，涂刷時涂量不穩(wěn)定會引起摩擦因數(shù)波動。②當(dāng)停止涂覆后，經(jīng)過約100 s，摩擦因數(shù)降低至干摩擦下的0.11附近。

圖9 軸重15 t時摩擦因數(shù)變化曲線

將試驗參數(shù)帶入到模型中驗證。由于脂類材料未標(biāo)明該參數(shù)，須通過試驗測定其取值。在滾動接觸試驗中，最小試驗力為1.2 kN，所產(chǎn)生的接觸應(yīng)力等效于現(xiàn)場的8.6 t軸重，而所用試樣的粗糙度僅有0.6 μm。在1.2 kN試驗力、10 r∕min的工況下進行了5組試驗，5組試驗的摩擦因數(shù)最大值的均值為0.22，故取脂介質(zhì)全膜彈流狀態(tài)時的摩擦因數(shù)為0.22。取涂覆穩(wěn)定階段的摩擦因數(shù)均值作為該軸重下摩擦因數(shù)的試驗值。得到試驗結(jié)果和數(shù)值計算結(jié)果對比，見表3。

表3 試驗結(jié)果與數(shù)值計算結(jié)果對比

由表3可知，試驗結(jié)果的最大相對誤差在4%以內(nèi)，平均相對誤差為2.30%，說明理論計算結(jié)果與試驗結(jié)果較吻合。

5 結(jié)論

1）相比于干摩擦下的Hertz分布，在脂介質(zhì)條件下，脂介質(zhì)承擔(dān)的壓力在接觸區(qū)出口附近會出現(xiàn)一個二次尖峰，粗糙峰承擔(dān)的壓力在接觸區(qū)內(nèi)變化不大。液膜厚度在出口處有頸縮現(xiàn)象。

2）在幾種影響因素中，脂介質(zhì)的黏度、輪軌滾動速度上升都會引起輪軌間摩擦因數(shù)上升，而軸重上升會引起摩擦因數(shù)下降。

3）脂介質(zhì)條件下，若處于相同的黏度、軸重、速度情況時，輪軌表面越光滑輪軌間摩擦因數(shù)越大。原因是脂介質(zhì)的摩擦因數(shù)高于干摩擦或水介質(zhì)，當(dāng)輪軌表面越光滑時，脂介質(zhì)形成的液膜厚度增加，將承擔(dān)更多的荷載，從而使得接觸區(qū)的摩擦因數(shù)的數(shù)值向脂介質(zhì)靠近，即發(fā)揮更好的摩擦控制效果。

4）涂覆脂介質(zhì)后，隨著試件上施加的試驗力增加，試驗測得的摩擦因數(shù)逐漸降低，模擬軸重由10.0 t增加至22.5 t時，摩擦因數(shù)大約由0.199降低至0.184。將實測數(shù)據(jù)與數(shù)值結(jié)果進行對比，最大誤差在4%以內(nèi)，數(shù)值計算結(jié)果與試驗結(jié)果較吻合。