金玉峰, 段鸞芳, 吳偉偉

(1.江蘇省常州技師學(xué)院 智能制造學(xué)院, 江蘇 常州 213032; 2. 南京工業(yè)職業(yè)技術(shù)大學(xué) 電氣工程學(xué)院, 江蘇 南京 210023; 3.揚(yáng)州大學(xué) 機(jī)械工程學(xué)院, 江蘇 揚(yáng)州 225127)

引言

液壓油是數(shù)控機(jī)床加工過程中不可或缺的材料,主要起到潤滑與冷卻的作用[1-2]。當(dāng)液壓油用于導(dǎo)軌、絲杠及軸承(如靜壓導(dǎo)軌、靜壓絲杠及靜壓軸承)中時,還需要具有一定的承載能力[3-7]。靜壓轉(zhuǎn)臺是數(shù)控機(jī)床中的常見部件,其以靜壓軸承為支承,液壓油作為工作介質(zhì),在壓力作用下形成壓力油膜,具有較好的承載剛性[8-9]。

靜壓轉(zhuǎn)臺的結(jié)構(gòu)是影響剛度和承載能力等的主要因素。張瑞濤等[10]利用Fluent軟件對節(jié)流器孔徑及軸承間隙對剛度的影響進(jìn)行了分析,基于分析結(jié)果設(shè)計了具有較高剛度的靜壓轉(zhuǎn)臺。仝志偉等[11]分析了支承布局對平面度誤差的影響,并對支承布局進(jìn)行重構(gòu),可減小平面度誤差,提高轉(zhuǎn)臺的承載能力。王京等[12]發(fā)現(xiàn)雙環(huán)形油腔較傳統(tǒng)油腔具有更好的承載能力。王永柱等[13]對扇形轉(zhuǎn)臺中油腔形狀對承載能力的影響進(jìn)行了分析對比,結(jié)果表明口字型油腔具有更好的承載能力。張慶鋒[14]分析研究了位移率和間隙對油膜承載能力的影響,根據(jù)結(jié)果優(yōu)化了轉(zhuǎn)臺結(jié)構(gòu)。申峰等[15]具體分析了封油邊形狀特征參數(shù)對轉(zhuǎn)臺承載能力的影響,利用分析結(jié)果優(yōu)化封油邊參數(shù)可顯著改善承載能力。MICHALEC M等[16]利用CFD模型和MATLAB提出一種雙參數(shù)優(yōu)化方法,對靜壓軸承的油腔形狀進(jìn)行了優(yōu)化,可有效降低能量損失。ZHANG Yanqin等[17]針對雙矩形腔的靜壓轉(zhuǎn)臺中載荷對油膜性能的影響展開了研究,研究結(jié)果表明,載荷大小與油膜的平均溫升之間存在非線性關(guān)系。YADAV S K等[18]利用有限元分析對比分析了圓形、橢圓形、方形及矩形油腔對油膜剛度系數(shù)的影響,結(jié)果表明油腔結(jié)構(gòu)對油膜承載能力具有重要影響。

在靜壓轉(zhuǎn)臺的實(shí)際運(yùn)行工作過程中,通過實(shí)驗(yàn)發(fā)現(xiàn)液壓油的油液品質(zhì)(油液的黏度)對機(jī)床零部件運(yùn)作時的性能也具有一定影響[19-21]。液壓油在使用過程中,隨著時間的推移,油液的黏度會發(fā)生變化,需要了解該變化對承載能力的影響。本研究將以承載能力較好的雙環(huán)形油腔結(jié)構(gòu)為研究對象,利用格子玻爾茲曼方法(LBM)分析討論封油邊厚度對油液黏度的時變特性的影響,為后續(xù)選擇合適的油液黏度和封油邊厚度提供理論依據(jù)。

1 流變測試

表1 兩種油品流變數(shù)據(jù)擬合結(jié)果

圖1 兩種油品實(shí)驗(yàn)數(shù)據(jù)與擬合結(jié)果對比

根據(jù)表1結(jié)果,將其與流變實(shí)驗(yàn)的數(shù)據(jù)相比較,結(jié)果顯示出較好的一致性,且表1中的相關(guān)系數(shù)無限接近于或等于1,因此擬合獲得的流變模型是可信的。由表1結(jié)果可知,1#和2#油品均呈現(xiàn)出非牛頓特性,具體為帶屈服特征的Herschel-Bulkley型流體,僅當(dāng)剪切應(yīng)力大于初始屈服應(yīng)力時,流體才會發(fā)生剪切運(yùn)動。

2 適應(yīng)Herschel-Bulkley流體的LBM方法

2.1 虛外力項LBM方法

LBM的演化方程為:

(1)

為了解決由非牛頓特性引起的仿真離散及精度差的問題,引入“虛外力項”用以描述Herschel-Bulkley特征,對于一般含外力項的LBM演化方程為[22]:

feq(r,t)]+δtF′

(2)

式(2)中包含了LBM的兩個主要過程——粒子遷移和粒子碰撞,其中粒子遷移具體可描述為:

f(r+eiδt,t+δt)=f+(r,t)

(3)

粒子碰撞具體可描述為:

(4)

式中,r—— 位置矢量

t—— 某時刻

f—— 粒子分布函數(shù)

feq—— 平衡態(tài)分布函數(shù)

其具體表達(dá)式為:

(5)

式中,ρ—— 格子密度,通常取1

u—— 速度矢量

ωi—— 權(quán)重系數(shù)

對于LBM的D2Q9模型,其具體取值如下:

(6)

其中,ei為D2Q9模型的速度配置,具體取值為:

(7)

(8)

黏度一般可根據(jù)流變模型計算獲得,對于牛頓流體,黏度不隨剪切率的變化而變化,其為恒定常數(shù),則松弛時間也為恒定值;而對于非牛頓流體,黏度隨剪切率的變化而變化,因此松弛時間為一變化值。對于液壓油所屬的Herschel-Bulkley流體,其流變模型為:

(9)

式中,τh—— Herschel-Bulkley流體的剪切應(yīng)力

τh0—— 初始屈服應(yīng)力

μh0—— 黏度系數(shù)

n—— 冪律指數(shù)

為了方便計算,前人針對Herschel-Bulkley流體提出了修正黏度-剪切率方程為:

(10)

(11)

其中,DII為應(yīng)變率張量第二不變量,其可通過下式獲得:

(12)

式中,l—— 計算時的維數(shù),此處l=2

Sαβ—— 應(yīng)變率張量

其具體計算公式為:

(13)

對于式(4)中的虛外力項,具體計算方式如下:

(14)

(15)

結(jié)合式(14)和式(15),可知離散虛外力項具體可以描述為:

(16)

2.2 理論驗(yàn)證

為了驗(yàn)證上述提出的虛外力項方法的可靠性,利用Herschel-Bulkley流體在泊肅葉流中的理論解進(jìn)行驗(yàn)證,其理論解為[23]:

(17)

式中,H為泊肅葉流中兩平行平板之間的距離,yτ為一臨界點(diǎn),表征了剪切應(yīng)力是否達(dá)到屈服應(yīng)力,具體值為yτ=τh0/(?p/?x),?p/?x為在x方向上的壓力梯度。具體參數(shù)設(shè)置如下:H=1,平板長度L=1,壓力梯度?p/?x=-4×10-3,黏度系數(shù)μh0=0.1,冪律指數(shù)n=0.9,初始屈服應(yīng)力τh0=6×10-4。在數(shù)值模擬計算時,格子數(shù)取為200×200,此處僅為了驗(yàn)證模型的可靠性,因此直接使用了無量綱數(shù)進(jìn)行模擬與對比。經(jīng)過仿真,對比結(jié)果如圖2所示,對比結(jié)果表明提出的虛外力項LBM具有較好的可靠性,可用于Herschel-Bulkley流體的仿真中。

圖2 泊肅葉流LBM數(shù)值解與理論解的對比

2.3 靜壓轉(zhuǎn)臺仿真步驟

將上述提出的虛外力項LBM用于液壓油在靜壓轉(zhuǎn)臺中的流動分析,其具體過程如下。

(1) 初始化參數(shù)。設(shè)定初始參數(shù),主要包括格子數(shù)(由轉(zhuǎn)臺的尺寸決定)、入口壓力或入口速度、格子密度、初始粒子分布函數(shù)、初始黏度、初始松弛時間等。在此步驟中,需進(jìn)行單位轉(zhuǎn)換,LBM是一種無量綱方法,此處以雷諾數(shù)為準(zhǔn)則數(shù),由實(shí)際尺寸和格子大小可確定長度比例,自定義一較小格子速度,由實(shí)際速度值和格子速度可確定速度比例,基于式(18),保證實(shí)際雷諾數(shù)與格子雷諾數(shù)相等,即可獲得運(yùn)動黏度比例,以此類推可獲得其他比例,進(jìn)而計算出其他格子物理量,完成實(shí)際物理量與仿真格子參數(shù)之間的轉(zhuǎn)換[24]:

(18)

其中,v為特征速度,d為特征長度,υ為運(yùn)動黏度。

(2) 初始化平衡態(tài)分布函數(shù)。在初始階段可以自定義一個函數(shù)作為初始平衡態(tài)分布函數(shù),后續(xù)根據(jù)式(5)實(shí)際計算迭代獲得新的平衡態(tài)分布函數(shù)。

(3) 執(zhí)行碰撞步式(4)和遷移步式(5)。在執(zhí)行碰撞步時,需根據(jù)式(13)～式(16)計算虛外力項,用以表征Herschel-Bulkley流體的非牛頓特性。

