非同心型擠壓油膜阻尼器空化流場特性數(shù)值模擬

崔穎， 王清智， 王永亮， 李婷

(大連海事大學 船舶與海洋工程學院，遼寧 大連 116026)

擠壓油膜阻尼器(squeeze film damper，SFD)作為軸系的減振元件，具有體積小、結構簡單、減振效果明顯等優(yōu)點，廣泛應用于航空和艦船燃氣輪機軸系減振結構中。SFD結構分為同心型和非同心型，由于非同心型阻尼器結構中內環(huán)重力載荷沒有彈支結構來平衡，阻尼器的內環(huán)與外環(huán)存在靜態(tài)偏心距。相比之下，非同心型阻尼器內的流動與減振特性較同心型擠壓油膜阻尼器更為復雜，掌握其空化流場特性可為對其減振特性分析奠定基礎。

擠壓油膜阻尼器普遍存在氣穴現(xiàn)象，氣泡的生成、發(fā)展和潰滅會對SFD的油膜壓力分布產(chǎn)生較大程度的影響。國內外學者對阻尼器的氣穴現(xiàn)象做了大量的研究，傳統(tǒng)的求解方法是基于單相流的J-F-O模型，認為正壓區(qū)均為液相，負壓區(qū)均為氣相，采用Elrod算法確定油膜的邊界[1-3]。White[4]通過實驗發(fā)現(xiàn)當內環(huán)的進動偏心距較大時，高壓區(qū)仍會出現(xiàn)氣泡，這表明傳統(tǒng)的求解方法存在不足。Zeidan等[5]和Vance等[6]應用高速攝像技術觀察到阻尼器在高速工作過程中所形成的氣穴為吸入空氣和滑油蒸汽的混合物，因而有必要在擠壓油膜阻尼器建模中對不同的空化類型和機制進行區(qū)分。針對兩端開口型SFD，Diaz等[7]最先把氣泡動力學(Rayleigh-Plesset，R-P)方程的引入到SFD的空化問題研究中，不考慮R-P方程的慣性項、阻尼項和表面張力項，提出了一種考慮從外界吸入空氣的空化模型。Someya[8]在R-P方程中引入了氣泡表面的膨脹黏度，補充并完善了Diaz 在研究時所采用的簡化形式的R-P方程。Gehannin等[9]提出將完全R-P方程應用于油蒸汽空化現(xiàn)象預測中，并分析了R-P方程中各個因素對SFD流場特性的影響。以上研究中所采用的空化模型，都是建立在雷諾方程的基礎上。若所求解的SFD的物理模型稍加復雜，則會存在一定程度的不足。

近年，計算流體力學(computational fluid dynamics, CFD)軟件的快速發(fā)展使得一些學者開始利用它來對SFD的流場特性深入計算分析[10-13]，CFD軟件所包含的豐富計算模型和方法可對空化流場進行精確求解，并且建立的物理模型不受幾何結構限制。本文針對文獻[14]中所采用的非同心型SFD結構，利用Mixture多相流模型、Z-G-B空化模型和動網(wǎng)格技術，通過Fluent軟件中的用戶自定義函數(shù)UDF(user-defined function)編輯SFD內環(huán)進動頻率和半徑，從而對流場不同時刻壓力分布和空化特性進行數(shù)值研究，利用實驗結果對數(shù)值模擬結果進行了驗證。

1 阻尼器流場的數(shù)值建模

1.1 阻尼器三維建模

采用Gambit軟件對文獻[14]中的非同心型SFD結構進行三維建模，SFD的結構示意如圖1所示。

圖1 SFD結構示意[14]Fig.1 Schematic diagram of SFD[14]

阻尼器的結構參數(shù)為：外環(huán)直徑D=140 mm，徑向間隙C=0.385 mm，靜偏心距d=0.016 mm，軸向長度L=30 mm。進油和出油均是軸向方向，基于文獻[14]中的第2組實驗(Ⅱ run)的工況進行邊界條件設置，并將后續(xù)的計算結果與實驗結果作對比。設置出口壓力是環(huán)境大氣壓，進油壓力選取不同的值，潤滑油的動力黏度μ=0.066 Pa·s。

