基于BB算法的左截斷威布爾分布可靠性分析

孔祥芬 張凱奇 張 俊 張 飛

1.中國民航大學(xué)航空工程學(xué)院，天津，3003002.深圳航空有限責(zé)任公司維修工程部，深圳，518000

基于BB算法的左截斷威布爾分布可靠性分析

孔祥芬1張凱奇1張 俊2張 飛2

1.中國民航大學(xué)航空工程學(xué)院，天津，3003002.深圳航空有限責(zé)任公司維修工程部，深圳，518000

針對壽命分布類型和參數(shù)估計方法問題，對國內(nèi)某航空公司現(xiàn)役機型B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)進行可靠性分析。結(jié)合數(shù)據(jù)類型特點，初步預(yù)判其壽命分布類型可能為左截斷兩參數(shù)威布爾分布、對數(shù)正態(tài)分布、指數(shù)分布。基于極大似然估計法分別估計三種分布的未知參數(shù)，并利用根據(jù)BB算法改進的柯爾莫哥洛夫-斯米爾諾夫假設(shè)檢驗(K-S檢驗)法判斷其壽命分布類型。實證分析結(jié)果表明，采用左截斷兩參數(shù)威布爾分布來擬合是合適的，并提出了具體維修建議。

起落架；可靠性分析；極大似然估計；K-S檢驗；左截斷威布爾分布；BB算法

0 引言

飛機起落架減振支柱是飛機安全飛行的重要前提和基礎(chǔ)[1]，其工作的可靠性直接影響到飛機起飛和著陸性能的實現(xiàn)及飛行安全[2]。針對各飛機系統(tǒng)或部件的不同故障模式和故障后果，擬定不同的維修方式和維修制度，對提高飛行安全和控制維修成本具有重要意義[3]。目前，學(xué)者們對飛機部件的可靠性分析問題已進行了大量的研究。在失效數(shù)據(jù)分布類型研究方面，周英彪等[4]對比分析了兩參數(shù)威布爾分布與三參數(shù)威布爾分布的優(yōu)點和缺點，并應(yīng)用球磨機缺陷數(shù)據(jù)進行了實證分析。在失效數(shù)據(jù)精度分析上，辛龍等[5]以最小二乘法為威布爾分布參數(shù)估計方法，從經(jīng)驗分布函數(shù)角度對比了平均秩次法和中位秩法，并證明平均秩次法精度更高。伍建軍等[6]進一步結(jié)合平均秩次法和加權(quán)最小二乘法的優(yōu)點，對兩參數(shù)威布爾分布進行了參數(shù)估計，對比了水平殘差平方和垂直殘差平方，得到形狀參數(shù)取值范圍不同應(yīng)該采取不同的參數(shù)估計方法的結(jié)論。秦明等[7]發(fā)現(xiàn)，威布爾分布參數(shù)的最小二乘估計并不是最優(yōu)方法，并且指出極大似然估計法(maximum likelihood estimation,MLE)對于較復(fù)雜樣本情況仍然具有優(yōu)良的適用性和精確的估計。吳江[8]應(yīng)用粒子群算法(particle swarm optimization,PSO)對混合威布爾分布的參數(shù)估計值進行了數(shù)值逼近，對比證明了該方法收斂速度相對較快，不易陷入局部最優(yōu)解，但該方法并沒有考慮K-S法擬合檢驗中參考統(tǒng)計量D的選取問題，僅對參數(shù)值的數(shù)值逼近方法進行了優(yōu)化。

研究和實踐結(jié)果表明，機械、電氣、氣動、液壓本身或其所包含零件和子系統(tǒng)的故障模式適合用威布爾分布來擬合[9]。在參數(shù)估計的數(shù)值分析方法中，MLE方法不僅可以得到理想的參數(shù)估計結(jié)果，而且適用于完全樣本以外的截尾數(shù)據(jù)和截斷數(shù)據(jù)情形，甚至對于小樣本情形也有較好的適用性。國內(nèi)某航空公司自2010年9月以來，引進了大量B737-800飛機，其數(shù)量約占公司機隊規(guī)模的一半左右。然而該機型飛機主起落架減振支柱滲油故障時有發(fā)生。截止到2016年3月底，在統(tǒng)計到的35架B737-800飛機起落架系統(tǒng)故障記錄中共發(fā)生70多次主起落架減振支柱漏油故障。因而從安全角度以及經(jīng)濟效益影響方面綜合考量，對B737-800飛機主起落架減振支柱漏油故障的壽命分布及其失效原因加以分析和研究是非常必要和有現(xiàn)實意義的。

