基于Workbench的高壓圓盤氣體軸承三維流場分析*

郭良斌 吳永良

(1.武漢科技大學冶金裝備及其控制教育部重點實驗室 湖北武漢 430081；2.武漢科技大學機械傳動與制造工程湖北省重點實驗室 湖北武漢 430081)

氣體軸承因其高精度、長壽命和低功耗等特點，已被廣泛運用于紡織機械、儀表機床等現代工業領域[1]。傳統的靜壓氣體軸承剛度低、承載力小、供氣壓力低，只適用于輕載精密支撐[2]。為此，研究人員對氣體軸承的發展與應用進行了深入研究[3-4]。徐凡[5]參考風洞實驗中有關收縮段的設計理論，設計出氣體軸承新型雙對稱收縮段壁面型線，使氣體軸承內部流場得到有效改善，消除了供氣壓力的限制。劉劍和郭良斌[6]對高壓圓盤氣體軸承可壓縮邊界層形態進行了詳細分析，研究了氣膜厚度、軸承半徑以及供氣壓力對氣膜邊界層發展特征的影響。郭良斌和劉劍[7]對不同供氣總壓時高壓圓盤氣體軸承邊界層的對流換熱系數進行了進一步研究，通過對氣膜內的流場進行數值模擬，分析了速度邊界層的特征，并預測了氣膜內邊界層的轉捩點。然而，因為該高壓圓盤氣體軸承具有雙對稱性，上述研究均是在雙對稱條件下取模型的1/2或1/4進行簡化，從而得到具有對稱面和對稱軸的二維簡化模型及其流場特性。目前，國內外已有不少學者對氣體軸承的三維流場特性進行了研究[8-14]，對未來研究氣體軸承特性具有很好的參考價值。由于高壓圓盤氣體軸承模型的特殊性，以及氣膜間隙內湍流的不確定性，對該模型進行三維流場分析顯得尤為重要。

為了研究高壓圓盤氣體軸承完整模型的三維流場特性，本文作者采用Ansys17.0仿真軟件內的Workbench工作平臺，利用Solidworks專業建模軟件對高壓圓盤氣體軸承進行三維建模，采用Workbench自帶的網格劃分工具Mesh對幾何模型進行網格劃分，并利用Fluent流體分析模塊進行計算，最終得到了高壓圓盤氣體軸承的三維流場，并與二維簡化模型流場計算結果進行了對比。

1 軸承結構原理與內流域抽取

1.1 軸承結構原理

由于文中主要研究的是高壓圓盤氣體軸承內流場狀況，故對軸承本身建模作簡化處理，忽略穩流管與工作圓盤之間的密封墊以及工作圓盤上、下端面的倒角。高壓圓盤氣體軸承的三維模型結構如圖1所示。

由圖1可知，完整的高壓圓盤氣體軸承由上、下圓盤組件構成；上圓盤組件包括上穩流管和上工作圓盤，二者通過螺紋連接緊固在一起，上圓盤組件與下圓盤組件完全對稱。穩流管長72 mm，工作圓盤厚30 mm，兩工作圓盤之間的工作氣膜間隙為0.2 mm，空氣由上、下穩流管的供氣孔進入，經氣膜間隙從兩工作圓盤間的外邊緣射流流出。由于采用Workbench中的Fill功能進行內流域抽取，該功能需保持入口和出口均為封閉狀態。為了使整個計算域包含氣膜出口縫隙射流區域，將圖1的原始三維模型修改為如圖2所示的結構。

圖1 高壓圓盤氣體軸承三維結構

圖2 修改后的高壓圓盤氣體軸承三維結構和剖面

即在軸承氣膜出口添加軸向長度20.2 mm、徑向長度15 mm的三維域作為氣體出口的遠場，用薄殼包裹，殼厚2 mm。

1.2 內流域抽取

將圖2模型導入Workbench中的Fluid Flow(Fluent)分析系統，利用Workbench自帶的Design Modeler建模模塊的Fill功能對該軸承內流場進行流體域抽取，再通過壓縮軸承實體得到軸承內流道的幾何模型，如圖3所示。

圖3 軸承內流道幾何模型

2 網格劃分與邊界條件設置

2.1 網格劃分

如圖4所示，Workbench中的Fluid Flow(Fluent)分析系統會自動保存Design Modeler模塊中的建模，并將其導入Mesh網格劃分模塊。采用的網格劃分方式為自動網格劃分，通過全局網格參數設置，將網格精度relevance、網格細化程度relevance center以及跨度中心角span angle center調至最優，網格平滑度smoothing調至中等；選擇尺寸函數類型為curvature，系統會根據曲率法向角度自動確定細化邊與曲面處的網格單元尺寸大小，設置曲率法向角度curvature normal angle為30°，最小網格尺寸為0.05 mm，最大面部尺寸為1.2 mm，最大四面體尺寸為1.2 mm，網格增長率為1.05，其余保持默認設置。最終劃分網格平均質量為0.83，平均Aspect ratio為1.9，網格節點數約30萬個，網格數量約160萬個。

