楊寶鋒，賈少鋒，李 斌，陳 暉

(1.液體火箭發(fā)動機技術重點實驗室，陜西 西安 710100；2.西安航天動力研究所，陜西 西安 710100；3.航天推進技術研究院，陜西 西安 710100)

0 引言

液體火箭發(fā)動機渦輪泵通常采用環(huán)形密封來控制流體介質(zhì)的泄漏，以提高渦輪泵的工作效率[1-2]。然而密封在減小泄漏量的同時，還會產(chǎn)生流體激振力。該激振力一方面可引入一定的剛度阻尼，提高轉(zhuǎn)子系統(tǒng)的臨界轉(zhuǎn)速，控制轉(zhuǎn)子振動響應，甚至可取代軸承作為支承[3]；另一方面，其引入的交叉剛度對轉(zhuǎn)子系統(tǒng)穩(wěn)定性構(gòu)成威脅。

對于密封轉(zhuǎn)子動力特性的研究，早期普遍采用基于整體流動理論(Bulk-Flow Theory)的求解方法，其忽略了密封徑向速度變化情況，通過試驗結(jié)果建立了壁面剪切力與平均流速之間的關系，極大地簡化了求解過程[4]。該方法因其求解速度快，一定條件下可獲得較為可靠的結(jié)果，在工程上得到了廣泛的應用[5-10]。但由于整體流動模型是一種線性模型，只適用于小擾動情況，超出線性擾動范圍之外將不再適用。因此，隨著計算流體力學(CFD)技術的發(fā)展，利用三維CFD求解密封流場受到關注。相比傳統(tǒng)Bulk-Flow模型，CFD方法無需做過多假設，理論上是更為精確的求解方法，并且可用于研究轉(zhuǎn)子大擾動下的密封動力特性。密封的CFD求解方法主要有兩種：CFD準穩(wěn)態(tài)法和CFD瞬態(tài)法。CFD準穩(wěn)態(tài)法最早由Tam等[11]提出，該方法假定轉(zhuǎn)子繞密封中心軸線做圓形渦動，通過引入固定在轉(zhuǎn)子上的動參考系，將瞬態(tài)問題轉(zhuǎn)化為穩(wěn)態(tài)問題，從而避免了瞬態(tài)模擬的復雜耗時，是近年來在密封動力學研究中使用最為廣泛的求解方法[12-18]；CFD瞬態(tài)方法則直接在靜止坐標系下對密封流場進行模擬，對密封結(jié)構(gòu)及渦動軌跡等無限制，是研究密封的一種普適方法。但由于瞬態(tài)方法復雜耗時，并且在處理小間隙網(wǎng)格變形等方面存在諸多難點，因此在工程應用方面受到限制，目前的研究還集中在提高瞬態(tài)求解效率對結(jié)構(gòu)較為復雜的阻尼密封開展相關研究[19-22]。

對于實際工作的渦輪泵小間隙密封，由于安裝精度以及轉(zhuǎn)子彎曲等原因，通常會在偏心狀態(tài)下運行；此外，渦輪泵內(nèi)復雜流動引起的激振現(xiàn)象往往會使轉(zhuǎn)子擾動超出線性小擾動范圍，基于線性小擾動的密封動特性計算方法不再適用。然而目前對于渦輪泵轉(zhuǎn)子大偏心尤其是大擾動下密封動力學特性的研究還比較匱乏。

本文分別采用修正的Bulk-Flow模型和CFD準穩(wěn)態(tài)法對我國某型液體火箭發(fā)動機渦輪泵浮動環(huán)密封在大偏心以及大擾動下的動力學特性進行了研究。并且通過試驗數(shù)據(jù)對兩種方法的求解精度及可靠性進行了驗證。

1 Bulk-Flow理論模型

渦輪泵密封兩端通常承受較高的壓差，由于介質(zhì)粘度較小，導致其內(nèi)部流動表現(xiàn)為高度的湍流流動。圖1給出了密封幾何結(jié)構(gòu)及相應坐標系。

