多孔介質空氣軸承靜態特性仿真分析

劉 浩， 須 穎， 王本明， 姚建華， 陶文彬

(廣東工業大學 機電工程學院， 廣東 廣州 510006)

引言

半導體檢測設備在芯片制造領域有重要作用，這些設備對運動導軌的精度和速度有著很高的要求。普通機械剛性接觸導軌很難實現高速高精度運動。氣浮導軌以其非接觸、低阻尼的優勢，可以實現亞微米甚至納米級的運動精度，因而成為實現高速高精密運動的主要導軌形式[1-4]。作為氣懸浮導軌的核心零件，氣浮軸承的合理設計能夠顯著提高運動臺的承載力和精度。氣浮軸承的靜態特性與氣膜厚度有著很大關系，氣膜太薄對導軌及軸承面的加工要求太高，氣膜太厚則其承載力和靜剛度又會顯著下降。有效預測氣浮軸承在不同氣膜厚度下的靜態特性有重要工程價值。

氣浮軸承的節流方式有小孔節流、表面節流及狹縫節流等。小孔節流氣浮軸承具有加工簡單、綜合性能良好等優勢；同時，小孔節流氣浮軸承的研究也比較成熟，這類氣浮軸承的設計準則已基本建立[5-7]。帶腔小孔節流氣浮軸承具有較高的承載力和剛度，但是容易發生振動[8-10]。帶均壓槽結構的氣浮支承軸承通過設計均壓槽形狀可以改善氣膜內壓力分布[11-12]。狹縫節流氣浮軸承具有更好的穩定性[13-14]。

多孔介質氣浮軸承有別于以上幾種節流方式，其特點在于采用多孔介質材料節流，具有更高的承載力、剛度和較好的穩定性[15-16]。其次，多孔材料表面上有成千上萬個出氣孔，因而多孔介質氣浮軸承對氣浮面刮傷不敏感。多孔介質氣浮軸承的另一個重要研究對象是多孔介質材料孔隙結構。通常采用滲透率綜合評價孔隙結構阻礙介質流動的程度[17-18]。而貫通孔隙度和滲透率有很大的聯系，其他孔隙對氣體流動沒有影響。

但是，上述研究在分析多孔介質滲透率時，多通過二維成像手段獲取多孔介質孔隙度，沒有考慮三維孔隙連通性的影響。本研究引入micro-CT斷層掃描方法獲得多孔介質孔隙結構并分離出孤立孔隙和貫通孔隙。理論計算了多孔介質的黏性系數和慣性系數，并進行了試驗驗證。采用CFD方法進行仿真，研究了多孔介質止推軸承在不同形狀、不同氣膜厚度下的承載力和剛度，為設計多孔介質氣浮軸承提供參考。

1 物性參數分析

1.1 滲透系數

多孔介質前后的壓差與介質內部微通道特征有很大的關系，比如孔隙大小、孔隙形狀、貫通孔隙度等。氣體通過多孔介質流動物理模型示意圖如圖1所示，淺色為固相，深色為孔隙相。

圖1 多孔介質流動示意圖

氣流從進口流入，沿著貫通孔隙流出，貫通區域壁面拖拽氣流。當多孔介質內部氣速較低時，流動為層流狀態，多孔介質流動方向壓差與速度成線性關系，方程為：

(1)

式中， Δp—— 壓差

L—— 多孔介質流動方向的長度

μ—— 流動介質的黏度

K—— 滲透率

v—— 多孔介質前端表面氣速

(2)

式中,Q—— 多孔介質前端工況體積流量

A—— 多孔介質截面積

當多孔介質壓差與速度遵循式(1)時，流動狀態為Darcy流；隨著流速增大，壓差與速度不再遵循式(1)，流動將脫離Darcy狀態進入Forchheimer狀態，壓差與速度的關系可用式(3)描述：

(3)

式中，ρ—— 流動介質密度

C—— 慣性系數

(4)

假設流動是穩態充分發展的，流動主要沿著水平方向。從孔隙尺度上，內部的速度場、壓力場比較復雜。Forchheimer方程的推導是基于管道流，因此在推導過程中采用了體積平均的方法。假設在微小控制容積長度dx內，壓力和速度是都是均勻分布的。

平均速度可以表示為質量的函數，其中G為質量流量，A為多孔介質截面積。

(5)

考慮到空氣的壓縮性，引入理想狀態方程，壓縮空氣的密度ρ可以表示為：

(6)

式中，p—— 壓縮氣體的絕對壓力

R —— 特定氣體常數

T—— 溫度

