吳 垚 郗文君 楊利花 張彩麗 曹巨江 王 哲

(1.陜西科技大學(xué)機(jī)電工程學(xué)院 陜西西安 710021；2.西安交通大學(xué)航天航空學(xué)院 陜西西安 710049；3.施耐德(西安)創(chuàng)新技術(shù)有限公司低壓事業(yè)部 陜西西安 710075；4.西安工程大學(xué)材料工程學(xué)院 陜西西安 710048)

隨著微型燃?xì)廨啓C(jī)、微電機(jī)、超精密機(jī)床主軸、諧振濾波器和慣性導(dǎo)航系統(tǒng)等具有高能量密度和高功率密度的微動(dòng)力機(jī)電系統(tǒng)(micro-electromechanical systems，MEMS)和微流體器件的快速發(fā)展，微氣體滑動(dòng)軸承因其結(jié)構(gòu)簡(jiǎn)單、高轉(zhuǎn)速、高精度、低功耗、低摩擦特性以及工作溫度范圍寬等優(yōu)點(diǎn)在微旋轉(zhuǎn)機(jī)械領(lǐng)域中被廣泛應(yīng)用[1-3]。傳統(tǒng)軸承設(shè)計(jì)中通常假設(shè)滑動(dòng)軸承軸瓦為剛性且潤(rùn)滑劑黏度保持不變。而微型軸承由于氣膜厚度和特征尺度較小，軸瓦彈性變形會(huì)顯著改變軸承表面幾何形狀甚至與初始軸承間隙達(dá)到相同數(shù)量級(jí)；此外，根據(jù)VEIJOLA和TUROWSKI[4]的有效黏度模型，可壓縮潤(rùn)滑劑黏度是Knudsen數(shù)的函數(shù),而Knudsen數(shù)是表征氣體稀薄程度的參數(shù)，定義為氣體分子平均自由程與微通道內(nèi)氣體流動(dòng)特征長(zhǎng)度的比值[5]。因此，剛性體和常黏度潤(rùn)滑劑假設(shè)在超薄氣膜潤(rùn)滑計(jì)算中將不再適用，非常有必要分析軸瓦彈性變形和黏度變化對(duì)稀薄氣體潤(rùn)滑性能的影響，為提高微氣體軸承-轉(zhuǎn)子系統(tǒng)的穩(wěn)定性提供理論依據(jù)。

