王 民，王樹華，孫鐵偉，石 龍

(1.北京工業(yè)大學(xué)機(jī)械工程與應(yīng)用電子技術(shù)學(xué)院，先進(jìn)制造技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室，北京 100124；2.電火花加工技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室，北京 100191)

現(xiàn)代制造業(yè)水平的飛速提升導(dǎo)致了人們對(duì)數(shù)控機(jī)床整體性能的要求越來(lái)越高，滾珠絲杠副是數(shù)控機(jī)床的重要傳動(dòng)部件，與數(shù)控機(jī)床整體的可靠性緊密關(guān)聯(lián)，在數(shù)控機(jī)床、自動(dòng)化裝備、車輛交通及航空航天等領(lǐng)域常常用到[1].因此，對(duì)滾珠絲杠副的可靠性分析工作得到越來(lái)越多的科研機(jī)構(gòu)和研究人員的重視.

近年來(lái)，隨著可靠性工程的不斷發(fā)展，相關(guān)的數(shù)據(jù)樣本量較小.對(duì)于滾珠絲杠副產(chǎn)品，由于價(jià)格昂貴、可靠性水平較高，而且試驗(yàn)周期長(zhǎng)，很難進(jìn)行大樣本的統(tǒng)計(jì)試驗(yàn)，甚至有的專項(xiàng)計(jì)劃只生產(chǎn)1至2臺(tái)產(chǎn)品.綜合考慮上述情況，本文通過引入小子樣、加速壽命試驗(yàn)的方法完成對(duì)單樣本滾珠絲杠副的壽命評(píng)估.

在實(shí)際工況下，數(shù)控機(jī)床進(jìn)給系統(tǒng)的滾珠絲杠在加減速運(yùn)行過程中，潤(rùn)滑會(huì)變差，加劇接觸部件的磨損，降低其精度和壽命[2]，進(jìn)而導(dǎo)致數(shù)控機(jī)床的可靠性水平下降、加工質(zhì)量變差.另外，根據(jù)Archard[3]磨損理論知道外界負(fù)載應(yīng)力是影響磨損量的主要因素之一.因此本文將綜合考慮轉(zhuǎn)速和負(fù)載雙應(yīng)力對(duì)滾珠絲杠副可靠性壽命的影響進(jìn)行加速壽命試驗(yàn)，然后通過對(duì)試驗(yàn)數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析得到在加速應(yīng)力條件下的壽命分布模型，最后利用加速方程外推出在正常應(yīng)力水平下的壽命特征，完成對(duì)滾珠絲杠副的壽命評(píng)估[4-5].

為了提高小子樣單樣本滾珠絲杠副的可靠性評(píng)估精度，本文使用虛擬增廣法對(duì)加速壽命試驗(yàn)得到的故障數(shù)據(jù)進(jìn)行等價(jià)擴(kuò)充，根據(jù)虛擬增廣法擴(kuò)充得到的數(shù)組和原試驗(yàn)數(shù)據(jù)的分布形式相同.然后對(duì)增廣后的數(shù)組進(jìn)行統(tǒng)計(jì)分析，從而達(dá)到提高試驗(yàn)評(píng)估精度的目的.

常用的可靠性壽命及故障分布模型有正態(tài)分布、指數(shù)分布、威布爾(Weibull)分布[6].其中Weibull分布是一類廣泛應(yīng)用于機(jī)械、電子、航天航空等領(lǐng)域的壽命分布模型，是目前較為常用的壽命模型之一.Weibull分布有多種參數(shù)估計(jì)方法[7]：其中圖估計(jì)法易受人為因素影響，精度偏低；極大似然估計(jì)法求解煩瑣，常常需要借助計(jì)算機(jī)；線性無(wú)偏估計(jì)法不適用于小子樣情況，在大樣本條件下還需借助查閱專用表格，雖然步驟簡(jiǎn)單，但表格量大且在無(wú)表格處無(wú)法使用此方法；最小二乘法求解簡(jiǎn)便，具有唯一最優(yōu)解和好的解析性質(zhì)，是Weibull分布參數(shù)估計(jì)的一種常用方法.

