蔣仁言, 張碧雯

(長(zhǎng)沙理工大學(xué) 汽車(chē)與機(jī)械工程學(xué)院,湖南 長(zhǎng)沙 410114)

0 引言

威布爾分布是建模非負(fù)、連續(xù)隨機(jī)變量使用最廣泛的分布模型[1]；威布爾分布的更新函數(shù)(簡(jiǎn)稱(chēng)威布爾更新函數(shù))在可靠性和可用度分析、維修決策優(yōu)化、質(zhì)保費(fèi)用分析、備件庫(kù)存管理、排隊(duì)系統(tǒng)、交通流分析等領(lǐng)域有廣泛的應(yīng)用[2～8]。威布爾更新函數(shù)沒(méi)有解析表達(dá)式，這給解各種涉及更新函數(shù)的優(yōu)化問(wèn)題帶來(lái)極大不便。因此，開(kāi)發(fā)威布爾更新函數(shù)的近似式已經(jīng)吸引了持續(xù)的注意。

根據(jù)近似式在時(shí)間軸t的適用范圍，可將各種威布爾更新函數(shù)的近似式分為以下兩大類(lèi)。第一類(lèi)近似式適合于中偏小的t值[9]，另一類(lèi)近似式適合于整個(gè)t值范圍(即0到∞)。本文著重討論第二類(lèi)近似式。

根據(jù)模型的構(gòu)建方式，第二類(lèi)近似式可進(jìn)一步分為三個(gè)子類(lèi)：分段函數(shù)近似[10～13]，級(jí)數(shù)近似[14～16]和混合模型形式的近似[9,17]。一個(gè)好的近似式應(yīng)同時(shí)滿(mǎn)足精確性和簡(jiǎn)單性?xún)蓚€(gè)要求。分段模型在結(jié)構(gòu)上簡(jiǎn)單，但其精度通常偏低。對(duì)于某些分布函數(shù)(如正態(tài)分布和伽瑪分布)，級(jí)數(shù)模型的精度隨項(xiàng)數(shù)的增加而增大，但過(guò)多的項(xiàng)數(shù)導(dǎo)致近似更復(fù)雜[18]。混合形式的近似采用一個(gè)權(quán)函數(shù)將兩個(gè)極限關(guān)系平滑地連接成一個(gè)函數(shù)；在結(jié)構(gòu)上相對(duì)地簡(jiǎn)單，且其精度比分段函數(shù)近似要高，即此類(lèi)近似能較好地滿(mǎn)足精確性和簡(jiǎn)單性要求。但是，現(xiàn)有混合形式的近似當(dāng)威布爾形狀參數(shù)較大時(shí)，其精度仍有一定的改進(jìn)空間[17]。

為了改進(jìn)現(xiàn)有近似式在形狀參數(shù)較大時(shí)精度欠高的問(wèn)題，本文提出一個(gè)新的混合形式的威布爾更新函數(shù)近似式，可供威布爾形狀參數(shù)是大的(>3.65)情況時(shí)使用。此外，通過(guò)對(duì)比現(xiàn)有各近似式在不同形狀參數(shù)取值范圍內(nèi)的精度，識(shí)別好的近似式，確定其適用范圍，從而確立一個(gè)基于形狀參數(shù)的近似式選擇方法。

本文結(jié)構(gòu)如下。第1節(jié)介紹現(xiàn)有威布爾更新函數(shù)的近似式，分析其精度和適用范圍。第2節(jié)提出新的近似式，分析其精度和適用范圍。第3節(jié)通過(guò)一個(gè)維修政策優(yōu)化的數(shù)例例證提出的近似式的精確性和有用性。第4節(jié)總結(jié)全文。

1 現(xiàn)有近似式及其精度

1.1 更新函數(shù)精確值的計(jì)算

計(jì)算更新函數(shù)精確值是評(píng)價(jià)一個(gè)近似式精確性的先決條件。令F(t)記非負(fù)隨機(jī)變量T的分布函數(shù)。本文著重考慮F(t)是威布爾分布的情況，其表達(dá)式為:

F(t)=1-exp[-(t/η)β]

(1)

這里，β是形狀參數(shù)，η是尺度參數(shù)。在可靠性和維修應(yīng)用領(lǐng)域，β一般大于1，且很少超過(guò)4[12]。平均(μ)、標(biāo)準(zhǔn)偏差(σ)和變異系數(shù)(cv)分別為

μ=ηΓ(1+1/β),

σ=η[Γ(1+2/β)-Γ2(1+1/β)]0.5,cv=σ/μ

