劉雨薇， 葉福浩， 李和言， 吳健鵬， 孫能祥

(1. 中國礦業(yè)大學(北京) 機電與信息工程學院，北京 100083；2. 華為機器有限公司，東莞 523808；3. 深圳技術大學 城市交通與物流學院，深圳 518118；4. 北京理工大學 機械與車輛學院，北京 100081)

機械系統(tǒng)的各零部件結合面之間存在各種各樣的接觸熱導，接觸熱導是其熱分析過程中最重要的不確定因素，對傳熱效率及熱力耦合起著至關重要的作用，在極端情況下會影響熱控系統(tǒng)和機器運行的可靠性。對接觸熱導的預測研究已經成為低溫冷卻、機械制造、微電子、高溫超導、航空航天等諸多工程領域的關鍵環(huán)節(jié)[1]。

工程應用中，固體材料的真實表面是不光滑的，表面上分布有大量微小的凸起，使得實際的固體接觸面積只占名義接觸面積很小的一部分，并且絕大部分熱流通過接觸點傳遞出去[2]。因此，熱流流經接觸面時產生收縮，從而形成了額外的換熱阻力，即接觸熱阻(Thermal Contact Resistance，TCR)。接觸熱阻的倒數(shù)形式為接觸熱導(Thermal Contact Conductance，TCC)，定義為熱流密度與固體界面間的溫差之比。接觸熱導可用于表征接觸界面間的傳熱特性。

半個多世紀以來，學者們從理論推導、數(shù)值分析和實驗測試方面對固體界面間的接觸熱導開展了大量的研究工作。影響接觸熱導的因素大多為非線性因素，包括固體材料的種類、接觸壓強、表面粗糙度、表面氧化、接觸面溫度、接觸點的變形情況等[3-5]。早期的研究學者Cooper等[6]提出的CMY(Cooper-Mikic-Yovanovich)熱傳導模型為之后開展的許多研究奠定了基礎。Whitehouse等[7]基于隨機過程理論，考慮了微凸體曲率半徑和峰頂高度的相關性，采用統(tǒng)計學參數(shù)對固體表面間的彈性接觸問題進行研究，建立了彈性接觸的表面接觸熱導模型(W-A模型)。Mikic[8]首次明確給出了微凸體分別發(fā)生純塑性變形和純彈性變形時兩種形變條件下的接觸熱導模型。Yovanovich等[9-10]通過對涂層表面的接觸熱導進行研究，在表面輪廓高度遵循高斯分布的條件下對CMY模型進行了擴展，建立了CMY塑性接觸熱阻模型。徐瑞萍等[11]對真空低溫環(huán)境下的接觸熱阻進行了研究，采用Cantor集分形理論描述粗糙接觸界面的拓撲形貌，考慮了真實接觸面的分形維度和變形對接觸熱阻的影響，推導出基于分形理論的接觸熱阻理論模型，該理論模型計算結果與實驗數(shù)據(jù)吻合良好。Liu等[12-13]以碳碳復合材料和超耐熱不銹鋼GH600為基體材料，對高溫度條件下的接觸熱導進行了試驗研究，結果表明界面溫度對接觸熱導值具有很大影響，接觸熱導隨著表面溫度的升高而減小。Fieberg等[14]基于瞬態(tài)的溫度測量方法，通過紅外攝像儀記錄的接觸表面附近的瞬態(tài)溫度反推接觸熱阻值，獲得了高溫、高壓條件下的接觸熱導表達式，其中接觸壓強的變化范圍為7.5～85 MPa，溫度范圍為60 ℃～280 ℃。研究發(fā)現(xiàn)接觸熱導與表面接觸壓強間呈線性關系，并且揭示了溫度對接觸熱導的影響不可忽略。Bi等[15]對界面溫度、接觸壓強以及基體材料對接觸熱導的影響進行了研究，溫度范圍從70～290 K，壓強從0.2～0.7 MPa，結果表明接觸熱導與界面溫度之間具有冪指數(shù)關系，與接觸壓強之間具有線性關系。Dou等[16]對表面粗糙度、接觸壓強以及界面溫度對接觸熱導的影響開展了實驗研究，指出相對較為光滑的表面，具有較高表面粗糙度的界面接觸熱導隨著基體溫度的升高增長較快。

