干式摩擦元件的熱彈不穩定性分析*

郭 超 張 昱 楊 凱

(1.中國航空發動機研究院 北京 101304；2.北京航空航天大學機械工程及自動化學院 北京 100191；3.北京衛星制造廠 北京 100094)

摩擦元件的熱彈不穩定性是指由于其制造、安裝、結構、控制等方面的不足，制動器內部的物理場存在一定的初始非均勻性，當制動速度超過臨界值時，系統將出現熱彈性不穩定(Thermoelastic Instability)狀態，簡稱TEI[1]。制動器出現熱彈性不穩定(TEI)現象，加速了材料的磨損及熱疲勞破壞，另外還引入了“熱顫振”(Hot Judder)導致系統產生摩擦振動與噪聲的可能性加大，因此，TEI是影響制動器穩定性能的關鍵因素之一，必須深入研究。

PEI等[2]分析了摩擦系統中摩擦熱導致的熱機械失穩機制，提出一種熱邊界控制方法，并研究了轉盤在氣流下的氣動彈性不穩定性問題。ABDULLAH等[3]用Perturbation方法研究了2個摩擦面接觸的穩定性問題。對于熱彈性失穩的研究，需要分析臨界轉速，超過臨界轉速是熱彈性失穩的主要原因。TANG等[4]采用有限元法建立了離合器熱彈性不穩定性的瞬態模型。馬彪、于亮等人[5- 6]建立離合器接合過程的二維理論模型，分別求解了采用不同摩擦材料的離合器發生TEI時的臨界速度等。SHPENEV[7]通過有限差分方法建模分析各向異性圓盤樣品在非定常摩擦下引發熱彈性失穩的問題。HAN和GU[8]采用數值方法研究TEI問題，將求解臨界速度問題轉化為求解系統矩陣特征值問題。YU等[9]首次采用擾動法研究了系統的TEI問題。WALIA等[10]研究發現系統相對滑動速度較高時，擾動會向表層集中。擾動的來源往往是多方面的，與系統結構、溫度場、磨損量、壓力分布以及應力狀態有關。

拋開系統狀態說系統穩定性問題是不全面的，WANG等[11-13]從定性和定量2個方面對制動器工作過程中接觸摩擦生熱、瞬態溫度場、對偶件的變形和應力及接觸壓力分布等制動器熱彈性耦合分析的基本問題，進行了較為廣泛的研究。

本文作者在上述研究的基礎上，進一步探討干式摩擦元件熱彈性失穩的原因，將理論分析和有限元仿真相結合，對干式摩擦元件進行熱彈不穩定分析。

1 干式摩擦元件熱彈不穩定分析模型

當兩個相互接觸的物體產生相對滑動時，由于表面不規則粗糙，真實接觸面積是由名義接觸區域內許多微凸體接觸產生的。因此界面上壓力分布本質上是不均勻的，這種微觀的不均勻是導致系統不穩定性的根本原因。相對摩擦所產生的熱將引起接觸體的熱彈性變形，進而改變接觸面之間的壓力分布。反之，接觸面間壓力分布不均勻，又使局部產熱更加不均勻，當系統能量輸入超過其能量耗散時，該過程將不可逆轉地惡化下去，使得微觀尺度下的不均勻性會隨著溫度累計形成宏觀的局部高溫點。因此，對于給定的摩擦因數，將存在著某個滑動速度，當高于此數值時，系統將變得不穩定，即產生TEI現象[14]。

所以對摩擦元件熱彈性不穩定性問題的研究，主要集中在求解摩擦系統進入熱彈性不穩定狀態的臨界速度，分析摩擦元件在熱彈性不穩定狀態時摩擦表面的溫度場分布，熱點分布規律，系統參數對熱彈性不穩定性的影響等方面。在理論方面對于干式摩擦系統的TEI問題研究一般采用擾動法，在Barber和Burton模型的基礎上，改進模型參數，利用伽遼金法[15]，分別建立摩擦副系統二維有限元模型，通過求解TEI矩陣的特征值得到摩擦副在不同滑摩速度下的擾動增長系數，并確定TEI臨界速度。

