基于勒讓德-掃描法的機械結合部法向接觸剛度區(qū)間分析

劉 勇, 渠 杰, 張 青, 王大偉, 閆方超

(1. 中國民航大學 航空工程學院,天津 300300;2. 天津布爾科技有限公司,天津 200392)

機械零件間緊密接觸的結合部分被稱為“機械結合部”,結合面間的接觸實際上是兩個粗糙表面間大小不一的微凸體之間的接觸[1]。航空發(fā)動機是由上萬個零件組成的復雜機械系統(tǒng),其內部存在大量的機械結合部,影響航空發(fā)動機的整機動態(tài)特性。因此,在研究航空發(fā)動機的整機動態(tài)特性時,將結合部的法向接觸剛度考慮在內就顯得尤為必要。

目前關于結合部法向接觸剛度模型的研究大多基于確定性假設,但在實際工程中,受材料表面的不均勻性、加工誤差、測量誤差等因素的影響,粗糙表面微觀形貌參數會存在較大的不確定性,導致結合部法向接觸剛度也存在不確定性。Fukuoka[2]的研究指出結合部法向接觸剛度的不確定性是連接結構剛度離散性的主要誘因,而傳統(tǒng)的確定性接觸剛度模型無法準確預測結合部法向接觸剛度。因此,考慮粗糙表面形貌參數的不確定,建立結合部法向接觸剛度的區(qū)間模型,對于研究結合部法向接觸剛度以及確保航空發(fā)動機的穩(wěn)定運行都具有至關重要的作用。

粗糙表面接觸建模是法向接觸剛度研究的基礎。Greenwood等[3]以Hertz接觸理論為基礎,建立了經典彈性接觸(Greenwood-Williamson,GW)模型,為粗糙表面統(tǒng)計學接觸研究提供了可靠的支撐。GW模型存在諸多缺點,為進一步完善GW模型,許多學者針對于GW模型的不足進行了修正[4-6]。基于GW模型及其改進模型,一些學者針對于結合部法向接觸剛度進行了研究。莊艷等[7]研究了微凸體粗糙峰之間的側向接觸,建立了考慮微凸體實際大小和空間分布的法向接觸剛度計算模型,并將實測的表面形貌參數代入所建立模型中,通過與試驗對比驗證了所提出模型的有效性。田小龍等[8]建立了考慮粗糙表面微凸體相互作用的法向接觸剛度統(tǒng)計學模型,在此基礎上研究了不同塑性指數下微凸體之間相互作用對法向接觸剛度的影響。Yin等[9]將微凸體變形階段重新劃分,建立了考慮硬度變化的法向接觸剛度統(tǒng)計學模型。考慮到傳統(tǒng)法向接觸剛度模型表達式復雜,難以獲得接觸載荷與法向接觸剛度之間的之間關系,Xiao等[10]基于不同微觀接觸模型的數值模擬結果,建立了關于法向載荷的法向接觸剛度顯式近似表達式。Bai等[11]考慮加工表面的實際形貌并基于高斯分布建立了更符合真實情況的虛擬表面,通過計算微凸體變形與接觸力的關系,建立了粗糙表面法向接觸剛度的統(tǒng)計學模型。綜上所述,關于結合部法向接觸剛度統(tǒng)計學模型的研究大多基于確定性假設,忽略了表面微觀形貌參數的不確定性。因此,為了更深入研究結合部法向接觸剛度,需要考慮表面形貌參數不確定性的影響。