本文以該航空公司35架B737-800飛機主起落架油氣式減振支柱漏油故障數(shù)據(jù)(平均無故障間隔時間)為研究對象進行壽命分布對比研究，以極大似然估計法為核心進行參數(shù)估計，并對參數(shù)估計結(jié)果進行了基于BB算法改進的K-S擬合檢驗。基于所獲得的壽命分布結(jié)果，進一步開展主起落架減振支柱漏油故障數(shù)據(jù)可靠性分析，并提出了相應(yīng)的維護措施。

1 相關(guān)理論

1.1 左截斷兩參數(shù)威布爾分布

兩參數(shù)威布爾分布的概率分布函數(shù)[10]如下：

(1)

其中，τ為隨機變量；非負參數(shù)α為尺度參數(shù)；β為形狀參數(shù)。形狀參數(shù)β取值不同，概率密度曲線的形狀不同，恰好能擬合浴盆形狀的壽命分布曲線的不同時期。如果假定設(shè)備的壽命符合兩參數(shù)威布爾分布，則可以通過形狀參數(shù)β的不同值來建議維修方式[11]。

如果一個樣本容量為n的獨立同分布左截斷數(shù)據(jù)集合(τ1，τ2，…τn)具有τL<τi(i=1,2,…,n)的性質(zhì)，其中τL為非負參數(shù)，表示數(shù)據(jù)的左截斷點，則左截斷威布爾分布的概率分布函數(shù)和密度函數(shù)分別為

(2)

(3)

τ>τL

本文假設(shè)非負參數(shù)左截斷點τL是已知或可以設(shè)定的，則利用極大似然估計法可以得到確定分布參數(shù)α、β似然估計值的方程：

(4)

(5)

(6)

如果不滿足式(6)要求，則會求出無效解α=β=0。

1.2 左截斷兩參威布爾分布的Kolmogorov-Smirnov檢驗

單樣本的Kolmogorov-Smirnov檢驗(K-S檢驗)用來檢驗一組樣本數(shù)據(jù)的觀測經(jīng)驗分布是否與某一已知的理論分布相符合，當(dāng)兩者的差距很小時，推斷該樣本取自已知的理論分布。

(7)

根據(jù)Brownian Bridge和Donsker’s定理以及數(shù)值分析法(Brownian Bridge-Donsker’s， BB)[12-14]，可以得到臨界值DL的數(shù)值方法逼近值。限于篇幅，僅給出相應(yīng)理論推導(dǎo)結(jié)論(以下結(jié)論公式均假定顯著性水平取值為0.05)。

當(dāng)參數(shù)α、β已知時，臨界值DL可看作關(guān)于樣本容量n的函數(shù)：

(8)

(9)

η≡(τL/α)β

其中，η為截斷參數(shù)。當(dāng)把DL當(dāng)作關(guān)于截斷參數(shù)η的函數(shù)時，有

(10)

其中，A(n)、B(n)、C(n)、D(n)、E(n)均為BB算法給出的式(10)中涉及的具體參數(shù)值。結(jié)合式(8)和式(10)，可以得到修正的臨界值DL函數(shù)形式：

(11)

表1給出了顯著水平為0.05時參數(shù)A、B、C、D、E、F的取值。

表1 威布爾分布BB算法K-S檢驗擬合參數(shù)對照表

2 主起落架減振支柱漏油故障數(shù)據(jù)壽命分布

2.1 故障數(shù)據(jù)的采集與處理

