考慮滑移效應的可傾瓦動壓氣體軸承承載性能分析

李青松，張珂，燕震雷

(上海應用技術(shù)大學 機械工程學院，上海 201418)

以氣體代替油脂作潤滑的動壓氣體軸承依靠其收斂間隙內(nèi)氣體流動產(chǎn)生的動壓效應來承載[1]?？蓛A瓦動壓氣體軸承由多個瓦塊構(gòu)成，各瓦塊可以自由傾斜擺動進行自我調(diào)節(jié)，具有轉(zhuǎn)速高，穩(wěn)定性好，能耗低，摩擦小等優(yōu)良性能， 在低溫、高速旋轉(zhuǎn)等機械中得到應用[1-5]；但其結(jié)構(gòu)復雜，制造困難[5-6]。

可傾瓦氣體軸承在國外廣泛研究運用而國內(nèi)研究并不成熟：文獻[7]提出了氣體軸承的理念；文獻[1]和文獻[8]首次研制出了以空氣作潤滑的氣體軸承；文獻[9-10]在等溫、考慮可壓縮性的條件下計算分析了無限長氣體軸承的性能；文獻[11]以可傾瓦動壓氣體軸承為研究對象計算分析了其靜態(tài)性能；文獻[12]基于矢量疊加原理求解分析氣體軸承靜態(tài)特性參數(shù)；文獻[13-14]提出以裝配法計算分析軸承的靜態(tài)性能；文獻[15]采用線性方法對可傾瓦氣體軸承的靜、動態(tài)特性進行了分析研究；文獻[16-17]利用間隙耦合法分析了可傾瓦氣體軸承靜、動態(tài)性能。上述研究大部分是先計算單一瓦塊承載力，再進行矢量合成的數(shù)值計算方法分析可傾瓦動壓氣體軸承的性能參數(shù)，但對氣體軸承性能分析大部分不考慮滑移邊界，僅利用連續(xù)介質(zhì)模型進行分析。

軸承高速旋轉(zhuǎn)時，轉(zhuǎn)、靜子兩者較小的氣膜間隙可導致兩者與氣膜潤滑作用面產(chǎn)生稀薄氣體滑移效應[18-19]。文獻[20]最早利用線性玻爾茲曼方程推導出克努森數(shù)Kn變化下的潤滑方程，證實滑移邊界影響動壓氣體軸承的性能。滑移邊界模型是通過引入速度滑移邊界條件修正雷諾方程的方法分析滑移效應下氣體軸承的各項性能，最常用的是1階滑移模型。文獻[21]考慮氣體滑移效應提出2階滑移修正邊界模型，文獻[22]建立了1.5 階滑移速度邊界模型。當0.0010.1，連續(xù)型模型和1階滑移模型對于間隙較小且Kn>0.1的氣流過渡區(qū)域的靜態(tài)性能預測結(jié)果準確性低，但當Kn>1時，1階和1.5階滑移模型的預測計算結(jié)果開始迅速偏離線性化玻爾茲曼方程解[23-25]。Frenkel-Kontorova(FK)模型被認為是最準確的模型，因為此模型比其他滑動模型更接近第一性原理，但是FK模型數(shù)據(jù)庫依賴于幾何形狀和流動參數(shù)，因而不便在實踐中應用[23,26]。

因此，文獻[27]從分子運動理論出發(fā)詳細分析了微小間隙下分子與分子壁面之間的碰撞運動，將碰撞分子劃分為臨界角θ和小于臨界角θ碰撞壁面，建立了基于任意克努森數(shù)Kn下的新滑移修正模型，該模型與傳統(tǒng)FK模型的計算結(jié)果相近；改進的滑移模型不僅便于應用且在高克努森數(shù)下所提供的預測結(jié)果優(yōu)于其他模型；文獻[28] 對微型動壓氣體軸承和三瓦可傾瓦動壓氣體軸承研究分析，比較了連續(xù)介質(zhì)模型、1階滑移修正模型和WU新滑移修正模型下軸承的性能參數(shù)，WU新滑移模型下計算的軸承承載力數(shù)值更接近玻爾茲曼方程解。

四瓦可傾瓦軸承主要安裝于第四代核電機組高溫氣冷堆氦風機中，其安全穩(wěn)定性能對冷卻堆芯系統(tǒng)正常運行，防止核電站發(fā)生事故極為重要。本文以四瓦可傾瓦動壓氣體軸承為研究對象，工作介質(zhì)為氦氣，考慮滑移影響，基于牛頓-拉弗森迭代法，利用有限差分離散方程求解連續(xù)介質(zhì)、1階滑移、WU新滑移模型修正后的主控雷諾方程；分析不同模型下軸承承載力變化以及滑移邊界對軸承性能的影響。