將阻尼器外環(huán)圓心設為坐標軸原點，阻尼器軸向方向為Z軸方向，Z軸負方向進油，正方向出油，該阻尼器的三維模型關于XY面對稱，如圖2所示。阻尼器的外環(huán)定義為靜止壁面，利用Fluent軟件中的UDF對內環(huán)進動頻率和半徑進行編輯。

圖2 SFD三維模型及邊界定義Fig.2 Three-dimensional model and boundary definition of SFD

1.2 多相流模型

Fluent軟件中Mixture模型是一種簡化的多相流模型，通過對連續(xù)方程和動量方程的求解，再結合相對速度方程以及氣體體積分數(shù)方程來模擬多相流動，具有較高的收斂性，故選取該模型來進行模擬。為研究方便，假設液相和氣相位均勻混合狀態(tài)，均不可壓縮，在靠近壁面處流體無滑移，并且不考慮相間的相對速度和滑油黏度的變化。其控制方程為：

(1)

(2)

(3)

?(ρmvm)/?t+·(ρmvm)=-p+

(4)

式中：ρv是氣相的密度；ρl是液相的密度；ρm是混合相的密度；V是氣相的速度；Re和Rc分別是氣體生成和凝結的速率；αv是流場中氣相體積分數(shù)；vm是質量平均速度；F為體積力；μm是混合相的動力黏度；p是壓力。

1.3 空化模型

采用Z-G-B空化模型來求解氣泡的蒸發(fā)和凝結，其忽略流體中的不可凝結氣體，假定全部氣泡尺寸相同，利用氣泡的密度數(shù)(n)來求解單位體積內總的相間質量傳遞率(R)，其表達式為：

(5)

根據(jù)廣義Rayleigh-Plesset方程推導出氣泡動力學方程：

(6)

式中：RB為氣泡半徑；PB為氣泡表面壓力；P為當?shù)剡h場壓力；ρl為液相密度；S為液體表面張力系數(shù)。

若忽略二階項與表面張力項，上式可簡化為：

(7)

將方程(5)與方程(7)聯(lián)立并結合氣相體積分數(shù)，可得：

(8)

若要將其應用于氣泡凝結過程，引入下式：

(9)

式中：F為經(jīng)驗校準系數(shù)，考慮當氣泡體積分數(shù)增加，成核位置的密度減小，用(1-αv)αnuc來替換原有的氣相體積分數(shù)，可得：

當P≤Pv時：

(10)

當P>Pv時：

(11)

式中：氣泡的半徑RB=1 μm，在成核位置處的體積分數(shù)αnuc=5×10-4，蒸發(fā)系數(shù)Fvap=50，凝結系數(shù)Fcon=0.01。在計算過程中本文將空化壓力取為絕對壓力0。

1.4 動網(wǎng)格模型

選取動網(wǎng)格模型中彈簧光順法(smoothing methods)，當邊界處的網(wǎng)格發(fā)生變形時，內部的網(wǎng)格結點也隨之變形，而結點數(shù)目和連通性維持不變，結點位置的更新方程為：

xnew=xold+uΔt

(12)

本文選擇Diffusing模型控制網(wǎng)格的變化，該模型能使網(wǎng)格在運動過程中保持較高的質量，控制網(wǎng)格形變的擴散方程為：

·(γu)=0

(13)

式中：γ是擴散系數(shù)；u表示網(wǎng)格移動的速度。

1.5 無關性驗證

兼顧計算機計算精度、計算速度以及計算的收斂性，需對網(wǎng)格數(shù)量進行無關性驗證，從而選取一個最佳平衡點。采用六面體單元對流場進行網(wǎng)格劃分，將網(wǎng)格數(shù)量分別設置為10萬、25萬、50萬、75萬和100萬，對一個周期內內環(huán)所受的油膜力最大值進行比較，如圖3所示。

圖3 不同網(wǎng)格數(shù)量下的油膜力Fig.3 Oil film force under different grid numbers

可以看出，當網(wǎng)格數(shù)量大于50萬，隨著網(wǎng)格數(shù)量上升油膜力升高幅度變緩，因此后續(xù)計算均采用50萬的網(wǎng)格數(shù)量。其中，在阻尼器的周向布置2 000個結點，在阻尼器的軸向布置50個結點，在阻尼器的油膜厚度方向上布置5個結點。