本文采集到的原始故障數(shù)據(jù)來自國內(nèi)某航空公司35架B737-800飛機起落架系統(tǒng)的電子版故障記錄，記錄時間從2010年9月14日到2016年3月29日。本文考慮的B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)主要包括以下兩種類型：①在故障描述中說明主起落架減振支柱漏油嚴重，并在排故措施中采取更換封嚴等措施的故障；②在故障描述中體現(xiàn)主起落架減振支柱發(fā)生漏油現(xiàn)象，并因后續(xù)漏油決定采取持續(xù)監(jiān)控或鏡面加熱措施，且在排故措施中不符合勤務(wù)范圍而采取更換封嚴等措施的故障。

利用單機同一主起落架減振支柱相鄰兩次漏油故障的飛行小時數(shù)(NF)記錄，得到其可靠性計算指標故障間隔時間，并按照故障間隔時間從小到大升序排列，具體統(tǒng)計見表2。

對表2的故障間隔時間以NF=800為間隔進行分組，繪制頻數(shù)分布直方圖(圖1)。

頻數(shù)分布直方圖可以直觀地體現(xiàn)數(shù)據(jù)的趨中性和分散性，可預(yù)測數(shù)據(jù)的分布類型。觀察圖1的數(shù)據(jù)分布趨勢，初步預(yù)測B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)分布符合指數(shù)分布、對數(shù)正態(tài)分布和左截斷兩參數(shù)威布爾分布，且由數(shù)據(jù)走勢預(yù)判左截斷兩參數(shù)威布爾分布的形狀參數(shù)β>1。

表2 B737-800飛機主起落架減振支柱漏油故障間隔時間統(tǒng)計表

圖1 B737-800飛機主起落減振支柱漏油故障頻數(shù)分布直方圖Fig.1 Frequency distribution histogram of main landing gear shock absorber pillar oil spill for B737-800

2.2 壽命分布類型的參數(shù)估計與檢驗

2.2.1 參數(shù)估計與檢驗步驟

基于極大似然估計法，利用MATLAB軟件編程計算估計分布函數(shù)對應(yīng)未知參數(shù)。具體參數(shù)為：指數(shù)分布的參數(shù)λ、對數(shù)正態(tài)分布的對數(shù)均值μ和對數(shù)標準差σ，左截斷兩參數(shù)威布爾分布的尺度參數(shù)α和形狀參數(shù)β。運用K-S檢驗法對各分布函數(shù)進行擬合優(yōu)度檢驗，比較各個統(tǒng)計觀察值Dn，若Dn相差較大，則Dn較小的分布擬合效果比較好，更符合實際情況；若Dn相差不大，則進一步利用Minitab軟件進行Personχ2擬合優(yōu)度交叉檢驗。計算出相關(guān)系數(shù)ρxy并繪制相應(yīng)的概率圖，最終驗證數(shù)據(jù)的分布類型。

2.2.2 參數(shù)估計與檢驗實例過程

首先對提供的B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)進行過濾，剔除無用數(shù)據(jù)。由于本文主要側(cè)重對左截斷類型數(shù)據(jù)的研究，故將截尾數(shù)據(jù)等同于無用數(shù)據(jù)處理。得到截斷水平p≈0.3，并根據(jù)式(9)計算得到相應(yīng)的截斷參數(shù)η≈0.35。則取定n=19、η=0.35的情況下，結(jié)合表1參數(shù)值代入式(11)，得到DL(0.35,19)≈0.9525。

運用MATLAB編程，計算得到B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)威布爾分布參數(shù)的估計值α=3445.7，β=3.0，再結(jié)合式(9)計算得到截斷點τL≈2428.29。則B737-800飛機主起落架減振支柱漏油故障的壽命分布函數(shù)

(12)

由表2數(shù)據(jù)和式(7)，計算得到Dn=0.1681。顯然Dn≈0.2185，所以接受原假設(shè)，即該故障數(shù)據(jù)分布符合左截斷兩參數(shù)威布爾分布。

2.2.3 參數(shù)估計與檢驗結(jié)果對比分析

同理，利用MATLAB編程分別計算另外兩種分布的參數(shù)，然后按K-S檢驗步驟進行擬合檢驗。三種分布的對比結(jié)果列于表3。