上述研究均集中于靜態(tài)條件下，外界加載情況、界面初始溫度、表面粗糙程度等因素對接觸熱阻的影響。但是，當固體接觸面之間存在相對運動時，運動速度對接觸熱導也具有非常重要的影響。例如活塞式發(fā)動機中的活塞與氣缸之間，在相對運動的過程中有很大一部分熱量需要通過接觸導熱傳遞出去。本文將在已有研究成果的基礎上，對滑動雙粗糙表面間的有效接觸熱導進行分析和預測。

1 滑動粗糙表面間接觸熱導分析

1.1 接觸熱導的理論模型

文獻[17]中，本文作者已經對滑動粗糙表面間的瞬態(tài)導熱問題進行了研究，假設微凸體的變形為彈性變形，峰頂高度服從高斯分布，并且忽略平行于滑動方向的熱傳導，即只考慮垂直于滑動方向的導熱。在Greenwood和Williamson(GW)統(tǒng)計接觸模型的基礎上，建立了一個高Peclet數(shù)(Pe=Va/k)下，粗糙表面間的滑動熱接觸模型。該理論模型對接觸面上的平均閃點溫度進行了預測，并且給出了平均名義熱流密度與體積溫差之間的關系式，即

(1)

(2)

(3a)

(3b)

對于兩個完全相同的接觸表面，即具有相同的材料屬性和表面粗糙參數(shù)，式(1)和式(2)可以化簡為

(4)

(5)

由此得到高Pe數(shù)下，滑動粗糙接觸面間單位面積上的接觸熱導的表達式為

(6)

1.2 接觸熱導的數(shù)值模型

作者本人對單對微凸體的接觸導熱過程進行了數(shù)值分析，并且對理論模型的研究范圍進行了擴展，對不同Pe數(shù)下單對微凸體的導熱量進行仿真計算，擬合得到由體積溫差ΔT引起的兩個微凸體間的導熱量公式為[18]

(7)

式中：a0為最大接觸半徑。根據(jù)Hertz接觸理論，最大接觸半徑為[19]

(8)

式中：R*為物體1和物體2的粗糙接觸面上相互作用的兩個微凸體峰頂?shù)膹秃习霃剑琑*=R1R2/(R1+R2)，d0由Dou等研究中的式(5)定義。

采用統(tǒng)計方法對單對微凸體間的導熱量進行積分，得到單位滑動距離下、單位名義接觸面積上的總導熱量為

(9)

式中：Ni為平面i上單位面積上的微凸體峰頂個數(shù)；φi(hi)為微凸體峰頂高度的概率密度函數(shù)。將式(7)代入式(9)，得到

(10)

其中，

(11a)

(11b)

(11c)

根據(jù)牛頓冷卻定律qc=hcΔT，兩滑動接觸面間的換熱系數(shù)，即單位面積上的接觸熱導為

(12)

對于完全相同的兩個滑動接觸表面，結合文獻[20]中的式(18)、式(19)、式(24)，式(12)可以化簡為

其中，

(14a)

F1(If)=0.987 7(If)0.908 1

(14b)

F2(If)=0.598 2(If)0.932 3

(14c)

2 接觸熱導比較驗證

2.1 接觸熱導模型主要參數(shù)確定

Nayak[21]提出任意各向同性，并且服從高斯分布的隨機表面都可以用三個功率譜密度矩來描述，分別為m0,m2和m4。這些譜密度矩可以通過表面輪廓方程z(x)來定義，即有

m0=2；m2=2；m4=2

(15)

式中：z(x)為輪廓表面相對于距離為x處的任意參考平面的高度偏差，可以通過表面輪廓測量儀或者原子力顯微鏡(Atomic Force Microscope，AFM)來提??；零階譜矩m0為輪廓的高度分布，m0的平方根等于表面的RMS(Root-Mean-Square)粗糙度；二階譜矩m2為輪廓的斜率分布，等于斜率的均方差；m4為輪廓的曲率分布，等于曲率的均方差。

接觸熱導模型參數(shù)Ni和Ri可以通過表面輪廓的譜矩獲得[22]，即

(16)

Bush等[23]給出的微凸體峰頂高度的標準偏差σi和RMS粗糙度σ′的表達式分布為

(17a)

(17b)

其中，

(18)

