基于威布爾分布的加速試驗(yàn)剖面設(shè)計(jì)方法

王海東,蔣 剛,陳志偉,賈 斐

(1.上海航天精密機(jī)械研究所,上海 201600;2.西安電子科技大學(xué) 機(jī)電工程學(xué)院,陜西 西安 710071)

0 引言

產(chǎn)品壽命是評(píng)價(jià)航空航天裝備使用性能的重要指標(biāo)。為準(zhǔn)確、快速地評(píng)價(jià)航空航天產(chǎn)品的使用壽命,目前通常采用加速壽命試驗(yàn)方法進(jìn)行產(chǎn)品定壽。加速壽命試驗(yàn)(ALT)是在合理的工程及統(tǒng)計(jì)假設(shè)的基礎(chǔ)上,采用提高工作應(yīng)力或環(huán)境應(yīng)力的方法來使產(chǎn)品快速地暴露故障,再利用加速應(yīng)力下的試驗(yàn)信息(失效數(shù)據(jù)、失效原因等)來外推產(chǎn)品在正常應(yīng)力水平下的可靠性指標(biāo)的一種壽命試驗(yàn)方法[1]。加速壽命試驗(yàn)采用加速應(yīng)力進(jìn)行試件的壽命試驗(yàn),縮短了試驗(yàn)時(shí)間,提高了試驗(yàn)效率,降低了試驗(yàn)成本,其研究使高可靠長(zhǎng)壽命產(chǎn)品的可靠性評(píng)定成為可能[2]。加速壽命試驗(yàn)受到國內(nèi)外可靠性工作人員和工程研究人員的高度重視,且已被廣泛應(yīng)用于機(jī)電產(chǎn)品、半導(dǎo)體材料、武器裝備及航空航天等實(shí)際問題中[3-7]。

加速壽命試驗(yàn)的類型很多,按照應(yīng)力施加方式的不同,一般分為恒定應(yīng)力加速壽命試驗(yàn)、步進(jìn)應(yīng)力加速壽命試驗(yàn)、序進(jìn)應(yīng)力加速壽命試驗(yàn)3種基本類型[8]。在NELSON[9]提出的累積損傷模型基礎(chǔ)上,我國學(xué)者研究了不同分布類型下步進(jìn)加速試驗(yàn)數(shù)據(jù)處理方法。序進(jìn)加速試驗(yàn)數(shù)據(jù)處理方法難度較大,當(dāng)前取得的研究成果相對(duì)較少,還未能成熟運(yùn)用。恒定加速試驗(yàn)和步進(jìn)加速試驗(yàn)應(yīng)用較廣泛,但步進(jìn)加速試驗(yàn)對(duì)試驗(yàn)時(shí)間和樣品數(shù)量需求更低,試驗(yàn)效率較高,已成為產(chǎn)品定壽試驗(yàn)的應(yīng)用趨勢(shì)。而二參數(shù)威布爾分布作為一種連續(xù)分布,可以很好地?cái)M合各種類型的試驗(yàn)數(shù)據(jù),因此被廣泛應(yīng)用于各個(gè)領(lǐng)域,特別是可靠性領(lǐng)域中的壽命數(shù)據(jù)分析,因此基于二參數(shù)威布爾分布的研究一直是可靠性領(lǐng)域研究的熱點(diǎn)之一。鎂合金是實(shí)際應(yīng)用中最輕的金屬結(jié)構(gòu)材料,在航天航空產(chǎn)品領(lǐng)域具有廣闊的應(yīng)用前景[10-11],所制造出的產(chǎn)品也具有更高的可靠性,同樣符合二參數(shù)威布爾分布。

1 試驗(yàn)剖面優(yōu)化設(shè)計(jì)

1.1 設(shè)計(jì)變量

1.1.1 確定試驗(yàn)剖面

為了便于分析,將對(duì)簡(jiǎn)單三步進(jìn)加速壽命試驗(yàn)剖面進(jìn)行研究。三步進(jìn)加速壽命試驗(yàn)剖面如圖1所示。圖中：粗黑線表示試驗(yàn)載荷施加方式;Tm為不改變失效機(jī)理的最大溫度應(yīng)力,可根據(jù)少量摸底試驗(yàn)或仿真確定;tc為試驗(yàn)截尾時(shí)間,可根據(jù)工程經(jīng)驗(yàn)或者少量摸底試驗(yàn)確定;T1和T2為加速溫度應(yīng)力;t1,t2為應(yīng)力轉(zhuǎn)換時(shí)間。除此之外,產(chǎn)品正常工作溫度應(yīng)力為T0。

圖1 三步進(jìn)加速壽命試驗(yàn)剖面Fig.1 Three-step accelerated life test profile

1.1.2 模型歸一化