圖4 軸承整體網格與局部網格

2.2 網格無關性驗證

一般而言，模型劃分網格數量越多，其計算結果精度越高，但過高的網格數量，會延長計算時間，降低工作效率，甚至普通的計算機難以進行正常工況計算。因此劃分合適的網格數量顯得尤為重要。

如圖5所示為83萬～207萬網格數量下的軸承氣膜對稱線上馬赫數的分布狀況。可知，不同網格數量下在軸承入口和收縮段部分(72.39～132.39 mm)以及平行氣膜區(42.39～72.39 mm、132.39～162.39 mm)其馬赫數曲線基本重合，僅在收縮段出口處(72.39、132.39 mm)以及縫隙射流區(24～42.39 mm、162.39～178 mm)等處出現小幅變化。但馬赫數總體變化趨勢是一致的，且160萬網格數量時在軸承氣膜出口處能更好反映出馬赫數最大值?？紤]到計算機配置及程序運行等情況，文中采用160萬網格數量的三維劃分模型進行工況計算。

圖5 網格無關性驗證結果

2.3 邊界條件設置

如圖6所示，軸承圓盤直徑為120.54 mm，穩流管氣流入口直徑為20 mm，供氣孔直徑為10 mm。將上、下穩流管入口面設置為壓力入口IN，計算域與遠場外壁的接觸面設置為壓力出口OUT，壓力出口OUT軸向寬度20.2 mm，其余面WALL-WL、WALL-BEARING、WALL-ZZCM、WALL-FARFIELD均設為壁面。入口總壓設置為2 MPa(20 atm)，初始表壓設置為1.998 MPa(19.98 atm)。湍流項選擇湍流強度和水力直徑，其中湍流強度設為7%，入口水力直徑為20 mm，入口總溫為300.15 K。出口表壓設置為0.1 MPa(1 atm)，湍流強度設為7%，出口水力直徑為150.54 mm，總溫為300.15 K。其中，完全發展的管流的核心湍流強度用以下經驗公式計算：

圖6 軸承各邊界命名

(1)

式中：I為湍流強度；u′為湍流脈動速度；uavg為平均速度；ReDH為按水力直徑計算得到的雷諾數，圓管水力直徑即為其自身直徑D：

DH=D

(2)

所有壁面均為固定無滑移恒溫壁面，壁溫300.15 K。

3 控制方程與湍流模型

3.1 控制方程

流體處于穩態流動時，其三維可壓縮的、笛卡爾坐標系下的控制方程[7]為

(3)

(4)

(5)

(6)

3.2 湍流模型

完整的控制方程需要引入湍流模型將其封閉，建立雷諾應力方程組。文中數值計算使用Realizablek-ε模型。 Realizablek-ε是近期出現的帶旋流修正的k-ε模型，對比標準的k-ε模型，其增加了一個與湍流黏性有關的公式，為湍流耗散率增加了新的傳輸方程，該模型直接的好處就是對于平板和圓柱射流的發散比率有著更精確的預測。Realizablek-ε模型的方程如下：

Gb-ρε-YM+Sk

(7)

(8)

4 數值計算與結果分析

4.1 計算參數設置

采用Workbench中的Fluent模塊，對模型進行求解。選取基于壓力的絕對穩態求解器，操作壓力設為0，開啟能量方程與湍流方程，設置流入氣體為空氣，并假設其為連續的理想氣體，其密度采用理想氣體狀態方程，黏度采用薩瑟蘭公式：

(9)

式中：T0=273.15 K；μ0=1.716×10-5Pa·s；Ts=110.56 K。空氣導熱系數、比定壓熱容以及分子質量設為默認常數，湍流模型使用標準壁面函數。求解方法采用SIMPLEC算法，壓力、密度以及湍流動能等差分格式均保持默認設置，各項欠松弛因子、各項殘差標準也保持默認設置。

4.2 結果及分析

4.2.1 三維氣膜流場靜壓分析