2 阻尼器空化流場的數(shù)值分析

數(shù)值模擬了不同時刻下隨著內環(huán)進動阻尼器流場油膜壓力與氣相體積分數(shù)分布，氣相體積分數(shù)為流場單位體積中氣相所占的體積比，并分析了進油壓力對流場油膜壓力與氣相體積分數(shù)的影響。

2.1 空化流場壓力與氣相體積分數(shù)的瞬態(tài)分布

將阻尼器靜偏心的位置與內環(huán)進動的起始位置均設置在X軸負方向，通過UDF功能設置內環(huán)沿逆時針方向進動，進動頻率為50 Hz，進動半徑為0.192 5 mm，進油壓力為相對壓力0.02 MPa，截取內環(huán)分別運動到1/4周期、1/2周期、3/4周期和1個周期時刻下的壓力與氣相體積分數(shù)分布云圖，如圖4所示。

圖4 不同時刻油膜壓力與氣相體積分數(shù)分布Fig.4 Distributions of oil film pressure and gas-phase volume fraction at different times

可以看出，隨著內環(huán)進動，在油膜厚度最小處的前方產(chǎn)生了高壓區(qū)，后方出現(xiàn)低壓區(qū)且有氣泡生成，在低壓區(qū)的中間位置存在條狀氣相帶，其寬度沿內環(huán)進動方向逐漸減小，高壓區(qū)域的范圍和油膜壓力峰值以及壓力梯度均在瞬時變化。

由于靜偏心的存在，在不同時刻內環(huán)對油膜的擠壓程度不同，油膜壓力和氣相體積分數(shù)隨之變化。如在t=T時刻，內環(huán)運動到X軸負方向，此時油膜厚度最小，高壓區(qū)油膜壓力達到最大；在t=T/2時刻，內環(huán)運動到X軸正方向，此時油膜厚度最大，高壓區(qū)油膜壓力達到最小。氣穴比為氣相體積占整個阻尼器內流場體積的百分比，其在一個周期內變化情況如圖5，由于油膜壓力的不斷變化，不同時刻整個流場中油膜壓力低于空化壓力的部分隨之變化，因此，滑油空化所產(chǎn)生的氣體占整個流場的氣穴比周期性變化。

圖5 氣穴比變化Fig.5 Variation of cavitation ratio

2.2 進油壓力對空化流場的影響

對阻尼器的進油壓力Pi分別為相對壓力0.2、0.5、1.0和2.0 bar條件下該阻尼器的壓力與氣相體積分數(shù)分布進行了計算模擬。為對比方便，將t=3T/4時刻的三維云圖在Y軸負方向處裁開并延展為二維云圖，向上方向為Z軸正方向，即下方為進油邊、上方為出油邊。圖6給出了不同進油壓力下阻尼器油膜壓力與氣相體積分數(shù)分布。可見，隨著進油壓力不斷升高，高壓區(qū)范圍的壓力值逐漸增大，低壓區(qū)范圍不斷縮小。并且條狀氣相帶的范圍和氣相體積分數(shù)隨進油壓力增大而不斷減少。另外，當進油壓力較小時，高壓區(qū)和低壓區(qū)處于軸向的中間位置，隨著進油壓力不斷升高，阻尼器進出口壓差對空化流場的影響越來越大，導致高壓區(qū)逐漸向進油邊移動，低壓區(qū)逐漸向出油邊移動。

圖6 不同進油壓力下的油膜壓力與氣相體積分數(shù)分布Fig.6 Distributions of oil film pressure and gas-phase volume fraction under different supply pressure

3 數(shù)值模擬與試驗結果對比分析

根據(jù)文獻[14]中的實驗監(jiān)測點布置，下面對這3個監(jiān)測點處的油膜壓力與氣相體積分數(shù)的變化規(guī)律進行分析。選取軸向的中間位置即Z=0平面，在外環(huán)上布置3個監(jiān)測點TP1、TP2和TP3，其相對位置如圖7所示。

圖7 監(jiān)測位置示意[9]Fig.7 Diagram of monitoring location[9]

