動壓箔片軸承氣膜厚度分布特點

王文升，鄭越青,2，徐剛，徐慶元

(1.中國工程物理研究院 機械制造工藝研究所，四川 綿陽 621900；2.西安交通大學 能源與動力工程學院，西安 710049)

動壓箔片軸承是一種由平箔和波箔組合的空氣動壓軸承，它以平箔作為軸承表面，以波箔作為彈性支承[1]，采用空氣潤滑。由于在運行過程中采用氣體潤滑，箔片軸承具有摩擦力低、損耗小及壽命長等特點；箔片的彈性可以平衡轉子的不平衡性，具有很強的環境適應性；選擇性能優良的軸承材料和涂層材料的箔片軸承可以耐受高溫。箔片軸承的這些特點使其在諸多性能上超越了傳統支承[1-3]，可以作為高速和超高速轉子系統的支承。

箔片軸承特性仿真涉及流體力學與彈性力學，并且存在二者的耦合。目前最常用的方法是簡單彈性基礎模型[2-3]，該模型將平箔與波箔的組合簡化為彈性均勻的柔性支承體，并且認為柔性支承體的剛度完全取決于波箔。該模型假設支承體表面任一位置的變形僅與支承剛度及該處壓強有關，而與周圍結構變形無關，即柔性支承體表面各點是相對獨立的。該假設顯然不符合平箔是連續介質這一基本事實。

為解決該模型存在的缺陷，在簡單彈性基礎模型基礎上，細化了對平箔與波箔的數學描述，采用Euler梁單元描述平箔，簡單彈簧單元描述波箔。以此模型為基礎，對箔片軸承的氣膜厚度與壓力分布進行了分析。

1 數學模型

箔片軸承結構示意圖如圖1所示，其由構成軸承表面的平箔和提供彈性的波箔組成。平箔的一端與波箔固定，另一端自由。正常工作時，轉子沿自由端至固定端的方向旋轉。

1.1 壓力場控制方程

箔片軸承氣隙中的氣體壓力場采用Reynolds方程描述，為減小舍入誤差，Reynolds方程通常采用無量綱形式，即

(1)

式中：X為無量綱周向坐標；Y為無量綱軸向坐標；P,H分別為無量綱氣膜壓力和無量綱氣膜厚度；Λ為軸承數。無量綱數與軸承數的定義可參考文獻[2]。

(1)式為非線性偏微分方程，需采用Newton-Raphson方法進行求解。將連續的Reynolds方程以中心差分形式展開，得到各差分節點處的離散Reynolds方程(角標規則如圖2所示)

圖2 網格與節點角標規則

(2)

構造除邊界點外的各節點上的離散Reynolds方程組，并采用Newton-Raphson迭代法進行求解，獲得各節點壓力值。(2)式中P值在首次計算時需給定一組初始值，然后以迭代過程中獲得的P分布不斷更新。

1.2 氣膜厚度數學模型

氣膜厚度等于初始氣膜厚度與平箔變形量的疊加[4]

h0=c+ecos(θ-θ0)+u，

(3)

式中：c為軸承名義間隙；e為偏心量；θ為周向角度；θ0為軸承偏位角；u為平箔變形量。

當轉子旋轉時，空氣在黏性力作用下隨轉子在氣隙中運動，當進入收斂區時氣體被壓縮而使壓力提高。氣體壓力作用于平箔表面，一方面波箔發生彈性變形，另一方面平箔本身也會發生變形。為描述以上現象，對平箔與波箔建立不同的有限元模型，平箔采用Euler梁單元描述，而波箔采用簡單彈簧描述(圖3)。

圖3 波箔與平箔的有限元模型

平箔的Euler梁單元剛度矩陣為

(4)

式中：E為彈性模量；Iz為平箔截面矩；l為單元長度。

波箔的簡單彈簧單元剛度矩陣為

(5)