在不確定性研究中,由于區(qū)間算法只需定義不確定變量的上下邊界,無需獲取不確定變量的精確概率分布,并且對試驗數據依賴性小,區(qū)間算法被廣泛用于不確定性研究中[12-15]?；贚egendre多項式所建立的區(qū)間代理模型可避免復雜運算并且計算效率較高,在多體動力學和轉子動力學領域取得了應用[16-17]。與傳統(tǒng)的切比雪夫包絡函數法相比,基于Legendre多項式的掃描法求解精度更高[18]。目前利用區(qū)間算法對于連接結構的研究大多集中于宏觀動力學方面[19-20],采用區(qū)間算法對于結合部接觸剛度進行研究的報道較少,Zhao等[21]基于測量所得的粗糙表面分形參數,并引入切比雪夫區(qū)間法,建立了粗糙表面法向和切向接觸剛度的分形區(qū)間模型。在此基礎上,李玲等[22]引入矩譜法計算了表面微觀形貌參數區(qū)間,并分析了微觀形貌參數對于粗糙表面法向接觸剛度區(qū)間的影響規(guī)律。上述關于粗糙表面接觸剛度的研究,作者重點關注結合部法向接觸剛度的區(qū)間表征與分析,所建立的接觸剛度區(qū)間范圍不夠精確,導致模型的實際參考價值不足。

有鑒于此,本文基于勒讓德-掃描法(Legendre scanning method, LSM)建立機械結合部法向接觸剛度區(qū)間模型并進行求解,即基于勒讓德多項式和粗糙表面統(tǒng)計學接觸剛度模型,建立考慮粗糙表面形貌參數不確定的法向接觸剛度區(qū)間模型,將微凸體曲率半徑、分布密度、微凸體高度方差視為區(qū)間變量,并采用掃描法求解區(qū)間模型。以45號鋼為算例,通過與莊艷等研究中的法向接觸剛度試驗結果以及傳統(tǒng)區(qū)間算法求解結果進行對比,驗證了本文建立接觸剛度區(qū)間模型的有效性和準確性?；谠撃Ｐ?分析了不同區(qū)間變量對結合部法向接觸剛度區(qū)間的影響規(guī)律。

1 粗糙表面法向接觸剛度模型

兩個粗糙表面的接觸如圖1(a)所示,可利用一個等效粗糙面與和一個光滑剛性平面的接觸來代替[23]。在圖1(b)中:β為微凸體曲率半徑;z為微凸體高度,且有z=ω+d;ω為微凸體變形量;d為光滑剛性平面與微凸體平均高度之間的距離。一般認為,隨著接觸壓力的增加,微凸體變形經歷三個階段:彈性、彈塑性和塑性變形。

圖1 粗糙表面接觸示意圖Fig.1 Rough surface contact

當載荷較低或者接觸變形量ω較小時,微凸體發(fā)生彈性變形。接觸載荷pe、接觸面積Ae、平均接觸壓力Ne與ω之間的關系可利用Hertz彈性接觸理論描述

(1)

Ae=πβω

(2)

(3)

式中:E*為兩個表面接觸的等效彈性模量,E*=E/(1-μ2);μ為材料的泊松比。

高載荷下,微凸體發(fā)生完全塑性變形,該階段平均接觸壓力pp、接觸面積Ap和接觸載荷Np可由AF模型[24]表示

pp=H

(4)

Ap=2πβω

(5)

Np=2πβHω

(6)

式中,H為材料硬度,H與屈服強度δy相關,滿足H=2.8δy。

金屬材料表面受到的接觸應力超過屈服強度時,微凸體存在彈塑性變形狀態(tài)。KE彈-塑性接觸模型考慮了微凸體的彈塑性變形,更貼近真實接觸狀態(tài)?；贙E接觸模型得到的接觸載荷F(hn)和法向接觸剛度K(hn)表達式為

(7)

(8)

式中:hn為粗糙面高度均線與剛性平面之間的距離,下標“n”為對變量進行無量綱化;zn為無量綱參數;N為粗糙表面參與接觸的微凸體總數,N=ηAn(An為名義接觸面積,η為微凸體分布密度);σ為微凸體表面高度方差,σ與微凸體高度分布標準差σs、粗糙表面高度分布函數φ(z)之間滿足以下關系

(9)

2 考慮形貌參數不確定性的結合部法向接觸剛度區(qū)間模型

由式(8)可看出,影響結合部法向接觸剛度的參數較多,傳統(tǒng)的統(tǒng)計學接觸模型為簡化計算,假設表面微凸體的曲率半徑相同,該假設與實際情況不符。在機加工表面上,微觀形貌參數的變化會使結合部法向接觸剛度產生較大的不確定性,不利于結合部法向接觸剛度的分析和優(yōu)化。為此,本文基于勒讓德多項式和KE彈塑性接觸剛度模型,建立考慮表面形貌參數不確定性的接觸剛度區(qū)間模型。