圖1 環(huán)形密封幾何結(jié)構(gòu)Fig.1 Geometric structure of annular seal

基于整體流動理論，建立密封內(nèi)的流體控制方程組[23-24]

(1)

其中

h=Cr(1+ecos(x/R))

式中:ρ為介質(zhì)密度；h為密封間隙；u為密封周向平均速度；w為軸向平均速度；p為流體壓力；us,ur分別為靜子、轉(zhuǎn)子表面速度；fs,fr分別為靜子、轉(zhuǎn)子表面阻力系數(shù)；γ為壁面摩擦系數(shù)；t為時間；R為轉(zhuǎn)子半徑；Ω為轉(zhuǎn)子表面角速度；Cr為密封半徑間隙；e為偏心率。

在獲得密封流場控制方程組后，基于線性攝動理論獲得零階以及一階攝動方程組，具體過程詳見文獻[25]。求解零階攝動方程組，并對密封出口軸向平均速度進行積分獲得密封的泄漏量

(2)

求解一階攝動方程組，對密封間隙流域壓力進行積分，基于最小二乘法線性擬合，即可獲得密封剛度、阻尼以及慣性系數(shù)

(3)

式中：ω為渦動轉(zhuǎn)速；K為剛度系數(shù)；C為阻尼系數(shù)；M為慣性系數(shù)；α,β=X,Y；ΓX=cosθ，ΓY=sinθ；當β=X時,ν=a，ζ=b；當β=Y時,ν=b，ζ=a。

本文采取文獻[26]提出的有限元方法對攝動方程進行求解，并對模型中幾個關鍵參數(shù)(如入口損失系數(shù)、預旋系數(shù)、湍流摩擦系數(shù)等)進行修正，以提高模型預測的準確性。

2 CFD準穩(wěn)態(tài)方法

2.1 求解思想

當轉(zhuǎn)子在繞靜偏心做小位移擾動時，密封流體激勵力可用下式表示

(4)

式中：Fx，F(xiàn)y分別為密封x，y方向的流體激勵力；Kxx,Kyy為主剛度系數(shù)；Kxy，Kyx為交叉剛度系數(shù)；Cxx，Cyy為主阻尼系數(shù)；Cxy，Cyx為交叉阻尼系數(shù)；Mxx，Myy為主慣性系數(shù)。

當轉(zhuǎn)子進行同心渦動時，密封動力特性系數(shù)中直接項相等，交叉項互為負值，此時，上式可轉(zhuǎn)換為

(5)

對于CFD準穩(wěn)態(tài)法，其假設轉(zhuǎn)子做同心圓形渦動，因此在轉(zhuǎn)子中心建立轉(zhuǎn)速等于渦動轉(zhuǎn)速的旋轉(zhuǎn)坐標系，即可將復雜的瞬態(tài)問題轉(zhuǎn)換為穩(wěn)態(tài)問題。此時作用在轉(zhuǎn)子上流體徑向力及切向力分別為

(6)

基于CFD仿真獲得不同渦動轉(zhuǎn)速下的Fr和Ft，代入式(6)進行擬合即可獲得相應的密封動特性系數(shù)。

2.2 物理模型及網(wǎng)格劃分

本文對兩種密封模型進行仿真。模型一為試驗條件下的密封結(jié)構(gòu)用于對兩種求解方法進行驗證，其幾何尺寸及運行參數(shù)如表1所示；模型二為某型火箭發(fā)動機渦輪泵浮動環(huán)密封，其幾何結(jié)構(gòu)及受力情況如圖2所示，密封半徑間隙為0.15 mm，進出口壓力分別為pin=23 MPa，pout=5 MPa；轉(zhuǎn)子轉(zhuǎn)速為18 000 rpm。

圖2 浮動環(huán)密封結(jié)構(gòu)示意圖Fig.2 Structural diagram of the floating ring seal

