溫玉全，張利敏，洪東跑

(北京理工大學(xué) 爆炸科學(xué)與技術(shù)國家重點實驗室，北京100081)

火工品發(fā)火可靠性試驗及其數(shù)據(jù)分析的主要目的是驗證或評估產(chǎn)品的可靠性［1］。由于火工品屬于一次性作用的敏感性產(chǎn)品，一般假設(shè)每個產(chǎn)品均存在一個臨界刺激量，當(dāng)外界施加的刺激量大于臨界刺激量時，該產(chǎn)品響應(yīng);否則不響應(yīng)。一般把產(chǎn)品的臨界刺激量的分布稱為感度分布。

火工品的臨界刺激量是不能直接測量的，其發(fā)火可靠性試驗是在若干刺激量下做刺激試驗，試驗數(shù)據(jù)是這些刺激量和相應(yīng)的響應(yīng)或不響應(yīng)數(shù)。為了對火工品的發(fā)火可靠性進行評估或驗證，目前國內(nèi)外最常用的試驗及其數(shù)據(jù)分析方法為:1)進行升降法試驗［2］，利用試驗數(shù)據(jù)，估計感度分布中的未知參數(shù)，進而求得分布函數(shù)某分位點估計或某刺激量處的可靠度估計。由于感度分布參數(shù)的估計值可能有較大的偏差，從而導(dǎo)致分位點估計或可靠度估計偏差較大，特別是對高可靠性要求的產(chǎn)品，這種偏差有時會達到難以容忍的程度。2)在工作刺激量處做成敗型試驗，應(yīng)用二項分布模型驗證或評估產(chǎn)品的可靠性。該方法需要的樣本量太大，如航空航天類火工品，在置信度0.90 或0.95 下，要求可靠度達到0.999 以上，為了鑒定產(chǎn)品的可靠性需要在工作刺激量處試驗2 303 發(fā)或2 996 發(fā)，且必須全部響應(yīng)。不僅試驗工作量大而且試驗費用昂貴。故本文結(jié)合工程應(yīng)用背景，綜合利用升降法試驗數(shù)據(jù)和固定刺激量處的成敗型試驗數(shù)據(jù)，給出了基于感度試驗的火工品發(fā)火可靠性試驗數(shù)據(jù)的分析方法。

1 感度分布參數(shù)估計

假設(shè)感度分布為F(μ，σ)，本文利用升降法試驗估計感度分布參數(shù)(μ，σ).升降法試驗包括3 個因素:試驗樣本量N、初始刺激量x0和步長d.x0和d 選定后，用x0作第1 次刺激-響應(yīng)試驗;第2 次及以后每次試驗所用刺激量的取法如下:如前一次試探的反應(yīng)結(jié)果為響應(yīng)，則本次試探用刺激量為xi+1=xi-d;如為不響應(yīng)，則為xi+1=xi+d.如此循環(huán)試驗，至完成預(yù)定試驗樣本量N 為止。由升降法試驗數(shù)據(jù)可以得到火工品感度分布參數(shù)(μ，σ)的極大似然估計，當(dāng)數(shù)據(jù)存在“混合區(qū)”(最大不響應(yīng)刺激量要大于最小響應(yīng)刺激量)時，極大似然估計(MLE)唯一［3］。在實際應(yīng)用中，只有數(shù)據(jù)存在混合區(qū)參數(shù)才有唯一的MLE，否則認(rèn)為試驗無效，需要重新進行升降法試驗。

有研究［3-4］表明是μ 的無偏估計，不是σ的無偏估計，且試驗方案(x0，n，d)對的影響較大，其中d 的影響最大。為了進一步研究試驗方案對參數(shù)估計的影響，利用蒙特卡羅方法模擬升降法試驗。模擬結(jié)果表明，試驗方案對的影響較小，對的影響較大，其中d 對試驗的有效性(是否出現(xiàn)數(shù)據(jù)“混合區(qū)”)及參數(shù)估計的精確性有較大的影響。通過對不同的試驗方案模擬比較表明，當(dāng)試驗樣本量相同時試驗方案x0=μ、d =σ 是最理想的，故在確定試驗方案應(yīng)該盡可能獲得(μ，σ)比較精確的預(yù)估。