綜合以上特點(diǎn)，本文針對(duì)小子樣單樣本滾珠絲杠副提出一種雙應(yīng)力三水平的加速壽命試驗(yàn)設(shè)計(jì)和統(tǒng)計(jì)分析方法，利用虛擬增廣法和最小二乘法完成在3個(gè)步進(jìn)加速雙應(yīng)力水平下的參數(shù)估計(jì)和可靠性評(píng)估，然后聯(lián)立求解出加速方程，進(jìn)而計(jì)算出滾珠絲杠副在正常應(yīng)力條件下的壽命指標(biāo)，完成對(duì)單樣本滾珠絲杠副的可靠度壽命分析.

1 加速壽命試驗(yàn)?zāi)Ｐ?/h2>

1.1 基本假定

滾珠絲杠副具有加速性的3個(gè)條件為：失效分布同一性、失效機(jī)理一致性和失效過程規(guī)律性[8].根據(jù)以上具有加速性的3個(gè)條件，本文提出3個(gè)基本假定.

假定1：樣本在任意應(yīng)力下的壽命T均服從Weibull分布，即T～W(m,η)，累積分布函數(shù)為

(1)

式中：η>0為尺度參數(shù)，表示特征壽命；m>0為形狀參數(shù).

假定2：樣本在任意應(yīng)力下的失效機(jī)理不變，由于Weibull分布中的形狀參數(shù)反映失效機(jī)理，因此有

m1=m2=…≡m

(2)

假定3：樣本特征壽命與試驗(yàn)雙應(yīng)力之間的關(guān)系模型選用Nelson廣義對(duì)數(shù)線性模型，其數(shù)學(xué)表達(dá)式為

lnη=β0+β1φ1(V)+β2φ2(P)

(3)

式中：η為特征壽命；φ1、φ2常為應(yīng)力的已知函數(shù)；β0、β1、β2為待估參數(shù)；V、P分別為轉(zhuǎn)速和負(fù)載.

考慮到本文的2個(gè)應(yīng)力均為非熱應(yīng)力，而且每個(gè)應(yīng)力的加速效果均呈逆冪律形式，因此應(yīng)力Vi、Pi與相對(duì)應(yīng)的特征壽命ηi可以通過式(3)轉(zhuǎn)化為雙應(yīng)力加速的逆冪律模型[9]：

(4)

對(duì)式(4)兩邊取對(duì)數(shù)，得到加速方程的對(duì)數(shù)形式為

lnηi=lnA-αlnVi-βlnPi

(5)

式中A、α、β均為常數(shù).

1.2 加速系數(shù)

加速系數(shù)是指滾珠絲杠副在正常應(yīng)力組合(V0,P0)與加速應(yīng)力組合(Vi,Pi)下達(dá)到特征壽命需要時(shí)間T0和Ti的比值，記為τ.即

(6)

由式(4)(6)可得加速系數(shù)理論公式為

(7)

1.3 可靠性指標(biāo)

滾珠絲杠副是可修復(fù)、可更換系統(tǒng)，出現(xiàn)故障后，可對(duì)零部件進(jìn)行修復(fù)或更換，并且記錄故障時(shí)間和修復(fù)時(shí)間，然后繼續(xù)進(jìn)行試驗(yàn).

將平均故障間隔時(shí)間作為可靠性定量指標(biāo).其表達(dá)式為

(8)

式中Γ為伽馬函數(shù)，其值可以通過查閱伽馬函數(shù)標(biāo)準(zhǔn)表[10]得到.

在特定工作條件下，將滾珠絲杠副零部件出現(xiàn)故障的概率稱為壽命可靠性.由假定1可知滾珠絲杠副可靠度的數(shù)學(xué)表達(dá)式為

(9)

2 加速壽命試驗(yàn)

2.1 加速壽命試驗(yàn)方案設(shè)計(jì)

