雙粗糙滑動(dòng)表面間接觸熱導(dǎo)研究

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

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

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

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

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

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

1 滑動(dòng)粗糙表面間接觸熱導(dǎo)分析

1.1 接觸熱導(dǎo)的理論模型

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

(1)

(2)

(3a)

(3b)

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

(4)

(5)

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

(6)

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

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

(7)

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

(8)

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

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

(9)

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

(10)

其中，

(11a)

(11b)

(11c)

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

(12)

對(duì)于完全相同的兩個(gè)滑動(dòng)接觸表面，結(jié)合文獻(xiàn)[20]中的式(18)、式(19)、式(24)，式(12)可以化簡(jiǎn)為

其中，

(14a)

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

(14b)

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

(14c)

2 接觸熱導(dǎo)比較驗(yàn)證

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

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

m0=2；m2=2；m4=2

(15)

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

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

(16)

Bush等[23]給出的微凸體峰頂高度的標(biāo)準(zhǔn)偏差σi和RMS粗糙度σ′的表達(dá)式分布為

(17a)

(17b)

其中，

(18)

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

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

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

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

(19)

(20)

(21)

(22)

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

Ks=2K1K2/(K1+K2)

(23)

(24)

(25)

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

(26)

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

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

2.3 高速滑動(dòng)條件下接觸熱導(dǎo)比較

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

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

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

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

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

圖4 不同滑動(dòng)速度下，接觸熱導(dǎo)的理論解、數(shù)值解 以及Mikic模型接觸熱導(dǎo)值比較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 接觸熱導(dǎo)在宏觀熱傳導(dǎo)分析中應(yīng)用

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

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

4 結(jié) 論

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

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

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