1)綜合利用產(chǎn)品的可靠性信息確定升降法試驗方案(x0，n，d)，隨機抽取產(chǎn)品進行升降法試驗，由試驗數(shù)據(jù)得感度分布參數(shù)的MLE和

2 假設(shè)檢驗

假設(shè)某火工品的感度分布為正態(tài)分布N(μ，σ)，則在刺激量x 處的響應(yīng)概率為

假設(shè)某火工品的可靠性指標(biāo)為:在置信水平γ下，在刺激量xB處，可靠度為RB.則由(1)式可知對于給定的μ，必存在唯一的σb，滿足

由(2)式可得

式中uRB為標(biāo)準(zhǔn)正態(tài)分布的RB分位點。

在實際中μ 通常是未知的，由于升降法能較好的估計參數(shù)μ，一般通過升降法試驗獲得其估計值。為了使估計更精確，一般通過3 組樣本量為50 的升降法試驗，再求參數(shù)μ 估計值的平均值得.由于火工品在設(shè)計時要求μ ＞3σ，而實際的火工品往往有μ?3σ.模擬試驗結(jié)果表明當(dāng)μ ＞3σ 且n =50時，其中var (μ)為μ 的方差。取3 組的平均值，有由此可知|(-μ)/μ|較小，特別對于機械類火工品，|-μ|小于試驗誤差。

對于上述火工品，為了判定其可靠性是否達到指標(biāo)要求，要做如下的假設(shè)檢驗(顯著水平為1 -γ)

定理1 ?x ＞μ 有

σ≤σb?F(x，μ，σ)≥F(x，μ，σb).

證明 F(x)對σ 求導(dǎo)有

證畢

由定理1 可知，當(dāng)刺激量xB＞μ 時，(4)式等價于如下假設(shè):

根據(jù)定理1，通過在某一較低的刺激量xC下進行試驗(要求xC＞μ)，對(5)式的假設(shè)檢驗進行判斷，從而對(4)式的假設(shè)檢驗進行判斷，就可以實現(xiàn)對產(chǎn)品響應(yīng)可靠性的驗證。由(1)式可知R(xC)=F(xC，μ，σ)為產(chǎn)品在刺激量xC的響應(yīng)概率，RC=F(xC，μ，σb)為結(jié)合指標(biāo)產(chǎn)品在刺激量xC響應(yīng)概率的估計。再由定理1 可知(5)式等價于

因此在工作刺激量xB處檢驗R(xB)≥RB的問題，就可轉(zhuǎn)變?yōu)樵诖碳ち縳C處檢驗R(xC)≥RC的問題。

隨機抽取樣本量為nc的產(chǎn)品，在xC處進行試驗，記不響應(yīng)數(shù)為f.由于不響應(yīng)數(shù)服從二項分布即f～B(nc，1 -R(xC))，在顯著水平1 -γ 下，若不能拒絕(6)式中的H0，則有

由此可以確定檢驗的拒絕域，根據(jù)β 分布與F分布分位數(shù)的關(guān)系，有

當(dāng)nc和RC已知時由(7)式可確定f，隨機抽取樣本量為nc的產(chǎn)品，在xC處進行試驗，記不響應(yīng)數(shù)為f'，如果f'≤f，則不能拒絕(6)式的H0，從而不能拒絕(4)式的H0;否則，拒絕(6)式的H0，從而拒絕(4)式的H0.

3 抽樣方案

在刺激量點xC處進行樣本量為nc的試驗，取f=0，如果產(chǎn)品全部響應(yīng)，則不能拒絕(4)式的H0;否則接受(4)式的H1.對應(yīng)抽樣方案(nc，0)則有

