郭 嘉，吳 飛，何朝輝

（浙江省機電設計研究院水泵所，浙江杭州310021）

0 引言

高壓多級離心泵機組作為流體輸送的關鍵設備廣泛應用于石油、化工、電力和鋼鐵等行業(yè)，如用于高壓鍋爐給水，大型化肥工業(yè)用液氨泵及甲氨泵等。這類多級離心泵中，由于實現平衡壓力及密封的需要，大量存在如密封口環(huán)、級間密封、平衡鼓及平衡盤等環(huán)形密封，高壓多級離心泵環(huán)形密封示意圖如圖1所示。該類環(huán)形密封結構一方面影響泵內的級間泄漏量從而直接影響離心泵的效率，另一方面在密封間隙流道內的高壓流體在高速轉子的作用下，產生了較大的間隙流體力，對機組的動力學性能及動力學行為具有明顯影響，這一現象在上世紀60年代首先由Loma?kin 提出，被稱為“洛馬金效應”［1］。目前，準確計算由Lomakin效應產生的口環(huán)密封等效動力學性能并預估考慮該效應的“濕轉子”臨界轉速已成為高壓多級離心泵設計中的重要內容，在API610等國際泵行業(yè)標準中明確均有明確規(guī)定。

20世紀60年代，環(huán)形密封動力學特性的相關研究開始逐漸興起。Black 和Jenssen 等人［2-4］針對“短密封”首次提出了主剛度的數值求解方法，但這一方法在求解中對流動方程進行了較大簡化，在口環(huán)軸向長度與半徑之比（簡稱長徑比）大于0.3的口環(huán)計算結果誤差較大，適用范圍較小。1983年，Childs［5-6］在Black研究的基礎上，基于Hirs［7］的整體流動理論補充了“短密封”求解方法并提出了“有限長”求解方法，將口環(huán)的求解范圍擴大到長徑比小于等于0.6。“有限長”求解方法考慮了慣性項及入口旋流等因素，求得結果與實驗結果對比較好，在工程設計領域具有廣泛的應用。在隨后的十幾年間，Childs 所領導的美國德州農工大學渦輪機械實驗室［8-11］繼續(xù)對離心泵密封口環(huán)、平衡鼓的動力學求解進行了理論模型的完善，試驗對比研究并將“有限長”求解方法應用于錐形及迷宮密封等多種型式的液體密封中。

本研究對Childs提出的“有限長”方法進行編程改進，在原有收斂條件的基礎上增加出口壓力恢復項，并采用“打靶法”對高壓多級離心泵口環(huán)密封內流體控制方程組進行Matlab編程求解，從而求得口環(huán)密封的動特性系數。同時，用VB編寫程序界面，通過動態(tài)鏈接庫的方式將界面與內核算法結合，使計算程序脫離Matlab 運行環(huán)境，程序計算結果與Childs 理論計算結果一致。

1 出口壓力恢復系數的添加

1967年，Domm 等提出離心泵口環(huán)出口位置同樣存在壓力的不連續(xù)，即壓力突變的現象。流體即將流出出口面的壓力比流體剛流出出口面時要低，即壓力的恢復效應，其壓力關系可定義為：

對壓力分布及軸向速度進行攝動量簡化，可表示為：

本研究將以上攝動量代入無量綱處理后的壓力關系式，可得密封口環(huán)內壓力分布在密封入口及出口位置的壓力邊界條件：

入口壓力損失系數與出口壓力恢復系數的無量綱化及攝動法求解對比流程圖如圖1所示。

2 Matlab 程序的編制

考慮到Matlab 對矩陣及微分方程求解運算的便捷性，筆者采用Matlab編寫環(huán)形間隙動特性求解的內核程序，求解中采用打靶法對“有限長”理論推導出的由周向、軸向動量方程及連續(xù)性方程組成的流體一階攝動控制方程組進行求解。

圖1 入口壓力損失系數與出口壓力恢復系數添加流程圖對比

