橡塑往復(fù)密封性能參數(shù)數(shù)值求解流程與物理模型探討

殷圖源 魏大盛 索雙富

(1.北京航空航天大學(xué)能源動力工程學(xué)院 北京100191；2.清華大學(xué)機(jī)械工程系 北京100048)

往復(fù)橡塑密封應(yīng)用范圍廣泛，主要應(yīng)用于液壓缸、柱塞泵、換向閥等液壓元件。物理參數(shù)如液體壓力為0.5 MPa～1 GPa[1]，往復(fù)速度為0.05～5 m/s，液體黏度為0.001～0.1 Pa·s，流量為0.1～10 L/min不等[2]。密封壽命是反映產(chǎn)品可靠性的重要因素之一，也是密封性能的核心問題。影響密封性能有結(jié)構(gòu)設(shè)計合理性、初始裝配預(yù)應(yīng)力、承載力、材料特性、材料配副、表面微織構(gòu)等多種因素[3-4]。

為獲取影響密封性能的參數(shù)，通常采用有限元仿真以與數(shù)值求解方法。仿真即設(shè)定密封初始工作參數(shù)并借用有限元軟件以獲取密封性能參數(shù)；數(shù)值計算即建立密封物理模型，基于已知設(shè)定參數(shù)量，采取數(shù)值方法求解目標(biāo)參數(shù)。盡管仿真方法理論上可獲得影響密封性能的參數(shù)，但缺少靈活性，尤其計算微米級尺度處理如：膜厚計算、表面粗糙度、表面織構(gòu)等，有限元軟件尺度精確性相對較低[4](一般長度單位最小尺度0.01 mm)。另外其處理非穩(wěn)態(tài)流程操作較繁瑣，雖可模擬多個靜態(tài)過程，但求解周期隨迭代次數(shù)非線性增加。而最大缺陷是仿真商業(yè)軟件缺少優(yōu)化密封設(shè)計參數(shù)算法，無法較快求出最佳參數(shù)解集域，只能通過枚舉法逐一設(shè)置，效率低且精度差[5-6]。數(shù)值求解不僅可以揭示密封微觀參數(shù)，還可以優(yōu)化密封性能與工作設(shè)計參數(shù)解集域，為密封具體設(shè)計提供量化參考。

盡管往復(fù)橡塑密封數(shù)值求解方法可在外文文獻(xiàn)中檢索到，但具體難點(diǎn)分析并不全面，大部分都是針對具體問題的研究闡述，文獻(xiàn)多雜不且易理解。因此本文作者分析了目前往復(fù)密封常用數(shù)值求解流程，以提高與實現(xiàn)測試結(jié)果吻合度；同時對求解流程做了進(jìn)一步解釋、總結(jié)、歸納，方便相關(guān)研究人員開展研究。

1 往復(fù)密封性能數(shù)值求解

目前橡塑往復(fù)密封求解方法是2010年佐治亞大學(xué)SALANT等[7]提出的，而在此之前一直沿用德國學(xué)者米勒和納烏[8]提出的定性分析，即采取實驗手段測得密封壓力分布，根據(jù)壓力分布特性提出密封泄漏判據(jù)。該方法未說明對密封結(jié)構(gòu)適用性，且在壓力分布做了較強(qiáng)假設(shè)，僅考慮了極值點(diǎn)與拐點(diǎn)處壓力與其對應(yīng)膜厚數(shù)組；并在力系平衡上做了液體作用力pf與接觸應(yīng)力pc等價假設(shè)[9]，而假設(shè)僅在膜厚不大于1 μm時成立，因此有極大局限性。但SALANT在此基礎(chǔ)上系統(tǒng)地給出了往復(fù)密封數(shù)值求解流程，對密封結(jié)構(gòu)型式、數(shù)值求解流程、Reynolds方程邊界設(shè)置、方程數(shù)值求解、接觸模型建立、表面微織構(gòu)耦合效應(yīng)、潤滑充足與否等方面作了量化分析。