1 四瓦可傾瓦動壓氣體軸承結(jié)構(gòu)

四瓦可傾瓦動壓氣體軸承結(jié)構(gòu)示意圖如圖1所示，主要由1個軸頸和4個弧形瓦塊組成，各瓦塊在圓周方向均勻布置且可繞其支點自由擺動以適應工作位置。可傾瓦動壓氣體軸承的瓦塊支點位置角和放置位置因瓦塊數(shù)量的不同而不同；支點位置角度和擺放位置影響軸承氣膜間隙和承載能力[30-31]。圖中：R為軸頸半徑，θ為軸承偏位角，γ為支點位置角，φ為周向角坐標，e為軸頸中心相對支點圓圓心偏移距離，β為瓦塊的支點的位置角，z為軸向坐標，h為氣膜厚度。

圖1 四瓦可傾瓦動壓氣體軸承結(jié)構(gòu)示意圖

2 滑移邊界下修正理論模型

2.1 修正雷諾方程

微觀尺度上存在著一些不同于宏觀尺度的物理現(xiàn)象，此時，宏觀尺度的分析理論將不再適用；微觀尺度上發(fā)生的現(xiàn)象因其在宏觀尺度上的影響較小而被忽略，但微觀尺度下，這些現(xiàn)象成為影響力學性能的重要因素[32-33]。

可傾瓦動壓氣體軸承支點位置角、擺放位置因瓦塊數(shù)目不同而有差異，并對軸承氣膜間隙、承載力造成影響[16,29]。在軸承的分析計算中，很多參數(shù)通常以量綱一化的形式存在；通過參數(shù)的量綱一化，可以將物理問題轉(zhuǎn)化為數(shù)學問題且避免計算過程中不必要的因素，這樣可有效提高計算精度。

(1)

式中：P為量綱一的壓力；p為壓力；p0為環(huán)境壓力；H為量綱一的氣膜厚度；Cb為瓦塊內(nèi)圓弧半徑與軸頸半徑之差；φ，λ為量綱一的周向坐標。

滑移速度為氣體流動速度方向與壁面切向速度之差，其大小正比于分子平均自由程?？伺瓟?shù)Kn是描述微觀尺度物理現(xiàn)象時常用的重要參數(shù)[20]，恒溫下：Kn=C1/(ph)，即分子自由程σ=C1/p與h的比值，其中C1取6.5×10-4。若氣膜間隙較小時，Kn可能增大至氣流過渡區(qū)域(0.1

表1 3種模型的修正系數(shù)

連續(xù)介質(zhì)下不考慮滑移效應，此時無邊界滑移的氣體速度為

(2)

式中：δ為氣膜邊界。

笛卡爾坐標系下氣體潤滑軸承的N-S方程為

(3)

(4)

(5)

式中：u，v，w分別為x，y，z方向的氣體速度；μ為氣體動力黏度系數(shù)；ρ為氣體密度；t為時間。

簡化(3)，(4)和(5)式可得

(6)

根據(jù)(2)，(6)式得連續(xù)介質(zhì)(無滑移)、1階滑移、WU新滑移模型的修正氣體速度為

(7)

(8)

(9)

式中：uR為速度項系數(shù)。

從氣體速度的表達式可以看出，不同模型的速度滑移方向相同。與連續(xù)介質(zhì)模型相比，滑移模型引入邊界速度條件使計算結(jié)果向真實結(jié)果接近[33-34]，其中WU新滑移模型系數(shù)的值隨著Kn的變化而發(fā)生改變，比1階滑移模型能更好的適應任意Kn的真實情況，也更好地與實際物理模型吻合。

等溫條件下主控雷諾方程[35-36]量綱一的形式為

(10)

聯(lián)立(7)—(10)式并化簡整理，得到經(jīng)滑移模型修正的可壓縮氣體的定常量綱一的雷諾方程為

Λ=6μω/(p0Ψ2)，

Ψ=Cb/R，

(11)

式中：Λ為量綱一的軸承數(shù)；ω為軸頸旋轉(zhuǎn)角速度。

2.2 氣膜厚度方程

對軸承而言，各瓦塊可以繞其支點自由擺動，擺動角度也會影響瓦塊內(nèi)表面與軸頸之間的間隙[16,29]。因此，氣膜厚度和偏心距與瓦塊自由運動有關(guān)[36-37]。單瓦塊量綱一的氣膜厚度為

H=1-(m+u0)cos(β-φ)+εcos(φ-θ)+

χ0sin(β-φ)，

(12)

χ0=χ/Ψ，m=1-Cp/Cb，

式中：ε為偏心率(軸頸中心距支點圓中心的距離)；χ0為量綱一的瓦塊擺角；m為預載系數(shù)；Cp為支點間隙；u0為量綱一的靜態(tài)徑向位移，取u0=0。

2.3 邊界條件

軸承運行時，瓦塊邊界與氣體相通[38]，氣膜壓力在瓦塊內(nèi)表面與軸頸表面間的量綱一化邊界條件為

P(φ=φ0,λ) =1,P(φ=φ1,λ) =1，

P(φ,λ=l/(2R)) =1,P(φ,λ=-l/(2R)) =1，

(13)

式中：l為瓦塊長度；φ0為瓦塊前緣進氣端的角度；φ1為瓦塊后緣出氣端的角度。

3 數(shù)值求解

3.1 牛頓迭代法轉(zhuǎn)化雷諾方程

數(shù)值算法的選擇影響求解精度和計算速度[39]。進行數(shù)值模擬計算前對單一瓦塊展開并進行離散網(wǎng)格，因為氣膜厚度的方向數(shù)量級非常大，相對來說其他方向上的數(shù)量級可以忽略不記，因此即使瓦塊展開成平面，也不會對最后的計算結(jié)果產(chǎn)生明顯影響?；谝陨峡紤]，忽略瓦塊表面曲率對計算結(jié)果的影響，將各瓦塊近似展開成矩形平面，并離散成40×40的網(wǎng)格。采用牛頓迭代法和有限差分法求解非線性雷諾方程。構(gòu)造量綱一的氣膜壓力P的函數(shù)F(P)為