由于參數(shù)較多,不便于分析。為了簡(jiǎn)便計(jì)算,且使各參數(shù)的取值與優(yōu)化模型結(jié)果之間的關(guān)系模型更具一般性,需對(duì)模型進(jìn)行標(biāo)準(zhǔn)化處理。

1)應(yīng)力水平標(biāo)準(zhǔn)化

已知正常工作應(yīng)力水平為T0,最高加速應(yīng)力水平為Tm,令則加速應(yīng)力水平T1和T2轉(zhuǎn)化為標(biāo)準(zhǔn)應(yīng)力水平u1和u2。

2)應(yīng)力轉(zhuǎn)換時(shí)間標(biāo)準(zhǔn)化

已知試驗(yàn)截尾時(shí)間tc,令

則應(yīng)力轉(zhuǎn)換時(shí)間t1和t2轉(zhuǎn)化為標(biāo)準(zhǔn)應(yīng)力轉(zhuǎn)換時(shí)間v1和v2。

1.2 約束條件

針對(duì)三步進(jìn)加速試驗(yàn)的設(shè)計(jì)變量有4個(gè),分別為加速應(yīng)力水平u1和u2,以及應(yīng)力轉(zhuǎn)換時(shí)間v1和v2。對(duì)應(yīng)的約束條件有

1.3 目標(biāo)函數(shù)

1.3.1 確定分布函數(shù)

產(chǎn)品在不同應(yīng)力水平下的壽命服從同一分布函數(shù),且試驗(yàn)中所有隨機(jī)試樣都具有統(tǒng)計(jì)獨(dú)立性。在常見產(chǎn)品壽命模型中,主要有指數(shù)分布、正態(tài)分布和威布爾分布,其中威布爾分布模型的應(yīng)用較為廣泛,因?yàn)橥紶柗植寄Ｐ途哂卸嘧兊暮瘮?shù)形態(tài),能夠很好地?cái)M合各類數(shù)據(jù)。威布爾分布作為一種連續(xù)的概率分布模型,有近半個(gè)世紀(jì)的應(yīng)用和研究歷史,在機(jī)械可靠性工程中得到了一定的應(yīng)用,二參數(shù)和三參數(shù)分布形式是威布爾分布的基本形式。常用的二參數(shù)威布爾分布的累積失效分布函數(shù)和概率密度函數(shù)分別為

式中：x為失效時(shí)間;η為尺度參數(shù),η＞0;m為形狀參數(shù)(根據(jù)工程經(jīng)驗(yàn),0.5≤m≤5)。

為了使產(chǎn)品在不同應(yīng)力水平下的失效機(jī)理保持一致,需假定產(chǎn)品在不同應(yīng)力水平下分布函數(shù)的形狀參數(shù)保持不變,即二參數(shù)威布爾分布的形狀參數(shù)m保持不變。

1.3.2 確定加速模型

加速壽命試驗(yàn)的實(shí)質(zhì)：當(dāng)產(chǎn)品壽命特征處于高應(yīng)力水平下時(shí)去外推產(chǎn)品處于正常應(yīng)力水平下的壽命特征,其推導(dǎo)關(guān)鍵在于建立壽命特征與應(yīng)力水平之間的數(shù)學(xué)關(guān)系,這種關(guān)系被稱為加速模型或加速方程。常用的單應(yīng)力加速模型有Arrhenius模型、逆冪律模型、Eyring模型和廣義Eyring模型。本文考慮的加速應(yīng)力只有溫度,故采用Arrhenius加速模型,即

式中：T為環(huán)境溫度,℃;Ea為激活能,e V;k為玻爾茲曼常數(shù),值為8.617×10-5e V/K;A為未知常數(shù)。

將式(6)兩邊分別取對(duì)數(shù),有

可將式(7)轉(zhuǎn)換成對(duì)數(shù)線性關(guān)系式,即

式中：γ1,γ2為未知常數(shù);φ(T)為跟溫度應(yīng)力有關(guān)的函數(shù)。

1.3.3 構(gòu)建累積失效模型

NELSON累積失效理論(CE模型),假設(shè)產(chǎn)品的剩余壽命僅依賴于當(dāng)時(shí)已累積失效的部分和當(dāng)時(shí)的應(yīng)力水平,與累積方式無關(guān)。這一假定是NELSON根據(jù)物理原理提出的。如果產(chǎn)品的壽命分布為F(t),則該假定的數(shù)學(xué)含義是：某產(chǎn)品在應(yīng)力水平Si下的工作時(shí)間為t i,則t i內(nèi)產(chǎn)品的累積失效概率為F i(t i),相當(dāng)于在應(yīng)力水平Sj下工作時(shí)間為tj時(shí)產(chǎn)品的累積失效概率F j(t j),即F i(t i)=F j(t j)。根據(jù)這一假定,可對(duì)不同應(yīng)力水平下的壽命數(shù)據(jù)進(jìn)行折算。