為考慮密封預旋效應以及入口節(jié)流降壓過程，對模型一入口遠場流域以及模型二密封入口前泄漏流域也進行了建模。采用ANSYS Workbench對兩種模型進行三維結(jié)構(gòu)化網(wǎng)格劃分，并對間隙壁面區(qū)域進行網(wǎng)格加密，兩種模型的網(wǎng)格劃分結(jié)果如圖3～圖4所示。

表1 試驗密封幾何尺寸及運行參數(shù)

圖3 密封模型一網(wǎng)格Fig.3 Grids of seal model 1

圖4 密封模型二網(wǎng)格Fig.4 Grids of seal model 2

采用了7套網(wǎng)格方案對兩種模型進行網(wǎng)格無關性驗證，間隙內(nèi)網(wǎng)格層數(shù)分別為4, 6, 8, 10, 12, 14, 16層，模型一對應網(wǎng)格數(shù)在15×104～200×104之間，模型二網(wǎng)格數(shù)在52×104～86×105之間。網(wǎng)格無關性計算結(jié)果表明當網(wǎng)格層數(shù)超過8層時，泄漏量結(jié)果基本穩(wěn)定，變化率不超過0.2%；而密封流體激勵力在網(wǎng)格層數(shù)超過12層時才逐漸穩(wěn)定。因此為了保證仿真結(jié)果的可靠性，選擇間隙網(wǎng)格層數(shù)為14層的網(wǎng)格方案進行研究，相應的模型一網(wǎng)格數(shù)為1 505×103，模型二網(wǎng)格數(shù)為5 604×103。

2.3 數(shù)值方法及邊界條件

采取ANSYS CFX對密封流場進行仿真，對于密封動特性的研究使用較多的湍流模型為標準k-ε模型或其改進模型，近年來SSTk-ω模型的應用逐漸增多[14-15]，均獲得了較好的計算結(jié)果。本文對標準k-ε，RNGk-ε以及SSTk-ω三種模型仿真結(jié)果與試驗進行對比，最終確定選取SSTk-ω模型進行后續(xù)仿真研究。利用本文數(shù)值方法與文獻[27]中基于Fluent軟件，Reliazablek-ε模型的三種仿真方法獲得的結(jié)果進行對比如表2所示?？梢钥闯觯疚幕贑FX軟件、SSTk-ω的準穩(wěn)態(tài)方法與文獻中的兩種準穩(wěn)態(tài)方法相比，獲得的結(jié)果更佳，其計算精度與瞬態(tài)方法精度一致。其中主剛度系數(shù)誤差較大，文中給出的原因是主剛度系數(shù)對密封上游狀態(tài)敏感性較大，而模型中未考慮上游狀態(tài)的影響。

表2 不同仿真結(jié)果對比

邊界條件的選取：模型一根據(jù)試驗測量結(jié)果給定入口靜壓以及出口靜壓邊界條件；模型二為了更好考慮葉輪旋轉(zhuǎn)引起的預旋效應，根據(jù)文獻[28]中渦輪氧泵全流場仿真結(jié)果給定速度入口以及壓力出口邊界條件。轉(zhuǎn)子壁面角速度設為(Ω-ω)，密封壁面角速度設為-ω。收斂殘差為1×10-7，同時監(jiān)控流量曲線以及流體激勵力作為收斂判定準則，當所有結(jié)果穩(wěn)定后即認為計算結(jié)果達到收斂。最后選取5種渦動頻率比ω/Ω= -1,-0.5,0,0.5,1計算密封動力學特性系數(shù)及泄漏量。

3 仿真方法驗證

密封動力學特性試驗在哈爾濱工業(yè)大學的高轉(zhuǎn)速間隙密封試驗臺[29]上進行。利用所獲得的試驗數(shù)據(jù)對修正的Bulk-Flow模型以及CFD準穩(wěn)態(tài)法計算結(jié)果進行驗證。

在Bulk-Flow模型求解中，利用CFD仿真結(jié)果對入口損失系數(shù)以及預旋系數(shù)進行修正；通過試驗結(jié)果并基于泄漏量相等的原則，對湍流摩擦系數(shù)進行修正。修正后的入口損失系數(shù)及預旋系數(shù)分別為0.5和0.22，湍流摩擦系數(shù)為3。仿真與試驗結(jié)果對比如表3所示。