圖7(a)所示為氣膜流場靜壓云圖?？芍S承氣流入口壓力最高，氣流到達收縮段前壓力下降并不明顯；從收縮段開始，壓力以穩流管軸線為中心，軸對稱向四周逐漸降低，但整個軸承氣膜內壓力仍維持在0.5 MPa(4.93 atm)以上，僅在氣膜出口處很窄的范圍內壓力快速下降，達到了良好的高壓潤滑效果；氣流經氣膜間隙出口流出后，氣流逐漸降至略大于環境壓力。利用Fluent軟件的Plane Tool工具新建氣膜對稱面，再將軸承圓盤及穩流管的軸對稱面與其相交，截得氣膜對稱線。圖7(b)所示為氣膜對稱線上靜壓分布曲線?？芍?，靜壓分布曲線以R=102.39 mm的鉛垂線為對稱軸左右對稱，以右半邊曲線為例，氣流壓力在供氣孔區域a(102.39～110 mm)下降緩慢，而在收縮段區域b(110～138.9 mm)氣流壓力快速下降至0.89 MPa(8.9 atm)，而后氣流到達平行氣膜區c(138.9～62.39 mm)，壓力下降至0.5 MPa(4.93 atm)。在氣膜出口處很窄范圍內壓力迅速降至0.11 MPa(1.1 atm)，最后在縫隙射流區d(162.39～178 mm)壓力維持在略大于環境壓力。

圖7 氣膜流場靜壓云圖及氣膜對稱線上靜壓分布曲線

4.2.2 三維氣膜流場馬赫數分析

圖8(a)所示為氣膜流場馬赫數云圖。可知，穩流管氣流入口處馬赫數最低，氣流到達收縮段前馬赫數上升不明顯，馬赫數不超過；從收縮段開始，馬赫數以穩流管軸線為中心，軸對稱向四周逐漸升高，可以觀察到，位于收縮段終點的馬赫數已升至高亞音速，進而在平行氣膜區有小幅度降低；而后在氣膜出口處因氣流壓力遠大于環境背壓，使氣流速度上升至音速，隨后，在氣膜出口外進一步形成超音速圓盤縫隙射流，最終因膨脹波和壓縮波的交替作用馬赫數逐漸趨于0，氣流逐漸靜止。圖8(b)所示為氣膜對稱線上的馬赫數分布曲線?？梢?，馬赫數曲線以R=102.39 mm的鉛垂線為對稱軸左右對稱，取右半邊曲線為例，氣流在供氣孔區域a(102.39～110 mm)的馬赫數為0.2，經過收縮段b(110～138.9 mm)后，馬赫數上升至0.92，且在平行氣膜區c(138.9～162.39 mm)的絕大部分區域，氣流都維持在馬赫數0.7～1的高亞音速狀態。在氣膜間隙出口處，由于流通截面的突然增大，氣流因膨脹波的作用升至超音速馬赫數1.36，經過縫隙射流區d(162.39～178 mm)后，馬赫數逐漸趨于0。

圖8 氣膜流場馬赫數云圖及氣膜對稱線上馬赫數分布曲線

4.2.3 三維模型與二維簡化模型的計算結果對比

4.2.3.1 二維簡化模型

二維模型及其邊界條件如圖9(a)所示，即OUT為壓力出口，IN為壓力入口，SYM為對稱面，AXIS為對稱軸，其余WALL-WL、WALL-BEARING、WALL-ZZCM、WALL-FARFIELD均為壁面。采用ICEM對二維模型進行網格劃分，并對邊界層進行局部加密，控制y+值在2.0以內，最終劃分網格平均質量為0.96，平均Aspect ratio為1.98，網格數量約90萬個，收縮段入口網格劃分如圖9(b)所示。利用Fluent求解二維流場。入口總壓為2 MPa(20 atm)，初始表壓為1.998 MPa(19.98 atm)，湍流項選擇湍流強度和水力直徑，其中湍流強度為7%，入口水力直徑為20 mm，入口總溫為300.15 K；流體出口設置為壓力出口，出口表壓為0.1 MPa(1 atm)，湍流強度為7%，出口水力直徑為120 mm，總溫為300.15 K；采用SIMPLE算法，其余設置同三維模型。

圖9 二維模型邊界條件及收縮段入口網格劃分

4.2.3.2 二維模型網格無關性驗證

如圖10所示為網格數量在65萬～97萬時二維模型對稱面上的馬赫數分布狀況。

圖10 二維網格無關性驗證結果

由圖10可知，隨著網格數量的增加，軸承氣膜出口處(60～62 mm)最大馬赫數會發生小幅變化，而不同網格數量時馬赫數整體變化趨勢是一樣的，其馬赫數曲線基本重合?？紤]到計算機配置及程序運行等情況，文中采用90萬網格數量的二維劃分模型進行工況計算。

4.2.3.3 2種模型的計算結果比較