根據(jù)NELSON 累積失效理論,由F1(t1)=F2(τ1),得

同理,由F2(t2-t1+τ1)=F3(τ2)得

式中：τ1為在溫度應(yīng)力水平T1下試驗(yàn)t1時(shí)間換算到溫度應(yīng)力水平T2下的累積等效試驗(yàn)時(shí)間;τ2為在溫度應(yīng)力水平T2下試驗(yàn)(t2-t1+τ1)時(shí)間換算到溫度應(yīng)力水平Tm下的累積等效試驗(yàn)時(shí)間。

因此,試件失效時(shí)間X的累積分布函數(shù)為

對(duì)應(yīng)的概率密度度函數(shù)為

式中：F1(x),F2(x),F3(x)分別為在溫度應(yīng)力水平T1,T2,Tm下的分布函數(shù);f1(x),f2(x),f3(x)分別表示在溫度應(yīng)力水平T1,T2,Tm下的概率密度函數(shù);η1,η2,η3分別為在溫度應(yīng)力水平T1,T2,Tm下的尺度參數(shù)。同時(shí)滿足

1.3.4 模型參數(shù)估計(jì)

由于加速試驗(yàn)普遍存在試驗(yàn)截尾、失效模式不確定等問題,因此突發(fā)型失效模式數(shù)據(jù)通常是非完整數(shù)據(jù)。極大似然估計(jì)(MLE)不僅適用于完整數(shù)據(jù),同時(shí)適用于非完整數(shù)據(jù),因此,選擇MLE 對(duì)加速試驗(yàn)突發(fā)型失效模式數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析。

假設(shè)第i個(gè)樣本在時(shí)間xi失效,定義3個(gè)指示函數(shù)

則第i個(gè)樣本的對(duì)數(shù)似然函數(shù)L i為

當(dāng)投入試驗(yàn)的樣本量為n時(shí),則所有樣本的的對(duì)數(shù)似然函數(shù)L為

解L的一階偏導(dǎo)數(shù)方程組,由

1.3.5 分位數(shù)壽命估計(jì)

根據(jù)加速模型公式,可以計(jì)算得到在正常溫度應(yīng)力水平T0下的尺度參數(shù)η0的估計(jì)值為

從而得到在正常溫度應(yīng)力水平T0下的累積分布函數(shù)為

故產(chǎn)品的P階分位數(shù)壽命x P的MLE估計(jì)值為

1.3.6 優(yōu)化準(zhǔn)則與漸進(jìn)方差

1)優(yōu)化準(zhǔn)則

對(duì)于具有可靠性高且壽命長(zhǎng)的產(chǎn)品,做加速壽命試驗(yàn)時(shí),通常選取產(chǎn)品在正常應(yīng)力下P階分位壽命MLE估計(jì)值的漸進(jìn)方差作為其可靠性壽命特征量。為提高該特征量估計(jì)精度,并減少試樣數(shù)量、縮短試驗(yàn)時(shí)間,需要優(yōu)化加速壽命試驗(yàn)方案。在正常應(yīng)力下,以產(chǎn)品P階分位壽命MLE 估計(jì)值的漸進(jìn)方差最小為準(zhǔn)則,對(duì)加速壽命試驗(yàn)方案進(jìn)行優(yōu)化。

2)壽命的漸進(jìn)方差

一般協(xié)方差矩陣難以求得,根據(jù)MLE 估計(jì)理論,協(xié)方差矩陣Σ與信息矩陣F為互逆,即Σ=F-1,故通常采取對(duì)信息矩陣求逆來求得協(xié)方差矩陣。而信息矩陣可由n個(gè)樣本的對(duì)數(shù)似然函數(shù)對(duì)各模型參數(shù)的負(fù)二階偏導(dǎo)數(shù)矩陣的數(shù)學(xué)期望求得,有

通常,極大似然估計(jì)量具有漸進(jìn)無偏性和漸進(jìn)正態(tài)性,故可得到在正常溫度應(yīng)力水平T0下,產(chǎn)品P階分位數(shù)的對(duì)數(shù)壽命的漸進(jìn)方差為

地鐵供電系統(tǒng)中剛性接觸網(wǎng)常見故障和防范措施解析……………………………………………………… 鐘人正（7-72）

式中：H為一階偏導(dǎo)矩陣,即

由式(1),(2)得

再將式(27)代入式(25)中,可得到只與4個(gè)設(shè)計(jì)變量u1,u2,v1,v2有關(guān)的漸進(jìn)方差表達(dá)式。