近幾十年來(lái)，國(guó)內(nèi)外學(xué)者關(guān)于超薄氣膜潤(rùn)滑技術(shù)已開展了許多卓有成效的研究并取得了一些重要成果。其中，BURGDORFER[6]、HSIA和DOMOTO[7]、MITSUYA[8]、FUKUI和KANEKO[9-11]、BAHUKUDUMBI、BESKOK[12]分別提出了不同的泊肅葉流流量模型，推導(dǎo)了相應(yīng)的修正Reynolds方程，研究了氣體稀薄效應(yīng)對(duì)軸承性能的影響。由于在求解磁頭/磁盤的超薄氣膜潤(rùn)滑性能時(shí)軸承數(shù)很大且剪切流項(xiàng)含有壓力，黃平和牛榮軍等[13-14]提出了求解大軸承數(shù)氣體潤(rùn)滑Reynolds方程的有效有限差分算法。ZHANG等[15]利用有限體積法和原子力顯微鏡數(shù)據(jù)研究了隨機(jī)粗糙表面對(duì)氣體微型軸承流動(dòng)特性的影響，計(jì)算表明分形表面的粗糙指數(shù)越大，壁面氣體滑移速度變化越大。YANG等[16]通過(guò)13矩方程將滑移系數(shù)視為壁面適應(yīng)系數(shù)和Knudsen數(shù)的函數(shù)，推導(dǎo)出了新的滑移速度邊界條件。ZHANG等[17]耦合Greenwood-Williamson模型和修正Reynolds方程計(jì)算了不同槽深和粗糙峰高度標(biāo)準(zhǔn)差的氣體螺旋槽推力微型軸承的軸承力和粗糙面接觸力。關(guān)于軸承彈流潤(rùn)滑分析，在DOWSON和HIGGINSON[18]首次求解等溫條件下受重載彈性圓柱體彈流潤(rùn)滑問(wèn)題后，彈性流體動(dòng)力潤(rùn)滑理論得到了長(zhǎng)足的發(fā)展。近年來(lái)，YAGI和SUGIMURA[19]對(duì)一維固定滑動(dòng)面軸承的彈流動(dòng)壓潤(rùn)滑性能進(jìn)行了數(shù)值仿真，結(jié)果發(fā)現(xiàn)小于100 nm的彈性變形對(duì)薄膜動(dòng)壓潤(rùn)滑的承載力有明顯影響，潤(rùn)滑薄膜形狀有從入口出現(xiàn)的恒定薄膜厚度區(qū)、位于中心的略微增加區(qū)和出口周圍的收斂區(qū)3個(gè)區(qū)域。肖乾等人[20]應(yīng)用COMSOL Multiphysics軟件建立了考慮彈性變形和黏壓、黏溫效應(yīng)的止推軸承有限元模型，仿真結(jié)果表明載荷增大使承載瓦上產(chǎn)生較大變形且油膜壓力峰值向出口移動(dòng)。SUN等[21]研究了軸頸偏斜、表面粗糙度、黏壓效應(yīng)以及彈性變形對(duì)徑向滑動(dòng)軸承性能的影響，結(jié)果表明軸承表面彈性變形時(shí)，表面粗糙度對(duì)油潤(rùn)滑偏斜軸承有輕微影響且軸承偏心率越大，軸瓦變形對(duì)軸承特性影響越明顯。LINJAMAA等[22]考慮了軸承表面的彈性變形和熱變形，建立了含有柔性聚合物的動(dòng)靜壓徑向滑動(dòng)軸承參數(shù)化計(jì)算模型，該模型評(píng)價(jià)了包括軸頸偏斜等不同運(yùn)行條件下主要軸承性能參數(shù)。根據(jù)湍流潤(rùn)滑理論，CHETTI[23]推導(dǎo)了考慮湍流狀態(tài)、耦合應(yīng)力流體和彈性變形影響的修正Reynolds方程，得到了不同耦合應(yīng)力系數(shù)、彈性系數(shù)和Reynolds數(shù)下的承載力、偏位角、摩擦因數(shù)以及端泄量(side leakage)，結(jié)果表明耦合應(yīng)力流體提高了剛性和彈性滑動(dòng)軸承在層流和湍流狀態(tài)下的穩(wěn)態(tài)性能。

綜上所述，已有的文獻(xiàn)對(duì)微型氣體軸承彈性流體動(dòng)力潤(rùn)滑性能的研究較少。另外，由于微型動(dòng)壓軸承復(fù)雜的氣彈耦合行為，在求解微型軸承微觀彈流潤(rùn)滑問(wèn)題時(shí)，需要考慮氣膜厚度、氣膜壓力、稀薄氣體黏度、軸瓦彈性變形以及Knudsen數(shù)等參數(shù)間的相互作用。因此，本文作者研究了軸瓦彈性變形和氣體有效黏度對(duì)微氣體滑動(dòng)軸承靜動(dòng)態(tài)特性的影響；利用有限元法求解控制軸瓦位移場(chǎng)的三維彈性方程，推導(dǎo)了包含稀薄氣體泊肅葉流流量的廣義修正Reynolds方程；通過(guò)耦合迭代求解超薄氣膜潤(rùn)滑Reynolds方程和三維彈性方程中不斷修正黏度項(xiàng)和氣膜厚度項(xiàng)，得到不同軸承參數(shù)和彈性模量下的微型軸承承載力、摩擦因數(shù)以及動(dòng)態(tài)剛度和阻尼系數(shù)，并與不考慮氣體薄膜有效黏度的剛性軸承結(jié)果進(jìn)行對(duì)比。這對(duì)于保持高速微旋轉(zhuǎn)機(jī)械穩(wěn)定性具有重要意義。

1 薄膜氣體潤(rùn)滑控制方程

對(duì)圖1所示的微氣體滑動(dòng)軸承，忽略潤(rùn)滑膜中的慣性效應(yīng)和熱效應(yīng)且假設(shè)可壓縮氣體流動(dòng)為層流。由于軸頸材料的彈性模量遠(yuǎn)大于軸瓦材料，軸頸可看作剛性體，不考慮其彈性變形。則考慮稀薄效應(yīng)和軸瓦彈性變形的修正Reynolds方程可表示為