首先建立摩擦副熱二維TEI有限元模型(如圖1所示)，通過有限元研究結果與解析法研究結果的對比，證明有限元模型計算結果的準確性。材料1、2分別表示摩擦材料與對偶材料，面積分別為A1、A2，邊界分別為C1、C2。2種材料接觸邊界為C3，材料1以速度v、沿X軸正方向運動。

圖1 2D有限元TEI模型示意

首先假設摩擦副內的溫度場、應力場、熱流密度、熱彈性變形等物理場的擾動均可以表示成余弦形式，擾動的幅值隨時間呈指數變化，在笛卡爾坐標系下可以表示為

R{F(y)exp(bt+iλx)}

(1)

同樣的溫度擾動為

R{T(y)exp(bt+iλx)}

(2)

式中：R表示公式取實部；F(y)表示應力；T(y)表示溫度；b表示擾動增長系數；λ表示擾動空間頻率，i=-1的開方。

考慮對流，那么熱傳導方程存在對流項，表示成：

(3)

式中：K為材料的導熱系數；ρ為密度；cp表示比熱容。

將擾動代入熱傳導方程可得：

(4)

然后對表面進行積分，表面考慮權函數w(y)可得：

(5)

利用分步積分法，

(6)

將其化簡得到：

(7)

熱流密度對溫度場的擾動可以表示成：

(8)

將其代入熱傳導方程中，得：

(9)

將其離散化，引入有限元形函數，

(10)

式中：Nj表示材料的節點數；Tm表示第m個節點上的溫度擾動值；Wm(y)表示形函數。

在伽遼金方法中，權函數與形函數相等，所以可以得到：

(11)

式中：Wm和Wn分別表示形函數和權函數；WmWn表示一個Nj×Nj的矩陣，表示成矩陣為WmWn=WWT。至此文中就建立了基于擾動關系的二維有限元熱傳導模型。

形函數可以表示為

Wm=[W1，…，WNj]T

(12)

以上方程可以改為矩陣形式：

(K+G+bL)T+Q=0

(13)

式中：矩陣K為N維傳導矩陣；G為N維系數矩陣；L為N維質量矩陣；T表示溫度向量；Q為熱流向量。

各矩陣和向量的表達式為

(14)

G=?Aj(Kjλ2+ρjcpjiλvj)WWTdAj

(15)

L=?Aj(ρjcpj)WWTdAj

(16)

(17)

綜上，建立干式制動器摩擦副系統熱彈不穩定分析模型。模型的求解，關鍵是確定過渡矩陣H，它是一個Nc×N矩陣。利用有限元分析軟件ABAQUS，計算各節點溫度場和接觸壓力，從而推導過渡矩陣。

2 干式摩擦元件熱彈不穩定分析結果

系統矩陣與特征值解法，根據邊界特征可知，2個表面之間的熱流只存在于接觸表面。所以接觸表面節點數為Nc時，如果另一個方向也是Nc，則Nc=N，否則Nc=N1+N2-N。那么此時一定存在矩陣S滿足：

Q=SQc

(18)

(19)

式中：Qc、I分別為接觸熱流密度Nc維向量、Nc×Nc單位矩陣；矩陣S為N×Nc矩陣。

根據熱流計算公式：

q=μvp

(20)

式中：p為壓力；v是相對滑動速度；μ為摩擦因數。

將式(20)寫成矩陣代入方程(18)中可得：

Q=SQc=μvSHT

(21)

代入方程(13)中可得：

(K+G+bL+μvSH)T=0

(22)

構造矩陣U，使得：

UT=bT

(23)

則矩陣U的特征值，就是擾動增長系數b。從式(22)中可以看到，擾動增長與摩擦因數和滑動速度相關。當擾動增長系數恰好為0時，對應的摩擦副的相對滑動速度v即為熱彈不穩定的臨界速度。根據制動器摩擦副熱固耦合有限元仿真計算結果，導出接觸表面的溫度矩陣T和壓力矩陣pc，所以過渡矩陣可以寫成如下形式：

pc=HT

(24)

摩擦副部分的熱固耦合有限元計算結果如圖2所示，其中圖2(a)所示為一對摩擦副(對偶片和摩擦片各一組)相互貼合下的正應力分布結果；圖2(b)所示為溫度場分布；圖2(c)所示為應力場分布。仿真工況如表1所示。

圖2 摩擦副熱固耦合有限元分析結果

表1 仿真工況

