考慮結(jié)合面特性影響的機床整機有限元靜力學分析模型*

王 豪,黃為彬,方 兵*

(1.福建農(nóng)林大學 機電工程學院,福建 福州 350002;2.福建省農(nóng)業(yè)信息感知技術(shù)重點實驗室,福建 福州 350002)

0 引 言

目前,基于計算機輔助設(shè)計的理念已經(jīng)廣泛地應(yīng)用于精密數(shù)控機床產(chǎn)品的開發(fā)中。設(shè)計者期望在設(shè)計階段準確掌握產(chǎn)品的剛度、熱特性、動態(tài)響應(yīng)以及加工精度等性能。但由于機床零部件數(shù)量眾多、結(jié)構(gòu)復(fù)雜,存在多種類型的結(jié)合面;結(jié)合面不連續(xù)的特點對機床剛度、熱傳導(dǎo)特性的影響顯著,且結(jié)合面自身的影響因素較多、規(guī)律不明顯。因此,結(jié)合面的建模準確與否在較大程度上影響整體分析的結(jié)果[1-3]。

在處理結(jié)合面問題時,“彈簧-阻尼”法是常用的方法之一。該方法用若干個“彈簧-阻尼”單元來模擬結(jié)合面的靜、動態(tài)特性。如朱堅民等人[4]使用了“彈簧-阻尼”單元模擬主軸刀柄結(jié)合面,并進行了動力學參數(shù)識別,結(jié)果表明實驗值與計算值相對誤差在4%以內(nèi),具有較高精度。李院生等人[5]采用了4節(jié)點的“彈簧-阻尼”單元模型來等效結(jié)合面,解決了螺栓固定結(jié)合面參數(shù)識別問題。蘭國生等人[6]采用了“彈簧-阻尼”法,建立了固定結(jié)合面法向接觸阻尼模型及結(jié)合面間阻尼損耗因子模型。

然而“彈簧-阻尼”單元法存在一些缺點。該方法是一種簡化自由度的方法,需要根據(jù)具體的問題,選擇單元的數(shù)量,確定單元的位置,并且忽略了各單元之間的相互影響。

為了解決上述問題,近年來研究人員已經(jīng)開始采用虛擬材料法建立固定結(jié)合面模型。

虛擬材料法是一種根據(jù)結(jié)合面接觸部位微觀結(jié)構(gòu),采用具備彈性模量、泊松比以及密度的材料來模擬結(jié)合面靜力學和動力學特性的方法。比如,張學良等人[7]提出了一種利用宏觀各向同性虛擬材料等效的固定結(jié)合部動力學參數(shù)化建模方法,即采用虛擬材料的彈性模量、泊松比等來模擬結(jié)合面靜力學和動力學特性,結(jié)果顯示,計算值與實驗值的絕對誤差在10%以內(nèi)。孫偉等人[8-9]采用了一層虛擬材料來對含螺栓復(fù)合材料結(jié)構(gòu)進行動力學分析,均達到了較高的精度。

然而,上述研究并未考慮彈塑性變形對虛擬材料屬性的影響。

因此,筆者基于Greenwood-Williamson(G-W)接觸分析模型,考慮彈塑性變形的影響,依據(jù)粗糙表面形貌統(tǒng)計參數(shù),提出一種新的考慮摩擦力影響的剛度模型;根據(jù)結(jié)合面的剛度特性,將結(jié)合面等效成一層虛擬材料,建立固定結(jié)合面的靜力學分析模型,將其應(yīng)用于整機的有限元分析,并且開展機床靜剛度測試實驗,以驗證上述分析結(jié)果的正確性。

1 結(jié)合面虛擬材料模型

結(jié)合面實質(zhì)上是由兩個相互接觸的粗糙表面組成,可以將其簡化為一個粗糙表面與一個剛性平面接觸。