(14)

由牛頓迭代法得

F′(Pn)δ=-F(Pn)；n=0,1,2…，

δ=Pn+1-Pn，

(15)

式中：δ為2次迭代值之差；n為迭代次數(shù)。

使用牛頓-萊布尼茨迭代法將P的非線性方程轉(zhuǎn)化為壓力增量ΔP的線性方程

(16)

3.2 基于有限差分法離散方程

基于有限差分法，利用中心差分格式，將(16)式展開(i=1～41，j=1～41)為

(17)

利用(17)式將(16)式整理為有規(guī)律的系數(shù)矩陣方程組，即

Ai,jPi,j+Bi,jPi+1,j+Ci,jPi-1,j+Di,jPi,j+1+

Ei,jPi,j-1=-Hi,j。

(18)

有限差分法離散得到各系數(shù)項為

根據(jù)已知條件迭代求解矩陣方程組以獲得單個瓦塊的量綱一的壓力分布。迭代過程中，為保證計算過程的穩(wěn)定收斂和加快收斂速度，一般采用松弛迭代法計算

(19)

為了確保迭代過程可終止，必須保證每次迭代計算的壓力分布滿足計算精度。按照收斂準則，當計算結(jié)果達到精度要求后，程序?qū)⑻鲅h(huán)，從而終止迭代計算，允許相對誤差eeps=1×10-4，收斂準則為

(20)

3.3 求解流程

3.3.1 計算瓦塊擺角

瓦塊可以自由擺動，計算時需要判斷擺角的位置。當瓦塊處于靜平衡位置時，瓦塊不再擺動，此時的角度就是計算軸承靜態(tài)特性時的擺角位置，瓦塊的擺動角度與瓦塊的支點力矩M、氣膜壓力分布和量綱一的力矩滿足平衡關(guān)系，即

(21)

當計算結(jié)果不滿足瓦塊支點力矩平衡時，使用超大范圍牛頓公式迭代產(chǎn)生新的瓦塊擺角，迭代更新瓦塊擺角直至支點力矩M為0，迭代公式和收斂準則為

(22)

(23)

3.3.2 計算氣膜合力

所有瓦塊的合力滿足靜態(tài)平衡條件時，可得各單個瓦塊的氣膜壓力分布，通過矢量疊加的方法計算軸承的氣膜合力。軸承垂直和水平方向的量綱一的氣膜合力Fx，F(xiàn)y及量綱一的承載力FCL計算公式為

(24)

3.3.3 計算偏位角

壓力分布計算時需確定偏位角θ，偏位角迭代格式和收斂準則為

θn+1=θn-ζarctan(Fx/Fy) ，

(25)

Fx/Fy≤η。

4 結(jié)果分析

