叢 明，王冠雄，謝書文，潘大慶

（1.大連理工大學(xué)機(jī)械學(xué)院，遼守大連 116024；2.長(zhǎng)安汽車股份有限公司工藝部，重慶 404100）

加工中心早期故障的兩重威布爾模型研究*

叢 明1，王冠雄1，謝書文2，潘大慶2

（1.大連理工大學(xué)機(jī)械學(xué)院，遼守大連 116024；2.長(zhǎng)安汽車股份有限公司工藝部，重慶 404100）

為了對(duì)加工中心早期故障進(jìn)行有效的預(yù)測(cè)和控制，采集了26臺(tái)加工中心從試運(yùn)行開始的102組故障數(shù)據(jù)，并以此為基礎(chǔ)，列出并比較了3種較常用的兩重威布爾分布函數(shù)模型，從中選取分段威布爾函數(shù)作為擬合模型。為提高分段威布爾模型特征參數(shù)估計(jì)精度，提出迭代尋址求最優(yōu)解法，以擬合精度為目標(biāo)優(yōu)化量，有效避免了時(shí)間分界點(diǎn)選取錯(cuò)誤，提高了擬合精度，計(jì)算得出早期故障期大約持續(xù)4個(gè)月，與廠家反饋一致，并對(duì)分析結(jié)果進(jìn)行了擬合優(yōu)度檢驗(yàn)和Laplace檢驗(yàn)，檢驗(yàn)結(jié)果符合要求。

早期故障；浴盆曲線；加工中心

0 引言

加工中心生命周期過程一般呈現(xiàn)出“浴盆曲線”，包括早期故障期、偶然故障期和耗損故障期，如圖1所示［1］。由于加工中心早期故障期多為3個(gè)月到半年，時(shí)間較短，給數(shù)據(jù)獲取帶來(lái)了一定的難度，所以加工中心的可靠性研究多側(cè)重于偶然故障期，而早期故障期的相關(guān)研究較少。但是早期故障分布模型的建立，特別是分界時(shí)間點(diǎn)的確定，有助于加工中心早期可靠性的提高，對(duì)于減小加工中心早期故障危害有指導(dǎo)作用和實(shí)際意義。

賈亞洲教授［2］首先提出對(duì)數(shù)控機(jī)床早期故障進(jìn)行專項(xiàng)研究，他采用FMEA方法對(duì)9臺(tái)加工中心的早期故障進(jìn)行采集、分析，得出主要的故障模式與故障原因，文獻(xiàn)3深入分析了早期故障的產(chǎn)生原因，建立了早期故障消除技術(shù)體系，為提高加工中心早期故障期可靠性有一定的借鑒作用。文獻(xiàn)［4-7］，利用兩重威布爾分布對(duì)故障數(shù)據(jù)進(jìn)行擬合，對(duì)比一重威布爾分布，擬合效果更準(zhǔn)確，但二重威布爾分布模型類型較多，上述文獻(xiàn)并沒有對(duì)各種二重威布爾模型的應(yīng)用范圍進(jìn)行討論，而且都采用圖解法求解特征參數(shù)，雖然求解效率高，但不能準(zhǔn)確求出特征參數(shù)，并且易造成假設(shè)的分界時(shí)間點(diǎn)與實(shí)際分界時(shí)間點(diǎn)不統(tǒng)一等問題。

本文收集了18個(gè)月中加工中心早期故障數(shù)據(jù)，針對(duì)此數(shù)據(jù)，擬合出合適的威布爾分布模型，采用迭代尋址求最優(yōu)解方法，計(jì)算出早期故障期與偶然故障期的時(shí)間分界點(diǎn)，便于企業(yè)針對(duì)不同故障期制定合適的維修策略，對(duì)于提高生產(chǎn)效率、降低維護(hù)成本有很大的意義。

圖1 浴盆曲線

1 故障數(shù)據(jù)采集與分析

1.1 故障數(shù)據(jù)采集