(2)

這里Γ(.)是伽瑪函數(shù)。

令F(k)(t)記F(t)的k重卷積(F(1)(t)=F(t))，N(t)記時(shí)間區(qū)間(0,t]內(nèi)的更新事件數(shù)。更新函數(shù)(記為M(t))是離散隨機(jī)變量N(t)的數(shù)學(xué)期望，可表達(dá)為下面的級(jí)數(shù)形式或積分形式[15,19]:

(3)

對(duì)于大多數(shù)分布函數(shù)，包括威布爾分布，上述的卷積和積分沒(méi)有解析表達(dá)式。因此，開(kāi)發(fā)更新函數(shù)的近似式成為一個(gè)重要的研究方向。最簡(jiǎn)單、最著名的兩個(gè)近似關(guān)系為:

(4)

前者適合小t的情況，后者適合于大t的情況。

對(duì)于一個(gè)給定的近似式(記為M∞(t)),一個(gè)重要的問(wèn)題是確定它在時(shí)間軸t上的適用范圍。這個(gè)問(wèn)題等價(jià)于評(píng)價(jià)近似式的精度。為此，需要知道一些β和t的組合下的更新函數(shù)的精確值。對(duì)于一個(gè)給定的β和t的組合，其精確值可通過(guò)數(shù)值積分法解式(3)中的積分方程獲得[20,21]。文獻(xiàn)[17]的補(bǔ)充材料(Supplementary Material)顯示了η=1，β=1.0(0.5)4.5 (即從1.0到4.5，步長(zhǎng)為0.5)和t=0.05(0.05)3.00(即從0.05到3.00，步長(zhǎng)為0.05)條件下的威布爾更新函數(shù)精確值。這樣，近似式Ma(t)的精確性可由下式評(píng)價(jià):

ε(t;β)=|1-Ma(t;β)/M(t;β)|

(5)

對(duì)于一個(gè)給定的β，令εM記上述時(shí)間范圍內(nèi)ε(t;β)值中的最大值，用它作為評(píng)價(jià)近似式精確性的測(cè)度。

方程(4)的第一個(gè)關(guān)系的適用范圍可定義為(0,τ1)[17]，其中

τ1=sup{t;ε(t;β)<0.01}

(6)

方程(4)的第二個(gè)關(guān)系的適用范圍可定義為(τ2,∞ )[17]，其中

τ2=inf{t:ε(t;β)≤0.01}

(7)

這樣，“小的t”意味著t<τ1；“中等大小的t”意味著τ1τ2。

表1 近似式的εM值

1.2 現(xiàn)有近似式及其精度

以下主要介紹分段形式和混合形式的近似式；級(jí)數(shù)形式的近似式因缺乏簡(jiǎn)單性不作介紹。

最早的分段形式的近似式可追溯到Spearman[10]，其表達(dá)式為

Ma(t)=max[F(t),M2(t)]

(8)

(9)

對(duì)于威布爾分布，F(xiàn)(t)和M2(t)有一個(gè)唯一的交點(diǎn)[12]。因此，式(8)實(shí)際上是一個(gè)兩重分段模型。表1第2列顯示了式(8)的εM值。從中可以看出，該近似的精度偏低。

Parsa和Jin[13]提出一個(gè)可積的三重分段線(xiàn)性近似:

(10)

最近，Jiang[11]提出下面的兩重分段近似式

(11)

其中，ts是分界點(diǎn)，根據(jù)該點(diǎn)處的連續(xù)和平滑條件予以確定

M1(t)=H(t)/[1+αH(t)],α≥0

(12)

H(t)是累積風(fēng)險(xiǎn)函數(shù)，α是ts的一個(gè)已知函數(shù)。表1第4列給出了這個(gè)近似式的εM值。從中可以看出，這個(gè)近似式比前述兩個(gè)近似式更簡(jiǎn)單、更精確。但是，隨著β的增大，εM值快速增大。

Jiang[12]修改式(12)為

M1(t)-F(t)/[1-αF(t)],α≥0

(13)

分界點(diǎn)ts以類(lèi)似的方式確定，α是ts的已知函數(shù)。表1第5列給出了該近似式的εM值。從中可以看出，當(dāng)β較大時(shí)，該近似式比由式(11)和(12)所構(gòu)成的近似式更精確。然而，當(dāng)β較大時(shí)，εM值仍然偏大。

Jiang和Chen[17]提出下面的混合形式的近似式:

Ma(t)=w(t)M1(t)+[1-w(t)]M∞(t)

(14)

`其中，由下式確定[9]：

M1(t)=pF(t)+(1-p)H(t)
p=1-exp{-{(β-1)/0.873}0.9269}

(15)

權(quán)函數(shù)是β的函數(shù)：

w(t)=1-Φ(t;μw,σw)

(16)

這里，Φ(t;μw,σw)是參數(shù)為μw和σw的正態(tài)分布函數(shù)，μw=(0.9139+0.2020β)η,σw=(|0.6302β-2.0001|+0.1226)η/6。

2 提出的近似式及其精度

2.1 提出的近似式

提出的近似關(guān)系取式(14)的形式，具有不同的M1(t)；權(quán)函數(shù)w(t)有不同的參數(shù)。

Baker[22]提到，當(dāng)t>2μ時(shí)，式(4)中的M∞(t)一般地是精確的。這意味，M1(t)應(yīng)該恰當(dāng)?shù)卮_定使其在時(shí)間區(qū)間(0,2μ)內(nèi)是精確的。下式可能滿(mǎn)足這個(gè)要求：

(17)

為簡(jiǎn)單起見(jiàn)，使用一個(gè)分布函數(shù)G2(t)近似F(2)(t)；使用另一個(gè)分布函數(shù)G3(t)近似F(3)(t)。這樣，M1(t)被定義為

M1(t)=F(t)+G2(t)+G3(t)

(18)

考慮威布爾、伽瑪、正態(tài)和對(duì)數(shù)正態(tài)分布作為G2(t)和G3(t)的備選模型。對(duì)于一個(gè)給定的β和η的組合，使用式(2)計(jì)算F(t)的均值μ和方差σ2。則F(k)(t)(k=2,3)的平均為kμ，方差為kσ；Gk(t)的參數(shù)用矩法確定。例如，如果Gk(t)是形狀參數(shù)為uk、尺度參數(shù)為v的伽瑪分布函數(shù)，則

uk=kμ2/σ2,v=σ2/μ

(19)

表2 備選模型的εA值

為了對(duì)比備選模型的性能，類(lèi)似于計(jì)算威布爾更新函數(shù)的精確值，使用數(shù)值積分方法計(jì)算F(k)(t)的精確值。對(duì)于一個(gè)給定的備選模型和β值，計(jì)算在t/η=0.05(0.05)2.00處的相對(duì)誤差εg(t)，εA令記這些相對(duì)誤差的平均，用它作為評(píng)價(jià)備選模型性能的測(cè)度。采用相對(duì)誤差(而非絕對(duì)誤差)的目的是為了強(qiáng)調(diào)備選模型的左尾特征，因?yàn)镸1(t)適用于中偏小的t值。

表2顯示了備選模型的值εA。它清楚地表明，伽瑪分布是最好的。因此，G2(t)和G3(t)是伽瑪分布，其參數(shù)由式(19)計(jì)算?，F(xiàn)在確定權(quán)函數(shù)的參數(shù)：μw和σw。仔細(xì)分析M1(t)和M∞(t)之間的交點(diǎn)發(fā)現(xiàn)有以下三種情況：

·對(duì)于β=1, 它們相交在原點(diǎn)。

·當(dāng)β大于且接近于1時(shí)，有一個(gè)交點(diǎn)(參見(jiàn)圖1，ΔM(t)=M1(t)-M∞(t))。

·當(dāng)β較大時(shí)，M1(t)和M∞(t)之間有多個(gè)交點(diǎn)(參見(jiàn)圖2)。用t1記最小的交點(diǎn)，用t2記最大的交點(diǎn)。在這兩點(diǎn)之間，隨著t的增大，和越來(lái)越接近。

圖1 ΔM(t)=M1(t)-M∞(t)的圖形(β=1.5)

對(duì)于β=1.0(0.5)4.5，表3第2和第3行給出了對(duì)應(yīng)的t1/η和t2/η的值。平均參數(shù)μw取t1和t2的加權(quán)和，即

μw=w1t1+w2t2

(20)

為使μw更接近于t2，假定wk(k=1,2)正比于tk。這得wk=tk/(t1+t2)。式(20)成為

(21)

表3 t1/η、t2/η、μw/η和σw/η的值

標(biāo)準(zhǔn)偏差參數(shù)σw取為

σw=(t2-μw)/Φ-1(0.999;0,1)=(t2-μw)/3.0902

(22)

這里Φ-1(.)是正態(tài)分布的逆函數(shù)。當(dāng)t1=0時(shí)，μw=t2，取σw=0.01η。表3最后兩行給出了μw/η和σw/η的值。

2.2 本文近似的精度

表1倒數(shù)第2列給出了本文提出的近似式的εA值；最后一列給出最好的近似式。從中可以看出，在β=1.5和2.0之間，最好的近似式由文獻(xiàn)[11]的近似變?yōu)槲墨I(xiàn)[17]的近似。使用二次多項(xiàng)式插值法推斷兩者在β=1.69處有相同的εA值。類(lèi)似地，在β=3.5和4.0之間，最好的近似式由文獻(xiàn)[17]的近似變?yōu)楸疚奶岢龅慕?，并且推斷兩者在?3.65處有相同的εA值。因此，我們可以得出以下結(jié)論：

·對(duì)于β<1.69，由Jiang[11]所提出的近似式提供最好的精度；

·對(duì)于β=1.69-3.65，由Jiang和Chen[17]所提出的近似式有最高的精度；

·對(duì)于β>3.65，本文所提出的近似式是最精確的。

·按照以上β值范圍選擇模型，最大相對(duì)誤差將小于2%。這個(gè)精度滿(mǎn)足一般應(yīng)用的要求。

3 示例

當(dāng)一個(gè)產(chǎn)品含有幾個(gè)相同的元件(如軸承)，作為一個(gè)基于時(shí)間的預(yù)防維修政策，批更換政策(block replacement policy)可能適用于這些元件的預(yù)防更換[23]。在這個(gè)政策下，每隔時(shí)間周期Tp，元件被預(yù)防地更換；如果元件在此之前發(fā)生失效，則只更換失效的元件。

在時(shí)間區(qū)間((k-1)Tp,kTp)(k=1,2,…) 內(nèi)，每個(gè)元件的失效事件形成一個(gè)更新過(guò)程，平均失效數(shù)為更新函數(shù)M(t)。令cp記每個(gè)元件的平均預(yù)防更換費(fèi)用，cf記失效更換費(fèi)用。通常，cf比cp大得多。在上述時(shí)間區(qū)間內(nèi)的期望費(fèi)用率為[23]

J(Tp)=[cp+cfM(Tp)]/Tp

(23)

如果采用失效更換政策，費(fèi)用率將是J(∞)=cf/μ。批更換政策的有效性可用維修費(fèi)用節(jié)省率描述：

(24)

此值越大越好。

設(shè)想元件壽命服從β=4.0和平均壽命為700小時(shí)的威布爾分布；cf和cp的值分別為100和250費(fèi)用單位。使用精確的更新函數(shù)求得的最優(yōu)解顯示在表4第一列。如果采用失效更換政策，費(fèi)用率將是0.3571。這表明，批更換政策能帶來(lái)ηe=21.42%的維修費(fèi)用節(jié)省。因此，對(duì)于本文數(shù)例，批更換政策的有效性十分明顯。

表4 基于不同近似式的最優(yōu)解和相對(duì)誤差值

為進(jìn)一步例證本文所提出的近似式的精確性，表4最后五列顯示了從其它五個(gè)近似式獲得的最優(yōu)解和相對(duì)誤差值。從中可以看出，只有從文獻(xiàn)[17]的近似式得到的最優(yōu)解的精度是可以接受的。這例證了基于β值選擇近似式的必要性。

表5 本文近似式的參數(shù)

3 總結(jié)和結(jié)論

本文系統(tǒng)地總結(jié)了威布爾更新函數(shù)的近似計(jì)算式，針對(duì)現(xiàn)有近似式在β較大時(shí)精度不高的問(wèn)題，提出了一個(gè)新的近似式。精度分析發(fā)現(xiàn)，當(dāng)β<1.69時(shí)，由Jiang[11]所提出的近似式是最精確的；當(dāng)β∈[1.69,3.65]時(shí)，由Jiang和Chen[17]所提出的近似式是最精確的；當(dāng)β>3.65時(shí)，本文所提出的近似式是最精確的。若按此選擇近似式，最大相對(duì)誤差將小于2%。

本文數(shù)例例證了所提出的近似式的精確性，基于β值選擇近似式的必要性，以及維修決策優(yōu)化的有效性。

注意到用伽瑪分布G2(t)和G3(t)近似F(2)(t)和F(3)(t)的精度仍有改進(jìn)的空間(見(jiàn)表2)，故尋找更好的分布模型供近似F(2)(t)和F(3)(t)是一個(gè)值得進(jìn)一步研究的課題。