采用McCool[25]文章中的材料參數(shù)和表面矩函數(shù)，結合式(17)和式(18)可以得到單個粗糙表面的微凸體峰頂半徑R、峰頂高度偏差σ和單位名義接觸面積上的微凸峰個數(shù)N，所有參數(shù)見表1。

表1 接觸界面材料參數(shù)和表面粗糙度參數(shù)Tab.1 Material properties of the contact interface and surface roughness parameters

2.2 靜態(tài)條件下接觸熱導比較

Sridhar等[26-27]給出了靜止狀態(tài)下，微凸體發(fā)生完全彈性變形條件下的表面接觸模型

(19)

(20)

(21)

(22)

式中：“erfc”為余誤差函數(shù)；Ar為實際接觸面積；Aa為名義接觸面積；a為平均接觸半徑，P為接觸壓強；Ks為接觸表面的復合導熱系數(shù)；σ′為表面有效RMS粗糙度；m為表面輪廓的有效平均絕對斜率。

Ks=2K1K2/(K1+K2)

(23)

(24)

(25)

m可以根據(jù)Lambert等[28]給出的經驗獲得，即

(26)

將表1中的參數(shù)分別代入式(5)、式(19)和式(22)～式(26)可以確定Mikic模型的接觸熱導值。圖1顯示了接觸壓強在0.01～1.0 MPa內，Mikic模型的接觸熱導值hc_Mikic與由式(13a)計算得到的接觸熱導數(shù)值解hc隨接觸壓強的變化趨勢，其中“+”表示hc_Mikic，“×”表示hc。從圖1中可以看出，隨著接觸壓強的增大，hc_Mikic逐漸大于hc，二者之間的誤差始終小于10%，且與接觸壓強成近似線性關系。

圖1 Mikic模型接觸熱導值與本文模型接觸熱導值比較Fig. 1 Thermal contact conductances estimated from Mikic model and the model proposed in this paper

2.3 高速滑動條件下接觸熱導比較

將表1中的參數(shù)代入式(6)和式(13c)分別計算得到高Pe數(shù)下接觸熱導的理論解與數(shù)值解，如圖2所示，其中hc表示由式(13c)確定的數(shù)值解，hc_theo表示由式(6)確定的理論解。不同接觸壓強下，數(shù)值模型的接觸熱導值均大于理論模型的接觸熱導值，前者約為后者的2～3倍。此外，隨著接觸壓強的增大，二者之間的差距逐漸減小。

圖2 高Pe數(shù)下，接觸熱導理論解與數(shù)值解比較Fig. 2 Comparison of theoretical solution and numerical simulation to the high Peclet limit

微凸體在相互作用的過程中，實際接觸面積大小隨著接觸時間不斷發(fā)生變化，同時接觸位置也隨之發(fā)生移動，這樣的運動形式使得更多的表面參與到熱交換過程中。理論模型中假設接觸面積的大小隨著接觸時間的增大先增大后減小，而接觸面積的中心位置保持不變，如圖3所示。在該假設下，采用解析法計算得到的總導熱量將小于實際的導熱量，因此理論接觸熱導值小于數(shù)值接觸導熱值是合理的。

圖3 單對微凸體接觸面積變化趨勢Fig.3 Change tendency of contact area during the single asperity interaction

圖4為接觸壓強為0.1 MPa時，高速滑動條件下接觸熱導的理論解、數(shù)值解以及Mikic模型的接觸熱導值隨滑動速度的變化曲線。由圖4可知，不同滑動速度下，數(shù)值模型的接觸熱導值均為理論模型接觸熱導值的2.6倍，因此可將2.6作為理論模型的修正系數(shù)。此外，與理論模型和數(shù)值模型的接觸熱導值相比，Mikic模型的接觸熱導值始終為一條直線，無法反映接觸熱導隨滑動速度的變化規(guī)律。在滑動速度由50 m/s增大到150 m/s的過程中，數(shù)值模型的接觸熱導值為Mikic模型的接觸熱導值的3.5～6.2倍。由此可知，現(xiàn)有的靜態(tài)模型(如Mikic模型)對于預測高速滑動下的接觸熱導偏差很大，并且滑動速度越大，偏差越顯著，不能用于對高速滑動狀態(tài)下的接觸熱導進行預測。