本文針對(duì)某發(fā)動(dòng)機(jī)生產(chǎn)企業(yè)的26臺(tái)加工中心進(jìn)行故障數(shù)據(jù)采集，數(shù)據(jù)采集從加工中心試生產(chǎn)開始，持續(xù)了18個(gè)月，采集了102組故障數(shù)據(jù)，具體數(shù)據(jù)如表一所示，。因?yàn)榧庸ぶ行呐既还收掀跁r(shí)間較長(zhǎng)，進(jìn)入耗損故障期較緩慢，在短時(shí)間內(nèi)無(wú)法采集到耗損故障期的故障數(shù)據(jù)，所以本文只針對(duì)浴盆曲線的左側(cè)，即早期故障期和偶然故障期的前期做詳細(xì)分析。

表1 故障時(shí)間表

將102組故障間隙時(shí)間觀測(cè)值t∈［8，4319］分為10組，以每組故障間隙時(shí)間中值為橫坐標(biāo)，故障頻數(shù)為縱坐標(biāo)做出故障分布時(shí)間圖，如圖2所示。從圖2可以看出，隨著運(yùn)行時(shí)間的增加，發(fā)生故障的次數(shù)呈現(xiàn)出先減小后不變的趨勢(shì)，這與威布爾分布的趨勢(shì)相符合，為了更加清晰地判斷出故障數(shù)據(jù)的變化趨勢(shì)，還需對(duì)故障數(shù)據(jù)進(jìn)一步分析。

圖2 故障分布時(shí)間圖

1.2 WPP圖

判斷故障數(shù)據(jù)符合哪種威布爾分布，最常用的辦法是把故障數(shù)據(jù)經(jīng)過線性處理后，畫在威布爾頻率紙（WPP）上，通過分析處理后的圖形，與已知分布模型變化趨勢(shì)相比對(duì)，選擇適合的分布模型，進(jìn)而求出曲線的分布參數(shù)［3］。

把故障數(shù)據(jù)進(jìn)行線性變換，首先，通過中位秩估計(jì)故障分布函數(shù)F1（ti）。

式中，i表示故障數(shù)據(jù)的序數(shù)，n表示樣本量。對(duì)t和做變化F2（ti），令

將表1中的觀測(cè)值轉(zhuǎn)化為x、y，做出WPP圖，如圖3所示。

圖3 故障數(shù)據(jù)的WPP圖

已知兩參數(shù)威布爾分布的WPP圖為一條直線，三參數(shù)的WPP圖為一條曲線，都不和圖3吻合，所以用簡(jiǎn)單的一重二參數(shù)或三參數(shù)威布爾分布模型擬合曲線的效果不理想。從圖3得出變化趨勢(shì)為兩條直線，這與兩重威布爾分布模型的WPP圖類型，所以采用兩重威布爾分布對(duì)未知曲線進(jìn)行線性估計(jì)。

2 兩重威布爾分布類型與選擇

兩重威布爾分布類型很多，包括兩重分段型威布爾分布，兩重混合型威布爾分布，兩重競(jìng)爭(zhēng)型威布爾分布等［4-9］。每種威布爾分布的適用情況及WPP圖都不盡相同，下文分別對(duì)常用的三種兩重威布爾分布做詳細(xì)敘述與歸納。

2.1 兩重分段威布爾分布

分段威布爾分布多用于故障數(shù)據(jù)模型較復(fù)雜，并沒有明顯規(guī)律，不能用簡(jiǎn)單的分布擬合的情況。

失效分布函數(shù)F（t）由兩參數(shù)威布爾函數(shù)F1（t）、F2（t）組成。當(dāng)0≤t≤t0時(shí)，F(xiàn)（t）=F1（t）；當(dāng)t0＜t＜∞時(shí)，F(xiàn)（t）=k F2（t）。F（t）表示為

式中，μ1，μ2分別為這兩個(gè)分布的尺度參數(shù)；β1，β2為兩個(gè)分布的形狀參數(shù)；為了保證概率密度函數(shù)f（t）和分布函數(shù)F（t）在是連續(xù)的，引入調(diào)節(jié)參數(shù)k，得

假設(shè)β1＞β2，按式（1）和（2）對(duì)t和F（t）做線性變化，得

由式（7）得出，WPP圖由直線L1和曲線C組成。當(dāng)x＜ln t0時(shí)，y=β1［x-ln（μ2）］是一條直線L1；當(dāng)時(shí)x→∞，漸進(jìn)線L2為y（x）=β2［x-ln（μ2）］，進(jìn)而，繪出WPP圖，如表2所示。