由表3可以看出，兩種方法均能較好地預測密封動力特性系數(shù)，其中CFD方法具有更高的求解精度。修正的Bulk-Flow模型對密封各動特性系數(shù)的預測均偏高，其中主剛度系數(shù)預測誤差為10.32%，交叉剛度誤差為13.65%，主阻尼誤差為26.84%，交叉阻尼誤差為43.04%；CFD方法對主剛度系數(shù)預測偏低，而對其他系數(shù)預測偏高，其中主剛度與交叉剛度預測誤差僅為4.07%和3.47%，主阻尼誤差為19.63%，交叉阻尼誤差為44.30%。在泄漏量預測上，由于修正的Bulk-Flow模型是基于與試驗泄漏量相等原則來修正湍流摩擦系數(shù)，因此與試驗結(jié)果保持一致；此外，CFD準穩(wěn)態(tài)法能夠準確預測泄漏量，與試驗結(jié)果非常接近。因此后續(xù)研究可通過CFD求解泄漏量結(jié)果對湍流摩擦系數(shù)進行修正。

值得一提的是，本文試驗以及仿真獲得的交叉剛度顯著大于主剛度系數(shù)。對于低壓差浮動環(huán)密封，文獻[30]結(jié)果表明，當鎖定偏心率<0.35時，密封交叉剛度系數(shù)Kxy,Kyx均大于主剛度系數(shù)Kxx,Kyy，其中負交叉剛度系數(shù)Kyx遠大于主剛度系數(shù)Kxx,Kyy。而本文試驗臺密封鎖定于同心位置，因此出現(xiàn)交叉剛度大于主剛度這種情況。

表3 仿真與試驗結(jié)果對比

4 浮動環(huán)密封動力學特性

渦輪泵浮動環(huán)密封在穩(wěn)定運行時一般會處于鎖定狀態(tài)，其鎖定位置由于重力等影響往往會處于偏心位置。此外，由于浮動環(huán)密封間隙相對更小，在渦輪泵內(nèi)流體復雜激勵作用下轉(zhuǎn)子擾動往往會超出線性擾動范圍，傳統(tǒng)線性化方法預測的動特性系數(shù)將不再適用。因此，本章利用修正的Bulk-Flow模型對轉(zhuǎn)子偏心狀態(tài)下密封動特性進行研究，并通過CFD方法對轉(zhuǎn)子大擾動下密封所表現(xiàn)的非線性動力學特性進行研究。

4.1 浮動環(huán)密封動特性系數(shù)求解

假設浮動環(huán)鎖定位置與轉(zhuǎn)子同心，基于小擾動理論對密封動特性系數(shù)進行求解。在利用Bulk-Flow模型求解時，基于CFD仿真結(jié)果對模型經(jīng)驗參數(shù)進行修正，修正后的入口損失系數(shù)為0.58，預旋系數(shù)為0.4，湍流摩擦系數(shù)為3.5。

表4分別給出了利用Bulk-Flow以及CFD方法獲得的浮動環(huán)密封動特性系數(shù)?？梢钥闯?，兩種方法獲得的密封動特性系數(shù)十分接近，這一定程度上也證明了仿真求解的可靠性。對于高壓高轉(zhuǎn)速渦輪泵，采用小間隙浮動環(huán)密封引入的主剛度系數(shù)達到1×108N/m以上，與滾動軸承的支承剛度量級相當，因此在渦輪泵轉(zhuǎn)子動力學特性的研究中，密封影響不應忽略。此外，對于采用剛性轉(zhuǎn)子低溫氧泵來說，軸承阻尼較小可忽略，密封作為唯一可引入較大阻尼的部件，其主阻尼系數(shù)達到40 kN·s/m以上，通過合理的密封設計，可能會對渦輪泵轉(zhuǎn)子系統(tǒng)振動的抑制起到良好的效果。

表4 浮動環(huán)密封動特性系數(shù)