2.1 勒讓德多項式

一個定義域是[-1,1]的函數f(x)可以近似采用p階截斷勒讓德多項式表示[25]

(10)

式中:fi為截斷勒讓德多項式常系數;Li(x)為勒讓德多項式。勒讓德多項式滿足以下遞推關系

(11)

同理,對于具有多個不確定性輸入參數的多元函數f(x1,x2,…,xm),可以使用多維Legendre多項式進行近似表征

f(x1,x2,…,xm)≈fL(x1,x2,…,xm)=

(12)

計算上述高維問題所需插值點數量極多,為減少計算量,提升計算效率,對式(12)化簡可得

(13)

式中,ψi(x)為i個一維多項式的張量積,可由下式進行計算

ψi(x)=Li1,i2,…,im(x1,x2,…,xm)=
Li1(x1)Li2(x2)…Lim(xm)

(14)

式中,γi為勒讓德多項式系數矢量,表達式為

γi=(TTT)-1TTY=[γ0,γ1,…,γk-1]T

(15)

式中:Y為插值點處的模型輸出矢量;T為勒讓德轉換矩陣?；谏鲜龉?T可表示為

(16)

(17)

2.2 考慮形貌參數不確定性的法向接觸剛度區(qū)間模型

對粗糙表面進行法向接觸剛度分析時,首要任務是確定真實粗糙表面形貌參數,即微凸體等效曲率半徑β,微凸體分布密度η,微凸體表面高度方差σ。由于區(qū)間算法具有對輸入參數的依賴性較小、無需獲取輸入變量的具體概率分布等優(yōu)點,被廣泛用于處理不確定性問題。因此,本文把上述三個參數視作區(qū)間變量xj(j=1,2,3)

(18)

(19)

根據式(8)、式(9)和式(12),可建立考慮形貌參數不確定的結合部法向接觸剛度區(qū)間模型,表達式為

(20)

(21)

為簡化模型提高計算效率,結合部法向接觸剛度區(qū)間模型可由下式近似表示

(22)

式中:fL([ε])為基于勒讓德多項式計算得到的法向接觸剛度;ε=[ε1,ε2,ε3]為區(qū)間變量標準化后的矢量,εj(j=1,2,3)與區(qū)間變量xj(j=1,2,3)之間存在以下關系

(23)

輸出變量,即法向接觸剛度區(qū)間的不確定度λK可定義為

(24)

式中:Ku和Kl分別為法向接觸剛度區(qū)間的上界和下界;λK的取值范圍是[0,1],不確定度越小,說明所建立法向剛度區(qū)間的精度越高。

一般情況下可利用傳統(tǒng)勒讓德區(qū)間算法計算接觸剛度的上下邊界,但是傳統(tǒng)的區(qū)間算法精度有限,在求解多個不確定變量的問題時“包裹效應”會導致區(qū)間解過度放大,影響區(qū)間解的精度。為此,采用掃描法[26]對基于勒讓德多項式建立的粗糙表面法向接觸剛度區(qū)間模型進行求解。

2.3 求解結合部法向接觸剛度區(qū)間模型

對于式(22),在計算ψi([ε])時,由于“包裹效應”的存在,會產生估計值偏高的現象,導致計算結果偏大。因此可通過掃描法直接計算fL([ε])的最值以獲取法向接觸剛度的上界和下界。掃描法可以很好地控制“包裹效應”,并且當區(qū)間變量的個數m≤3時,掃描法可精確獲得區(qū)間模型的邊界?；趻呙璺ǖ膮^(qū)間上下界求解公式為

(25)

圖2為求解結合部法向剛度區(qū)間模型的算法流程。以下是模型求解步驟:

圖2 計算流程圖Fig.2 Flowchart of computation process

步驟1輸入區(qū)間向量的維數及Legendre多項式階數(j、p)。