GREENWOOD J A等人[10]早在20世紀60年代就提出了基于統(tǒng)計理論的G-W接觸分析模型。該模型假定凸峰的高度服從某種特定的分布(如高斯分布),凸峰具有相同的曲率半徑。當凸峰的變形量ω大于彈性極限值ωe時,凸峰開始呈現(xiàn)塑性變形的特點。然而,KOGUT通過有限元分析的方法發(fā)現(xiàn),凸峰由彈性變形階段開始進入塑性變形階段時,先要經(jīng)歷一個彈塑性混合變形階段;然后,直到其變形量大于等于110×ωe時,凸峰才呈現(xiàn)完全塑性變形的特點[11]。因此,單個凸峰的接觸過程可由彈性、彈塑性和塑性三個階段來描述。

1.1 彈性接觸變形

根據(jù)Hertz接觸理論,可建立彈性接觸面積ae、彈性接觸載荷pe,以及彈性變形量ω間的關(guān)系如下:

ae=πRω

(1)

(2)

當結(jié)合面間存在摩擦力時[12],由凸峰發(fā)生屈服的臨界平均接觸壓力pave=kυkμH可得到彈性變形極限,即:

(3)

式中:kυ為平均接觸壓力系數(shù),由材料泊松比υ確定[13],kυ=0.464 5+0.314 1υ+0.194 3υ2(若取υ=0.3,則kυ≈0.6);kμ為摩擦力影響系數(shù),其值與摩擦系數(shù)μ有關(guān)[9];H為材料硬度。

kμ的計算公式如下:

(4)

式中:μ為摩擦系數(shù)。

(5)

根據(jù)文獻[14]可知,單個凸峰的切向剛度為:

(6)

Qe和Pe分別表示該階段的總體切向載荷和法向載荷,其中:

(7)

式中:η為凸峰密度;A0為名義接觸面積;z為凸峰高度;d為剛性平面與粗糙表面中線的距離,且滿足關(guān)系式:ω=z-d;ωe為彈性極限值;φ(z)為高度分布的概率密度函數(shù)。

文獻[15]給出了單個凸峰所承受的切向載荷為:

(8)

式中:σy為與摩擦力方向相同的應(yīng)力;a為接觸點面積;p為接觸點上的載荷。

根據(jù)接觸面積與變形的關(guān)系,可得該階段總切向載荷為:

(9)

1.2 彈塑性接觸變形

當變形量在ωe和完全塑性臨界值ωp之間時,微凸峰處在彈塑性混合變形階段。

李玲等人[16]在接觸面積、平均接觸載荷滿足連續(xù)變化的基礎(chǔ)上,構(gòu)造了以變形量ω為自變量的多項式表達式,其接觸面積可表示為:

aep=πRκ1(ω)

(10)

其中:

接觸載荷可表示為:

(11)

類似的可定義剛度為接觸載荷對變形量的導(dǎo)數(shù),于是法向接觸剛度可以表示為:

(12)

其中:

κ3(ω)=-2g3(ω)+6(2-ω)g2(ω)+2(3ω-4)g(ω)+2。

則,切向接觸剛度可以表示為:

(13)

式中:Pep,Qep為該階段的總體切向載荷和法向載荷。

Pep和Qep表達式如下:

(14)

(15)

1.3 塑性接觸變形

KOGUT L等人[11]通過有限元計算發(fā)現(xiàn),當變形量ω大于110×ωe時,凸峰呈現(xiàn)完全塑性變形的特點。此時,接觸面積、載荷與變形量的關(guān)系可以表示為:

ap=2πRω

(16)

pp=Hap

(17)

類似的,可以求得該階段的法向接觸剛度和切向接觸剛度分別為:

(18)

(19)

Pp和Qp表示該階段的總體切向載荷和法向載荷,分別為:

(20)

(21)

綜合式(5)、式(12)和式(18),可得到總的法向接觸剛度為:

(22)

綜合式(6)、式(13)和式(19),可得到總的切向接觸剛度為:

(23)

1.4 虛擬材料屬性確定