試驗(yàn)選用國(guó)產(chǎn)雙螺母預(yù)緊式滾珠絲杠副為研究對(duì)象，其結(jié)構(gòu)參數(shù)見表1[11].本文綜合考慮轉(zhuǎn)速以及軸向負(fù)載作為加速壽命試驗(yàn)的2個(gè)應(yīng)力，分別將其記為V和P.額定轉(zhuǎn)速V0和額定負(fù)載P0分別為500 r/min、1 000 N.如果每個(gè)應(yīng)力選定3個(gè)加速水平，則共有3×3=9個(gè)應(yīng)力組合，如果將這9個(gè)應(yīng)力組合全部進(jìn)行試驗(yàn)，則會(huì)大大延長(zhǎng)試驗(yàn)周期，降低試驗(yàn)效率.因此按照低、中、高3個(gè)水平對(duì)加速應(yīng)力進(jìn)行合理等級(jí)劃分，根據(jù)對(duì)角線原則分別選取(V1,P1)、(V2,P2)、(V3,P3)作為3個(gè)步進(jìn)加速雙應(yīng)力進(jìn)行定時(shí)截尾試驗(yàn).由于試驗(yàn)過程中滾珠絲杠副的失效機(jī)理不能改變，綜合考慮試驗(yàn)方案和加載能力，結(jié)合工程經(jīng)驗(yàn)，本文選取軸向負(fù)載大小分別為5%Ca、10%Ca、15%Ca(Ca為額定動(dòng)載荷)，具體載荷大小如表2所示.

表1 滾珠絲杠副的參數(shù)Table 1 Parameters of the ball screw

表2 加速試驗(yàn)應(yīng)力組合方案Table 2 Stress combination scheme of accelerated test

考慮到試驗(yàn)對(duì)象為小子樣單樣本滾珠絲杠副，為了滿足虛擬增廣法的適用要求，提高其可靠度壽命的評(píng)估精度，提出在不同步進(jìn)加速應(yīng)力組合下的試驗(yàn)截尾時(shí)間設(shè)定準(zhǔn)則：在每個(gè)步進(jìn)加速雙應(yīng)力組合下出現(xiàn)至少2～3個(gè)非關(guān)聯(lián)性故障，并將最后一個(gè)故障時(shí)間記為該步進(jìn)加速雙應(yīng)力組合下的試驗(yàn)截尾時(shí)間，即當(dāng)在第1個(gè)步進(jìn)加速雙應(yīng)力下發(fā)生至少2～3個(gè)非關(guān)聯(lián)性故障時(shí)，提升應(yīng)力，繼續(xù)在第2個(gè)步進(jìn)加速雙應(yīng)力水平下進(jìn)行試驗(yàn)，以此類推，直至試驗(yàn)結(jié)束.

整個(gè)試驗(yàn)過程中環(huán)境溫度保持在(20±5) ℃，濕度最大為90%，要定期進(jìn)行潤(rùn)滑，潤(rùn)滑條件包括潤(rùn)滑劑種類、潤(rùn)滑間隔時(shí)間、潤(rùn)滑量等必須按照具體說明書執(zhí)行.首先對(duì)絲杠進(jìn)行加載跑合，排除由于設(shè)計(jì)等問題引起的早期故障，確保在各個(gè)加速應(yīng)力載荷階段無(wú)異常，然后在滾珠絲杠副可靠性試驗(yàn)儀器上進(jìn)行模擬加載，對(duì)其進(jìn)行較為全面的試驗(yàn)，在試驗(yàn)過程中對(duì)摩擦力矩、噪聲、振動(dòng)、溫升等性能參數(shù)進(jìn)行實(shí)時(shí)監(jiān)測(cè)和記錄，如果滾珠絲杠副出現(xiàn)超過允許值的噪聲和振動(dòng)，則判定為故障，修復(fù)或更換受損部件后繼續(xù)進(jìn)行試驗(yàn)并記錄時(shí)間作為一次壽命點(diǎn).如果經(jīng)過修復(fù)后仍無(wú)法正常工作即發(fā)生不可修復(fù)故障，則視為失效，停止試驗(yàn).

對(duì)可靠性進(jìn)行評(píng)估時(shí)，關(guān)鍵的步驟之一是弄清故障的判定和分析方法，這直接影響到評(píng)估結(jié)果的準(zhǔn)確性.對(duì)于滾珠絲杠副來(lái)說，故障是指產(chǎn)品在規(guī)定的工作運(yùn)行條件下，一個(gè)或某幾個(gè)功能喪失，或是性能參數(shù)超過允許界限，也包括滾珠絲杠副螺母、滾珠、回珠器等元部件的破損和斷裂.具體故障判據(jù)如下[12]：