4.1 程序驗證

為驗證本文數(shù)值計算模型的準確性，以文獻[10]和文獻[40]的計算結(jié)果與本文3種模型的計算結(jié)果進行對比，數(shù)值計算流程如圖2所示，后文若無特殊說明，所得結(jié)果皆為量綱一的值。

圖2 數(shù)值計算流程

文獻[10,40]與本文3種模型計算結(jié)果的對比如圖3所示：將文獻[10]中的可傾瓦氣體軸承參數(shù)代入數(shù)值計算流程，利用連續(xù)介質(zhì)模型計算不同偏心率下的軸承承載力，并與文獻[10]中的計算結(jié)果相比較，可以看出不同偏心率下的誤差值均小于5.8%；將文獻[40]中的可傾瓦氣體軸承參數(shù)代入數(shù)值計算流程，分別利用連續(xù)介質(zhì)模型、1階滑移模型和WU新滑移模型計算不同偏心率下的軸承承載力，并與文獻[40]中的計算結(jié)果相比較，可以看出連續(xù)介質(zhì)模型與文獻[40]結(jié)果之間的誤差最小，3種模型的誤差值均小于7.2%，考慮滑移邊界條件下所計算出的可傾瓦氣體軸承的承載力相比文獻[40]和連續(xù)模型顯著降低。由上述結(jié)果可知，圖2程序計算流程具有一定的準確性。

圖3 本文3種模型與文獻的計算結(jié)果對比

4.2 計算結(jié)果

給定預載系數(shù)m=0.2，偏心率ε=0.5下的軸承參數(shù)見表2。

表2 四瓦可傾瓦動壓氣體軸承參數(shù)[41-42]

4.2.1 氣膜厚度計算結(jié)果

m=0.2，ε=0.5時，瓦塊的氣膜厚度如圖4所示：氣膜厚度沿軸頸旋轉(zhuǎn)方向逐漸減小。瓦塊氣膜厚度隨周向節(jié)點的變化如圖5所示：任意位置瓦塊2，3周向節(jié)點方向氣膜厚度小于瓦塊1，4；瓦塊1氣膜厚度最大，瓦塊3氣膜厚度最小。

(a)瓦塊1

圖5 瓦塊量綱一的氣膜厚度隨周向節(jié)點的變化

4.2.2 氣膜壓力和克努森數(shù)的計算結(jié)果

以瓦塊1為研究對象，考慮滑移效應，取m=0.2，ε=0.5，3種模型下瓦塊1氣膜壓力分布(左圖)和等高線(右圖)如圖6所示，連續(xù)介質(zhì)模型下瓦塊1的克努森數(shù)Kn分布如圖7所示?；菩雇邏K1的氣膜壓力明顯降低，且WU新滑移模型計算所得氣膜壓力分布小于1階滑移模型；在瓦塊1的邊界處，氣膜壓力等于大氣壓力，3種模型在邊界處計算的氣膜壓力都相同，接近出口位置時的氣膜壓力達到最大值；除了邊界處，連續(xù)介質(zhì)模型計算出的氣膜壓力大于滑移模型，WU新滑移模型計算的氣膜壓力值小于1階滑移模型。另外，氣膜壓力沿軸向方向先增大后減小，分析原因：瓦塊1擺角在氣體沖擊下沿軸承旋轉(zhuǎn)方向擺動，瓦塊內(nèi)表面與軸頸表面間隙逐漸減小，氣膜壓力逐漸增大，并在瓦塊后段接近后緣部分的最小間隙處達到最大；瓦塊后緣出氣端與氣體環(huán)境相同，后緣部分壓力急劇減小至和外界壓力相同；由于可傾瓦塊沿軸向兩端都與大氣環(huán)境相同，所以氣膜壓力沿軸向表現(xiàn)為先增大后減小。圖7中一部分潤滑氣膜處于氣流過渡區(qū)，另一部分潤滑氣膜則處于滑移區(qū)；WU新滑移模型在處理氣流過渡區(qū)的計算性能優(yōu)于連續(xù)介質(zhì)模型和1階滑移模型。

4.2.3 偏心率對量綱一的氣膜厚度、克努森數(shù)和承載力的影響