經(jīng)過推導(dǎo)計(jì)算,信息矩陣F里面的元素過于復(fù)雜,無法求出具體的定積分表達(dá)式,故采用數(shù)值積分的方法進(jìn)行近似求解。考慮到使用數(shù)值積分需要已知積分區(qū)間,且整個(gè)推導(dǎo)過程中涉及到諸多變量替換,同時(shí)存在約束條件,故考慮使用非線性規(guī)劃算法進(jìn)行優(yōu)化求解。而經(jīng)典非線性規(guī)劃算法大多采用梯度下降的方法求解,需知道顯式的函數(shù)表達(dá)式,顯然不適用,故綜合考慮下決定使用隨機(jī)Monte-Carlo模擬的方法進(jìn)行優(yōu)化求解,求使?jié)u進(jìn)方差取得最小值所對(duì)應(yīng)的優(yōu)化變量

2 仿真算例

結(jié)合部分工程經(jīng)驗(yàn),假設(shè)參數(shù)真值為m=2,γ1=10,γ2=2 300;溫度應(yīng)力水平分別為T0=20 ℃,T1=50 ℃,T2=80 ℃,Tm=110 ℃;試 驗(yàn) 時(shí)間為t1=40 min,t2=80 min和tc=120 min;試驗(yàn)樣本量為n=40。

現(xiàn)以產(chǎn)品在正常溫度應(yīng)力水平T0下的對(duì)數(shù)中位數(shù)壽命估計(jì)值的漸進(jìn)方差為優(yōu)化準(zhǔn)則,以加速應(yīng)力水平u1,u2和應(yīng)力轉(zhuǎn)換時(shí)間v1,v2為設(shè)計(jì)變量,進(jìn)行三步進(jìn)加速試驗(yàn)方案的優(yōu)化設(shè)計(jì)。

采用Monte-Carlo模擬的方法進(jìn)行1 000次隨機(jī)尋優(yōu)求解的結(jié)果為=0.33,=0.82,此時(shí)所對(duì)應(yīng)的漸進(jìn)方差最小值為0.11。若采用傳統(tǒng)均勻設(shè)計(jì)的加速試驗(yàn)方案：=1/3,,此時(shí)所求得的漸進(jìn)方差為0.34。

顯然,經(jīng)過優(yōu)化后的壽命估計(jì)值的漸進(jìn)方差比傳統(tǒng)均勻設(shè)計(jì)的漸進(jìn)方差要小,這表明經(jīng)過優(yōu)化的加速試驗(yàn)方案比傳統(tǒng)均勻設(shè)計(jì)的試驗(yàn)方案具有更高的估計(jì)精度。

3 結(jié)束語

為實(shí)現(xiàn)對(duì)高可靠性、長(zhǎng)壽命產(chǎn)品在正常環(huán)境應(yīng)力水平下壽命與可靠性的快速評(píng)估,本文在單應(yīng)力情形下,研究了二參數(shù)威布爾分布的三步進(jìn)加速試驗(yàn)的優(yōu)化設(shè)計(jì)理論與試驗(yàn)數(shù)據(jù)的統(tǒng)計(jì)分析方法。通過建立Arrhenius與二參數(shù)威布爾分布試驗(yàn)的可靠性統(tǒng)計(jì)模型,以產(chǎn)品在正常溫度應(yīng)力水平下中位壽命估計(jì)的漸近方差最小為優(yōu)化準(zhǔn)則,給出了二參數(shù)威布爾分布的三步進(jìn)加速試驗(yàn)應(yīng)力水平與應(yīng)力轉(zhuǎn)換時(shí)間的優(yōu)化設(shè)計(jì)方法,并通過與傳統(tǒng)均勻設(shè)計(jì)加速試驗(yàn)方案的對(duì)比,表明經(jīng)優(yōu)化后的加速試驗(yàn)方案具有更高的估計(jì)精度。采用隨機(jī)Monte-Carlo模擬的方法代替?zhèn)鹘y(tǒng)非線性規(guī)劃算法進(jìn)行加速試驗(yàn)方案設(shè)計(jì)變量的優(yōu)化求解,以解決函數(shù)表達(dá)式過于冗長(zhǎng)、計(jì)算復(fù)雜等問題。今后將考慮基于二參數(shù)威布爾分布進(jìn)行多應(yīng)力(如溫度、濕度、振動(dòng)等)綜合加載方式下的步進(jìn)加速試驗(yàn)統(tǒng)計(jì)分析與試驗(yàn)設(shè)計(jì)優(yōu)化方法研究,并考慮基于二參數(shù)威布爾分布進(jìn)行多步進(jìn)(即加速應(yīng)力水平數(shù)大于3)情況下的步進(jìn)加速試驗(yàn)統(tǒng)計(jì)分析與試驗(yàn)設(shè)計(jì)優(yōu)化方法研究。