1.1 密封結(jié)構(gòu)型式研究

YANG和SALANT[10]研究的結(jié)構(gòu)型式有斯特封、雙唇U形密封、串聯(lián)密封等線接觸，如圖1(a)、1(b)所示。具體密封結(jié)構(gòu)工作參數(shù)見表1。

圖1(a) 斯特封結(jié)構(gòu) 圖1(b) 雙唇U形密封結(jié)構(gòu)

表1 密封具體結(jié)構(gòu)與工作參數(shù)

SALANT提出的密封結(jié)構(gòu)和工作設(shè)計參數(shù)包含了結(jié)構(gòu)尺寸、壓力載荷、柱塞往復(fù)速度、材料模型、介質(zhì)黏度、零件粗糙度[11-12]，尤其在結(jié)構(gòu)方面由單一O形圈擴(kuò)展到廣義線接觸密封型式，并且考慮了潤滑充足與否等不同工作狀態(tài)，與實際工作狀況一致。但缺陷是未分析密封裝配初應(yīng)力，初始間隙配合缺少量化，未分析唇口傾角位置對膜厚分布影響等，而這些參數(shù)也影響密封工作性能如液體壓力分布、摩擦力、磨損率等。

1.2 數(shù)值求解流程

具體求解流程見圖2，步驟包含限元仿真與數(shù)值模擬2個模塊。有限元仿真步驟：建立密封圈有限元模型，設(shè)定單元類型與材料模量，劃分網(wǎng)格，對密封模型邊界加載壓力，提取與柱塞之間靜態(tài)接觸應(yīng)力psc。

數(shù)值模擬流程：首先通過有限元計算出密封圈節(jié)點(diǎn)靜接觸應(yīng)力psc。計算時接觸模型采用GW模型，假設(shè)粗糙度形貌，將膜厚數(shù)值代入接觸應(yīng)力公式求出接觸應(yīng)力pc。設(shè)置初始截斷膜厚ht，將其代入Reynolds方程求解液體作用壓力pf。采取彈流潤滑求解膜厚方法：根據(jù)力系平衡與壓力收斂準(zhǔn)則，不斷調(diào)節(jié)與修正膜厚以得到最終膜厚。也可以選用擬合影響因子與壓力變量，即影響系數(shù)法求得出膜厚分布[13-15]，總體分為4個過程：提取有限元分析，建立接觸模型，載荷平衡約束求解壓力，影響系數(shù)法求解膜厚。但需要補(bǔ)充的是截斷膜厚與實際膜厚沒有必然聯(lián)系，原因是截斷膜厚未考慮粗糙度波峰幅值分布。

圖2 往復(fù)密封數(shù)值求解流程

1.3 有限元分析

首先建立密封圈與柱塞幾何模型并劃分網(wǎng)格。定義密封圈接觸邊線，將壓力載荷加載到密封圈，求解出靜接觸應(yīng)力psc;其次設(shè)定密封圈與柱塞接觸容差ε，將節(jié)點(diǎn)與柱塞徑向距離小于ε時設(shè)為接觸，提取該接觸數(shù)據(jù)以作為靜態(tài)壓力接觸應(yīng)力邊界，壓力節(jié)點(diǎn)提取方法亦同[16-17]。具體分析及流程見圖3(a)、3(b)，將小于ε的節(jié)點(diǎn)如781、782等提取。但經(jīng)過壓力閾值與接觸容差篩選，節(jié)點(diǎn)序號可能不連續(xù)，需補(bǔ)充相應(yīng)節(jié)點(diǎn)使其連續(xù)，否則無法作為靜態(tài)求解壓力邊界[18]。

圖3(a) 有限元接觸節(jié)點(diǎn)設(shè)置

1.4 Reynolds方程求解