表3 B737-800飛機主起落架減振支柱漏油故障壽命分布參數(shù)估計及K-S檢驗

由表3可見，B737-800飛機主起落架減振支柱漏油故障間隔時間符合左截斷兩參數(shù)威布爾分布和對數(shù)正態(tài)分布，并且擬合效果相近，兩參數(shù)威布爾分布擬合效果相對更好。進一步利用Minitab軟件對對數(shù)正態(tài)分布和左截斷兩參數(shù)威布爾分布進行Personχ2擬合優(yōu)度檢驗，擬合結(jié)果如圖2所示。

(a)對數(shù)正態(tài)分布 (b)威布爾分布 (95%置信區(qū)間) (95%置信區(qū)間)圖2 對數(shù)正態(tài)分布、威布爾分布擬合優(yōu)度檢驗Fig.2 Goodness of fit test of lognormal distribution and Weibull distribution

圖2表明，兩種檢驗對應(yīng)的P值均大于顯著性0.05，表明兩種分布假設(shè)均成立，但對B737-800飛機主起落架減振支柱漏油故障數(shù)據(jù)擬合效果較好的為左截斷兩參數(shù)威布爾分布。此時運用MATLAB程序求得其相關(guān)系數(shù)ρxy=0.9654，棄真概率小于0.05，有理由相信用左截斷兩參數(shù)威布爾分布擬合B737-800飛機主起落架減振支柱漏油故障的壽命分布更好。運用MATLAB程序繪制其威布爾概率圖，如圖3所示。

圖3 威布爾概率圖 Fig.3 Weibull Probability Chart

圖3中，數(shù)據(jù)點(xi,yi)是均勻且趨向于一條直線分布的，這也充分說明B737-800飛機主起落架減振支柱漏油故障的壽命分布用左截斷兩參數(shù)威布爾分布來擬合是準確的。

3 主起落架減振支柱漏油故障數(shù)據(jù)可靠性分析

B737-800飛機主起落架減振支柱漏油故障較好壽命分布為左截斷兩參數(shù)威布爾分布，其中，β=3.0，1<β<4，其失效率隨時間的延長而增大，故障類型為早期損耗，可以通過優(yōu)化相應(yīng)的預(yù)定維修計劃來提高可靠性。具體可以采用硬時限的控制方式，即確定一個由合理的可靠度計算出的使用時限，當(dāng)使用時間達到該時限時進行維護或翻修。

應(yīng)用MATLAB軟件分別繪出B737-800飛機主起落架減振支柱漏油故障左截斷兩參數(shù)威布爾分布(α=3445.7，β=3)的失效率函數(shù)曲線和可靠度函數(shù)曲線，如圖4、圖5所示。

圖4 失效率函數(shù)曲線Fig.4 Failure rate function curve

圖5 可靠度函數(shù)曲線Fig.5 Reliability function curve

由失效率函數(shù)曲線可知，B737-800飛機主起落架減振支柱漏油故障大約在NF=1500時失效率明顯上升，同時由可靠度函數(shù)曲線可知其可靠度在NF=1500左右急劇下降。假設(shè)該減振支柱已工作NF=1500，則當(dāng)前的失效率λ(1500)=1.65×10-4/Fh，可靠度R(1500)=0.92，這表明B737-800飛機主起落架減振支柱在工作時間達NF=1500后，單位時間內(nèi)的失效率為1.65×10-4，也表明每10個這樣的減振支柱在工作時間達NF=1000后就有約1～2個產(chǎn)生漏油故障。假設(shè)可靠性維修大綱要求B737-800飛機主起落架減振支柱可靠度為92%才能正常工作，則該作動器工作時間大約在NF=1500左右就必須進行檢修，這顯然符合工程實際中起落架部件的可靠度往往大于90%的要求。