由于入口與出口壓力分布均與軸向速度有關，本研究在求解中假設軸向速度為基礎變量，將入口與出口位置壓力值用基礎變量表示，并采用壓力出口大小為收斂邊界條件，設定收斂準則為相鄰兩時間步長內三組未知數求解殘差小于10-8。介于收斂入口與出口位置壓力分布均與軸向速度有關，故?。?/p>

式中：uz1，uθ1—軸向及周向速度一階攝動量；p1—軸向壓力分布一階攝動量；b—僅與口環(huán)密封幾何參數與運行工況相關的系數。

同時，設：

則：

原方程組各式，對γk求偏導數，可得：

其中：

由以上推導可知，邊界條件為M1(0)=1,M2(0)=0,M3(0)=k。

設：

采用牛頓法對初值進行后續(xù)改進加速收斂：

故原方程組的求解可化為對γk初值進行改進，并驗證原邊界條件是否滿足的迭代求解過程。最終迭代結束，便得到環(huán)形間隙內流體壓力沿軸向所在位置的分布情況：

此時，對環(huán)形液體環(huán)的動特性系數進行定義并做相應分析。作用于轉子上面的反作用力可表示為：

根據液體環(huán)內壓力分布情況，對反作用力進行徑向與周向的分解分析，并進行無量綱化處理，如下：

由以上周向與軸向力的分析可知，在任意渦動頻率下，均可通過所求得的軸向與周向無量綱壓力分布函數f3c(z),f3s(z)沿軸向積分求得。在求解過程中，6個動力特性系數組成唯一的一組由兩個方程組成的六元一次方程組。對于某一固定工作轉速N，可取渦動頻率為0，0.5，1.0，1.5，2.0 倍的工作轉速，組成5 組六元一次方程組，每3組方程可求解出一組動特性系數，5組方程排列組合共求解10組動特性系數，求其平均值，并通過Matlab輸出其計算結果。

本研究在Matlab 運行環(huán)境下輸入高壓多級離心泵葉輪口環(huán)各幾何尺寸、操作工況，介質屬性等參數即可實現其當前工況下動特性系數的自動化求解，包括主剛度、主阻尼、交叉剛度、交叉阻尼、主附加質量及交叉附加質量。Matlab環(huán)境下簡易的輸入界面如圖2所示。

圖2 Matlab環(huán)境下程序自動化輸入界面

3 操作界面與數據傳輸的實現

Matlab 作為一種高性能數值計算軟件，提供了強大的數值分析、矩陣運算等功能，作為該軟件的核心求解程序，編程便捷，計算效率高，但其界面開發(fā)能力較差，與其他后處理軟件、功能開發(fā)等的兼容性較差，且不能脫離Matlab 環(huán)境運行。相比之下，VisualBasic（VB）在設計開發(fā)Windows應用程序方面界面友好，方便快捷。本研究集成VB 和Matlab 各自的優(yōu)點，內核計算由Matlab 獨立完成，VB 編寫程序界面，完成數據傳輸及用戶交互，極大地提高程序的開發(fā)效率，并使其脫離于Matlab獨立運行。

本研究對中間文件傳遞法、ActiveX自動化服務技術、DDE技術進行比較，發(fā)現以上3種混合編程在調用的過程中比較復雜繁瑣，且無法脫離Matlab的編譯環(huán)境，這直接影響了所開發(fā)軟件的獨立運行和發(fā)布。為實現脫離Matlab 運行環(huán)境，方便以后程序的推廣與用戶使用，本次研究采用動態(tài)鏈接庫的方式完成數據傳輸與程序間的調用。

通過Matlab 調用C++編譯器，將核心計算程序編譯成可供VB程序調用的后綴名為dll的動態(tài)鏈接庫文件，并將此動態(tài)鏈接庫文件在VB中進行注冊，即可在VB 編譯過程中對此核心程序完成調用，且完全脫離Matlab運行環(huán)境。程序用戶操作界面如圖3所示。