式中:α 為生產(chǎn)方風(fēng)險;β 為使用方風(fēng)險。當(dāng)產(chǎn)品的響應(yīng)概率R≥R1時，則判斷產(chǎn)品是合格的，以高概率(大于或等于1 -α)接收;當(dāng)產(chǎn)品的響應(yīng)概率R≤R0時，則判斷產(chǎn)品不合格的應(yīng)該以低概率(小于或等于β)接收。結(jié)合產(chǎn)品可靠性指標(biāo)取R0=F(xC，，σb)，R1=F(xC，，)，則當(dāng)給定α 和β 時，(8)式無解或存在唯一解(xC，nc)，即當(dāng)(8)式有解的時候能唯一確定試驗刺激量點xC和抽樣方案(nc，0).一般取使用方風(fēng)險為β=1 -γ，代入(8)式可得

由于實際應(yīng)用中當(dāng)α 取值較小時，nc往往過大，導(dǎo)致試驗成本過高，在此把抽樣方案由一次抽樣(nc，0)，改為二次抽樣(nc2/nc1，0，2，1)，即隨機抽取樣本量為nc1產(chǎn)品在xC處試驗，如果全部響應(yīng)則不能拒絕(4)式H0;如果不響應(yīng)數(shù)為1，則再抽取樣本量為nc2產(chǎn)品在xC處試驗，如果全部響應(yīng)則不能拒絕(4)式的H0，否則接受(4)式的H1.其中nc1和nc2滿足

當(dāng)F(xC，，σb)已知時，解(10)式可確定唯一的抽樣方案(nc2/nc1，0，2，1).

與抽樣方案(nc，0)同理定義R0=F(xC，，σb)，R1=F(xC，，)，有下式成立

同理，解方程組(10)和(11)，如果方程無解，可適當(dāng)提高α;否則可得(nc1，nc2，xC).

4 基于小樣本的發(fā)火可靠性驗證方法

基于升降法試驗數(shù)據(jù)估計感度分布參數(shù)，利用假設(shè)檢驗原理，在某一較低刺激量點進行小樣本量的試驗，實現(xiàn)對產(chǎn)品發(fā)火可靠性的驗證或評估，并通過對雙方風(fēng)險的分析，確定試驗刺激量點及抽樣方案。假設(shè)某火工品的感度分布為正態(tài)分布(或通過變換后為正態(tài)分布)N(μ，σ)，可靠性指標(biāo)為:在置信水平γ 下，在工作刺激量xB處，響應(yīng)概率為RB，其小樣本試驗方法為

1)綜合利用已有的產(chǎn)品可靠性信息確定升降法試驗方案(x0，d，n)，進行有效的升降法試驗，并對σ 的估計值進行修正，得參數(shù)估計值(，).

3)在工程應(yīng)用中，由生產(chǎn)方對產(chǎn)品性能的了解確定生產(chǎn)方風(fēng)險α.

4)根據(jù)(α，γ)，由(10)式和(11)式(當(dāng)有解時)可得唯一的試驗方案(nc1，nc2，xC)，如果無解，適當(dāng)?shù)靥岣擀?，代入重新計算?/p>

隨機抽取樣本量為nc1的產(chǎn)品在xC處試驗，如果全部響應(yīng)則認(rèn)為該批產(chǎn)品的發(fā)火可靠性水平達到指標(biāo)要求;如果不響應(yīng)數(shù)大于1，則認(rèn)為該批產(chǎn)品的發(fā)火可靠性水平達不到指標(biāo)要求;如果不響應(yīng)數(shù)為1，則再抽取樣本量為nc2的產(chǎn)品在xC處試驗，如果全部響應(yīng)則認(rèn)為該批產(chǎn)品的發(fā)火可靠性水平達到指標(biāo)要求，否則認(rèn)為該批產(chǎn)品的發(fā)火可靠性水平達不到指標(biāo)要求。

5 實例

