郭良斌 吳永良

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

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

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

1 軸承結構原理與內(nèi)流域抽取

1.1 軸承結構原理

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

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

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

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

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

1.2 內(nèi)流域抽取

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

圖3 軸承內(nèi)流道幾何模型

2 網(wǎng)格劃分與邊界條件設置

2.1 網(wǎng)格劃分

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

圖4 軸承整體網(wǎng)格與局部網(wǎng)格

2.2 網(wǎng)格無關性驗證

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

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

圖5 網(wǎng)格無關性驗證結果

2.3 邊界條件設置

如圖6所示，軸承圓盤直徑為120.54 mm，穩(wěn)流管氣流入口直徑為20 mm，供氣孔直徑為10 mm。將上、下穩(wěn)流管入口面設置為壓力入口IN，計算域與遠場外壁的接觸面設置為壓力出口OUT，壓力出口OUT軸向?qū)挾?0.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。其中，完全發(fā)展的管流的核心湍流強度用以下經(jīng)驗公式計算：

圖6 軸承各邊界命名

(1)

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

DH=D

(2)

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

3 控制方程與湍流模型

3.1 控制方程

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

(3)

(4)

(5)

(6)

3.2 湍流模型

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

Gb-ρε-YM+Sk

(7)

(8)

4 數(shù)值計算與結果分析

4.1 計算參數(shù)設置

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

(9)

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

4.2 結果及分析

4.2.1 三維氣膜流場靜壓分析

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

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

4.2.2 三維氣膜流場馬赫數(shù)分析

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

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

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

4.2.3.1 二維簡化模型

二維模型及其邊界條件如圖9(a)所示，即OUT為壓力出口，IN為壓力入口，SYM為對稱面，AXIS為對稱軸，其余WALL-WL、WALL-BEARING、WALL-ZZCM、WALL-FARFIELD均為壁面。采用ICEM對二維模型進行網(wǎng)格劃分，并對邊界層進行局部加密，控制y+值在2.0以內(nèi)，最終劃分網(wǎng)格平均質(zhì)量為0.96，平均Aspect ratio為1.98，網(wǎng)格數(shù)量約90萬個，收縮段入口網(wǎng)格劃分如圖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 二維模型邊界條件及收縮段入口網(wǎng)格劃分

4.2.3.2 二維模型網(wǎng)格無關性驗證

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

圖10 二維網(wǎng)格無關性驗證結果

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

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

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

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

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

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

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

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

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

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

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

5 結論

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

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

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