2.2 兩重混合型威布爾分布

Mendenhall［9］在1958年最早使用兩重混合型威布爾分布對(duì)故障數(shù)據(jù)進(jìn)行擬合，從那以后，兩重混合型威布爾分布模型被廣泛使用。

兩重威布爾分布模型的分布函數(shù)為

式中，p+q=1。經(jīng)過對(duì)上式按式（2），（3）線性轉(zhuǎn)化，可得該模型的WPP圖趨勢(shì)為：當(dāng)x=-∞，曲線C漸近線L1為y（x）=β1［x-ln（μ1）］+ln p；當(dāng)x→+∞，曲線C漸近線L1為y（x）=β1［x-ln（μ1）］，可見兩條漸近線相互平行，WPP圖如表2所示。

2.3 兩重相乘型威布爾分布

兩重相乘型威布爾分布由兩種類型組成：其一，為兩重競(jìng)爭(zhēng)型威布爾分布，故障分布函數(shù)如式（8），主要應(yīng)用于兩個(gè)同時(shí)符合標(biāo)準(zhǔn)威布爾分布的系統(tǒng)串聯(lián)之后的故障分布擬合；其二，為兩重并聯(lián)型威布爾分布，故障分布函數(shù)如式（9），主要應(yīng)用于兩個(gè)同時(shí)符合標(biāo)準(zhǔn)威布爾分布的系統(tǒng)并聯(lián)之后的故障分布擬合。

式（8），（9）中，F(xiàn)（t）為總系統(tǒng)故障分布函數(shù)，F(xiàn)1（t）和F2（t）分別為兩個(gè)分系統(tǒng)的故障分布函數(shù)。兩重競(jìng)爭(zhēng)型威布爾分布和兩重并聯(lián)型威布爾分布分析類似，由于篇幅有限，本文僅對(duì)兩重競(jìng)爭(zhēng)型威布爾分布做詳盡分析。

對(duì)式（8）按照式（2）、（3）進(jìn)行線性轉(zhuǎn)化得

當(dāng)x→∞，曲線C漸近線L1為y（x）=β1［x-ln（μ1）］；當(dāng)x→∞，曲線C漸近線L1為y（x）=β2［x-ln（μ2）］，繪出WPP圖如表2所示。

表2 兩重威布爾分布模型及WPP圖

2.4 選擇模型

經(jīng)過與本文故障數(shù)據(jù)的WPP圖比較，兩重分段型威布爾分布和兩重相乘型威布爾分布的WPP圖與本文故障數(shù)據(jù)的WPP圖形狀相似，又因?yàn)楣收蠒r(shí)間較短時(shí)，分段型威布爾分布的WPP圖為一條直線，相乘型威布爾分布為一條曲線，經(jīng)過與原故障數(shù)據(jù)相比較，發(fā)現(xiàn)兩重分段型威布爾分布與本文故障數(shù)據(jù)分布趨勢(shì)最為相似，所以采用兩重分?jǐn)嘈屯紶柗植紝?duì)本文故障數(shù)據(jù)進(jìn)行擬合。

3 兩重分段威布爾分布參數(shù)確定

常用求解兩重威布爾參數(shù)的方法為圖解法，圖解法優(yōu)點(diǎn)是簡(jiǎn)單、直觀，但參數(shù)求解不準(zhǔn)確，受主觀影響大，同時(shí)極有可能造成假設(shè)的時(shí)間分界點(diǎn)與計(jì)算得出的時(shí)間分界點(diǎn)不統(tǒng)一，前后矛盾等問題。

本文提出的迭代尋址求最優(yōu)解方法，有效克服了上述問題。

3.1 擬合右端漸近線與確定時(shí)間分界點(diǎn)苑圍

由式（7）得出，當(dāng)時(shí)x→∞，曲線C的漸近線為y（x）=β2［ x-ln（μ2）］，由故障數(shù)據(jù)的WPP圖，求出第二段的漸近線方程為

進(jìn)而，求出β2=1.0851，μ2=2558.1。