圖3 口環(huán)計算程序用戶交互界面

4 可執(zhí)行文件的生成與打包

為保證高壓多級離心泵動力學分析計算軟件在其他非編譯計算機上的安裝及計算的正常執(zhí)行，研究者需通過VB 6.0將工程文件生成可執(zhí)行文件，利用軟件Setup Factory 7.0 對可執(zhí)行文件及其相關的窗口文件、鏈接庫文件、注冊文件等進行封裝，通過Setup Factory 7.0 內部的代碼編寫實現鏈接庫的自動注冊，最終打包成安裝程序，并進行軟件授權和安裝界面美化等方面的處理，打包文件列表如圖4所示。

圖4 打包文件列表

本研究目前針對已調試程序進行了簡略的打包試驗，安裝過程完整順暢，軟件安裝界面如圖5所示。安裝完成后，通過單獨執(zhí)行注冊程序對動態(tài)鏈接庫進行注冊，安裝完成后軟件能夠正常運行。

圖5 軟件安裝界面

5 結束語

本研究在“有限長”理論求解方法的基礎上，將口環(huán)密封出口壓力恢復系數添加到邊界條件中，并利用“打靶法”對該邊界條件下的口環(huán)內流體控制方程組進行了Matlab 編程求解。同時，筆者利用VB 編寫程序界面，實現Matlab核心計算程序與整體界面間的數據交互與程序調用，使計算程序脫離Matlab 運行環(huán)境。此外，本研究還對可執(zhí)行文件及其相關的窗口文件、鏈接庫文件、注冊文件等進行封裝，實現了非編譯計算機上的安裝及計算的正常執(zhí)行。

下階段，本研究將通過補充具有較大長徑比的平衡鼓的計算程序，并嘗試實現動態(tài)鏈接庫自動注冊，在今后的研究中，通過與更多的試驗結果進行對比，完善口環(huán)動特性系數的計算程序。

（References）：

［1］LOMAKIN A A.Calculating of critical speed and securing of dynamic stability of high pressure pumps with reference to forces arising seal gaps［J］.energomashinostroenic，1958，4（1）：1158-1162.

［2］BLACK H F.Effects of high-pressure ring seals on pump rotor vibrations［J］.Fluid Engineering Division.1971，184（3）：92-100.

［3］JENSSEN D N.Dynamics of rotors systems embodying high pressure ring seals［D］.Edinburgh：Heriot-Watt University，1970.

［4］BLACK H F.Calculation of forced whirling and stability of pump rotor vibrations［J］.Journal of Manufacturing Sci?ence and Engineering.1974，96（3）：1076-1081.

［5］CHILDS D W.Dynamic analysis of turbulent annular seals based on hirs' lubrication equation［J］.ASME Journal of lubrication technology，1983，105（3）：429-436.

［6］CHILDS D W.Finite-length solutions for rotordynamic co?efficients of turbulent annular seals［J］.Journal of Tribol?ogy，1983，105（3）：437-444.

［7］HIRS G G.Fundamentals of a bulk flow theory for turbu?lent lubricant films［D］.Pelft：PhD thesis，Delft University，1970.

［8］CHILDS D W.Turbomachinery rotordynamics：phenomena，modeling，and analysis［M］.Hoboken：John Wiley & Sons，1993.

［9］CHILDS D W.Dressman J B.Convergent-Tapered Annular Seals：Analysis and Testing for Rotordynamic Coefficients［J］.Journal of Tribology.1985，107（3）：307-316.

［10］KIM C H，CHILDS D W.Analysis for rotordynamic coeffi?cients of Helically- Grooved turbulent annular seals［J］.Journal of tribology，1987，109（1）：136-143.

［11］曹樹謙，陳予恕.現代密封轉子動力學研究綜述［J］.工程力學，2009，26（Ⅱ）：68-79.