圖5給出CFD仿真獲得浮動環(huán)密封軸向截面壓力分布。可以看出，密封入口處由于流道截面迅速減小，流速增大，壓力迅速降低(見局部放大圖)，導致一部分能量損失。密封間隙內(nèi)部壓力隨著密封長度增加逐漸減小，在出口處達到最低值，符合密封的實際工作狀態(tài)。

圖5 密封模型壓力分布Fig.5 Pressure distribution of the seal model

4.2 大偏心下的密封動力學特性

實際運行過程中，密封通常會處于不同程度的靜偏心狀態(tài)。大偏心下，其對密封動力學影響將不可忽視。利用修正的Bulk-Flow模型對浮動環(huán)在靜偏心率0～0.8范圍下的動力學特性進行研究。圖6給出密封剛度系數(shù)隨靜偏心率的變化曲線?？梢钥闯?，當偏心率小于0.2時，剛度系數(shù)基本保持不變；當偏心率大于0.2時，剛度系數(shù)隨偏心率增加而增大，大偏心下增大幅度加快；主剛度系數(shù)Kxx較Kyy增加幅度更大，偏心率0.8時，Kxx增加了56.94%；交叉剛度kyx較kxy增加幅度更大，偏心率0.8時，kyx增加110.8%。

圖6 剛度系數(shù)隨偏心率變化Fig.6 Variation of stiffness coefficients with eccentricities

圖7給出密封阻尼系數(shù)隨靜偏心率的變化曲線。同樣在靜偏心率小于0.2時，阻尼系數(shù)變化較?。欢斊穆蚀笥?.2時，阻尼系數(shù)隨靜偏心率的增大顯著變化。其中主阻尼系數(shù)Cxx，Cyy均隨靜偏心率增大，且Cxx增加幅度更大，偏心率0.8時，增加了85.57%；而交叉阻尼系數(shù)cxy，cyx卻呈現(xiàn)不同的變化規(guī)律：cxy隨偏心率增大逐漸減小，偏心率0.8時減小48.90%；而cyx則隨偏心率增加而增大，且在大偏心時增大幅度迅速增加，當偏心率為0.8時，Cyx增大了190.29%；因此總的來說，密封的交叉阻尼隨著轉(zhuǎn)子靜偏心的增加體現(xiàn)出增大的趨勢。

圖7 阻尼系數(shù)隨偏心率變化Fig.7 Variation of damping coefficients with eccentricities

圖8給出慣性系數(shù)及泄漏量隨靜偏心率的變化曲線。同樣在靜偏心率小于0.2時，慣性系數(shù)及泄漏量基本保持不變；隨著靜偏心率的增大，泄漏量緩慢增加，而慣性系數(shù)則顯著增加。其中Mxx增大幅度更大，當偏心率為0.8時，Mxx增加了102.37%。但總體上由于浮動環(huán)長度較短，其引入的慣性系數(shù)并不是很大。

圖8 慣性系數(shù)及泄漏量隨偏心率變化Fig.8 Variation of inertia coefficient and leakage with eccentricities

4.3 大擾動下的密封動力學特性

在線性密封動力學研究中，基于小擾動理論，轉(zhuǎn)子的擾動半徑控制在0.1Cr以內(nèi)，以滿足線性假設，來獲取密封動特性系數(shù)。而對于高壓高轉(zhuǎn)速渦輪泵，流體誘發(fā)的結(jié)構(gòu)振動劇烈，加之浮動環(huán)密封間隙較小，因此轉(zhuǎn)子擾動將超出線性范圍，處于大擾動狀態(tài)。本節(jié)利用CFD方法對不同擾動偏心(偏心率范圍0.05～0.8)下密封間隙內(nèi)流場進行仿真，利用式(5)和式(6)擬合獲得大擾動下密封的等效動特性系數(shù)。仿真過程中，假設浮動環(huán)鎖定位置與轉(zhuǎn)子同心。