接下來(lái)，確定時(shí)間分界點(diǎn)t0范圍。根據(jù)WPP圖選取分界時(shí)間點(diǎn)t0的變化范圍，選取t0合適的變化氛圍可以減少計(jì)算量，初步去除掉不合適的點(diǎn)，由圖3可以得出，ln（t0）∈7.4，8］，t0∈［1652，2275］，得分界點(diǎn)序數(shù)n0的取值范圍n0∈［60，80］。

3.2 迭代尋址求最優(yōu)解

首先，通過迭代符合范圍的n0，采用最小二乘法或者極大釋然法對(duì)序數(shù)從1到n0的點(diǎn)進(jìn)行擬合，擬合出直線L1，分別求出不同n0對(duì)應(yīng)的擬合曲線特征參數(shù)μ2和β1，如表3所示。

表3 選取不同分界點(diǎn)序數(shù)對(duì)應(yīng)的特征參數(shù)

確定約束條件，選取合理解。根據(jù)式（12），通過表3求出的不同n0對(duì)應(yīng)的μ2和β1，計(jì)算得出與之對(duì)應(yīng)的t0，t0為分界時(shí)間點(diǎn)，需滿足

式中，tn0

為序數(shù)n0的故障時(shí)間值，tn0+1為序數(shù)n0+1的故障時(shí)間值。如果沒有符合要求的t0，稍微修改μ2和β1后，從第一步重新計(jì)算。

本文的故障數(shù)據(jù)，經(jīng)過計(jì)算后，求出符合要求的為62和68，具體參數(shù)如表4所示。

表4 特征參數(shù)表

一般會(huì)有多個(gè)解滿足約束條件，同時(shí)，各解的t0間距較大，可能會(huì)間距數(shù)百小時(shí)，如何選取合適的t0尤為關(guān)鍵。

本文以曲線擬合質(zhì)量為目標(biāo)優(yōu)化量，尋求與故障數(shù)據(jù)分布最相符的函數(shù)為最優(yōu)解。擬合質(zhì)量指標(biāo)e可以比較不同威布爾分布的擬合數(shù)據(jù)的良好程度。e值越小，說明擬合質(zhì)量越好。擬合質(zhì)量e的計(jì)算公式為

式中，F(xiàn)（tj）為擬合曲線的故障分布函數(shù)，即式（4）；F1（ti）為中位秩估計(jì)的故障分布函數(shù)，即式（1）。

當(dāng)時(shí)n0=60，e=0.1248；當(dāng)時(shí)n0=68，e=0.1233。所以，n0=68是最優(yōu)解，擬合形狀和故障分布最相似。

綜上所述，對(duì)分段威布爾分布的特征參數(shù)確定進(jìn)行了詳細(xì)的敘述，借助于matlab等計(jì)算工具，采用迭代尋址求最優(yōu)解方法可以有效提高計(jì)算效率，避免了求得的分界時(shí)間點(diǎn)不在假設(shè)的分界范圍等不合理情況的發(fā)生，同時(shí)，提出以曲線擬合質(zhì)量e為衡量指標(biāo)，取得和原故障數(shù)據(jù)最相似的擬合曲線。

3.3 故障擬合曲線與故障率函數(shù)

把求得特征參數(shù)代入式（4），分別求得分布函數(shù)如下式所示。

繪出擬合曲線L1和L2，如圖4所示。

圖4 分段威布爾模型的WPP擬合曲線

由圖4得出，擬合曲線基本與故障分布曲線相似，表明該分段模型是一個(gè)恰當(dāng)模型。

故障率函數(shù)λ（t）=f（t）/F（t），式中f（t）為分布密度函數(shù)，可求得

把求得的β1、μ1、β2、μ2代入到式，求得

加工中心的故障率曲線見圖5，從圖中可以得出，該故障率曲線包含早期故障期和偶然故障期。當(dāng)運(yùn)行時(shí)間小于1840h時(shí)，故障率呈現(xiàn)下降趨勢(shì)；當(dāng)運(yùn)行時(shí)間大于1840h時(shí)，故障率基本不變，與浴盆曲線變化趨勢(shì)基本相符。以每周7天，兩班工作制計(jì)算，早期故障期為115天，大約為4個(gè)月，與工廠反應(yīng)的故障發(fā)生情況基本相符。