利用Fluent自帶截面工具截取氣膜對稱線，并獲取線上的馬赫數與靜壓數據。因為二維簡化模型流場的氣膜對稱線只有整個氣膜的1/2，為了方便對比，使用了MatLab軟件編程，對二維計算數據進行對稱擴展，得到整個氣膜的對稱線上的數據，在origin中繪圖進行對比，如圖11所示。

由圖11(a)可知，二維簡化模型與三維模型的靜壓變化趨勢基本一致，且變化的位置也基本一致，但在平行氣膜區開始部分(130～140 mm)，二維簡化模型靜壓曲線壓力明顯高于三維模型靜壓曲線的壓力。其原因在于三維模型使用Workbench自帶的mesh模塊進行網格劃分，其網格劃分能力稍差，難以添加高質量邊界層網格，對于高壓圓盤氣體軸承的微小氣膜間隙(0.2 mm)難以劃出較好的結構化網格，整體網格質量一般(0.83左右)；而二維簡化模型使用ICEM專業網格劃分軟件，網格劃分能力強，可添加邊界層并合理控制y+值，網格質量較好(0.96以上)，因而造成了計算結果上的小幅偏差，以及出現了不一樣的極值點。因此，二維模型網格質量比三維模型網格質量更高，其最大值顯示更為精確。

圖11 2種模型氣膜對稱線上靜壓和馬赫數對比

由圖11(b)可知，二維簡化模型與三維模型的馬赫數變化趨勢基本一致，且變化的位置也基本一致，但在氣膜間隙出口處，二維簡化模型馬赫數曲線的最大值明顯高于三維模型馬赫數曲線的最大值，其原因仍是二維網格質量更高使其最大值更為精確；二維模型縫隙射流區較窄(70 mm-58.9 mm=11.1 mm)，三維模型縫隙射流區較寬(178 mm-165 mm=13 mm)，且縫隙射流區(165～178 mm)三維模型馬赫數曲線明顯高于二維簡化模型馬赫數曲線，其原因在于三維模型中氣體流出氣膜間隙后，在縫隙射流區不僅有徑向的流動，還存在周向發散的流動，這是二維簡化模型難以預測到的現象。

4.2.4 不同湍流強度對高壓圓盤氣體軸承靜壓和馬赫數的影響

雖然湍流強度理論上按公式(1)進行計算，但實際上通過管道供氣，其湍流強度可能會發生改變。研究不同湍流強度下高壓圓盤氣體軸承靜壓和馬赫數的變化情況，對后續實驗有著重要指導作用。不同湍流強度下高壓圓盤氣體軸承靜壓和馬赫數分布情況如圖12所示。

圖12 不同湍流強度下氣體軸承靜壓和馬赫數分布情況

由圖12(a)可知，不同湍流強度時供氣孔區域a、收縮段區域b、平行氣膜區c以及縫隙射流區d的靜壓分布曲線基本重合，可見不同湍流強度對高壓圓盤氣體軸承內氣膜壓力無明顯影響。其原因在于湍流強度和軸承氣膜壓力之間并無直接關聯。

由圖12(b)可知，不同湍流強度時供氣孔區域a、收縮段區域b以及平行氣膜區c的馬赫數分布曲線基本重合，而在縫隙射流區d及氣流遠場區域(12.39～26 mm、178～192.39 mm)馬赫數存在較明顯的變化。其原因在于氣流經平行氣膜出口流出時加速為超音速，且湍流流動存在不確定性，湍流脈動速度在不斷發生變化。

5 結論

(1)高壓圓盤氣體軸承內，氣膜壓力隨著半徑方向R的增大而逐漸減小。平行氣膜區絕大部分區域氣膜壓力在0.5 MPa(4.93 atm)以上，僅在氣膜出口處很窄的范圍內壓力快速下降。經超音速圓盤縫隙射流流出后，氣流壓力逐漸降至略高于環境壓力。由于超音速流動時外部環境低壓不影響上游壓力分布，故文中的高壓圓盤氣體軸承可以實現良好的高壓潤滑作用。

(2)就整體變化趨勢而言，二維簡化模型與三維模型計算結果基本保持一致，出現小范圍差別的主要原因在于，三維模型使用Workbench自帶的mesh模塊進行網格劃分，其網格劃分能力稍差，難以添加高質量邊界層網格。但二者馬赫數、靜壓曲線變化趨勢基本一致，較好地反映了高壓圓盤氣體軸承內部流道的流動規律。

(3)不同湍流強度對高壓圓盤氣體軸承內流場的靜壓無明顯影響，但對縫隙射流區及氣流遠場區域的馬赫數有著小幅影響。因此后續實驗中，測量氣膜內流動參數時，可以排除供氣管道布置方式改變造成入口湍流強度變化這一影響因素。