采用虛擬材料模擬結(jié)合面的靜力學特性,需要知道虛擬材料的密度ρ(取兩接觸材料的平均密度)、彈性模量Ev、泊松比μv、虛擬材料的面積A和厚度t等參數(shù)。

根據(jù)機械加工表面形貌特征以及有限元分析的需求,參考文獻[17],可選取虛擬材料的厚度t為1 mm。

等效虛擬材料模型如圖1所示。

圖1 等效虛擬材料模型

若單位面積上的法向剛度和切向剛度分別為(量綱為[N·m-3]):

(24)

(25)

取一面積為ds的微元,在法向載荷Fn的作用下,其法向應(yīng)力和應(yīng)變分別為σ和ε,則有:

(26)

式中:δn為法向變形量。

由彈性模量定義及化簡式(14),得到:

(27)

式中:Ev為虛擬材料的等效彈性模量。

同理,可得虛擬材料的等效切變模量Gv為:

(28)

式中:τ為切變應(yīng)力;γ為切變應(yīng)變。

泊松比μv為:

(29)

2 整機靜剛度仿真分析

某立式加工中心的床身與立柱材料均為HT200,接觸表面采用銑削方式加工。在機床X、Y、Z三個方向各布置一根絲杠分別用于驅(qū)動主軸箱、工作臺和滑板。主軸前端安裝錐度為7:24的BT50刀柄,軸承采用4+2的布置方式,即前端采用4列并列,后端采用2列并列的方式配置,型號分別為7020C和7018A。

筆者對模型進行必要的簡化,并在CATIA V5r16中建立床身、立柱、滑板、工作臺、主軸箱以及導(dǎo)軌等主要零部件模型,按設(shè)計尺寸裝配成整機三維模型,再將CATIA V5r16中的整機三維模型導(dǎo)入ANSYS Workbench中,進行網(wǎng)格劃分、參數(shù)設(shè)定以及邊界條件設(shè)置等,以建立有限元分析模型[18]。

得到結(jié)合面的虛擬材料屬性是建立整機靜剛度分析模型的關(guān)鍵。此處以床身立柱結(jié)合面為列,該結(jié)合面由16只螺栓固定,正壓力約為56 000 N,接觸面積為0.238 m2,平均接觸壓力為2.35×105Pa。采用輪廓儀(泰勒霍普森Form Talysurf i120)測量所研究表面的形貌,并根據(jù)文獻[19]所提供的方法,經(jīng)計算得到服從高斯分布的表面凸峰的高度分布函數(shù)φ(z),平均半徑β為42.68 μm以及單位面積凸峰個數(shù)N為8.45×103個/mm2。

平均接觸壓力Pav與彈性模量Ev的關(guān)系如圖2所示。

圖2 平均面壓與彈性模量間的關(guān)系

由式(7)、式(14)、式(20)和高度分布函數(shù)φ(z)可得到總載荷P與偏差d間的關(guān)系,由式(24)和式(25)可分別得到彈性模量Ev、μv和偏差d間的關(guān)系。

綜合上述兩個關(guān)系,便可得到平均接觸壓力Pav與彈性模量Ev。隨著平均接觸壓力Pav的增加,虛擬材料的彈性模量也在增加,但趨勢趨于緩慢。

所建立的有限元分析模型主要考慮了床身立柱結(jié)合面、軸承內(nèi)圈與主軸的結(jié)合面、軸承外圈與箱體的結(jié)合面,以及主軸前后兩個端面;軸承、絲杠以及導(dǎo)軌則采用彈簧單元法進行等效。

根據(jù)實際工況,可得到主要結(jié)合面的虛擬材料參數(shù)如表1所示。

表1 主要結(jié)合面的虛擬材料屬性

仿真計算時,在主軸刀柄下端面沿著X、Y、Z三個方向施加大約1 960 N的載荷,反向載荷施加在工作臺面相應(yīng)位置上,并將床身底座與地面的接觸面固定約束。