經(jīng)計算得到B737-800飛機主起落架減振支柱在NF=500時R(500)=0.9969；在NF=800時R(800)=0.9876；在NF=1000時R(1000)=0.9759；在NF=1200時R(1200)=0.9586。假設(shè)工程實際中要求起落架部件的可靠度為99%以上，則建議在NF=500以后就應(yīng)密切關(guān)注其實際運行狀態(tài)。此外，狀態(tài)檢修方面，如點檢管理等，需要給出適當(dāng)?shù)臓顟B(tài)檢測時間間隔，即需要計算出平均無故障工作時間tMBF。通常情況下，狀態(tài)檢測周期Tm=0.10tMBF，此處給出根據(jù)威布爾分布兩參數(shù)計算得到的tMBF計算公式：

tMBF=ηΓ(1+1/β)

(13)

代入數(shù)據(jù)α=3445.7，β=3，得tMBF=3077.01，則狀態(tài)檢測周期Tm=307.70。因此，建議每隔NF=308左右，對主起落架減振支柱進行漏油方面的狀態(tài)檢測。

4 結(jié)語

工程實際中要求起落架部件的可靠度為99%以上，建議在NF=500以后就應(yīng)密切關(guān)注其實際運行狀態(tài)。在點檢方面，建議每隔NF=308左右，對主起落架減振支柱進行漏油方面的狀態(tài)檢測。同時航空公司的航材部門可以根據(jù)本文結(jié)論預(yù)測需要儲備航材的種類和數(shù)量，為進行合理備件提供決策依據(jù)，進而節(jié)約公司維修成本，提高飛機的使用效率。研究結(jié)果可為飛機其他系統(tǒng)、子系統(tǒng)及其部件的壽命分布研究和可靠性分析提供參考。

[1] 王紅，凌云，何超．對172R型飛機前起落架減振支柱滲油故障診斷及解決方法的研究[J].機床與液壓，2016，44(4):142-144． WANG Hong, LING Yun, HE Chao. The Study of the Leakage Fault Diagnosis and Solution of 172R Front Landing Gear Shock Strut Oil[J]. Machine Tool & Hydraulics,2016,44(4):142-144．

[2] LU Lixian, LEI Jianhe. Design and Reliability Prediction of a Distributed Landing Gear Control System [J].Aircraft Engineering and Aerospace Technology, 2010,82(1):15-22.

[3] 孫見忠,左洪福,梁坤. 基于民航發(fā)動機狀態(tài)數(shù)據(jù)的渦輪葉片剩余壽命評估[J]. 機械工程學(xué)報,2015,51(23):53-59. SUN Jianzhong, ZUO Hongfu, LIANG Kun. Remaining Useful Life Estimation Method for the Turbine Blade of a Civil Aircraft Engine Based on the QAR and Field Failure Data[J]. Journal of Mechanical Engineering, 2015,51(23):53-59.

[4] 周英彪,段權(quán)鵬,范杜平. 利用威布爾分布模型對球磨機可靠性分析[J]. 熱能動力工程,2009,24(1):95-99. ZHOU Yingbiao, DUAN Quanpeng, FAN Duping. The Reliability Analysis of Ball Mill by Weibull Distribution Model[J].Journal of Engineering for Thermal Energy and Power, 2009,24(1):95-99.

[5] 辛龍,周越文,翟穎燁,等. 基于Weibull分布的航空裝備部件壽命預(yù)測研究[J].電光與控制,2014,21(12):102-105. XIN Long, ZHOU Yuewen, ZHAI Yingye, et al. Weibull Distribution Based Lifetime Prediction for Components of Aerial Equipments[J]. Electronics Optics & Control, 2014,21(12):102-105.

[6] 伍建軍,游雄雄,吳事浪,等.典型礦冶廢舊零部件剩余壽命預(yù)測模型與可再制造性評估決策方法[J]. 機械科學(xué)與技術(shù),2014,33(12):1859-1864. WU Jianjun, YOU Xiongxiong, WU Shilang, et al. Remanufacturing Assessment of Typical Mining Waste Parts Based on Remaining Life Prediction Model[J]. Mechanical Science and Technology for Aerospace Engineering, 2014,33(12):1859-1864.

[7] 秦明,巫世晶,彭瀟, 等. 一種服從威布爾分布裝備的可靠性評估方法[J]. 武漢大學(xué)學(xué)報(工學(xué)版),2008,41(6):100-102. QIN Ming, WU Shijing, PENG Xiao, et al. A Reliability Evaluation Method for Equipment with a Weibull Failure Distribution[J]. Engineering Journal of Wuhan University,2008,41(6):100-102.