早期米勒和納烏[8]通過實驗測試密封內(nèi)外行程壓力分布并做出定性分析，針對壓力極值點(diǎn)、拐點(diǎn)與邊界約束代入Reynolds方程近似地求出密封膜厚分布，同時根據(jù)壓力分布給出了密封泄漏判據(jù)。SALANT對初始Reynolds方程做了較大修正，建立的方程考慮了粗糙度、空化、流量因子等[10]，如式(1)所示。該方程與潤滑Reynolds方程有實質(zhì)不同。

(1)

式中：φ為流體壓力/密度函數(shù)；H為量綱一油膜厚度；F為空化指數(shù)；HT為量綱一截面平均油膜厚度；φxx、φs,c,x為密封與柱塞綜合表面粗糙度的流量因數(shù)；ζ為量綱一柱塞往復(fù)速度；σ為量綱一平均粗糙度。

具體研究：(1)將柱塞內(nèi)外行程2種狀態(tài)Reynolds方程歸一化，從而避免繁瑣分類討論并且可計算單周期膜厚與壓力分布；(2)采用有限差分法與迎風(fēng)格式求解思路，將實際物理問題轉(zhuǎn)化數(shù)學(xué)PDE問題求解，以求解密封壓力，如圖4(a)、4(b)所示。具體原理為：當(dāng)柱塞處于內(nèi)行程工作狀態(tài)時，左邊靠近液體區(qū)單元信息已知，依次求出右邊靠近空氣區(qū)未知單元信息；當(dāng)柱塞處于外行程工作狀態(tài)時，右邊單元信息已知，依次求出左邊未知單元信息。

由于需要一次性求解所有密封節(jié)點(diǎn)壓力數(shù)值pf，因此需要將差分離散后Reynolds方程未知量pf處理成隱格式，再經(jīng)矩陣分解求出液體壓力。該方法不僅求解時間少于顯格式，而且克服PDE顯格式數(shù)值不穩(wěn)定缺陷[19]。

1.5 方程邊界處理

Reynolds方程邊界需要根據(jù)物理意義確定，為分析空化影響，在Reynolds方程中用邊界變量φ描述潤滑區(qū)的量綱一壓力，用變量F描述方程內(nèi)外行程不同狀態(tài)壓力邊界[20-21]。

內(nèi)行程φ≥0時，F(xiàn)=1和P=φ；P為液體壓力；外行程空化區(qū)：φ<0時，F(xiàn)=0和P=0。

產(chǎn)生空穴的原因是，當(dāng)柱塞處于外行程過程中，由于臨近空氣位置X處液體壓力低于外界壓力，且缺少潤滑介質(zhì)補(bǔ)充[22]，因此產(chǎn)生空穴現(xiàn)象，如圖5所示。

圖5 空化現(xiàn)象產(chǎn)生原因

1.6 仿真與數(shù)值求解區(qū)別

上文論述了數(shù)值求解流程，為更好地理解數(shù)值求解與仿真分析的區(qū)別，表2給出了2種方法的比較。可見，采用數(shù)值求解方法可以更好地量化密封工作設(shè)計參數(shù)，揭示更多密封工作性能。

表2 往復(fù)密封有限元仿真與數(shù)值模擬區(qū)別

2 計算結(jié)果及分析

文獻(xiàn)展示較多采用SALANT方法的計算結(jié)果分析，文中遴選有代表性結(jié)果進(jìn)行分析，如圖6所示。

圖6(a) 速度0.5 m/s下內(nèi)行程壓力分布 圖6(b) 速度0.5 m/s下外行程壓力分布

如圖6(a)所示，柱塞內(nèi)行程運(yùn)行過程中，液體壓力pf在密封位置0～0.5 mm處為0，主要原因是柱塞在下止點(diǎn)處速度為0，成膜較難，接觸應(yīng)力pc較大，因此pf為0。如圖6(b)所示，在柱塞外行程中，在密封位置1.7～2.5 mm處發(fā)生空穴，因此pf為0。此外pf與柱塞往復(fù)速度成正比，與粗糙度成反比。