仿真結(jié)果如圖3所示。

圖3 施加1 960 N載荷的整機變形

可見,主軸在承受載荷時,變形主要體現(xiàn)在主軸及其周邊的零部件,其中主軸前端影響最為明顯。在相同的載荷作用下,X、Y和Z方向的位移分別為68.6 μm、44.5 μm和47.3 μm,三個方向?qū)?yīng)的剛度分別為28.57 N/μm、44.04 N/μm和41.43 N/μm?？傮w來說,Y方向的剛度最大,Z方向上的剛度次之,X向剛度最小。

3 實驗與比較

為驗證整機靜剛度有限元分析結(jié)果的正確性,筆者開展了機床靜剛度測試實驗。

實驗測試原理如圖4所示。

圖4 機床靜剛度測量

筆者將機床各零部件放置在與仿真分析時對應(yīng)的位置上,以98 N為間隔逐步增加載荷至1 960 N,采用精度為0.98 N的BK-2型壓力傳感器測量載荷,由二次儀表顯示測量結(jié)果,采用分辨率為2 μm的千分表測量相對變形。

實驗結(jié)果如圖5所示。

圖5 載荷-變形曲線

隨著外加載荷的增加,主軸末端的相對變形也在逐漸增加,測得的“載荷-變形”數(shù)據(jù)基本呈線性規(guī)律?？蓴M合得到X、Y、Z三個方向的剛度,分別為26.60 N/μm、41.87 N/μm和40.17 N/μm。

為研究結(jié)合面對整體剛度特性的影響,筆者對未考慮結(jié)合面特性的模型進行分析計算。所謂不考慮結(jié)合面是指忽略由連接部位(也就是結(jié)合面)引起的剛度的變化。

在有限元模型中,筆者直接將連接部位黏接成一體,使得結(jié)合部與被連接部件的材料屬性一致。

實驗結(jié)果與仿真結(jié)果如表2所示。

表2 整機靜剛度仿真與實驗對比

其中,仿真的結(jié)果分為兩個模型,一個是考慮了結(jié)合面影響的模型,另一個是沒有考慮結(jié)合面影響的模型。

對比結(jié)果表明:不考慮結(jié)合面的影響時,機床工藝系統(tǒng)的靜剛度與實驗值的相對誤差達到17%左右;考慮結(jié)合面影響之后,其相對誤差基本在5%以內(nèi)。以上結(jié)果驗證了所提出模型的準確性。

4 結(jié)束語

針對機床數(shù)字化設(shè)計中典型結(jié)合面建模精度誤差較大的問題,筆者提出了一種采用虛擬材料等效結(jié)合面靜力學特性的方法,建立了考慮結(jié)合面特性影響的機床整機有限元靜力學分析模型,并完成了相應(yīng)的靜剛度測試實驗。

筆者基于粗糙表面接觸理論,推導(dǎo)了考慮摩擦力影響的固定結(jié)合面的法向剛度和切向剛度解析表達式,并確定了用于模擬結(jié)合面靜力學特性的虛擬材料的等效彈性模量和泊松比。

研究結(jié)果表明:

1)機床典型結(jié)合面等效后的虛擬材料彈性模量較基體材料小1到2個數(shù)量級。因此,對整機分析時,必須考慮結(jié)合面的影響;

2)與實驗測試結(jié)果對比,不考慮結(jié)合面影響的機床整機有限元靜力學分析模型的誤差在13%～17%之間,與實際偏差較大;

3)考慮了結(jié)合面特性影響后,得到機床整機有限元靜力學分析模型的誤差基本在5%以內(nèi),精度滿足工程設(shè)計的要求。

目前,筆者主要采用虛擬材料法研究了結(jié)合面的靜力學特性。后續(xù)工作中,筆者將在本研究的基礎(chǔ)上,將該模型應(yīng)用于結(jié)合面的動態(tài)特性以及熱-結(jié)構(gòu)耦合分析等方面。