以瓦塊1為研究對象，在m=0.2，ε取值為0，0.2，0.4，0.5，0.6，0.8時，探討偏心率的變化對氣膜厚度、克努森數(shù)、不同修正模型下軸承承載力的影響。瓦塊1氣膜厚度隨偏心率的變化如圖8所示：隨偏心率增大，瓦塊氣膜厚度呈增大趨勢。瓦塊1的Kn值沿軸向節(jié)點隨偏心率的變化如圖9所示：偏心率越大，沿周向節(jié)點Kn值越大。不同模型下軸承承載力隨偏心率的變化如圖10所示：隨偏心率增大，承載力FCL增大；WU新滑移修正模型計算所得承載力結(jié)果最小，1階滑移模型次之；隨偏心率增大，3種模型計算結(jié)果數(shù)值偏差也增大。

圖8 瓦塊1量綱一的氣膜厚度隨偏心率的變化規(guī)律

圖9 瓦塊1克努森數(shù)隨偏心率的變化

圖10 軸承承載力隨偏心率的變化

因瓦塊關(guān)于支點對稱布置，所以瓦塊1，2承載力相同，瓦塊3，4承載力相同，計算不同模型下單個瓦塊承載力隨偏心率的變化規(guī)律如圖11所示：對于單個瓦塊，隨偏心率增大，瓦塊1，2承載力減小，瓦塊3，4承載力增大?？紤]滑移修正下計算所得承載力顯著下降，且WU新滑移模型計算所得軸承承載力小于1階滑移模型。

(a)瓦塊1,2

4.2.4 預載系數(shù)和量綱一的軸承數(shù)對軸承承載 力的影響

以瓦塊1為研究對象，探討預載系數(shù)m和量綱一的軸承數(shù)對不同模型下軸承承載力的影響。ε=0.5時，不同模型軸承承載力隨預載系數(shù)m的變化如圖12所示：對同一個模型，預載系數(shù)m增大，軸承承載力隨之增大；對相同預載系數(shù)，WU新滑移模型計算所得的承載力最小。

圖12 軸承承載力隨預載系數(shù)的變化

m=0.2，ε=0.5時，3種模型下承載力隨量綱一的軸承數(shù)的變化如圖13所示：隨軸承數(shù)增加，3個模型的軸承承載力均增加；相同軸承數(shù)時，考慮滑移效應的模型軸承承載力較??；隨著軸承數(shù)增加，3種模型的軸承承載力差值減小。這是因為量綱一的軸承數(shù)只與角速度ω有關(guān)，因此隨著角速度速的增加，軸承數(shù)增加，瓦塊內(nèi)表面與軸頸表面間隙產(chǎn)生的動壓效果增強，滑移效應對軸承的影響增大，氣膜壓力增大，軸承承載力隨之增大；量綱一的軸承數(shù)一定時，考慮滑移效應的影響，瓦塊內(nèi)表面與軸頸表面間的氣膜壓力減小，導致承載力減小。

圖13 軸承承載力隨量綱一的軸承數(shù)的變化

m=0.2，ε=0.5時，3種模型單個瓦塊承載力隨量綱一的軸承數(shù)的變化如圖14所示：隨軸承數(shù)增大，瓦塊1，2的承載力呈先增大后減小；瓦塊3，4的承載力隨軸承數(shù)的增大而增大?？紤]滑移效應的修正模型所計算承載力結(jié)果明顯較小，即WU新滑移模型在承載力計算上有優(yōu)勢。

(a)瓦塊1，2

5 結(jié)論

以四瓦可傾瓦動壓氣體軸承為研究對象，考慮滑移效應的影響并分析滑移模型對量綱一的氣膜壓力分布和承載力的影響。從氣體潤滑基本方程出發(fā)，利用1階滑移模型和WU新滑移模型的氣體速度推導出不同模型的定常量綱一化靜態(tài)雷諾方程，討論連續(xù)模型、1階滑移模型、WU新滑移模型在不同軸承氣膜厚度、壓力等條件下的承載力。結(jié)論如下：

1)單個瓦塊量綱一的氣膜厚度均沿軸頸旋轉(zhuǎn)方向逐漸減??；任意位置瓦塊2，3周向節(jié)點量綱一的氣膜厚度小于瓦塊1，4。

2)可傾瓦動壓氣體軸承的承載力與偏心率、預載系數(shù)和量綱一的軸承數(shù)都成正比；滑移修正使軸承承載力顯著下降，且WU新滑移模型計算的承載力小于1階滑移模型；隨著偏心率、預載系數(shù)的增大，滑移效應對承載力結(jié)果的影響也增大。

3)WU新滑移模型計算的承載力、氣膜壓力結(jié)果最小，1階滑移模型次之，連續(xù)模型最大；基于滑移效應客觀存在的事實，使用WU新滑移模型對微小間隙下可傾瓦動壓氣體軸承的靜態(tài)性能計算更為準確。