撓性支承可傾瓦軸承動力特性研究

楊期江， 李偉光， 趙學智， 郭明軍

(1. 廣州航海學院 輪機工程學院，廣州 510725; 2. 華南理工大學 機械與汽車工程學院，廣州 510640)

傳統機械支點可傾瓦軸承通過軸瓦搖動或滑動，實現軸承剛度與阻尼系數低交叉耦合。撓性支承可傾瓦軸承(Flexure Pivot Tilt Pad Bearing，FPB)通過撓性支承(如圖1所示)，實現軸瓦支承點的彎扭耦合，達到相同的低交叉耦合和高穩定性，同時消除了傳統機械支點的磨損和高接觸應力。

圖1 傳統機械支點與撓性支點Fig.1 Traditional mechanical pivot and flexible pivot

軸承間隙是滑動軸承設計的重要參數，影響軸承轉子穩定性，臨界轉速和動態響應。撓性支承可傾瓦軸承采用電火花技術加工軸瓦，瓦塊與殼體之間通過撓性梁連接，使其成為一個整體。這種設計消除了瓦塊與殼體的裝配誤差(如圖2所示)，支點疲勞，瓦塊卸載的顫振問題等傳統可傾瓦軸承的缺陷。

圖2 傳統機械支承與撓性支承可傾瓦軸承公差帶Fig.2 Traditional mechanical pivot and flexible pivottilting pad bearing tolerance

國外Armentrout等[1]推出了全新的一體式可傾瓦軸承的設計，即撓性支點可傾瓦軸承(FPB)，相對與傳統可傾瓦滑動軸承(Tilting Pad Journal Bearing, TPJB)，它在轉子動力學方面更具有優勢。利用電火花加工(Electrical Discharge Machining, EDM)工藝整體線切割出 FPB撓性支點、軸瓦基體和瓦座。Kepple等[2]報道第一個FPB已安裝于1994年，討論了可傾瓦軸承的缺點，如支點磨損，設計，生產，安裝和維護的復雜性，并將傳統可傾瓦軸承與新的軸承設計“撓性支承可傾瓦軸承(FPB)”進行比較，提供了FPB應用在高速渦輪壓縮機上所測得的軸承-轉子穩定性試驗數據，發現撓性支承與傳統可傾瓦軸承在穩定性方面基本一致。Zeidan[3]討論了存在于所有的旋轉機械中油膜軸承交叉耦合特性及其對轉子系統穩定性的影響。接著介紹了之前所廣泛使用可傾瓦軸承的一些問題和缺點，并討論用FPB來解決這些問題。Chen等[4]在高速一體齒輪壓縮機將其與傳統的軸承進行比較，并闡述了優化設計之后FPB的效果，即降低了摩擦損失并具有更寬的工作溫度范圍。De Choudhury等[5]研究了在高速離心壓縮機上，具有一個新的徑向和推力軸承撓性支承5瓦可傾瓦軸承和推力軸承系統中的效果。他們的測試顯示出較低的油溫上升，較少的摩擦功率損失。目前該領域發展的方向：①與擠壓油膜阻尼器串聯組成全新的阻尼減振滑動軸承[6-7]；②將該類型結構的軸承應用于無油透平機械-如微型燃氣輪機，改用空氣作為潤滑劑，大大拓寬了其工作轉速與溫度[8]；③引入流體靜壓，發展成為動靜壓混合潤滑軸承[9-10]。

國內Feng等[11]提出了一種撓性支撐可傾瓦與金屬橡膠串聯的空氣軸承，建立了理論計算模型，并對其阻尼特性進行了研究。但目前國內在這一領域處于起步時期，本文從撓性支承可傾瓦軸承的瓦塊與軸頸的動力學特性出發，研究撓性支承可傾瓦軸承頻率縮減動力學模型，提出一種基于PDE工具箱快速求解軸承靜態以及頻率擾動壓力Reynolds方程，考慮可傾瓦軸承的瓦塊自由度，采用Newton-Raphson 迭代法可計算得到給定載荷和轉速工況下的軸承的靜平衡位置的數值計算方法。并與國外的最新試驗數據進行對比分析，驗證本文的動力學模型及仿真計算方法。為該類型軸承在透平渦輪上的應用做好理論工作準備。