圖5 分段威布爾模型的故障率擬合曲線

4 擬合優(yōu)度與趨勢(shì)檢驗(yàn)

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

常用的解析檢驗(yàn)法有皮爾遜檢驗(yàn)法、d檢驗(yàn)法等［11］。d檢驗(yàn)法具有更加精細(xì)，適用于小樣本等優(yōu)點(diǎn)。本文采用d檢驗(yàn)法對(duì)擬合的優(yōu)度檢驗(yàn)。

按照上式，計(jì)算得出Dn=0.093。取顯著性水平位α=0.10，查表的Dn，α=0.1554。因?yàn)镈n＜Dn，α，因此擬合函數(shù)可以被接受。由圖6可以看出68＜n＜90，F(xiàn)（ti）和F1（ti）偏差較大，這與WPP圖的擬合結(jié)果相同，可以看出WPP圖可以直觀反應(yīng)故障分布情況。

圖6 分段威布爾模型的故障分布函數(shù)擬合曲線

4.2 趨勢(shì)檢驗(yàn)

Laplace檢驗(yàn)也稱為U檢驗(yàn)，可以判明產(chǎn)品在實(shí)驗(yàn)中是否出現(xiàn)可靠性增長(zhǎng)、下降或者不變［12］。本文采用Laplace檢驗(yàn)法對(duì)早期故障期，和偶然故障期故障數(shù)據(jù)進(jìn)行檢驗(yàn)，驗(yàn)證模型以及時(shí)間分界點(diǎn)選取是否正確。

計(jì)算多臺(tái)定時(shí)截尾試驗(yàn)的統(tǒng)計(jì)量u為

式（10）中，T為試驗(yàn)總時(shí)間；M為樣本總數(shù)；t1為故障發(fā)生時(shí)間。

當(dāng)時(shí)u≤uα/2，表明在顯著水平α下，有顯著的可靠性增長(zhǎng)；當(dāng)時(shí)u≥u1-α/2，表明在顯著水平α下，有顯著的可靠性下降；若uα/2＜u＜u1-α/2，表明沒有明顯的變化趨勢(shì)，u值越接近0，越?jīng)]有變化趨勢(shì)。

4.2.1 早期故障期趨勢(shì)檢驗(yàn)

早期故障期為運(yùn)行時(shí)間小于1840h，發(fā)生故障總數(shù)M=68，T=1840h，代入式（10），求得u=-2.1593。取顯著程度α=0.05，查表得u0.025=-1.96，u1-0.025=1.96。由于u＜u0.025，得出在顯著水平0.05下，可靠性有顯著性增長(zhǎng)，故障率下降，與本文分析相符合。

4.2.2 偶然故障期趨勢(shì)檢驗(yàn)

偶期故障期的截尾時(shí)間為7300h，共發(fā)生了42次故障，把M=42，T=7300h，代入式（10），求得u=-0.9087。還是取顯著程度α=0.05，由于u0.025＜u＜u1-0.025，得出在顯著水平0.05下，可靠性沒有顯著性變化，處于偶然故障期，與本文分析相符合。

綜上所述，通過采用d檢驗(yàn)法和Laplace檢驗(yàn)法分別對(duì)擬合曲線的擬合優(yōu)度和故障數(shù)據(jù)可靠性變化趨勢(shì)進(jìn)行檢驗(yàn)，驗(yàn)證結(jié)果與原分析結(jié)果完全相同，擬合曲線及時(shí)間分界點(diǎn)可以被接受。

5 結(jié)論

本文以18個(gè)月采集的26臺(tái)加工中心的故障數(shù)據(jù)為基礎(chǔ)，列出了3種常用的兩重威布爾函數(shù)的WPP圖，從中選擇了兩重分段威布爾分布作為擬合模型，并以擬合質(zhì)量為優(yōu)化目標(biāo)，求出分布模型的特征參數(shù)和時(shí)間分界點(diǎn)，使用d檢驗(yàn)法和Laplace檢驗(yàn)法對(duì)結(jié)果進(jìn)行擬合優(yōu)度和趨勢(shì)檢驗(yàn)，檢驗(yàn)結(jié)果滿足要求，與原分析結(jié)果相一致，證明擬合函數(shù)可以被接受，對(duì)早期故障的預(yù)測(cè)和控制有一定指導(dǎo)意義。