1) 滾珠絲杠副的任一零部件發(fā)生松動(dòng)、損傷、破裂或局部斷裂；

2) 檢測(cè)系統(tǒng)采集到的摩擦力矩、振動(dòng)、噪聲等數(shù)值在某一段監(jiān)測(cè)時(shí)間間隔內(nèi)偏離正常范圍；

3) 絲杠傳動(dòng)不平穩(wěn)，無(wú)法正常運(yùn)行.

根據(jù)以上判據(jù)對(duì)可靠性試驗(yàn)過程中滾珠絲杠副出現(xiàn)的故障進(jìn)行統(tǒng)計(jì)分析，并將其按照文獻(xiàn)[13]中表5.2滾珠絲杠副的關(guān)聯(lián)故障和非關(guān)聯(lián)故障完成分類.

2.2 加速壽命試驗(yàn)統(tǒng)計(jì)分析

本文在定時(shí)截尾雙應(yīng)力加速壽命試驗(yàn)的基礎(chǔ)上，利用虛擬增廣法對(duì)試驗(yàn)數(shù)據(jù)進(jìn)行等價(jià)擴(kuò)充以滿足參數(shù)估計(jì)方法的要求，提高在單樣本情況下的壽命評(píng)估精度.虛擬增廣法的立論依據(jù)：增廣后子樣均值應(yīng)與原來(lái)的子樣均值相等；增廣后的子樣標(biāo)準(zhǔn)差應(yīng)與類似件的子樣標(biāo)準(zhǔn)差相等[14].為使虛擬增廣得到的故障時(shí)間數(shù)據(jù)樣本更為合理，使用近似公式

(10)

對(duì)試驗(yàn)得到的故障數(shù)據(jù)進(jìn)行增廣至m.式中：T0為樣本均值；c為待定系數(shù)；σ為樣本標(biāo)準(zhǔn)差；T為增廣后的樣本值.

利用近似公式(10)將試驗(yàn)樣本數(shù)據(jù)增廣至m=10個(gè)，那么增廣后的樣本數(shù)據(jù)T1～T10可以表示為

(11)

根據(jù)虛擬增廣法的2個(gè)立論依據(jù)可以列出方程：

(12)

(13)

根據(jù)式(12)(13)可求解得c=0.135 0，綜合試驗(yàn)數(shù)據(jù)及式(12)可以求出增廣后的樣本數(shù)據(jù)分別為

(T0-1.831 0σ)，(T0-1.089 0σ)，
(T0-0.559 0σ)，(T0-0.241 0σ)，
(T0-0.135 0σ)，(T0+0.135 0σ)，
(T0+0.241 0σ)，(T0+0.559 0σ)，
(T0+1.089 0σ)，(T0+1.831 0σ)

然后，將虛擬增廣后的數(shù)據(jù)利用最小二乘法進(jìn)行統(tǒng)計(jì)分析，具體做法如下.

將累積故障分布函數(shù)做如下線性變換：

(14)

其中，令

(15)

則式(14)可轉(zhuǎn)換為Weibull分布線性模型

y=bx+a

(16)

在每個(gè)步進(jìn)加速雙應(yīng)力情況下，通過試驗(yàn)得到的故障時(shí)間對(duì)應(yīng)的累積故障概率采用近似中位秩公式表示，其經(jīng)驗(yàn)分布函數(shù)為

(17)

式中n為樣本數(shù)據(jù)總數(shù).那么對(duì)于一組增廣后的試驗(yàn)數(shù)據(jù)

(ti,Fn(ti)),i=1,2,…,n

(18)

根據(jù)式(15)可轉(zhuǎn)換為

(19)

將線性模型y=bx+a通過最小二乘法進(jìn)行擬合，可以求出待估參數(shù)b、a的表達(dá)式分別為

(20)

(21)

綜上，求出Weibull分布中兩參數(shù)的點(diǎn)估計(jì)為

(22)

求出模型參數(shù)后進(jìn)行模型驗(yàn)證，判定故障間隔時(shí)間是否符合Weibull分布模型，常用的驗(yàn)證方法有擬合優(yōu)度檢驗(yàn)和假設(shè)檢驗(yàn).