圖4 不同滑動速度下，接觸熱導的理論解、數(shù)值解 以及Mikic模型接觸熱導值比較Fig.4 Comparison of thermal contact conductance predicted by the theoretical solution, by the numerical simulation and by the Mikic model with different sliding velocity

3 接觸熱導在宏觀熱傳導分析中應用

在對宏觀粗糙表面間的導熱問題進行數(shù)值分析時，通常情況下有限元網格難以準確地描述表面的粗糙度情況，即網格數(shù)量和細化程度都難以滿足要求。而且數(shù)量龐大的網格有可能導致計算不收斂，或者計算時間被大大延長。現(xiàn)有的大多數(shù)理論是用一個假想的具有復合粗糙度的表面與一個光滑表面來代替實際的雙粗糙接觸表面，并且假設接觸面積可以維持足夠長的時間以保證穩(wěn)態(tài)或者準穩(wěn)態(tài)導熱的建立。在這樣的假設下，粗糙表面上的熱源為固定熱源，光滑表面上的熱源是移動熱源，對于摩擦學應用中的典型情況(Pe≥1)，摩擦產生的絕大部分熱量將流入“移動的物體”，僅有非常少的熱量流入靜止的物體，除非移動的物體具有很低的熱導率。此外，如果兩個滑動粗糙表面間的Pe數(shù)(Peclet number)較大，尤其是兩個接觸表面的粗糙度和硬度相差不多的情況下，典型的微凸體間的相互作用將是瞬態(tài)的，持續(xù)的時間極短，無法滿足建立穩(wěn)定導熱狀態(tài)的條件，基于穩(wěn)態(tài)或者準穩(wěn)態(tài)導熱的分析方法將不再適用?；诒疚氖?13)中的接觸熱導表達式，可以通過在兩個光滑表面間加入有效接觸熱導的有限元模型來反映真實雙粗糙表面間的熱傳導過程，該有效接觸熱導包含了粗糙表面的所有微觀信息，采用該方法可以極大地提高計算效率。

此外，本文的研究結果為準確計算摩擦副表面的熱流分配系數(shù)提供了新的方法和思路。目前熱流分配系數(shù)主要通過理論公式[29]計算獲得，一旦材料確定，熱流分配系數(shù)也隨之確定。實際上，熱量在傳遞過程中涉及到摩擦表面接觸熱阻及接觸不充分等問題，并且由于摩擦副雙方的導熱性及幾何形狀不同，摩擦過程中熱量的傳遞速度及散熱性能也存在較大差異，這些均會對摩擦熱流分配產生較大影響。根據(jù)熱傳導的本質特性可知，兩個完全相同的表面上微凸體間的瞬態(tài)導熱過程是對稱的，摩擦生熱量會均等地流入兩個接觸面。然而，工程實際中摩擦副間的尺寸和材料特性往往是不對稱的，將會導致兩個接觸面間出現(xiàn)體積溫差，并且使得摩擦熱流量不再均等分配。對于此類問題，在對摩擦熱流分配系數(shù)進行求解時，采用疊加原理，將導熱問題等效于對稱的瞬態(tài)導熱問題與由于體積溫差產生的附加熱流量的疊加。其中，附加的熱流量取決于接觸面間的有效接觸熱導值，而附加熱流量又決定了最終的摩擦熱流分配系數(shù)。

4 結 論

粗糙表面之間的不完全接觸造成熱流線收縮而形成接觸熱阻，接觸熱阻會對摩擦熱流的分配產生影響，對機械設備的可靠性和使用壽命具有非常重要的影響。本文從微觀角度對具有相對滑動的兩個粗糙界面間的接觸熱導進行了分析研究，得到了適用于不同Pe數(shù)下的接觸熱導表達式，既可用于預測滑動接觸面間的有效接觸熱導，也可預估靜態(tài)下的有效接觸熱導，從而為探究接觸表面熱流分配系數(shù)和強化接觸換熱設計提供依據(jù)。主要結論如下：

(1) 高Pe數(shù)下，滑動粗糙表面間的有效接觸熱導不僅與接觸壓強、材料熱參數(shù)以及表面粗糙度有關，而且與滑動速度平方根成正比。

(2) 低Pe數(shù)或者靜態(tài)下，材料參數(shù)、微凸體峰頂半徑一定時，有效接觸熱導與接觸壓強成近似線性關系。