基于有限元法的水輪發電機組臨界轉速研究

姬聯濤，袁浩強，荊岫巖，王 璞，劉 昊，楊雍杰，李超順

（1.中國電力科學研究院有限公司，江蘇 南京 210003；2.華中科技大學 土木與水利工程學院，湖北 武漢 430074；3.國家電網有限公司，北京 100031；4.國網湖南省電力有限公司，湖南 長沙 410004）

0 引言

水電作為一種可再生能源，由于其諸多天然優勢，在我國使用極為廣泛，各流域大中小型水電站建設數量眾多。隨著國內更多中高水頭水電站的建設與增容改造，機組運行工況更加惡劣、穩定運行要求更高，在實際運行過程中，當機組轉速上升到一定值時，激振力頻率與軸系固有頻率相近，且作用方向與振型方向一致時，軸系會產生強烈振動，該轉速被稱為臨界轉速［1］。若機組在此轉速下長時間運行，可能造成結構損壞、噪聲加劇，嚴重危及廠房及工作人員的生命財產安全，因而對水電機組軸系振動特性尤其是臨界轉速的研究愈發重要。

目前針對此類轉子系統振動特性分析，國內外學者已經做了大量的研究：榮吉利等［2］基于有限元法自主編制程序ROTOR，對紅石和隔河巖水電站軸系的自振特性進行計算與分析，并與試驗結果進行對比；HUA 等［3］利用Muszynska 非線性密封力，構建轉子與密封的動力模型，并通過分岔圖和龐加萊圖分析了密封如何干擾系統的穩定性；楊曉明等［4］采用有限元法建立水電機組軸系橫向振動力學模型，并在考慮軸承簡單耦合關系下求解了固有頻率與臨界轉速；LUNENO 等［5］針對組合軸承對水電機組轉子不同方向的耦合作用，分析其改變機組固有頻率與振動模式的機理；翟黎明等［6］利用SAMCEF ROTOR 程序建立了抽水蓄能發電機組轉子-軸承-電磁系統三維有限元模型并分析了不平衡磁拉力與導軸承剛度變化對臨界轉速的影響；馬晨原等［7］采用有限元分析軟件探究了發電機質量偏心及導軸承剛度變化對機組振動特性的影響，并與現場測量數據進行對比；徐瑞紅等［8］利用ANSYS 軟件對某水電站軸系進行模態分析與臨界轉速計算，得出其發生共振可能性很小的結論。

以上研究為水電機組軸系的動態特性分析與結構優化提供了借鑒，但多數模型過于簡化，各類軸承與密封裝置的耦合作用未充分考慮，不能夠真實地反映出機組的邊界條件與工作狀態。鑒于此，本文針對某水電站的混流式機組，考慮上導、水導和推力軸承的耦合關系，分別構建無機械故障、轉子或轉輪質量偏心、聯軸法蘭不對中的有限元分析模型，進行轉子動力學仿真求解臨界轉速，從而分析軸系的動態性能。

1 轉子動力學仿真原理

轉子動力學分析主要包括穩定性分析、臨界轉速計算和不平衡響應三部分［9］，起初科研工作者們多數采用傳遞矩陣法，但隨著計算機運算性能的提高，研究者們更加青睞于精度更高、直觀性更強的有限元法。基于有限元動力學仿真推導出來的完整轉子動力學微分方程為：

求解轉子系統臨界轉速，就是在給定轉子系統轉速Ω的前提下，求解系統渦動轉速ω。當所得渦動轉速ω與給定轉速Ω相等時，便稱此轉速為臨界轉速，該渦動轉速下的振型就是臨界轉速的振型。

2 物理模型

2.1 幾何參數

由水輪機與發電機圖冊可知該機組為二導懸式混流式水輪發電機組，即存在兩個導軸承與一個推力軸承，利用三維可視化建模軟件UG，采用傳統的自底向上建模方法，即先對機組的轉輪、水輪機主軸、發電機大軸、轉子機架、磁極、磁軛、推力頭等部件進行精細化建模隨后裝配，最終得到軸系計算簡圖與三維裝配模型如圖1所示。

圖1 軸系計算簡圖與三維裝配模型Fig.1 Shafting calculation diagram and three-dimensional assembly model

2.2 材料參數

根據圖紙結合各部分的物理特性，軸系各部件材料參數統計如表1所示，在ANSYS中建立材料庫，并根據表1對機組各部分進行材料賦值。

表1 軸系各部件材料參數統計Tab.1 Material parameter statistics of shafting components

3 有限元模型