1 流體動力潤滑的基本方程

潤滑油膜的慣性力，以及熱效應影響油膜動力及減振特性。本文引入帶慣性項的不可壓縮非定常工況Reynolds 方程[12]

(1)

滑動軸承油膜二維偏微分能量方程為[13]

(2)

式中：x為軸承周向坐標；y為徑向坐標(膜厚方向)；z為軸承間隙軸向坐標；h為油膜厚度；μ為滑油動力黏度， Pa·s；P為油膜壓力，Pa；t為時間，s；Θ為周向位置角。

2 撓性支承可傾瓦軸承動力學模型

可傾瓦軸承幾何及坐標系如圖3所示，x方向為軸瓦圓周方向，z方向為軸線方向。 全局坐標系(X,Y,Z)以軸承中心為原點，eX,eY為軸心在X,Y方向上的位移。 瓦塊位置角θ=x/R。

圖3 撓性支承瓦塊坐標關系及支點運動位移Fig.3 Flexure pivot tilting pad geometry and pivot displacement

圖3中Θl為軸瓦起始邊角度，Θt為軸瓦終止邊角度，Θp為軸瓦支點位置角； 撓性支承可傾瓦軸承支承位置具有柔性，軸頸在旋轉過程中，瓦塊在支點處具有擺角(δk)， 徑向(ξk)以及切向(ηk)位移。

2.1 瓦塊運動方程

所有瓦塊油膜力的矢量和必須與軸頸外部載荷(WX,WY)平衡：

(3)

第kth塊瓦運動方程為

(4)

(5)

2.2 油膜力與力矩

油膜力與力矩可分為靜態與動態部分

(6)

其中靜態平衡位置(0階)油膜力與力矩為

(7)

式中： {Zαβ}為油膜力動剛度系數， 油膜力動剛度系數是可以根據軸心平衡位置(eX0,eY0)得到

(8)

式中：hx=cosθ，hy=sinθ。

2.3 擾動分析

當軸頸中心在靜態平衡(eX0,eY0)處做頻率為ω， 幅值為(ΔeX,ΔeY)的小擾動。瓦塊也做同頻小擾動。因此軸頸、瓦塊的擾動方程可表示為：

eX(t)=eX0+ΔeXeiωt，eY(t)=eY0+ΔeYeiωt，

式中：k表示瓦塊數， 即為：k=1, 2, …,Npad軸瓦油膜厚度與油膜壓力分布都由靜態域(0階)與擾動域(1階)組成，即為

(9)

式中： Δhk=ΔeXcosθ+ΔeYsinθ+Δξkcos(θ-Θp)+(Δηk-RΔδk)sin(θ-Θ), ΔPk=PXΔeX+PYΔeY+PξΔξk+PηΔηk+PδΔδk。

2.4 縮減動力系數

可傾瓦軸承動力系數取決于軸頸與軸瓦平衡位置，以及特定的激振頻率(ω)，通常認為同步激振頻率(ω=Ω)，或者不同步激振頻率(ω<Ω)。軸承的總自由度為2+3Npads即軸頸中心位移沿{X,Y}軸線，為2自由度，另一方面，每塊瓦有1個轉動自由度和2個徑向及切向位移自由度。由此可知，可傾瓦軸承包含了軸頸及軸瓦的耦合運動。

假設軸瓦與軸頸做頻率為(ω)的簡諧運動，將式(6)代入式(3)、式(4)中得到頻率縮減的動剛度系數。

(10)

3 支點剛度