式中：k為波箔的剛度系數[4]；b為軸承寬度；l0為波拱寬度的一半；ν為泊松比；t為波箔片厚度。

2 計算方法

2.1 邊界條件

計算過程采用有限元與有限差分2種計算方法，因此計算過程中使用兩類邊界條件：一類是壓力邊界條件；另一類是位移邊界條件。

氣膜間隙區域的4個邊界，即軸承兩側以及箔片固定端與自由端，均與大氣接觸，因此環境壓力邊界條件為

(6)

簡單彈簧單元的固定端施加位移邊界為

u(n+1：n+imax)=0。

(7)

2.2 收斂標準

求解量包括氣膜壓力P、氣膜厚度H與偏位角θ0，它們需要滿足的收斂標準為

(8)

2.3 迭代過程

整個計算由3層迭代過程組成，最底層是Newton-Raphson迭代，用于求解Reynolds方程；中間層是氣彈耦合迭代，用于獲得滿足耦合條件的氣膜壓力分布與氣膜厚度；最上層為偏位角修正迭代，用于校正轉子姿態，以滿足受力平衡要求。詳細計算過程如圖4所示。

圖4 程序流程框圖

給定初始的軸承名義間隙與偏心量，假定偏位角為0，且平箔沒有變形，根據(3)式計算出初始氣膜厚度；然后以該初始氣膜厚度求解Reynolds方程，計算壓力分布；由壓力分布獲得平箔受力情況，將其帶入有限元模型計算平箔變形，從而獲得新的氣膜厚度分布；再以該氣膜厚度計算壓力分布，如此迭代，直至收斂。

獲得收斂解后，得到穩態下的氣膜壓力，然后通過對x和y方向進行積分，得到軸承承載力以及偏位角。承載力與偏位角計算方法見文獻[4]。由于初始假定偏位角為0，該解并不滿足受力平衡要求，因此以計算獲得的偏位角替代之前的偏位角，重復上述過程。持續迭代，直至2次求解獲得的偏位角誤差小于收斂標準。

3 結果與分析

箔片軸承的氣膜壓力與氣膜厚度分布分別如圖5和圖6所示，其轉速為30 000 r/min，偏心率為0.9，其他軸承和氣體參數見表1。由圖5和圖6可知，隨著氣膜厚度的減小，氣體被壓縮使壓力迅速增大。由于重力作用氣膜厚度減小區域位于軸承下方，因此壓力增大區域也位于軸承下方，從而形成向上的承載力。通過氣膜厚度減小段后，氣體進入氣膜厚度增大段，氣壓逐漸恢復到大氣壓，在氣膜厚度較大區域(軸承上部分)氣膜壓力基本等于大氣壓[5]。

表1 軸承和氣體參數

圖5 無量綱氣膜厚度分布

圖6 無量綱氣膜壓力分布

以上結果與簡單彈性基礎模型的計算結果基本一致[6-7]，不同之處在于該模型的高壓區出現小尺度空間上的壓力波動，這一現象是由波箔下陷變形引起的。從圖5和圖6可見，在氣膜厚度較小處，也就是壓力較大區域，平箔出現局部下陷，該分析結果與文獻[8]中的試驗結果基本一致。

試驗后的平箔表面存在摩擦形成的類似斑馬線狀磨痕(圖7)，這是由運轉過程中轉子與定子碰撞摩擦所致。由于平箔的下陷效應，下陷區域沒有發生接觸摩擦，從而使磨痕呈現出斑馬線狀的條紋。

圖7 試驗后箔片的斑馬線磨痕

4 結束語

建立了考慮平箔片下垂效應的動壓箔片軸承氣彈耦合模型，通過迭代方法耦合求解了箔片軸承的氣膜壓力分布與氣膜厚度分布。基于仿真結果，對箔片軸承的氣膜厚度分布與壓力分布特點進行了分析，為平箔下陷這一試驗現象提供了合理解釋。