3.1 網格劃分

為了節約計算資源，并保證仿真結果的精度，采用映射劃分、自由劃分與掃掠劃分相結合的方法，并在網格不連續處設置綁定接觸，單元類型為六面體20 節點的SOLID186 與四面體10 節點的SOLID187 單元，經劃分的模型節點數為359 150，單元數為164 660，結果如圖2所示。

圖2 水輪發電機組網格劃分結果圖Fig.2 Grid division result of hydro-generator set

3.2 邊界條件

該二導懸式水輪發電機組的固有頻率很大程度上受到水導、上導、推力軸承等邊界條件的影響，故在保證高精度物理模型的基礎上合理簡化邊界條件尤為重要。對于實際的水電機組，軸承座、軸瓦和機墩等都具有一定彈性，在使用彈簧阻尼器模擬各軸承時，其彈簧剛度可視為油膜剛度、軸承座與機墩剛度三者的串聯，由于本研究著重分析軸系的振動特性，故假設軸承座與機墩為剛性，僅考慮油膜動力特性參數。

3.2.1 導軸承

為了準確模擬上導和水導軸承的油膜剛度與阻尼，本文采用Lund 提出的八參數模型［10］，其對應ANSYS 的彈簧阻尼單元COMBI214，包括線性化的4 個剛度系數與4 個阻尼系數，COMBI214單元幾何結構如圖3所示。

圖3 COMBI214單元幾何結構Fig.3 The geometry of COMBI214 unit

將支承視為各向同性，忽略交叉剛度與阻尼的影響，其主剛度與主阻尼可通過油膜力求得，但是油膜力表達式無法直接由流體潤滑的Reynolds 方程導出，故選擇無限短軸承理論［11］對軸承做近似處理，得到導軸承無量綱剛度系數和阻尼系數的解析表達式為：

式中：K為無量綱剛度系數；C為無量綱阻尼系數；B為導軸承長度；D為導軸承直徑；ε為軸頸中心的偏心率，等于偏心距離與軸承半徑間隙的比值。

導軸承剛度系數與阻尼系數的有量綱表達式為：

式中：k為有量綱剛度系數；c為有量綱阻尼系數；μ為潤滑油度；n為軸頸轉速；L為軸承軸向長度；φ為軸承間隙比，等于軸承半徑間隙與軸頸半徑的比值。

將涉及到的相關參數代入以上公式，軸頸中心偏心率取0.4，求得額定轉速下各導軸承剛度系數與阻尼系數：上導軸承剛度系數為8.55×106N/mm，阻尼系數為1.39×105（N·s）/mm；水導軸承剛度系數為7.09×107N/mm，阻尼系數為1.37×106（N·s）/mm。

3.2.2 推力軸承

考慮到推力軸承對軸系不僅有軸向的支撐作用，當鏡板隨軸的彎曲振動而傾斜時，油膜力還會提供恢復力矩，相當于一個扭轉彈簧。因而將推力軸承約束視為伸縮彈簧與扭轉彈簧的共同作用，其油膜剛度系數與阻尼系數可通過文獻［12，13］所提供的解析表達式并結合自身參數求得，軸向剛度系數取為8.50×108N/mm，阻尼為1.01×106（N·s）/mm；扭轉剛度系數取為3.18×109（N·s）/rad，阻尼為7.64×107（N·s）/mm。

最終得到軸系三維有限元模型如圖4所示。

圖4 軸系三維有限元模型Fig.4 Three dimensional finite element model of shafting

4 計算結果與分析

4.1 無機械故障的軸系振動特性分析

所謂無機械故障是指理想狀態下水電機組各部件的加工與設計階段完全一致，不存在材料的缺損、發電機轉子與水輪機轉輪的偏心與變形、聯軸法蘭不對中、軸系整體不對中等故障，且設備的安裝精度達到100%。利用軸系有限元模型考慮陀螺效應在八個轉速點下進行復模態分析，采用QR Damped 法提取各階模態，水電機組軸系前八階振型如圖5所示。

由圖5 可知，前兩階振型表現為以兩個導軸承為支點的擺動，轉輪下環處位移最大、推力頭次之；第三階為轉輪的扭轉變形，下環處位移最大且沿周向均勻分布；第四、五階振型的危險點位于推力頭上端，其次是轉輪，其余部位幾乎不發生振動；第六、七階主要表現為推力頭、轉子、轉輪的振擺，其中位移最大處仍位于推力頭上端與轉輪下沿；第八階振型下水輪機很安全，主要表現為發電機整體的擺動，推力頭處位移最大。隨后利用復模態分析結果繪制該轉速范圍內固有頻率與轉速相關激振頻率的關系曲線，即Campbell圖，結果如圖6所示。由于本研究主要針對軸系的徑向振動特性，故排除扭轉振型對應的模態分析結果，前四階固有頻率隨機組轉速的變化情況如表2所示。