對于撓性支承軸瓦，將撓性支承軸瓦視作一端固定，一端自由的懸臂梁，而軸瓦本身被視為集中慣量[14]，具體如圖4所示。

圖4 撓性支點簡化懸臂梁Fig.4 Cantilever bearm model of a tilting pad with flexural web

(11)

(12)

式中：mp是瓦塊質量;ro,p和rr,p是內、外瓦塊半徑；ρ為軸承材料的密度；l為瓦塊軸向長度；Larc,p為墊的弧長。

忽略撓性支承的圓角，同理可得撓性支點的壓縮與切向剛度為[14]

(13)

式中：E為材料彈性模量；w為撓性支承寬度；L為軸瓦寬度；h為撓性支承高度。

4 平衡位置迭代算法

可傾瓦軸承動態特性的傳統差分法仿真思路多是在給定偏心率的條件下進行編程計算，需經過插值方法，曲線擬合出實際軸承工作轉速與載荷下的油膜力以及平衡位置，再調用擾動壓力Reynolds方程計算程序，進而積分得到實際軸承工作轉速與載荷下的動力系數。但此方法過于繁瑣。本文提出采用Newton-Raphson 迭代法可直接計算得到給定載荷和轉速工況下的軸承的靜平衡位置。聯合Matlab以及Comsol PDE工具包，設計程序結構和循環迭代算法。最終可獲得平衡狀態下靜態特性、縮減動力學系數。其靜態平衡位置迭代計算的基本流程如圖5所示，具體過程如下：

圖5 撓性支承可傾瓦軸承的瓦塊及軸頸中心靜平衡位置迭代算法Fig.5 Flow chart of an algorithm to find the equilibriumposition of the the pads and journal in the flexible-pivottilting pad bearing

(14)

瓦塊擺角平衡位置計算迭代公式為

(15)

(2) 采用的Newton-Raphson方法對作用于軸心位置的加載力進行平衡迭代計算，即

F0-W0

(16)

軸心位置的迭代計算公式為

(17)

(18)

徑向位移迭代公式為

(19)

5 算例分析

撓性支承可傾瓦軸承動力特性仿真計算流程具體如圖6所示。首先根據撓性支承可傾瓦的動力學模型，采用PDE數值計算非定常工況Reynolds方程及能量方程，根據圖5所示的平衡位置計算及本文第5部分的縮減動剛度系數計算方法，獲得給定工況下的軸心及瓦塊平衡位置與頻率縮減動剛度系數，再根據最小二乘法曲線擬合得到剛度Kij、阻尼Cij和附加質量Mij。

圖6 撓性支承可傾瓦軸承動力特性仿真計算流程Fig.6 The flexuret ilting pad bearing support dynamiccharacteristics simulation process

5.1 撓性支承可傾瓦軸承參數

本章采用文獻[16]中的四瓦撓性支承可傾瓦軸承進行仿真計算對比，具體軸承工況條件如表1所示。

表1 軸承工況條件

其中四瓦撓性支承可傾瓦軸承的參數具體如表2所示。

表2 軸承參數表

撓性支承可傾瓦軸承潤滑油牌號為ISO VG32-透平油，其溫黏方程為

μ(T)=0.045 3e-0.030 07(T-21.11)

(20)

5.2 撓性支承可傾瓦軸承動力參數識別

軸承動力參數實驗采用倒置式試驗臺， 軸承安裝在軸承座內部。圖7描繪了軸承油膜八個線性化剛度和阻尼矩陣模型。

圖7 線性化剛度與阻尼系數Fig.7 Linearized bearing coefficients

軸承座(質量為Ms)運動方程可以寫為

(21)

(22)

式(22)是軸承油膜線性力-位移模型， 動力系數包括剛度Kij、 阻尼Cij和附加質量Mij。 Δx、 Δy表示軸頸和軸承座之間的相對運動位移。

將式(21)代入式(22)中，并進行傅里葉變換可得