圖7(a)、7(b)示出了U形密封結(jié)構(gòu)和斯特封結(jié)構(gòu)不同粗糙度下成膜臨界速度。可見線接觸斯特封比U形密封的成膜臨界速度小很多，主要原因是線接觸密封摩擦因數(shù)較低，摩擦接觸域短，而面接觸摩擦接觸域?qū)挘赡ぽ^難。但是粗糙度增加，成膜難度加在的原因，SALANT并未揭示。

圖7(a) U形密封結(jié)構(gòu)不同粗糙度下成膜臨界速度 圖7(b) 斯特封結(jié)構(gòu)不同粗糙度下成膜臨界速度

3 數(shù)值求解流程與數(shù)值模擬方法修正

下面分析數(shù)值求解流程與數(shù)值模擬方法修正，具體在密封結(jié)構(gòu)、數(shù)值求解流程、Reynolds方程建立、邊界確定與求解、接觸模型、載荷力系平衡等方面進(jìn)行探討。

3.1 不同工業(yè)背景應(yīng)用上

目前，橡塑往復(fù)密封主要應(yīng)用在工業(yè)柱塞泵與科研實驗室反應(yīng)裝置上。在泵工況下，當(dāng)處于外行程時，密封不受液膜作用力，因此可以忽略外行程對密封的影響；在反應(yīng)裝置工況下，當(dāng)處于外行程時，密封仍然受液膜作用力。具體如圖8(a)、8(b)所示。所以在不同工況下密封工作情形需要加以區(qū)分，不可混為一談。

圖8(a) 柱塞外行程時密封工作狀態(tài)(泵)

圖8(b) 柱塞外行程時密封工作狀態(tài)(實驗測試裝置)

3.2 密封結(jié)構(gòu)分析上

微觀結(jié)構(gòu)密封尺寸也需要考慮，如圖9(a)、(b)所示，如優(yōu)化密封唇口位置及高壓傾角等參數(shù)時，則以摩擦力為約束量，得出最佳壓力分布。

在結(jié)構(gòu)裝配上需考慮初始裝配預(yù)應(yīng)力，以及與主密封之間裝配公差。主密封與柱塞初始間隙等設(shè)置至關(guān)重要，在工作過程中溫度上升，密封圈徑向膨脹導(dǎo)致間隙減小，或裝配過程中出現(xiàn)偏差等，都易導(dǎo)致主密封與柱塞出現(xiàn)磨損卡滯。而初始間隙過大易出現(xiàn)泄漏，并且泄漏與膜厚呈3次方關(guān)系[1]，因此需要對密封微結(jié)構(gòu)參數(shù)開展深入探討[1]。

圖9(a) 往復(fù)密封系統(tǒng)模型 圖9(b) 往復(fù)密封傾角

3.3 數(shù)值求解流程方面

SALANT提出的方法中，液體膜厚采用影響因子法求解。但實際全彈流潤滑方式求解膜厚時，接觸應(yīng)力pc與液體作用力pf對膜厚有耦合作用，應(yīng)添加彼此約束，以防止計算膜厚時出現(xiàn)較大偏差。此外有限元仿真應(yīng)考慮動態(tài)接觸密封與柱塞間隙摩擦力、柱塞內(nèi)外行程工作狀況等，具體修正求解流程見圖10。

由于柱塞是往復(fù)運(yùn)動，因此運(yùn)動狀態(tài)等需要量化，可借助液壓仿真與流體商業(yè)軟件，共同模擬出液體對密封作用力。線接觸橡塑密封由于彈性模量相對較低，唇口處網(wǎng)格分布對其壓力分布有重要影響。同時柱塞往復(fù)速度、液體溫度等影響摩擦，因此需要調(diào)用子程序模塊建立動態(tài)摩擦因數(shù)[25]，具體步驟不再闡述。

3.4 Reynolds方程建立、邊界確定與求解