針對3個測點的油膜壓力變化規(guī)律進行分析，分別設置進油壓力為相對壓力0.2和2.0 bar，得到各監(jiān)測點處的油膜壓力關于內環(huán)進動角度β的變化規(guī)律，并與實驗結果作對比。通過幾何關系推導出各個監(jiān)測點處的油膜厚度關于內環(huán)進動角度的變化，進一步分析油膜壓力與油膜厚度之間的關系，如圖8所示。

由圖8可知，由于3個監(jiān)測點在布置上彼此相差120°的相位，所以其壓力變化曲線也各有120°的相位差，且由于徑向間隙不同，各個監(jiān)測點的油膜壓力的峰值也不相同。油膜壓力在上升階段為緩慢上升過程，下降階段為迅速下降過程。隨著進油壓力由0.2 bar增大到2.0 bar，監(jiān)測點油膜壓力隨之上升，但在下降階段曲線基本重合，且當監(jiān)測點處于阻尼器內壓力場的低壓區(qū)時，由于增大了進油壓力，低壓區(qū)的壓力得到補償，因此壓力曲線的低壓區(qū)域逐漸變窄。從油膜壓力曲線與油膜厚度曲線之間的關系可以看出，在油膜厚度達到最小之前壓力已達到最大值，當油膜厚度為最小值時壓力接近于空化壓力。本文采用的Z-G-B空化模型在油膜壓力低于給定空化壓力時產(chǎn)生空化現(xiàn)象，因而通過數(shù)值方法得到的壓力曲線在低壓區(qū)的壓力峰值以及持續(xù)時間與實驗結果略有不同，表明所采用的空化模型與實際的空化現(xiàn)象仍存在一定的偏差。但總體來看，油膜壓力的數(shù)值仿真結果與實驗結果基本一致。

圖8 不同監(jiān)測位置油膜壓力的計算結果與實驗數(shù)據(jù)對比Fig.8 Comparisons of numerical and experimental data of oil film pressure at different monitoring locations

在此基礎上，進一步數(shù)值模擬分析各個監(jiān)測點氣相體積分數(shù)的變化規(guī)律及其與油膜壓力之間的關系，在進油壓力為0.2 bar的條件下，對3個監(jiān)測點的油膜壓力與氣相體積分數(shù)進行分析，如圖9所示。

圖9 不同監(jiān)測位置油膜壓力與氣相體積分數(shù)的瞬時變化Fig.9 Instantaneous variations of oil film pressure and gas-phase volume fraction at different monitoring locations

由于3個監(jiān)測點處的最小油膜厚度不同，因而各處在產(chǎn)生空化時的氣相體積分數(shù)隨最小油膜厚度的增大而減小。隨著內環(huán)進動，對比3個監(jiān)測點處的油膜壓力變化曲線和氣相體積分數(shù)變化曲線可知：當油膜壓力低于空化壓力時出現(xiàn)氣泡，氣泡的生成過程是較為緩慢的過程；隨著監(jiān)測點處的油膜壓力不斷升高，當其高于空化壓力時氣泡開始減少直至全部消失，氣泡的潰滅過程相比于氣泡生成過程更為迅速。

4 結論

1)基于Mixture多相流模型、Z-G-B空化模型和動網(wǎng)格技術建立了非同心型擠壓油膜阻尼器三維非定?？栈鲌銮蠼饽Ｐ停瑪?shù)值模擬得到由于各個不同時刻內環(huán)對油膜擠壓程度不同，使得油膜壓力與氣相體積分數(shù)分布和阻尼器流場中的氣穴比隨內環(huán)進動而瞬時周期性變化。

2)數(shù)值模擬發(fā)現(xiàn)該型阻尼器的進油壓力對流場的壓力分布和氣穴比影響顯著，進油壓力增大使得油膜高壓區(qū)域范圍增大，油膜壓力峰值及壓力梯度增大，而低壓區(qū)縮小，氣穴比減少；并且隨著進油壓力的增大，流場內高壓區(qū)向進油邊移動，低壓區(qū)向出油邊移動。

3)通過對阻尼器外環(huán)監(jiān)測點瞬時油膜壓力計算結果與文獻中的實驗結果進行對比，表明各監(jiān)測點的模擬結果與實驗基本吻合，驗證了本文所采用的擠壓油膜阻尼器流場數(shù)值模擬方法具有較高的精確性，可用于更復雜結構的擠壓油膜阻尼器空化流場的預測分析。