將式(5)、式(6)帶入到式(4)中，然后沿著x方向進行積分，邊界條件為：x=0，p=p1；x=L，p=p2，最終得到如下可壓縮Forcheimer方程：

f=aG+bG2

(7)

其中，f表示壓力的平方在流動方向上的變化率。

(8)

(9)

(10)

實驗采用的多孔介質是多孔石墨，材料實驗參數如表1所示。

表1 多孔介質參數

實驗裝置及多孔介質夾具如圖2所示，多孔介質密封在金屬套管中，用速干膠水T1401密封周邊，防止漏氣。空氣壓縮機提供氣源，并經過了一系列干燥去油處理；減壓閥調節多孔介質前端供氣壓力；氣罐穩定多孔介質前端供氣壓力；差壓表檢測多孔介質前后的壓差；多孔介質前后放置2個溫度傳感器用來檢測溫度變化；流量計為差壓式質量流量計，其內部的多層柵格將氣流狀態調整為層流；為了避免雜質顆粒堵塞柵格，在流量計前端放置了2個雜質過濾器，過濾能力為5 μm。

圖2 實驗示意圖

圖3是試驗采集到的數據。空氣在多孔介質內部的流動主要受到貫通孔隙壁面黏性力和慣性力的作用，這兩種作用可以通過黏性滲透系數和慣性滲透系數來表征。多孔介質空氣軸承使用的材料滲透率一般介于3.13e-15～8.44e-14 m2，在這個數量級下，多孔介質內部的雷諾數很小，黏性力占主導地位。

圖3 質量流量及壓差的關系曲線

圖4為f和質量流量的關系曲線。可以看出，f和質量流量幾乎接近線性的關系，這說明多孔介質內部的氣體流動主要是黏性力主導，也就是滲透率系數起到了主要的作用，而慣性滲透率系數可以忽略。采用最小二乘法對實驗數據進行線性擬合得到方程為：

圖4 f與質量流量的關系曲線

f=9.936e19G2+1.64e17G

(11)

方差為0.999，得到多孔介質的黏性滲透系數為4.7447e-15 m2，慣性滲透系數為312.2390 m-1。

1.2 孔隙度

孔隙度的定義為多孔介質單位體積中孔隙所占的比例。孔隙度進一步細化可以分為貫通孔隙度和孤立孔隙度，如圖5所示。貫通孔隙在流動方向上是貫通的；孤立孔隙指的是多孔介質中閉死的孔隙，與外部大氣沒有任何連通通道。

圖5 孔隙結構細分示意圖

文獻[20-21]計算孔隙度多采用SEM圖像法、壓汞/壓水法、煮沸法等。SEM圖像分析對材料表面二維圖像進行閾值分割出孔隙，進而得到面孔隙度，將面孔隙度視為三維空間上的孔隙度。這種方法沒有考慮孔隙在三維結構上的連通性。

因此，本研究采用micro-CT計算機斷層掃描方法分離多孔介質的貫通孔隙度，其掃描原理如圖6所示。射線源發射的X光穿透樣品，部分射線被樣品所吸收，剩余量會打在探測板上被感應元件探測到。樣品旋轉360°，探測板共采集數百至數千幀數據，之后數據被重構，從而還原出樣品的三維結構。試驗數據采集使用的設備為天津某公司制造的nanoVoxel 3000型X射線三維顯微鏡。

圖6 micro-CT掃描原理示意圖

三維多孔介質圖像能夠顯示孔隙在三維方向上的連通關系，圖像后處理進行中值濾波降低噪聲，改善圖片質量。采用分水嶺算法分割孔隙，通過邏輯運算將孔隙結構細分為貫通孔隙、孤立孔隙。圖7a為貫通孔隙，圖7b為孤立孔隙。

圖7 孔隙結構圖

多孔介質孔隙度分析的關鍵之一在于找到其最小重復單元(Representative Elementary Volume,REV)。REV參考BEAR[22]提出的定義：不斷增大截取的子樣本體積并計算對應的孔隙度，當孔隙度不再隨著子樣本體積V增大而變化時，認為這個體積為REV。圖8顯示了REV尋找的曲線圖，當體積為0.25 mm3時孔隙度φ基本保持恒定。

圖8 樣本最小重復單元及孔隙度的關系曲線

2 靜態特性分析

多孔介質的阻尼作用以沉積項的形式附加在標準控制方程。因此，在計算流場時有如下假設：

(1) 多孔介質內部的堵塞效應不會被考慮；

(2) 多孔介質內部的流動假設為層流；