(23)

式(23)左邊為已知測量量， 式中Fx,Fy為輸入激勵力頻域值；Ax、Ay為軸承座加速度頻域值。因此4個軸承油膜動剛度未知量Hxx、Hxy、Hyx、Hyy可通過式(23)求出。

而軸承油膜動剛度與動力參數之間的關系為：

Hij=(Kij-Ω2Mij)+j(ΩCij)

Re(Hij)=Kij-Ω2Mij, Im(Hij)=ΩCij

式中： 下標i,j表示x,y方向，Ω是激振頻率。因此動力參數(剛度Kij、 阻尼Cij和附加質量Mij)可采用最小二乘法對軸承油膜動剛度的實部與虛部進行曲線擬合得出[17]。

5.3 油膜動剛度系數

根據圖6的撓性支承可傾瓦軸承動力特性仿真計算流程，采用考慮流體慣性的Reynolds方程，通過建立撓性支承可傾瓦軸承的頻率-縮減動力學模型，得到在不同轉速，載荷的工況條件下動剛度系數實部與虛部的動剛度系數計算結果，具體如圖8所示。

轉速為10 000 r/min，輕載(348.8 kPa)與重載(1 038.2 kPa)情況下，四瓦撓性支承可傾瓦軸承動剛度系數理論計算與試驗對比分析如圖7所示。隨著激振頻率的增加，動剛度實部直接項(Re(Hxx)與Re(Hyy))呈二次方下降趨勢，輕載時曲線下降斜率比重載時要小，說明重載時油膜剛性大。隨著激振頻率的增加，動剛度實部交叉項(Re(Hxy)與Re(Hyx))逐漸降低，但降幅較小。隨著激振頻率的增加，動剛度虛部直接項(Im(Hxx)與Im(Hyy))呈一次方上升趨勢。輕載時曲線上升斜率比重載時要大，說明重載時油膜阻尼大。隨著激振頻率的增加，動剛度虛部交叉項(Im(Hxy)與Im(Hyx))先上升后逐漸降低，但升降幅度十分微小。由此可知動剛度實部直接項(Re(Hxx)與Re(Hyy))理論計算與試驗分析比較接近，但動剛度實部交叉項(Re(Hxy)與Re(Hyx))理論計算與試驗分析結果有一定的差異。且動剛度虛部直接項(Im(Hxy)與Im(Hyx))理論計算與試驗分析理論計算與試驗分析結果的較吻合，但動剛度虛部交叉項(Im(Hxy)與Im(Hyx))理論計算與試驗分析結果的差異比較明顯。總體上10 000 r/min重載時，四瓦撓性支承可傾瓦軸承動剛度實部直接項(Re(Hxx)， Re(Hyy))與虛部直接項(Im(Hxx)， Im(Hyy))理論計算與試驗分析結果更加吻合，主要原因是重載時油膜厚度比輕載時油膜厚度小，油膜慣性力影響減小，因此理論計算更接近于試驗結果。

圖8 10 000 r/min輕載與重載情況下理論動剛度與試驗動剛度隨頻率變化曲線Fig.8 Theoretical and Experimental dynamic stiffnesses versus the excitation frequency at 10 000 r/min and light, high unit load

5.4 油膜剛度、阻尼與附加質量系數

應用本文算法(如圖6所示)，采用最小二乘法擬合得到在不同轉速，載荷的工況條件下油膜剛度、阻尼與附加質量系數。對不同轉速下剛度、阻尼、附加質量系數隨載荷變化的理論計算結果與文獻[16]中試驗數據曲線進行對比分析具體如圖9、圖10、圖11所示。

(1) 油膜剛度系數