［1］申桂香，賈亞洲，馬健.CNC機(jī)床故障分析及其可靠性［J］.中國(guó)機(jī)械工程，1996，7（6）：67-69.

［2］WANG Yi-qiang，JIA Ya-zhou，JIANG Wei-wei.Early Failure Analysis of Machining Centers：a Case Study［J］. Reliability Engineering and System Safety，2001，72：91-97.

［3］范秀君，許靜林，張根保.數(shù)控機(jī)床早期故障消除技術(shù)［J］.中國(guó)機(jī)械工程，2013，24（16）：2241-2247.

［4］鄒輝，黃筱調(diào)，洪榮晶.基于兩重威布爾分布的大型數(shù)控銑齒機(jī)可靠性評(píng)價(jià)［J］.機(jī)械設(shè)計(jì)與制造，2010（2）：198-200.

［5］陳其紅，闞樹林，秦臻.基于兩重威布爾分布的鑄造造型機(jī)可靠性分析［J］.機(jī)械制造，2012，50（2）：20-24.

［6］張根保，唐賢進(jìn)，廖小波.數(shù)控機(jī)床故障率分布曲線定量化建模及應(yīng)用［J］.重慶大學(xué)學(xué)報(bào)，2013，36（6）：119-136.

［7］陳殿生，王田苗，魏洪興.數(shù)控車床故障分布的兩重威布爾分段模型［J］.北京航空航天大學(xué)學(xué)報(bào)，2005，31（7）：767-769.

［8］顧蓓青，王蓉華，徐蹺嶺，等.威布爾兩重分段模型的統(tǒng)計(jì)分析［J］.山西大學(xué)學(xué)報(bào)（自然科學(xué)版），2008，31（3）：359-363.

［9］Murthy D.N.Weibull Models［M］.A JOHN WILEY& SONS，PUBLICTION，2004.

［10］張英芝，申桂香，吳甦.數(shù)控機(jī)床可用度分析［J］.制造技術(shù)與機(jī)床，2010（1）：124-126.

［11］Zequeira R，Valdes J.An approach for the Bayesian estimation in the case of ordered piecewise constant failure rates［J］.Reliability Engineering and System Safety，2001，72（3）：227-240.

［12］PulciniG.An exponential reliability-growth model in multicopy testing program［J］.IEEE Transactions on Reliability，2001，50（4）：365-373..

（編輯 趙蓉）

Analysis of Machining Center by Models Involving two Weibull Distributions

CONG Ming1，WANG Guan-xiong1，XIE Shu-wen2，PAN Da-qing2
（1.School of Mechanical of Engineering，Dalian University of Technology，Dalian Liaoning116023，China；2. Chongqing Changan Cars Co.，Ltd，Chongqing 404100，China）

In order to decrease and predict the damage of early failures，the failure datas of 26 machining centers were collected.Three models involving two weibull distributions were studied.Depending on comparing the plots on Weibull plotting paper（WPP），the sectional model is selected as the best one to model the failure data.A new method of calculating parameters is introduced，which is useful to avoid the error of demarcation time selection and improve the accuracy of fitting.It can be calculated that earlier failure period is about4 months，which is coincident with the consumer feedback.D test method and Laplace test method was used to examine the model，and the results validated the sectional model.

early failure；bathtub curve；machining centers

TH166；TG659

A

1001-2265（2015）04-0078-05 DOI：10.13462/j.cnki.mmtamt.2015.04.020

2014-07-09；

2014-10-17

國(guó)家＂高檔數(shù)控機(jī)床與基礎(chǔ)制造設(shè)備＂科技重大專項(xiàng)課題（2011ZX04015-021）

叢明（1963—），男，遼寧大連人，大連理工大學(xué)教授，博士，研究方向?yàn)闄C(jī)器人技術(shù)、可靠性技術(shù)，優(yōu)化設(shè)計(jì)，（E-mail）reliability1@163. com。