圖1 考慮軸瓦彈性變形的微氣體滑動(dòng)軸承結(jié)構(gòu)和三維模型

(1)

式中：量綱一氣膜壓力P=p/pa；量綱一氣膜厚度H=h/c；p為氣膜壓力；pa為環(huán)境氣體壓力；h為氣膜厚度；c為軸承半徑間隙；φ=x/R,λ=z/R分別是量綱一周向和軸向坐標(biāo)，R為軸頸半徑；Λeff=6μeffωR2/(pac2)為軸承數(shù)，μeff為氣體有效動(dòng)力黏度；ω為軸頸轉(zhuǎn)動(dòng)角速度；T為量綱一時(shí)間；Q為量綱一流量修正系數(shù)[24]。

Q=1+0.108 42Kn+9.359 3/Kn-1.174 68

(2)

式中：Knudsen數(shù)Kn=λ0/h，λ0為氣體分子平均自由程。

學(xué)者VEIJOLA和TUROWSKI[4]提出了氣體動(dòng)力黏度也受稀薄效應(yīng)的影響，即氣體有效黏度是Knudsen數(shù)的函數(shù)，氣體有效黏度的計(jì)算公式為

(3)

式中：μ為氣體潤(rùn)滑劑動(dòng)力黏度。

為了得到考慮軸瓦彈性變形和氣體有效黏度的微滑動(dòng)軸承靜動(dòng)態(tài)性能，通過(guò)數(shù)學(xué)變換PH=S,(PH)2=S2=Π，將方程(1)轉(zhuǎn)化為橢圓型偏微分方程形式求解氣膜壓力。假設(shè)軸頸在彈性微型軸承靜平衡位置(ε0,θ0)以擾動(dòng)頻率Ω作周期性小擾動(dòng)[25]，任意時(shí)刻的軸頸位置可表示為

(4)

式中：E0和Θ0是擾動(dòng)偏心率和擾動(dòng)偏位角定義在復(fù)數(shù)范圍內(nèi)的幅值；Ω=ν/ω為量綱一擾動(dòng)頻率，ν軸頸擾動(dòng)頻率；i為虛數(shù)單位。

在軸頸動(dòng)態(tài)小擾動(dòng)條件下，假設(shè)量綱一氣膜壓力P和氣膜厚度H包含靜動(dòng)態(tài)兩部分：

(5)

(6)

軸承在水平和豎直方向的氣膜合力為

(7)

微氣體軸承的偏位角、承載系數(shù)和摩擦因數(shù)的計(jì)算公式為

(8)

將式(5)代入式(1)可得到考慮稀薄效應(yīng)、軸瓦彈性變形以及氣體有效黏度的動(dòng)態(tài)修正Reynolds方程的一般形式：

(9)

為求解式(9)得到微氣體軸承的動(dòng)態(tài)剛度和阻尼系數(shù)，采用偏導(dǎo)數(shù)法定義：

(10)

式(6)和(9)分別對(duì)E0和Θ0求偏導(dǎo)，可得到關(guān)于PE、HE、Pθ和Hθ的偏微分方程，求得PE和Pθ后，按式(11)計(jì)算軸頸小擾動(dòng)時(shí)微氣體滑動(dòng)軸承的動(dòng)態(tài)剛度和阻尼系數(shù)。

(11)

將式(9)中的動(dòng)態(tài)系數(shù)進(jìn)行轉(zhuǎn)換，得到相應(yīng)笛卡爾坐標(biāo)系下的動(dòng)態(tài)剛度和阻尼系數(shù)為

(12)

通過(guò)MATLAB軟件中偏微分方程工具箱和有限元法編寫相關(guān)程序，采用松弛迭代法提高計(jì)算的精度和效率，收斂精度均為 0.000 1，微型氣體軸承彈流潤(rùn)滑的計(jì)算流程如圖2所示。

圖2 微型氣體軸承彈流潤(rùn)滑計(jì)算流程

2 軸瓦彈性變形和氣膜厚度變化