直接剛度項Kxx與Kyy理論計算結果是相等的，隨著載荷(>3.1 kN)增大， 直接剛度項Kxx與Kyy理論計算與試驗結果都呈線性增大趨勢；當高轉速(10 000 r/min與12 000 r/min)與低載荷(<6.16 kN)下理論計算數值比試驗結果要低，高轉速(10 000 r/min與12 000 r/min)與重載荷(≥6.16 kN)下理論計算數值與試驗結果更加接近，說明高速重載下的直接剛度項Kxx、Kyy本文理論計算數值與試驗結果較吻合； 交叉剛度項Kxy與Kyx在上述工況載荷下理論計算數值與試驗分析結果也都較吻合。

圖9 不同轉速下剛度系數隨載荷變化Fig.9 Stiffness coefficients vs. static loads fordifferent rotor speeds

(2) 油膜阻尼系數

直接阻尼項Cxx與Cyy理論計算結果是相等的。隨著載荷(>3.1 kN)增大， 直接阻尼項Cxx與Cyy理論計算呈線性緩慢增大趨勢。高轉速(10 000 r/min與12 000 r/min)與低載荷下(<6.16 kN)理論計算數值比試驗結果要低，高轉速(10 000 r/min與12 000 r/min)與重載荷下(≥6.16 kN)理論計算數值與試驗結果更接近，如說明高速時， 重載下的直接阻尼項Cxx與Cyy本文理論計算數值與試驗結果較為吻合。 交叉剛度項Cxy與Cyx理論計算結果變化不明顯，且在上述工況條件下理論計算數值與試驗分析結果也都較吻合。

(3) 油膜附加質量系數

直接附加質量項Mxx與Myy理論計算結果是相等的。 直接附加質量項Mxx與Myy理論值隨著載荷的增大呈線性增大趨勢，而試驗結果在高轉速(10 000 r/min與12 000 r/min)與低載荷(<6.16 kN)下理論計算數值比試驗結果要低，高轉速(10 000 r/min與12 000 r/min)與重載荷下(≥6.16 kN)理論計算數值與試驗結果更接近，如工況(轉速12 000 r/min，輕載3.1 kN)下， 直接質量Mxx理論值為21.24 kg， 試驗值為25.5 kg， 誤差為16.7%； 工況(轉速12 000 r/min，重載荷9.24 kN)下， 直接阻尼Mxx理論值為29.01 kN·s/m，試驗值為27.94 kN·s/m，誤差為3.82% ，說明高速時，重載下的直接質量項Mxx與Myy本文理論計算數值與試驗結果較為吻合。 交叉質量項Mxy與Myx理論計算結果變化不明顯，且在上述工況條件下理論計算與試驗分析結果也都較吻合。

圖10 不同轉速下阻尼系數隨載荷變化Fig.10 Damping coefficients vs. static loads in fordifferent rotor speeds

圖11 不同轉速下質量系數隨載荷變化Fig.11 Added-mass coefficients vs. static loads in fordifferent rotor speeds

6 結 論

(1) 引入帶慣性項的Reynolds方程，考慮動壓潤滑油膜溫黏效應；首先在軸瓦與軸頸作同頻振動的假設的基礎上，根據瓦塊支點受力的平衡方程，靜平衡位置瓦塊運動的微分方程等，建立了柔性阻尼支承可傾瓦軸承頻率-縮減動力學模型。

(2) 提出一種撓性支承可傾瓦軸承的靜平衡位置迭代計算方法，即采用PDE工具箱有限元數值方法可快速求解軸承非定常工況Reynolds 方程及二維能量方程，將軸承工作載荷和轉速作為初始條件，并考慮可傾瓦軸承的瓦塊擺角自由度，通過Newton-Raphson 迭代法可直接計算得到給定載荷和轉速工況下的軸承的靜平衡位置。最終可獲得平衡狀態動力學系數。

(3) 采用本文提出的動力學模型及平衡位置迭代計算方法對撓性支承可傾瓦軸承進行了動態特性參數分析研究，并與相關文獻的理論及試驗仿真結果進行了對比，驗證了本文提出的動力學模型，平衡位置迭代算法。