1) 擬合優(yōu)度檢驗(yàn)

可用相關(guān)系數(shù)r表明數(shù)據(jù)xi、yi的擬合情況，其表達(dá)式為

(23)

相關(guān)系數(shù)越趨近于1，說明壽命分布服從Weibull模型越好.顯著性水平為α?xí)r通過查表確定相關(guān)系數(shù)起碼值rα.當(dāng)r>rα?xí)r，認(rèn)為回歸效果顯著；反之，認(rèn)為回歸不顯著.可以據(jù)此來(lái)判定故障數(shù)據(jù)服從Weibull分布的好壞情況.

2) 假設(shè)檢驗(yàn)

基本思想：先提出兩個(gè)對(duì)立的假設(shè)，根據(jù)試驗(yàn)數(shù)據(jù)計(jì)算出檢驗(yàn)統(tǒng)計(jì)量，然后根據(jù)檢驗(yàn)統(tǒng)計(jì)量決定其中某一個(gè)假設(shè)成立，完成檢驗(yàn).常用的模型檢驗(yàn)方法有卡方檢驗(yàn)和柯爾莫戈洛夫(Kolmogorov-Smirnov,KS)檢驗(yàn)[15].其中，卡方檢驗(yàn)對(duì)小樣本的檢驗(yàn)精度不如柯爾莫戈洛夫檢驗(yàn)，因此選用KS檢驗(yàn)法對(duì)試驗(yàn)?zāi)Ｐ瓦M(jìn)行假設(shè)檢驗(yàn).

KS檢驗(yàn)法將經(jīng)驗(yàn)分布函數(shù)Fn(t)和統(tǒng)計(jì)分析得到的模型累積分布函數(shù)F0(t)進(jìn)行對(duì)比，原假設(shè)H0與備擇假設(shè)H1如下.

H0：故障間隔時(shí)間服從Weibull分布

H1：故障間隔時(shí)間不服從Weibull分布

其檢驗(yàn)統(tǒng)計(jì)量

Dn=sup|F0(ti)-Fn(ti)|=max{di}

(24)

通過查KS檢驗(yàn)臨界值表可以得到在顯著性水平為α?xí)rDnα的值.如果有Dn≤Dnα，則H0成立，故障間隔時(shí)間數(shù)據(jù)服從Weibull分布；反之，H1成立，Weibull壽命分布模型不適用.

3 加速壽命試驗(yàn)仿真算例

對(duì)樣本進(jìn)行加速壽命仿真試驗(yàn)，其結(jié)構(gòu)參數(shù)如表1所示，所施加的載荷應(yīng)力如表2所示.然后借助數(shù)學(xué)方法和計(jì)算機(jī)得到需要的模擬故障數(shù)據(jù)，其具體步驟如下.

對(duì)Weibull累積分布函數(shù)式(1)進(jìn)行逆變換得到其反函數(shù)為

(25)

其中隨機(jī)變量y服從均勻分布U(0,1)，y與1-y為同分布，所以可轉(zhuǎn)換為

(26)

可以通過計(jì)算機(jī)MATLAB軟件產(chǎn)生基于均勻分布U(0,1)的隨機(jī)數(shù)y，然后計(jì)算出F-1(y)的值.

結(jié)合上述方法以及截尾時(shí)間設(shè)定準(zhǔn)則，根據(jù)文獻(xiàn)[16]中的滾珠絲杠副的形狀參數(shù)m=1.65，特征壽命η分別取800、480、350 h得到在3個(gè)步進(jìn)加速雙應(yīng)力水平下的模擬故障間隔時(shí)間數(shù)據(jù)并按從小到大排列，如表3所示.

表3 加速應(yīng)力下的模擬數(shù)據(jù)Table 3 Simulation data under accelerated stress

對(duì)加速應(yīng)力下的故障間隔時(shí)間數(shù)據(jù)進(jìn)行分析，根據(jù)式(11)(12)(13)得到每一個(gè)步進(jìn)加速雙應(yīng)力下的虛擬增廣后的數(shù)據(jù)如表4所示.

表4 經(jīng)過虛擬增廣后的數(shù)據(jù)Table 4 Data after virtual augmentation