膜厚求解如圖11所示。膜厚求解可參考彈流潤滑方法，先假設(shè)初始膜厚與初始壓力，不斷調(diào)節(jié)膜厚計算出壓力；當(dāng)計算壓力與承載力誤差小于許用值時，則停止迭代。上述求解出的膜厚是截斷膜厚ht，如何轉(zhuǎn)化為實際膜厚待探討。而SALANT求解采用影響因子法求解，脫離了Reynolds方程求解膜厚，只適合低壓液體場合，具有局限性。該方法求解出的同樣是截斷膜厚ht，與實際膜厚h關(guān)系缺少闡述，因此需要對2個膜厚關(guān)系進(jìn)行關(guān)聯(lián)，以簡化參數(shù)。

改進(jìn)后的膜厚求解流程，避免了之前求解的封閉性，接觸應(yīng)力與液體作用力相互約束，最后求解迭代；并且將截斷膜厚ht與實際膜厚h關(guān)聯(lián)，減少了未知變量數(shù)量。

3.5 接觸模型方面

采用GW模型，將粗糙度表面中心簡化為圓弧，圓弧中心到柱塞徑向距離lc采取高斯分布，但由于粗糙度有較大的隨機(jī)性，因此lc1應(yīng)該按公式(2)積分平均。考慮到lc最值對實際膜厚h有顯著影響，因此應(yīng)該將lc最值考慮如圖12所示。

lcmax=max(z(x)|L(x)|)

lcmin=min(z(x)|L(x)|)

lc=(lc1+c1·lcmax+c2·lcmin)/3

(2)

式中：L為密封表面粗糙峰到柱塞徑向距離，也稱為實際膜厚；c1、c2為粗糙峰值影響系數(shù)，c1>>c2；z(x)為表面粗糙峰概率密度-高斯分布。

接觸域容差設(shè)置：GW模型為固體表面直接接觸，但實際密封圈之間存在液體，因此需要考慮微流體之間潤滑影響。由于柱塞往復(fù)運(yùn)動是動態(tài)過程，接觸容差ε設(shè)定應(yīng)與往復(fù)速度、液體黏度正相關(guān)，應(yīng)該對不同條件下接觸域給予動態(tài)設(shè)計，而非根據(jù)經(jīng)典膜厚比判據(jù)設(shè)定。

3.6 密封力學(xué)平衡

如圖13所示，根據(jù)密封徑向力平衡，考慮粗糙度影響[10]，則密封圈靜態(tài)接觸應(yīng)力psc等于液體作用力pf與接觸應(yīng)力pc之和,即:

psc=pf+pc

(3)

由于粗糙峰對接觸應(yīng)力與液體作用力均有影響，彼此間權(quán)重與潤滑充足性有關(guān)，因此需要修正公式(3)。引入系數(shù)并考慮速度、溫度、液體黏度等的影響，同時不同工況下柱塞外行程液體壓力不同，即對于柱塞泵外行程壓力較低并為常數(shù)，而實驗測試密封性能裝置外行程壓力為變量[25]，因此接觸容差不能直接歸一化為式(4)，如圖14所示。

圖13 密封圈徑向力載平衡

圖14 改進(jìn)密封圈平衡力系

psc=αpc+βpf+I1(pc,pf,γ,cl，wc)·I2(wc)

(4)

約束:

αpc≤psc,βpf≤psc；αpc≤βpf潤滑充足，βpf≤αpc潤滑不充足；

α:(U(t)m2,μ(T)m3,Fm1,Zm3,cl)；

β:(Um1(t),μm2(T),Fm4,Z-m3,cl,S)

式中：S為柱塞工作狀態(tài)；U為柱塞往復(fù)速度；μ為液體黏度；F為載荷；Z為粗糙峰均值；i為耦合項；m1～m4為影響敏度，且m1>m2>m3>m4>0；γ為耦合因子；cl為潤滑充足系數(shù)；wc為柱塞工作狀況。

在柱塞與密封圈裝配過程中存在偏心誤差，因此需要建立三維模型。但由于三維模型徑向膜厚分布不均勻，如圖15所示，因此需要同時分析徑向和軸向膜厚分布域。粗糙度由于具有隨機(jī)性，因此需要建立三維分布，運(yùn)用保角變換投影以等效建立模型，具體流程不再闡述。