步驟2將區(qū)間變量進行表示,同時根據式(17)確定插值點個數,并通過式(13)～式(16)利用一維多項式的張量積構造勒讓德轉換矩陣。

步驟3運用Matlab軟件設計循環(huán),在i+1個步驟中,根據式(20)～式(21)分別求解給定插值點處的結合部法向接觸剛度值,基于求解結果構造插值點處的模型輸出矢量Yi+1。將Yi+1代入式(15)中,結合步驟2得出的勒讓德轉換矩陣,可構造勒讓德多項式系數矢量γi+1。

步驟4基于上述參數,對于矩陣γi+1的每一列元素采用掃描法以及傳統(tǒng)區(qū)間算法求解法向接觸剛度區(qū)間模型,得到對應的法向接觸剛度邊界。

步驟5重復步驟3、步驟4直至循環(huán)結束,可得到結合部法向接觸剛度區(qū)間。

3 模型驗證

3.1 與試驗結果對比

為驗證基于LSM的接觸剛度區(qū)間模型的有效性,將模型結果與莊艷等研究中試驗測得的剛度值進行對比,結果如圖3所示。試驗所用材料為45號鋼,材料屬性參數為:彈性模量E=200 GPa,泊松比v=0.29,材料硬度H=2 025 MPa。莊艷等研究測量得到的45號鋼表面微觀形貌參數區(qū)間范圍,如表1所示。

表1 表面形貌參數區(qū)間分布表

圖3 法向接觸剛度理論區(qū)間與試驗結果對比Fig.3 Comparison of theoretical intervals of normal contact stiffness and experimental results

如圖3所示,隨著無量綱接觸距離ωn的減小,基于本文方法得出的法向接觸剛度區(qū)間與試驗數據基本吻合,證明了本文建立區(qū)間模型的有效性。當接觸間隙為0.4時,本文所建立結合部剛度區(qū)間模型的計算值偏離試驗值的主要原因在于:未充分考慮微凸體之間的相互作用;忽略了材料自身形變時的強化。

3.2 與現有不確定接觸模型對比

為驗證基于LSM的接觸剛度區(qū)間模型的準確性,將模型結果與基于切比雪夫包絡函數法所得結果進行對比。如圖4所示,LSM和切比雪夫包絡函數法均可獲得較可靠的結果?；贚SM所得法向接觸剛度區(qū)間不確定度是0.263,利用切比雪夫包絡函數法求解所得法向接觸剛度區(qū)間不確定度是0.468,說明LSM相較于切比雪夫包絡函數法在法向接觸剛度區(qū)間模型計算時具有更高的精度。

圖4 兩種區(qū)間算法計算結果對比Fig.4 Comparison of calculation results using different interval algorithms

3.3 與傳統(tǒng)區(qū)間算法對比

根據表1中的數據,由式(20)～(25)得到法向接觸剛度區(qū)間。為方便對比,采用傳統(tǒng)區(qū)間算法對式(22)進行求解。如圖5所示,對結合部法向接觸剛度區(qū)間模型進行求解時,掃描法和傳統(tǒng)區(qū)間算法均可獲得較可靠的結果。利用傳統(tǒng)區(qū)間算法求解所得法向接觸剛度區(qū)間的不確定度為0.633,利用掃描法求解所得法向接觸剛度區(qū)間的不確定度為0.263,求解所得區(qū)間的不確定度大幅降低,說明掃描法相較于傳統(tǒng)區(qū)間算法在求解法向接觸剛度區(qū)間模型時具有更高的精度。

圖5 不同算法下法向接觸剛度隨ωn變化關系Fig.5 Effect of contact gap on normal contact stiffness under different algorithms

法向接觸剛度區(qū)間上下界之差(即界差)與接觸間隙之間的關系,如圖6所示。隨著接觸間隙的減小,界差逐漸增大,掃描法計算得出的界差小于勒讓德法得出的界差。這是由于掃描法是基于密集網格樣本,計算不同參數區(qū)間內所有采樣點的響應,從而獲得較為精確的法向接觸剛度區(qū)間模型邊界。因此,掃描法可有效抑制區(qū)間擴張導致的“包裹效應”。