微型滑動(dòng)軸承在氣膜壓力作用下導(dǎo)致軸瓦表面發(fā)生彈性變形，而表面彈性變形直接改變氣膜厚度反過(guò)來(lái)也會(huì)影響氣膜壓力分布。由于軸瓦表面彈性變形量很小，因此可認(rèn)為變形和載荷呈線性關(guān)系，利用有限元法建立軸瓦柔度矩陣求解軸瓦表面彈性變形。

如圖3所示，將微氣體滑動(dòng)軸承沿圓周方向展開，由于軸瓦厚度遠(yuǎn)小于軸瓦展開后長(zhǎng)度，軸瓦模型近似為規(guī)則的長(zhǎng)方體，采用八節(jié)點(diǎn)六面體線性單元將模型離散；根據(jù)最小位能原理得到各單元體的節(jié)點(diǎn)位移和節(jié)點(diǎn)載荷間的關(guān)系，再將各單元的剛度矩陣進(jìn)行組裝，建立軸瓦模型結(jié)構(gòu)的整體剛度矩陣為

圖3 微氣體軸承的有限元模型

K=?VBTDBdxdydz

(13)

式中：

B=[B1B2B3……B8],

(i=1,2,3,...,8)

(14)

(15)

其中N為形函數(shù)。

對(duì)整體剛度矩陣求逆矩陣得到微型軸承的柔度矩陣C來(lái)確定氣膜壓力引起的軸瓦上表面各節(jié)點(diǎn)徑向變形。

δt=∑C·P=∑K-1·P

(16)

式中：δt為氣膜壓力矩陣P下微型軸承內(nèi)表面所有節(jié)點(diǎn)的徑向位移矩陣。

通過(guò)計(jì)算微型軸承表面的彈性變形對(duì)氣膜厚度進(jìn)行修正

H=H0+δt=1+εcos(φ-θ)+δt

(17)

在求解修正Reynolds方程得到氣膜壓力的每次循環(huán)迭代中計(jì)入有限元法求解軸瓦彈性變形和氣體有效黏度的影響，最終同時(shí)達(dá)到收斂，得到考慮氣體有效黏度的彈性微氣體滑動(dòng)軸承的靜動(dòng)態(tài)性能。

3 結(jié)果及討論

文中對(duì)考慮和不考慮氣體有效黏度的剛性和彈性微氣體滑動(dòng)軸承計(jì)算結(jié)果進(jìn)行了分析和討論。微型軸承的幾何參數(shù)：c=1 μm，R=1 mm，軸承寬度B=200 μm，pa=1.033×105Pa，長(zhǎng)徑比B/D=0.1。

為驗(yàn)證當(dāng)前模型的準(zhǔn)確性，如圖4所示，將文中計(jì)算結(jié)果與文獻(xiàn)[26]給出的軸承中間截面量綱一氣膜壓力分布進(jìn)行對(duì)比，可以看出2種方法得到的仿真結(jié)果非常接近，誤差小于3%。

圖4 文中數(shù)值結(jié)果與文獻(xiàn)[26]量綱一氣膜壓力分布對(duì)比(R=2.0 mm,B=0.4 mm,pa=1.033×105 Pa,Λ=2.4)

圖5所示為分別計(jì)入軸瓦彈性變形和同時(shí)計(jì)入氣體有效黏度條件下微型軸承量綱一承載力CL隨軸承數(shù)Λ的變化規(guī)律?？梢钥闯觯瑲怏w稀薄效應(yīng)顯著降低軸承承載力，微氣體軸承的承載力隨軸承數(shù)的增加近似線性增加，這與潤(rùn)滑表面的高軸頸轉(zhuǎn)速和氣體的可壓縮性有關(guān)。而CL隨軸瓦彈性模量降低而增加，這是因?yàn)闅怏w彈流潤(rùn)滑狀態(tài)在一定程度上削弱了氣體稀薄效應(yīng)的影響。不考慮氣體有效黏度的彈性微型軸承量綱一承載力大于考慮氣體有效黏度的情況。軸瓦彈性變形和軸承數(shù)越大CL增加越明顯。