圖15 考慮偏心密封端面模型

針對SALANT提出的往復(fù)密封流程與數(shù)值模擬方法，黃樂等人[14]在此基礎(chǔ)上提出需要考慮老化影響并對密封結(jié)構(gòu)做了參數(shù)優(yōu)化。但盡管如此，目前仍有不完備之處，在求解流程、仿真技術(shù)、Reynolds方程、數(shù)值模擬方法、模型建立等方面的修正探討見表3。

表3 GW模型修正

4 往復(fù)密封模型修正

4.1 界面科學(xué)的引入

SALANT提出了GW接觸模型，但并未從微尺度角度表征接觸模型，更未從機(jī)制方面深入揭示粗糙度對往復(fù)密封的影響，而近年來材料界面學(xué)科發(fā)展較快，因此借助其可解決該問題。建立往復(fù)密封模型需考慮的因素如圖16所示。為更好地量化目前GW模型不完備之處，提出修正需考慮的因素，具體見表3。

圖16 影響GW模型因素

清華大學(xué)胡松濤[26]針對密封表面粗糙度處理，采用統(tǒng)計學(xué)手段將柱塞與密封圈粗糙度分離解耦，如圖17所示，采用傅里葉快速變換表征形貌，將粗糙度形貌分解重構(gòu)，進(jìn)行軸向和徑向形貌表征，更符合實際加工情形，避免之前密封與柱塞粗糙度歸一化假設(shè)。但此研究僅僅基于機(jī)械密封并只適用于較低壓力高速場合，與往復(fù)密封工作場合不符，因此不能完全照搬，但該文獻(xiàn)提供了一種有價值的研究方法，可以參考。

圖17 粗糙表面分離解耦

如圖18所示，GW模型需考慮干摩擦情形、黏彈性材料、Reynolds方程、表面粗糙度、微觀、宏觀接觸應(yīng)力共同影響，而非僅和接觸應(yīng)力關(guān)聯(lián)。

圖18中，pcdry為干摩擦接觸應(yīng)力；u為柱塞往復(fù)速度；E為彈性體作用外力功；Et為應(yīng)變率；Ek為彈性應(yīng)變釋放變能。

4.2 修正GW模型算法流程

考慮黏彈性模型密封數(shù)值求解流程：輸入密封結(jié)構(gòu)參數(shù)，材料設(shè)置黏彈性，設(shè)置黏彈性時間步長，進(jìn)行有限元求解；將求解的非穩(wěn)態(tài)接觸應(yīng)力，代入數(shù)值計算Reynolds邊界，根據(jù)動態(tài)設(shè)定接觸域準(zhǔn)則，求解出非穩(wěn)態(tài)情形接觸應(yīng)力與液體作用力；然后依次推進(jìn)直至求解結(jié)束，以獲取考慮黏彈性時間效應(yīng)的密封圈膜厚、承載力、泄漏率等參數(shù)，從而更接近密封實際工作狀態(tài)，并能較好地揭示密封性能機(jī)制。

5 結(jié)論

(1)闡述SALANT提出的橡塑往復(fù)密封求解流程與數(shù)值模擬方法發(fā)展，分析流程特點(diǎn)、難點(diǎn)及不完備之處，指出該方法較系統(tǒng)地建立往復(fù)密封求解流程與數(shù)值求解流程，極大地完善了往復(fù)密封理論。

(2)對密封結(jié)構(gòu)、數(shù)值求解流程、Reynolds方程邊界及數(shù)值求解、密封載荷力系平衡等方面，對橡塑往復(fù)密封求解流程進(jìn)行了修正，使其更接近實際工作狀況。

(3)根據(jù)前沿界面科學(xué)，修正往復(fù)密封目前GW模型，分析往復(fù)密封接觸模型需要考慮的因素與求解流程，并對修正接觸GW模型建立數(shù)值求解流程，結(jié)果表明，考慮黏彈性GW模型修正更接近往復(fù)密封實際工作性能。