依據(jù)最小二乘法原理，通過MATLAB軟件分別對(duì)3個(gè)步進(jìn)加速雙應(yīng)力下增廣后的數(shù)據(jù)進(jìn)行直線擬合，擬合結(jié)果如圖1所示.

圖1 數(shù)據(jù)擬合曲線Fig.1 Data fitting curve

從圖1可以看出，在3個(gè)步進(jìn)加速雙應(yīng)力水平下擬合得到的直線基本保持平行，斜率近似相等，即形狀參數(shù)保持不變，從而驗(yàn)證了假定2的正確性，即在不同的加速應(yīng)力組合下，樣本的失效機(jī)理沒有發(fā)生改變.

由表5分析可知，3條擬合直線的相關(guān)系數(shù)都很接近于1，說明線性相關(guān)度很高，擬合效果較好.

表5 擬合直線數(shù)據(jù)Table 5 Data of the fitting line

結(jié)合式(22)可求得在3個(gè)步進(jìn)加速雙應(yīng)力下的形狀參數(shù)分別為

m1=1.611 0,m2=1.611 4,m3=1.609 8

為滿足假定2，形狀參數(shù)取各加速應(yīng)力組合下所求得的形狀參數(shù)的平均值，即

將m的值代入式(22)中的第2個(gè)式子，進(jìn)而求出在3個(gè)步進(jìn)加速雙應(yīng)力下的特征壽命分別為

η1=828.92,η2=498.06,η3=361.06

再將以上特征壽命數(shù)據(jù)代入式(5)，聯(lián)立求解可得

A=6.145 3×105,α=0.469 2,β=0.460 4

因此加速方程為

lnηi=13.328 6-0.469 2lnVi-0.460 4lnPi

(27)

將(V0,P0)代入式(27)求出正常應(yīng)力條件下的壽命為η0=1 383.0 h.

根據(jù)式(9)可知滾珠絲杠副在額定載荷下的可靠度函數(shù)為

由式(8)可知滾珠絲杠副的MTBF為

4 應(yīng)用討論

滾珠絲杠副在實(shí)際使用中，一般可分為高速輕載、低速重載等情況，以下結(jié)合加速方程對(duì)不同工況下的壽命可靠度進(jìn)行分析討論.

由加速方程式(27)可以求出在高速輕載(V3,P0)和低速重載(V0,P3)工況下的特征壽命分別為

η30=721.68,η03=691.93

根據(jù)式(9)可知在這兩種工況下的可靠度函數(shù)分別為

不同工況下的可靠度函數(shù)圖形如圖2所示.

圖2 滾珠絲杠副不同工況下的可靠度函數(shù)Fig.2 Reliability function diagram of ball screw under different working conditions

從圖2中可以看出，在相同可靠度的情況下，滾珠絲杠副在高速輕載下的壽命比在低速重載下的壽命略大，均遠(yuǎn)小于在額定載荷下的可靠度壽命，與實(shí)際相符.

5 結(jié)論

1) 為了在滾珠絲杠副可靠性壽命試驗(yàn)中綜合考慮速度和負(fù)載兩個(gè)應(yīng)力水平對(duì)可靠性壽命的影響，本文提出一種雙應(yīng)力三水平的加速壽命試驗(yàn)設(shè)計(jì)和可靠性分析方法，通過建立雙應(yīng)力加速壽命試驗(yàn)?zāi)Ｐ筒⒉捎米钚《朔▽?duì)故障數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析，建立滾珠絲杠副的加速壽命試驗(yàn)的雙應(yīng)力加速方程.為進(jìn)一步開展?jié)L珠絲杠副實(shí)際不同工況條件下壽命評(píng)估提供方法和依據(jù).

2) 本文提出利用虛擬增廣法對(duì)小子樣單樣本條件下雙應(yīng)力壽命試驗(yàn)數(shù)據(jù)進(jìn)行等價(jià)擴(kuò)充后得到不同應(yīng)力水平下的威布爾分布參數(shù)，并建立加速壽命試驗(yàn)的加速方程的方法.該方法可以滿足小子樣單樣本試驗(yàn)條件下滾珠絲桿副可靠性壽命評(píng)估的精確性要求.并且，本文給出了仿真算例驗(yàn)證了所提出方法的有效性.