圖5 軸承數(shù)、氣體有效黏度和軸瓦彈性模量對(duì)承載力的影響(ε=0.7,υ=0.3)

對(duì)于不同柔性的微型軸承，其軸頸表面摩擦因數(shù)隨軸承數(shù)的變化曲線如圖6所示。由于微氣體滑動(dòng)軸承的摩擦來(lái)自稀薄氣流的黏性剪切力，軸承數(shù)越高，即軸頸轉(zhuǎn)速越高，黏性剪切力則越大，因此軸承數(shù)對(duì)摩擦因數(shù)的影響就越明顯。隨著軸瓦彈性模量E降低，微型軸承柔性增加，軸頸表面摩擦因數(shù)增大。這是由于氣膜厚度隨彈性模量的減小而增大，更多的氣流進(jìn)入超薄潤(rùn)滑間隙導(dǎo)致軸頸表面摩擦因數(shù)增加。此外，考慮氣體潤(rùn)滑劑的有效黏度時(shí)，微型軸承摩擦因數(shù)低于氣體黏度不變時(shí)的結(jié)果。

圖6 軸承數(shù)、氣體有效黏度和軸瓦彈性模量對(duì)摩擦因數(shù)的影響(ε=0.7,υ=0.3)

圖7描述了計(jì)入不同軸瓦柔性和氣體有效黏度的量綱一承載力CL隨偏心率ε的變化情況。結(jié)果表明，隨偏心率的增大，微型軸承的承載能力越來(lái)越強(qiáng)。當(dāng)偏心率ε<0.7時(shí)，隨軸瓦彈性模量從5 GPa增加到200 GPa，量綱一承載力CL逐漸增大。而在ε較高時(shí)，CL隨E的增大呈現(xiàn)出相反的變化趨勢(shì)，這是因?yàn)檩S瓦彈性變形增加了氣膜承載區(qū)域從而彈性微型軸承的承載能力增強(qiáng)?？紤]稀薄氣體有效黏度的影響，隨著ε的增加進(jìn)一步降低了彈性微型軸承的承載力CL。

圖7 偏心率、氣體有效黏度和軸瓦彈性模量對(duì)承載力的影響(Λ=30,υ=0.3)

圖8和9分別示出了在不同軸瓦彈性模量E和量綱一擾動(dòng)頻率Ω條件下，考慮和不考慮稀薄氣體有效黏度的微型軸承動(dòng)態(tài)剛度和阻尼系數(shù)計(jì)算結(jié)果?？梢钥闯?，剛性和彈性2種微型軸承情況下，隨著Ω的增加，正交剛度系數(shù)Kxx和Kyy均增大，交叉剛度系數(shù)Kyx減小而Kxy先增大后減小。由于潤(rùn)滑薄膜主要在垂直方向支承轉(zhuǎn)子重力，因此Kyy略大于Kxx。與剛性微型軸承相比，彈性微型軸承中軸瓦彈性模量的降低使動(dòng)態(tài)剛度系數(shù)有所增加且隨著軸承柔性上升其剛度系數(shù)增加幅度更顯著。這是由于軸瓦彈性變形帶來(lái)氣膜厚度增加導(dǎo)致流過(guò)軸承間隙內(nèi)的氣流增加。微氣體滑動(dòng)軸承的動(dòng)態(tài)阻尼系數(shù)除了交叉阻尼Dyx均隨擾動(dòng)頻率的增加而減小。當(dāng)彈性模量減小，彈性微型軸承的直接阻尼系數(shù)小于剛性微型軸承。與恒定黏度氣體潤(rùn)滑相比，當(dāng)擾動(dòng)頻率Ω>2時(shí)，氣體有效黏度效應(yīng)降低了微型軸承的直接剛度，同時(shí)增加了其直接阻尼系數(shù)，剛性和彈性微型軸承的交叉耦合阻尼系數(shù)的絕對(duì)值均隨Ω的增大而趨于0。

圖8 擾動(dòng)頻率、氣體有效黏度和軸瓦彈性模量對(duì)剛度系數(shù)的影響(ε=0.7,Λ=80,υ=0.3)

圖9 擾動(dòng)頻率、氣體有效黏度和軸瓦彈性模量對(duì)阻尼系數(shù)的影響(ε=0.7,Λ=80,υ=0.3)