圖9給出剛度系數(shù)對擾動偏心率的變化曲線?？梢钥闯?，主剛度及交叉剛度系數(shù)均隨擾動偏心的增加呈現(xiàn)出復雜的非線性關系。在小偏心擾動下剛度系數(shù)保持一致，表明傳統(tǒng)線性化動力學系數(shù)假定的合理性。當擾動偏心增加時，主剛度系數(shù)首先迅速下降，在偏心率0.3時達到最低值，此后隨著偏心率繼續(xù)增加，主剛度系數(shù)一致處于增大趨勢；交叉剛度隨偏心率增加先緩慢下降，在偏心率大于0.3后開始增大，且在大擾動下增大趨勢迅速增加?？梢姶髷_動下密封引入的交叉剛度會對轉(zhuǎn)子系統(tǒng)的穩(wěn)定性構(gòu)成威脅。

圖9 剛度系數(shù)隨擾動量的變化Fig.9 Variation of stiffness coefficients with disturbances

圖10給出了阻尼系數(shù)隨擾動偏心率的變化曲線。同樣不同擾動偏心下阻尼系數(shù)呈現(xiàn)出復雜的非線性關系，且在小擾動下，阻尼系數(shù)基本保持不變。隨著擾動的增大，主阻尼呈現(xiàn)出先降低后增大的趨勢，當擾動偏心率為0.3時，主阻尼系數(shù)達到最低；而交叉阻尼隨擾動偏心的增加呈現(xiàn)出先增大后迅速降低的趨勢，當擾動偏心率為0.4時，交叉阻尼達到最大值。

圖10 阻尼系數(shù)隨擾動量的變化Fig.10 Variation of damping coefficients with disturbances

圖11給出不同擾動偏心率下密封慣性系數(shù)及泄漏量曲線?？梢钥闯觯瑪_動偏心對泄漏量沒有影響，但卻對慣性系數(shù)影響較為復雜。隨著擾動偏心的增加，慣性系數(shù)先降低，在偏心率0.3時達到最低，而后開始增大，在偏心率0.4時達到最大，最后又隨著擾動量的增大迅速減小。

圖11 慣性系數(shù)及泄漏量隨擾動量的變化Fig.11 Variation of inertia coefficient and leakage with disturbances

總的來說，轉(zhuǎn)子在大擾動下，密封動特性系數(shù)呈現(xiàn)出復雜的非線性關系。對于采用浮動環(huán)密封的渦輪泵，在其轉(zhuǎn)子系統(tǒng)振動特性的研究中密封與轉(zhuǎn)子振動位移之間的非線性耦合作用不應被忽視。

5 結(jié)論

本文采用修正的Bulk-Flow模型以及CFD準穩(wěn)態(tài)法對轉(zhuǎn)子大偏心以及大擾動下的渦輪泵浮動環(huán)密封的動力學特性進行了研究，得出結(jié)論如下：

1)本文采用的修正的Bulk-Flow模型及CFD準穩(wěn)態(tài)法均能較好地預測密封動特性系數(shù)，且CFD法具有更高地求解精度。其中Bulk-Flow模型對各動特性系數(shù)的預測均偏高，而CFD法對主剛度系數(shù)的預測偏低，其余系數(shù)的預測均偏高。

2)對于高壓高轉(zhuǎn)速渦輪泵，浮動環(huán)密封引入的主剛度系數(shù)達到108N/m量級，與滾動軸承支承剛度量級相當，因此在渦輪泵轉(zhuǎn)子系統(tǒng)動力學研究中，浮動環(huán)密封的作用不應被忽略。

3)轉(zhuǎn)子靜偏心量在一定程度上改變了密封動力學特性。總體來說，隨著靜偏心量增大，密封動力學系數(shù)均呈現(xiàn)出增大的趨勢，且越大偏心下增大的趨勢越快。

4)轉(zhuǎn)子大擾動對密封動力學特性影響顯著。不同擾動偏心下，密封各動特性系數(shù)呈現(xiàn)出復雜的非線性關系。在渦輪泵轉(zhuǎn)子系統(tǒng)振動特性研究中，應當考慮密封的非線性耦合作用。