(4) 計算松弛時間。由式(13)計算應(yīng)變率張量,結(jié)合式(11)、式(12)計算剪切率,由式(8)和式(10)計算下一次迭代所需的松弛時間。

(5) 邊界處理。不同于傳統(tǒng)有限元方法的邊界處理,LBM的邊界處理主要是對邊界粒子的狀態(tài)進(jìn)行定義,目前比較成熟的方法有反彈格式、非平衡外推等,此處選用非平衡外推進(jìn)行邊界處理。

(6) 判斷是否計算結(jié)束。設(shè)定一較小值為對照標(biāo)準(zhǔn),以某一物理量為基準(zhǔn),計算上一循環(huán)與本次循環(huán)的結(jié)果差值的絕對值,當(dāng)此值小于上述設(shè)定的對照標(biāo)準(zhǔn)時,即認(rèn)為計算已滿足條件,跳出迭代計算,否則進(jìn)行下一輪迭代循環(huán)計算。

(7) 計算宏觀物理量。參考式(19)進(jìn)行物理量計算,并再次進(jìn)行量綱轉(zhuǎn)換獲得宏觀物理量。

(19)

式中,ρ′ —— 仿真中油液密度

u—— 仿真中獲得的速度

3 不同靜壓轉(zhuǎn)臺結(jié)構(gòu)中的流動分析

3.1 靜壓轉(zhuǎn)臺結(jié)構(gòu)

基于雙環(huán)形靜壓轉(zhuǎn)臺較傳統(tǒng)靜壓轉(zhuǎn)臺具有更好的承載能力,采用如圖3所示的雙環(huán)形油腔結(jié)構(gòu)的靜壓轉(zhuǎn)臺為研究對象。

圖3 雙環(huán)形靜壓轉(zhuǎn)臺結(jié)構(gòu)

油液入口1的半徑r1設(shè)置為4 mm,油液入口2和3的半徑r2均設(shè)置為2 mm,內(nèi)圓環(huán)的內(nèi)圓半徑R2設(shè)置為20 mm,外圓半徑R3設(shè)置為24 mm,外圓環(huán)的內(nèi)圓半徑R4設(shè)置為46 mm,外圓半徑R5設(shè)置為50 mm,油腔的總深度為H,封油邊深度為h,為了研究油液黏度時變特性下深度對油液在腔體內(nèi)部流動的影響,其具體設(shè)置如表2所示。由于雙環(huán)形油腔是一典型的對稱結(jié)構(gòu),為提高運(yùn)算效率,在以下的仿真中均取一半結(jié)構(gòu)進(jìn)行分析,仿真時徑向的格子數(shù)設(shè)置為1000,軸向的格子數(shù)由具體尺寸等比例設(shè)定。

表2 油腔深度的設(shè)置

3.2 1#油品流動分析

根據(jù)2.1節(jié)及2.3節(jié)中的描述,結(jié)合MATLAB 2017b軟件將1#油品的流變方程代入LBM方法中進(jìn)行仿真分析,設(shè)定油液在3個入口的速度均為0.1 m/s,可以得到如圖4～圖7的結(jié)果,其中圖4為速度云圖,圖5為頂層油膜速度大小的分布情形,圖6為頂層油膜速度矢量的角度情形,圖7為頂層油膜的壓力分布結(jié)果。通過對比圖4a～圖4c可知,在保持油腔深度h2(h2=H-h)不變的情形下,不斷增加h,速度云圖的變化不明顯。根據(jù)圖4d～圖4f可以發(fā)現(xiàn),在保持h不變的情形下,不斷增加h2,速度云圖差異較為明顯,在入油口1處,由1個主渦逐漸發(fā)展成2個主渦相連,在入油口2處,主渦強(qiáng)度逐漸增強(qiáng)。在環(huán)形臺階上方,速度增強(qiáng)的現(xiàn)象逐漸減弱。

圖4 1#油品在6種結(jié)構(gòu)參數(shù)下的速度云圖

圖5 1#油品頂層油膜的速度大小分布

圖6 1#油品頂層油膜的速度矢量角度

圖7 1#油品頂層油膜的壓力分布

圖5a和圖5b分別描述了Case 1～Case 3和Case 4～Case 6情形下頂層油膜的速度大小分布情況,兩圖具有相似的變化趨勢,在兩個環(huán)形臺階上方的頂層油膜具有較大的速度,且入油口2輸入的油液增強(qiáng)了其左右兩側(cè)環(huán)形臺階上方頂層油膜的速度,造成同一臺階上方的速度大小有較明顯的不同。圖5b表明,相同的h下,較小的h2會導(dǎo)致環(huán)形臺階上方出現(xiàn)較大的速度,與圖4d～圖4f展現(xiàn)出的速度云圖相吻合。圖6中速度矢量與水平線的夾角展示了近乎一致的變化,該角度的變化與形成的渦流在頂層處的形狀相關(guān),由圖4可看出在頂層油膜渦流具有相似形狀,與圖6量化的角度大小基本吻合。