圖10和11分別給出了當(dāng)擾動(dòng)頻率Ω=3.5和偏心率ε=0.7時(shí)，計(jì)入不同軸瓦彈性模量和氣體有效黏度效應(yīng)的微氣體軸承動(dòng)態(tài)剛度和阻尼系數(shù)隨軸承數(shù)Λ的變化曲線?？梢钥闯觯械膭?dòng)態(tài)剛度系數(shù)均隨著軸承數(shù)的增加而增加，而動(dòng)態(tài)阻尼系數(shù)除了交叉阻尼Dyx均呈現(xiàn)出先增加后減小的趨勢(shì)。與剛性微氣體軸承相比，當(dāng)Λ>30時(shí)，剛度系數(shù)隨軸瓦彈性模量減小而顯著增大，相反直接阻尼系數(shù)變小。這是因?yàn)檩^高的軸承數(shù)，即高軸頸轉(zhuǎn)速增強(qiáng)了微氣體滑動(dòng)軸承中的動(dòng)壓效應(yīng)，有效減少了超薄軸承間隙中的能量耗散。另外，氣體潤(rùn)滑劑的有效黏度效應(yīng)顯著降低了微型軸承的直接剛度系數(shù)，并增加了微型軸承的直接阻尼系數(shù)。這表明如果應(yīng)用恒定氣體黏度計(jì)算超薄氣膜潤(rùn)滑特性，會(huì)導(dǎo)致計(jì)算的動(dòng)態(tài)阻尼系數(shù)偏小。

圖10 軸承數(shù)、氣體有效黏度和軸瓦彈性模量對(duì)剛度系數(shù)的影響(Ω=3.5,ε=0.7,υ=0.3)

圖11 軸承數(shù)、氣體有效黏度和軸瓦彈性模量對(duì)阻尼系數(shù)的影響(Ω=3.5,ε=0.7,υ=0.3)

4 結(jié)論

基于VEIJOLA提出的稀薄氣體有效黏度模型，利用偏導(dǎo)數(shù)法、松弛迭代法和有限元法對(duì)修正的Reynolds方程和軸瓦彈性方程進(jìn)行了耦合數(shù)值求解。結(jié)果表明稀薄效應(yīng)、軸瓦彈性變形以及氣體有效黏度對(duì)微氣體滑動(dòng)軸承潤(rùn)滑特性有顯著影響。一些重要規(guī)律總結(jié)如下：

(1)氣體稀薄效應(yīng)顯著降低了微型軸承的承載能力。量綱一承載力和摩擦因數(shù)隨軸承數(shù)的增加而增加。由于氣體彈流潤(rùn)滑狀態(tài)一定程度降低了稀薄效應(yīng)的影響，軸瓦柔性的增加可提高軸承承載力和摩擦因數(shù)。相比恒定氣體黏度潤(rùn)滑的剛性微型軸承，考慮氣體有效黏度的彈性微型軸承承載力和摩擦因數(shù)較小。

(2)微氣體軸承的量綱一承載力隨偏心率增加越來(lái)越高。當(dāng)偏心率較小時(shí)，隨著軸瓦彈性模量的增加，軸承承載力增大。而在大偏心率情況下，量綱一承載力隨微型軸承柔性增加略有增加。

(3)正交剛度系數(shù)Kxx和Kyy隨擾動(dòng)頻率的增加而增加，而動(dòng)態(tài)阻尼系數(shù)除交叉阻尼系數(shù)Dyx均隨擾動(dòng)頻率的增加而減小。微氣體滑動(dòng)軸承的動(dòng)態(tài)剛度系數(shù)均隨著軸承數(shù)的增加而增加，而相應(yīng)的動(dòng)態(tài)阻尼系數(shù)除了交叉阻尼Dyx均呈現(xiàn)出先增加后減小的趨勢(shì)。相比剛性微氣體軸承，軸瓦柔性的提高顯著增加了動(dòng)態(tài)剛度系數(shù)且降低了微型軸承的直接阻尼系數(shù)。稀薄氣體的有效黏度效應(yīng)降低了微型軸承的直接剛度，同時(shí)增加了其直接阻尼系數(shù)。