某撞擊火帽要求可靠性指標(biāo)γ =0.90，發(fā)火可靠度R≥0.999，響應(yīng)刺激量為:落錘質(zhì)量(388±1)g，落高100 mm.

根據(jù)火工品感度分布模型研究結(jié)果，其感度分布服從對數(shù)正態(tài)分布。因其感度分布參數(shù)未知，故先確定試驗升降法試驗方案，進行了1 組有效試驗，得參數(shù)的MLE.再綜合歷史信息和試驗獲得的參數(shù)的MLE，確定下一組試驗方案，由此一共進行3 組試驗，每組試驗樣本量為50，試驗數(shù)據(jù)列于表1中。

表1 3 組撞擊火帽升降法試驗數(shù)據(jù)Tab.1 Three groups of data in up-and-down test

對表1的試驗數(shù)據(jù)進行分析得參數(shù)的MLE，并對σ 的MLE 進行修正，結(jié)果如表2所示。

表2 3 組升降法分布參數(shù)估計值Tab.2 Sensitivity distribution parameter estimation

為了方便工程應(yīng)用，一般取生產(chǎn)方風(fēng)險為α =0.1，則對應(yīng)于試驗方案(nc2/nc1，0，2，1)，可由(11)式和(12)式得nc1=32，nc2=17，xC=6.88.

隨機抽取32 發(fā)產(chǎn)品在xC=6.88 處進行試驗，產(chǎn)品全部響應(yīng)，表明該產(chǎn)品的可靠度達到了指標(biāo)要求。

由于步進法試驗樣本量較大，其參數(shù)估計較為穩(wěn)定，根據(jù)步進法試驗數(shù)據(jù)估計分布參數(shù)，求出滿足可靠性指標(biāo)的刺激量上限，結(jié)合產(chǎn)品的技術(shù)指標(biāo)，就可以可判定產(chǎn)品是否達到了可靠性指標(biāo)要求。某撞擊火帽的步進法試驗數(shù)據(jù)如表3所示。

表3 撞擊火帽步進法試驗數(shù)據(jù)Tab.3 Test data of stepping method

6 結(jié)論

根據(jù)火工品感度分布函數(shù)的特性，結(jié)合假設(shè)檢驗原理，綜合利用了升降法試驗數(shù)據(jù)和固定刺激量處的成敗型試驗數(shù)據(jù)，給出了基于感度試驗的火工品發(fā)火可靠性小樣本評估與驗證方法。通過對某撞擊火帽進行大小樣本對比試驗，驗證了本方法的正確性和可行性，利用本方法，可以用200 以下的樣本量實現(xiàn)對高可靠度要求的火工品進行發(fā)火可靠性驗證。

References)

［1］田玉斌，李國英，房永飛.火工品可靠性試驗數(shù)據(jù)的綜合分析方法［J］.系統(tǒng)科學(xué)與數(shù)學(xué)，2006，26 (2):147 -158.TIAN Yu-bin，LI Guo-yin，F(xiàn)ANG Yong-fei.The synthetically analytical method for data sets on pyrotechnics reliability test［J］.Journal of Systems Science and Mathematical Science，2006，26(2):147 -158.(in Chinese)

［2］Dixon W J，Mood H M.A method for obtaining and analyzing sensitivity data［J］.Journal of the American Statistical Association，1948，43(241):109 -126.

［3］Jeff Wu.Efficient sequential designs with binary data［J］.Journal of the American Statistical Association，1985，80(392):974-984.

［4］Wetherill G B.On the existence of maximum likehood estimators for the binomial response models［J］.Journal of the Royal Statistical Society B，1981，43(3):310 -313.

［5］Efron B.Better bootstrap confidence interval［J］.Journal of the American Statistical Association，1987，82(397):171 -185.

［6］Chao M T，F(xiàn)uh C D.Bootstrap method for the up and down test on pyrotechnology sensitivity analysis［J］.Journal of the Royal Statistical Society B，2001，11(1):1 -21.