圖7的壓力分布結(jié)果間接展示了轉(zhuǎn)臺的承載能力,6種情形的結(jié)果均表明,中心凹槽的壓力遠(yuǎn)大于環(huán)形凹槽的壓力,間接說明了入油口2處油液的輸入形成了外圍渦流大大增加了中心凹槽的壓力,起到了增壓的作用,提升了轉(zhuǎn)臺的承載能力。此外,通過對比發(fā)現(xiàn)H值越小,頂層油膜所受的壓力越大,且呈現(xiàn)指數(shù)級增長的趨勢,中心凹槽的壓力與環(huán)形凹槽的壓力差值也越大。

3.3 2#油品流動分析

圖8～圖11為2#油品仿真分析的結(jié)果,圖8為速度云圖,其中圖8a～圖8c中,保證h2不變,逐漸增大封油邊h,速度云圖未有明顯差異。圖8d～圖8f中,保證h不變,h2較小時,中心凹槽內(nèi)的渦流較明顯,隨著h2逐漸增大,中心凹槽內(nèi)的渦流逐漸減弱,但均比圖8a～圖8c中的渦流明顯。

圖8 2#油品在6種結(jié)構(gòu)參數(shù)下的速度云圖

圖9 2#油品頂層油膜的速度大小分布

圖10 2#油品頂層油膜的速度矢量角度

圖11 2#油品頂層油膜的壓力分布

圖9中展示了頂層油膜速度大小的分布情況,圖9a表明隨著h的逐漸增大,環(huán)形臺階上方的速度逐漸減小,圖9b表明隨著h2的逐漸增大,頂層油膜的速度在徑向分布上均有一定程度的減小,在環(huán)形臺階上方的速度減幅較大。圖10中所展示的速度矢量角度表明在6種情形下頂層油膜速度與水平線夾角的變化趨勢幾乎一致,與圖8所示的速度云圖基本吻合。

圖10中展示了頂層油膜受壓情況,由圖可知,入油口2中油液的輸入同樣增加了中心凹槽內(nèi)的壓力,起到了增壓作用,且Case 4～Case 6的結(jié)果表明,H較小時,壓力的增幅較明顯。

通過綜合對比1#油品與2#油品的仿真結(jié)果可以發(fā)現(xiàn),隨著時間的推移,油液中混入的雜質(zhì)等改變了油液的流變特性,增加了油液的黏度,對比圖3與圖7、圖4與圖8可知,油液黏度的增加使油液的流速明顯減小,形成的主渦強(qiáng)度較弱。對比圖6與圖10的壓力分布可知,油膜的承壓能力也大大減弱。

4 結(jié)論

經(jīng)過實(shí)際加工過程觀察及流變實(shí)驗(yàn)分析發(fā)現(xiàn),液壓油在經(jīng)過一段時間的使用后,油液的流變特性發(fā)生了明顯變化,且剛使用與使用一段時間后的油液均呈現(xiàn)出典型的非牛頓特性,針對非牛頓油液流變隨時間變化的特點(diǎn),提出了一種虛外力項LBM方法,對油液在雙環(huán)形轉(zhuǎn)臺的內(nèi)部流動進(jìn)行分析,通過仿真結(jié)果對比后,可以得出以下結(jié)論:

(1) 雙環(huán)形轉(zhuǎn)臺的封油邊深度h變化時,速度云圖及頂部油膜速度變化主要集中在環(huán)形臺階上方,但壓力分布差異較明顯,隨著封油邊深度h的增加,入油口2油液的輸入對中心凹槽的增壓作用減弱;

(2) 雙環(huán)形轉(zhuǎn)臺的油腔深度h2不斷增大時,兩種油品的頂層油膜速度沿徑向分布減小,且在環(huán)形臺階上方的速度減幅較大。兩種油品中心凹槽內(nèi)的速度云圖有明顯不同,承載能力呈現(xiàn)“斷崖式”減小;

(3) 液壓油的流變特性隨時間的推移發(fā)生變化,黏度不斷增大,減弱了油液的流動,導(dǎo)致2#油品情形下,頂層油膜的速度減小,承載能力也大幅減弱。對于剛使用的1#油品,隨著油腔深度增加,中心凹槽的主渦區(qū)域不斷擴(kuò)大,而2#油品的中心凹槽主渦區(qū)域不斷減小,此現(xiàn)象也是由油液黏度的變化導(dǎo)致。為保證轉(zhuǎn)臺的正常工作,應(yīng)定期更換轉(zhuǎn)臺用液壓油。