根據圖3所示的有限元的計算結果，輸出二維接觸表面的溫度和壓力矩陣。文中構造了一個接觸表面15×15的過渡矩陣，如表2所示。按照矢量關系提取摩擦副的壓力與溫度云圖中的數據，得到每一條路徑上的值，其結果如圖4所示。

圖3 不同階段摩擦副表面應力與溫度

圖4 摩擦副表面接觸壓力提取

提取表面節點沿自定義向量方向的取值，然后根據節點的取值結果，計算過渡矩陣。根據溫度場分布結果計算出來的過渡矩陣如表2所示。

表2 過渡矩陣(×105)

將過渡矩陣代入式(22)中計算獲取臨界速度的值，其計算結果如圖5所示。

圖5 臨界速度典型曲線

臨界速度先減小，后增大。在臨界速度產生極小值之后，隨著擾動頻率的增加，臨界速度是增大的，當擾動減小時系統的臨界速度也會增大。當沒有擾動時，理論上就不會產生系統失穩，那么臨界速度理論上是無限大。也就是說無論轉速多大，都不會超過臨界速度，系統也就不會產生失穩。但是這在實際當中是不可能的，因為擾動必然存在。

針對典型工況，研究溫度造成熱應力、壓力變化等對臨界速度的影響。由于系統的輸入就是由溫度-壓力計算獲得的過渡矩陣，所以溫度對于熱應力的影響會通過壓力矩陣傳遞到過渡矩陣H上。根據摩擦元件的常用工作狀態，文中制定了6種典型工況參數(如表3所示)。這6種典型工況涵蓋了制動器工作過程的大部分的運作狀態。

表3 6種典型工況的參數

根據熱固耦合的有限元計算結果，結合6種工況獲得了不同工況下的臨界速度曲線，如圖6所示。可以看到，并不是熱流密度越高的工況造成系統失穩的可能越大。在散熱條件一定的情況下，熱流密度越大會造成制動器熱功率越大，那么在制動器使用過程中產生的溫度一般也越高。從圖中可以看到工況1和工況5的最小臨界速度較大，但是工況1的溫度最高，工況5的溫度最低(如圖7所示)，這說明臨界速度跟摩擦片整體的溫度高低沒有必然關系。6種工況中，工況6更容易產生系統失穩，這是因為工況6的制動時間更短，熱量的集中釋放比較明顯，造成了較為明顯的局部溫差。

圖6 不同工況下臨界速度變化

圖7顯示了6種工況下摩擦副表面溫度分布曲線和局部溫度的差異。圖7(a)所示為6種工況的最高與最低溫度曲線，可見工況1的溫度最高，工況5和工況6的溫度較低。

圖7(b)所示為不同工況下的摩擦片溫度差，左側的縱坐標代表絕對溫差，右側的縱坐標代表相對溫差，相對溫差為絕對溫差與最高溫度的比值。可見工況6的絕對溫差和相對溫差都比較大，工況2和工況4相對溫差高于工況1和工況5，這與臨界速度結果一致。這進一步說明了臨界速度跟摩擦片整體的溫度高低沒有必然關系，而是與局部溫差尤其是相對溫差相關。因此，引起失穩的重要原因之一是局部溫差，因為局部溫差增大了摩擦片熱彈性變形，所以系統也更加容易失穩。

圖7 不同工況下摩擦片表面溫度分布曲線和表面溫差

3 結論

提出一種理論和仿真相結合的計算方法，可以將系統擾動因素和系統結構因素相統一，計算結果可以直接用來指導摩擦片熱彈性失穩的優化設計。主要結論如下：

(1)擾動頻率通過臨界轉速影響系統失穩，造成系統失穩的臨界速度存在一個極小值點，在該點系統最容易產生失穩，高頻擾動和低頻擾動對系統失穩的影響較小。

(2)摩擦副局部熱應力過大是造成系統失穩的主要原因，摩擦副整體溫度的高低與系統失穩之間不存在必然的聯系。

(3)摩擦片結構對系統失穩也有明顯影響，對于文中所研究的摩擦片結構，其在通常摩擦片的線速度范圍內(一般不超過50 m/s)，容易引起系統失穩的擾動頻率區間為50～550 Hz。可以通過改變設計工況，或者改變摩擦片結構的方法進行調整，以避免系統失穩。