(3) 多孔介質的孔隙度是各向同性的；

(4) 不考慮多孔介質與激波的相互作用。

單相空氣流動的體積平均質量和動量守恒方程如下所示：

(12)

(13)

熱平衡方程為：

(14)

式中，U—— 速度向量

I—— 單位矩陣

S—— 多孔介質源項

γ—— 貫通孔隙度

Ef—— 總流體能量

Es—— 總固體能量

ρf—— 流體密度

ρs—— 固體密度

keff—— 有效熱傳導率

2.1 網格無關性驗證

采用全結構六面體網格對多孔介質空氣軸承流場特性進行數值計算。設計了3種相同的承載面積和介質厚度不同物理形狀的多孔介質空氣軸承，網格模型如圖9所示。在氣膜厚度為10 μm時，進行網格無關性驗證，驗證結果如圖10所示，其中橫坐標表示網格數量n，縱坐標表示承載力F。可以看出，多孔介質空氣軸承總網格數在8×104左右，繼續增加網格數量對結果幾乎沒有影響。

圖9 不同形狀多孔介質空氣軸承網格

圖10 網格數量及承載力的關系曲線

2.2 靜態特性分析

靜態特性包含承載力和靜剛度。靜剛度指在外部載荷變化過程中，不考慮氣膜的動態變化，單位氣膜厚度變化所需要的力。通常，靜剛度越高則運動部分抗載荷變化的能力越強，保持高精度的能力越強。

流體控制方程為復雜的偏微分方程組，獲得方程組的分析解是很困難的。數值解能夠在一定的精度下為工程應用提供有效參考價值。數值模擬通常需要模型網格化、方程組離散化等一系列操作。其中，網格的質量對迭代計算是否收斂影響很大，各個形狀網格的質量如表2所示。

表2 網格質量

基于正交模擬的原則，除了結構形狀外，其他物理邊界條件保持相同，詳細數值模擬參數如表3所示。

表3 CFD模擬條件

由圖11可以看出，隨著氣膜厚度l的減小多孔介質空氣軸承的承載力F一直增加。這是由于隨著氣膜厚度的減小，氣膜腔的表面積體積值之比增大，阻尼增大，氣膜內的平均壓力升高，空氣的彈性模量增大，進一步壓縮氣膜需要更大的力。

圖11 氣膜厚度及承載力的關系曲線

由圖12可以看出，隨著氣膜厚度的減小，靜剛度會出現1峰值，之后氣膜厚度的再次減小導致靜剛度下降。不同形狀的多孔介質空氣軸承靜剛度峰值及對應的氣膜厚度也有略有差異。其中，圓形多孔介質氣膜厚度為6 μm時靜剛度峰值為82 N/μm，方形多孔介質在氣膜厚度為5 μm時靜剛度峰值為82 N/μm，三角形多孔介質氣膜厚度為5 μm時靜剛度峰值為88 N/μm。相同的截面積情況下，三角形的靜態特性值略高一些。合理設計多孔介質截面的大小以及氣膜的厚度能夠很好均衡承載能力和靜剛度。

圖12 氣膜厚度及剛度的關系曲線

為了進一步理解相同面積不同形狀多孔介質氣浮軸承的承載差異，采用有限體積法計算了氣膜內的壓力場，CFD模擬參數如表3所示，氣膜厚度選定為5 μm，壓力云圖如圖13所示。

圖13 氣膜內壓力場力

沿著圖13中3條對稱線繪制壓力值曲線，如圖14所示,橫坐標表示對稱線長度x，縱坐標表示對應點靜壓力p。可以看出，圓形的高壓區壓力高且保持線最長，正方形高壓區壓力保持較高但保持線短，三角形的高壓區壓力相對不大。由此可見，高壓區壓力大小和面積決定了承載力的大小。

圖14 對稱線上壓力分布

3 結論

完整提出了多孔介質靜態特性預測方法，從多孔介質的物性參數探測、貫通孔隙度計算到基于有限體積的靜態特性數值模擬，得到以下結論：

(1) 采用可壓縮Forchheimer方程獲得多孔介質黏性系數和慣性系數；

(2) 采用micro-CT方法得到多孔介質的最小重復單元，體素在1.48 μm，樣本REV尺寸為0.25 mm3；

(3) 多孔介質空氣軸承的承載能力隨著氣膜的減小而增大，剛度有1個峰值，圓形多孔介質峰值剛度為82 N/μm在6 μm處，方形為82 N/μm在5 μm處，三角形為88 N/μm在5 μm處，相同的截面積情況下，三角形的剛度較高。