圖6 掃描法和區(qū)間算法計算法向接觸剛度界差對比Fig.6 Comparison of the boundary difference calculated by scanning method and interval algorithm

結合上述驗證,進一步說明了本文建立區(qū)間模型的有效性和準確性。最后基于LSM分析了表面形貌參數對結合部法向接觸剛度區(qū)間的影響規(guī)律。

4 分析與討論

微凸體分布密度(η)、微凸體曲率半徑(β)、微凸體表面高度方差(σ)分別為區(qū)間參數時,采用LSM所得法向接觸剛度變化范圍如圖7～圖9所示。由圖7～圖9可知,不同區(qū)間變量對法向接觸剛度的影響不同。如圖7(a)所示,結合部法向接觸剛度變化范圍受微凸體分布密度(η)的影響較小。圖7(b)表明隨著微凸體分布密度的增大,不同載荷下的結合部法向接觸剛度值也相應增加,但變化不明顯。造成上述現象的原因是,在相同接觸間隙下,當微凸體分布密度增大時,單位面積上實際參與接觸的微凸體數量增多,但由于微凸體曲率半徑較小,真實接觸面積變化不明顯,所以結合部法向接觸剛度范圍受微凸體分布密度(η)的影響較小。

圖7 η對法向接觸剛度區(qū)間影響Fig.7 Influence of η on normal contact stiffness interval

如圖8(a)所示,結合部法向接觸剛度變化范圍受微凸體曲率半徑(β)的影響較明顯。圖8(b)表明隨著微凸體曲率半徑的增大,不同載荷下的結合部法向接觸剛度值也相應增加且變化較明顯。這是由于當微凸體曲率半徑變化時,相同接觸間隙下真實接觸面積的變化較明顯,真實接觸面積的改變會影響接觸剛度,導致結合部接觸剛度呈區(qū)間變化。

圖8 β對法向接觸剛度區(qū)間影響Fig.8 Influence of β on normal contact stiffness interval

如圖9(a)所示,結合部法向接觸剛度變化范圍受微凸體表面高度方差σ的影響最顯著。圖9(b)表明隨著微凸體表面高度方差σ的增大,不同載荷下的結合部法向接觸剛度值也相應減小,并且變化較為明顯。這是因為σ影響粗糙表面微凸體高度分布,σ增大時,表面粗糙度增加,微凸體與接觸表面之間的真實接觸面積減小,從而導致接觸剛度減小。由圖7～圖9可以看出,相較于微凸體曲率半徑和微凸體分布密度,表面高度方差對于法向剛度區(qū)間的影響更為顯著。

圖9 σ對法向接觸剛度影響Fig.9 Influence of σ on normal contact stiffness interval

5 結 論

(1) 考慮粗糙表面微觀形貌參數的不確定性,基于勒讓德多項式和傳統(tǒng)KE彈塑性接觸剛度模型,建立了結合部法向接觸剛度區(qū)間模型。

(2) 采用基于勒讓德多項式的掃描法對上述區(qū)間模型進行求解,并與傳統(tǒng)區(qū)間算法求解結果進行了對比。結果表明,兩種方法對結合部法向接觸剛度值的估計均較為可靠,但是掃描法相較于傳統(tǒng)區(qū)間算法在求解法向接觸剛度區(qū)間模型的精度上有明顯優(yōu)勢,使不確定度由原來的0.633減小到0.263。

(3) 以45號鋼為例進行試驗驗證,可看出采用勒讓德掃描法計算所得法向接觸剛度區(qū)間結果與宏觀接觸剛度試驗結果具有相同變化趨勢且預測準確性較高。利用LSM建立并求解法向接觸剛度區(qū)間模型可為航空發(fā)動機的結合部設計提供參考。

(4) 將區(qū)間變化的表面微觀形貌參數代入粗糙表面法向接觸剛度區(qū)間中并利用掃描法求解。結果表明,微凸體表面高度方差對結合部法向接觸剛度區(qū)間的影響最大,微凸體曲率半徑次之,影響較小的是微凸體分布密度。