表2 前四階固有頻率隨機組轉速的變化情況Tab.2 Variation of the first four order natural frequencies with unit speed

圖6 無機械故障下軸系Campbell圖Fig.6 Campbell diagram of shafting without mechanical failure

結合圖6 與表2 分析可知，隨著機組轉速增大，正向渦動固有頻率上升，反向渦動固有頻率下降，前四階臨界轉速分別為：74.35、105.64、124.55、169.38 rad/s，但第一、三階為反向渦動，在實際運行過程中該共振并不會被激發，因而其并不能稱作真正意義上的臨界轉速。對于水輪發電機組，一般要求第一階臨界轉速為飛逸轉速的1.4 倍以上［14］，由資料可知該機組的飛逸轉速為300 r/min，即31.42 rad/s，遠小于正向渦動的第一階臨界轉速，故滿足設計要求，機組在正常運行時不會發生共振現象。

4.2 質量偏心對軸系振動特性的影響

由于水電機組在安裝階段的誤差或者長期運行過程中的磨損變形，在軸系的主要部位，即發電機轉子和轉輪往往會出現質量偏心問題，為探究其偏心率其對軸系振動的影響，首先計算出理想工況下兩者的質量、體積、質心位置、轉動慣量等物理參數，再以該發電機定轉子間的平均氣隙長度，即27 mm 為基準，由平行軸定理［15］計算出轉子與轉輪分別偏心20%、40%、60%后的參數變化，結果如表3所示。

表3 轉子與轉輪物理參數隨偏心率的變化情況Tab.3 Variation of rotor and runner physical parameters with eccentricity

在ANSYS 中采用mass21 單元分別模擬發電機轉子或水輪機轉輪偏心對軸系的影響，隨后以同樣方式繪制Campbell 圖并計算各階臨界轉速，得到臨界轉速隨轉子或轉輪偏心率的變化情況如表4所示。

表4 臨界轉速隨轉子或轉輪偏心率的變化情況Tab.4 The variation of critical speed with the eccentricity of the rotor or runner

由表4 可知，隨著發電機轉子或者水輪機轉輪偏心率的增大，前四階臨界轉速的變化微乎其微，變化幅度不超過0.01%，因而可以得出結論：至少在該偏心范圍內，轉子與轉輪的質量偏心不會對軸系臨界轉速產生影響，即臨界轉速與單純的質量偏心無關。

4.3 聯軸法蘭不對中對軸系振動特性的影響

由電站資料可知，水輪機軸與發電機軸通過法蘭結構連接，聯軸螺栓為20 條鉸制孔銷釘螺栓，其規格為M115×6 mm，長度為480 mm，在機組吊裝過程中不可避免地會造成一定的誤差，因而聯軸法蘭不對中現象普遍存在，為了探究偏心距離對臨界轉速的影響，本文分別建立了聯軸法蘭偏心10、50、100 mm的有限元模型并進行轉子動力學分析，得到臨界轉速隨聯軸法蘭偏心距離的變化情況如表5所示。

表5 臨界轉速隨聯軸法蘭偏心距離的變化情況Tab.5 Variation of critical speed with eccentric distance of coupling flange

由表5可知，隨著聯軸法蘭偏心距離的增大，低階臨界轉速呈減小態勢，高階臨界轉速略微增大，當聯軸法蘭偏心距離為100 mm時，前三階臨界轉速減小，第四階臨界轉速增大，與對中情況相比，增減幅度最大不超過0.33%，且實際運行過程中聯軸法蘭偏心距離遠達不到該數值，因而可以得出結論：水電機組臨界轉速受到聯軸法蘭不對中的影響不大。

5 結論

基于限元法建立了水輪發電機組軸系高精度動力學模型，考慮了上導、水導和推力軸承的耦合關系并計算出各向剛度與阻尼值，在無機械故障、轉子或轉輪質量偏心、聯軸法蘭不對中3種工況下進行了轉子動力學仿真，得到各階模態固有頻率、振型以及臨界轉速。通過分析可知：該水電機組軸系在上述3 種工況下的臨界轉速均滿足要求，且受質量偏心與聯軸法蘭不對中的影響很小，該研究對水輪發電機組結構優化設計具有一定的意義。