[8] 吳江.飛機可靠性分析中的混合威布爾分布參數(shù)估計方法[J]. 工程設(shè)計學(xué)報,2015,22(1):26-29. WU Jiang. Parameter Estimation Method of Mixed Weibull Distribution in Reliability Analysis of Aircraft[J]. Chinese Journal of Engineering Design, 2015,22(1):26-29.

[9] LI S M, DONG C L, HUANG Y X.Evaluate Reliableness of Aero-engine Turbine Blades Based on Weibull-distributing [J]. Journal of Civil Aviation University of China, 2008,26(4):14-17.

[10] BLOCKJ B. A Study of Aviation Field Data Fitting to a Two Parameter Weibull Distribution [C]//Risk, Reliability and Societal Safety: Proceedings of the European Safety and Reliability Conference 2007. Stavanger,2007:1339-1445.

[11] 徐微,胡偉明,孫鵬. 基于兩參數(shù)威布爾分布的設(shè)備可靠性預(yù)測研究[J]. 中國工程機械學(xué)報,2013,11(2):112-116. XU Wei, HU Weiming, SUN Peng. Equipment Reliability Prediction Based on Dual-parametric Weibull Distribution [J].Chinese Journal of Construction Machinery, 2013,11(2):112-116.

[12] PANC C, CHU L. Reliability Assessment for One-shot Product with Weibull Lifetime Components[J]. International Journal of Quality & Reliability Management,2010,27(5):596-610.

[13] KIZILERSU A, KREER M, THOMAS A W. Goodness-of-fit Testing for Left-truncated Two-parameter Weibull Distributions with Known Truncation Point[J]. Austrian Journal of Statistics,2016,45(3):15-42.

[14] O’Connor P T.Practical Reliability Engineering [M]. 4th ed. London:Heyden & Son Ltd.,2002.

(編輯 陳 勇)

Reliability Analysis of Left-truncated Weibull Distribution Modified by BB Algorithm

KONG Xiangfen1ZHANG Kaiqi1ZHANG Jun2ZHANG Fei2

1.College of Aeronautical Engineering, Cival Aviation University of China, Tianjin，300300 2.Maintenance & Engineering Department of Shenzhen Airlines, Shenzhen, Guangdong，518000

A reliability analysis was done based on the failure data of B737-800 main landing gear shock strut oil leakages herein, which came from a domestic airlines’ electronic records, in consideration of the life distribution type and parameter estimation method. The characteristics of data preliminarily demonstrated that the life distribution type might be left-truncated two parameter Weibull distribution, logarithmic normal distribution or exponential distribution. Next, the parameters corresponding to three types of distributions were estimated by using the MLE method. Then, the Kolmogorov-Smirnov(K-S) hypothesis test, which modified by BB algorithm and goodness-of-fit test , was used to identify and confirm the life distribution function model. The analysis results show that the distribution of the data is left-truncated two parameter Weibull distribution. In the end, reasonable maintenance plan and optimized maintenance scheme are provided.

landing gear；reliability analysis；maximum likelihood estimation(MLE)；Kolmogorov-Smirnov (K-S)hypothesis testing；left-truncated Weibull distribution；BB(Brownian Bridge-Donsker’s)algorithm

2016-10-11

航空科學(xué)基金資助項目(20130267001)

V226.3

10.3969/j.issn.1004-132X.2017.15.012

孔祥芬，女，1974年生。中國民航大學(xué)航空工程學(xué)院副教授。主要研究方向為工業(yè)工程、質(zhì)量管理與服務(wù)質(zhì)量。E-mail:759185381@qq.com。張凱奇，男，1989年生。中國民航大學(xué)航空工程學(xué)院碩士研究生。張 俊，男，1993年生。深圳航空有限責(zé)任公司維修工程部航線維護工程師。張 飛，男，1989年生。深圳航空有限責(zé)任公司維修工程部航線維護助理工程師。