基于最大熵和仿真的導彈裝備備件需求預測

趙建忠，徐廷學，李海軍，徐衡博

(1．海軍航空工程學院 兵器科學與技術系，山東 煙臺，264001； 2．中國人民解放軍91423部隊 裝備部，山東 萊陽 265200)

0 引言

備件是導彈裝備使用和維修的重要物質基礎，其保障水平對于保持導彈裝備的戰備完好性有著直接的影響。準確地預測備件的消耗數目才能在有限的經費條件下，最大限度地滿足導彈裝備保障的備件需求。目前，大部分的備件預測方法或模型將備件視為在冷儲備條件，即不考慮備件在儲備期間發生故障的情況，而事實上導彈裝備經常會出現零部件或備件在儲備期間損壞而導致短缺的問題，所以它們的預測精度往往不高。因此，建立間斷工作導彈裝備備件需求預測模型，更準確地確定備件需求量成為亟待解決的問題。

在建立模型之前，首先需要確定導彈裝備的壽命分布。導彈裝備故障機理十分復雜，而且其備件消耗大都呈現“短周期、小樣本”的特點，故采用通過大量數據樣本進行建模的方法計算導彈備件的壽命分布參數，效果往往不太理想。實踐表明，對于常見統計分布函數，壽命分布假定有時在0.05的水平被拒絕[1]。為了解決這一問題，本文將導彈壽命分布形式看作未知，提出用最大熵法確定壽命分布，并利用遺傳算法求解系統建模過程中的非線性規劃問題。

文獻[2]利用最大熵法來確定統計分布函數參數。由于傳統最大熵模型中的約束條件是基于經典矩構建的，因此，基于經典矩的最大熵法能夠有效擬合大樣本數據的分布函數，然而導彈裝備故障的數據往往不能滿足大樣本條件[3]。本文結合最大熵原理與概率加權矩，提出了基于最大熵原理和概率加權矩的導彈裝備壽命分布確定方法。它直接從樣本信息出發，不需要對待估隨機變量的統計分布類型作任何假定，從而為建立統計分布函數提供了便利。

1 基于最大熵原理與概率加權矩的部件導彈單元壽命分布確定

為了克服基于不同樣本量估計的模型對于系統統計特性的描述偏差相差極大的缺陷,本文引入概率加權矩(probability weighted moments, PWMs)方法。根據最大熵原理[4]，對于定義域R上的連續隨機變量x(x的概率密度函數為f(x))，假設其累積分布函數F(x)=P(X≤x)及其逆函數x=x(F)存在，則其概率加權矩可以表示為

(1)

式中：變量i,n,k為實數。

目前常用的概率加權矩有[5]

(2)

(3)

并設αk與βn對應的無偏估計為ak和bn。

由式(1)可得

(4)

在上述分析的基礎上，可得最大熵逆累積分布函數模型

(5)

求解得

(6)

式中：N為樣本矩的最高階數；λi,i=0,1,…,N為與第i階矩約束相對應的拉格朗日算子。

如果式(5)中約束等于各階經典矩，式(6)概率密度函數的約束等于各階概率加權矩，則式(6)為逆累積分布函數。

在上述分析的基礎上，bn可表示為[7-8]

(7)

式(7)是關于λi(i=0,1,…,N)的N個方程。

為便于數值求解,將式(7)改寫為如下形式：

(8)

確定如下優化問題:

(9)

通過遺傳算法求解式(9),可得λi,i=0,1,…,N的值。具體的算法步驟為：①對樣本數據進行排序；②計算概率加權矩；③建立殘差優化模型，確定優化初始值；④應用遺傳算法對模型進行求解；⑤判斷是否收斂，若收斂則轉⑦，否則轉⑥；⑥重新設定初始值，轉④；⑦輸出計算結果；⑧計算結束。

2 基于保障度的間斷工作導彈裝備備件需求模型的建立

為建立導彈裝備維修備件需求量模型作如下假設：

(1) 某導彈備件供應保障系統由S+1個單元組成。當導彈裝備處于工作狀態時，認為其中1個單元工作，其余單元作為溫儲備。當工作單發生故障時，由其中一個沒有發生故障的溫儲備單元替換；如果導彈裝備處于停機狀態，所有未發生故障的單元均作為溫儲備。

(2) 初始狀態，所有零部件及其備件均為新品，系統中1個零部件在工作，S個備件處于溫儲備用于替換故障件。

(3) 在每個階段狀態中，系統部件和備件只有正常/故障2種狀態，且相同零部件或備件壽命分布相同，零部件在工作期間和備件在儲備期間壽命服從同一分布。

(4) 部件只有2種狀態：“工作”或“不工作”，即維修保障周期為1年(365天)，裝備每周工作tw(h)，其余時間ts=168-tw(h)停機。導彈裝備工作時，其零部件處于工作狀態，相應的備件處于儲備狀態。導彈裝備停機時，零部件的故障率與備件的故障率相同。

(5) 零部件及其備件均不可修，導彈裝備發生故障時進行換件維修，且修復后性能如初，換件時間忽略不計。

(6) 導彈裝備如果在工作期間發生故障立即實施維修，如果在停機期間發生故障等到再次開機時實施維修。

(7) 裝備系統各組成器件的失效相互獨立，其失效不會發生在同一時刻。

(8) 任一器件在發現其不能工作之前總是完好的，即不能工作時間從故障發現時開始。

(9) 裝備數量及其所屬的零部件數量一定。

可推導出

(10)

(11)

式中:t1為第1個零部件發生故障前最后一次開機保持良好狀態的時間;i表示天數。

則第1個零部件的壽命可表示為

T1=24n1+t1.

(12)

如果處于儲備狀態的備件沒有發生故障，處于工作狀態的第1個零部件故障后直接用其替換，則第2個備件的累計故障分布函數可表示為

(13)

壽命T2=24(n1+n2)+t2。以此類推，第N個備件的累計故障分布函數為

(14)

通過迭代計算，可以推導出單元的壽命為TN=TN-1-(tw-tN-1)+24n+tN。其中tN為該零部件最后一次開機后良好運行的時間。

這樣依次循環直到所以備件都用完，就認為導彈裝備不能工作，則最后1個沒有在儲備期間發生故障的零部件壽命就是備件的可保障時間。

3 仿真思路規劃

根據MC直接抽樣法[11]，可得到每個單元的累計故障概率抽樣值ηj，即[0,1]區間均勻分布的隨機數(j=1,2,…,S+1)。由于處于溫儲備狀態的備件也可能發生故障，所以當第1個零部件在工作期間發生故障需要進行換件維修時，首先要判斷第2個單元是否已經發生故障。如果沒有故障，則立即進行換件維修；否則按以上方法對第3個單元進行判斷。以此類推，直到所以備件都用完，就認為導彈裝備不能工作。當最后一個零部件的壽命不小于備件供應保障周期時，就認為備件的數量滿足需要，否則就是不滿足需要。根據Monte-Carlo思想，將上述過程重復Nmn次，記備件滿足需要的次數為Nsuc，當Nmn足夠大時，Nsuc/Nmn即為備件的保障概率。針對故障概率和平均壽命的計算，進行1 000次仿真運行即可，本文取Nmn=1 000次。

3.1 實例分析

根據調研分析，某型導彈裝備及其零部件平時主要處于停機和使用2種狀態，即工作和儲備2種模式，裝備平均每周使用大約12 h，其余時間處于停機儲備狀態，且裝備和備件的儲備條件一樣。裝備發生故障后采用換件維修的方式。對于該型導彈裝備某電子模塊工作和儲備期間的故障情況，這里統計分析了有連續故障記錄的20個樣本，以該模塊的故障間隔時間表示，單位為h，具體如表1所示。要求備件保障度為0.9。

表1 故障統計表

由于該電子模塊故障數據較少，適合采用基于概率加權矩的最大熵方法確定壽命分布函數。其中，所用矩的階數m取為2，具體的m選取參見文獻[7-8]。在構造遺傳算法時，由于目標函數是求最小值，同時為增大適應度函數間的差距，以提高優秀個體被選中的概率，因此選擇目標函數值倒數的1010倍作為遺傳算法的適應度函數，求得其參數為：λ0=6.25×10-3，λ1=4.49×10-3，λ2=2.63×10-3，其算法尋優過程如圖1所示。經擬合優度檢驗知，該電子模塊工作壽命服從參數為λ=6.25×10-3的指數分布。在此基礎上，運用隨機生成服從相同參數分布的50組數據樣本，分別基于概率密度函數和逆累積分布函數來計算1階矩，并用在定義域內的積分對于理論總體均值的相對誤差來評判不同方法的優劣，具體計算結果如圖2所示。圖2表明，當樣本量少于10時，采用經典矩方法相對誤差達到50.033 7%，而采用概率加權矩方法僅為3.969 4%。因此，當樣本≤30時，基于經典矩的最大熵方法不太適用；若樣本量>30，兩種方法均可行，但概率加權矩方法精確度相對高些。同理，可求得儲備壽命為服從μ=1.851 8×10-4的指數分布。

圖1 適應度函數圖Fig.1 Fitness function figure

圖2 不同樣本下2種方法計算相對誤差Fig.2 Bias of two different methods obtained from various sample

3.2 需求仿真

利用Matlab編程進行1 000次仿真，該模塊1年的備件需求量部分計算結果如表2所示。其中，該導彈裝備每周工作20 h時，1年所需備件需求量如圖3所示。

3.3 結果分析與評價

表2中，tw=0表示備件一直處于儲備狀態，可作為溫儲備看待，即所有單元均以μ=1.851 8×10-4h的故障率進行儲備，可利用并聯系統公式[12]進行計算，結果為S=1，P=0.928 6。tw=24可視為連續工作的熱儲備系統，利用熱儲備系統公式[12]可解得S=6，P=0.928 1。由此可見，2種方法的計算結果是一致的，這說明了本文仿真方法的正確性。

表2 備件需求量計算結果

圖3 備件需求量仿真圖Fig.3 Simulation figure of spare parts demand

4 結束語

針對導彈裝備故障數據樣本量較小的情況，提出了基于最大熵原理與概率加權矩的分析方法，能夠在確定概率分布的同時估計分布參數。該方法的優勢在于適用范圍較廣，而且與基于經典矩的最大熵方法相比，對于樣本量小于30的情形擬合效果更為理想。備件需求模型既分析了導彈裝備工作期間發生情況，也考慮了導彈裝備在停機期間及備件在儲備期間也會發生故障的情況，更加貼近實際。

參考文獻：

[1] PANDEY M D. Direct Estimation of Quantile Functions Using the Maximum Entropy Principle [J]. Structural Safety, 2000,22(1): 61-79.

[2] 俞禮軍,嚴海,嚴寶杰. 最大熵原理在交通流統計分布模型中的應用[J].交通運輸工程學報, 2001, 1(3): 91-94.

YU Li-jun, YAN Hai, YAN Bao-jie. Maximum Entropy Principle and Its Application in Probability Density Function of Traffic Flow [J].Journal of Traffic and Transportation Engineering, 2001, 1 (3): 91-94.

[3] HOSKING R M, WALLIS J R. Regional Frequency Analysis[M].New York: Cambridge University Press, 1997.

[4] TAN L, TANIAR D. Adaptive Estimated Maximum Entropy Distribution Model[J].Information Sciences, 2007,15(177):3110-3128.

[5] 俞禮軍,徐建閩. 出行時間價值最大熵分布估計模型[J].交通運輸工程學報, 2008, 8(1):83-96.

YU Li-jun, XU Jian-min. Estimation Models of Distribution Functions for Travel Time Value with Maximum Entropy Principle[J]. Journal of Traffic and Transportation Engineering, 2008, 8(1):83-96.

[6] 俞禮軍, BenHeydecker. 基于概率加權矩與最大熵原理的交通流統計分布估計方法[J].公路交通科技, 2008, 25(2):113-133.

YU Li-jun, BenHeydecker. Probability Density method of Traffic Flow Based on Probability Weighted Moments and Maximum Entropy Principle[J]. Science and Technology of Road Transportation, 2008, 25(2):113-133.

[7] 朱堅民，郭冰菁，王中宇，等.基于最大熵方法的測量結果估計及測量不確定度評定[J].超重運輸機械，2005，42(8)：5-8.

ZHU Jian-min,GUO Bing-jing,WANG Zhong-yu. Study on Evaluation of Measurement Result and Uncertainty based on Maximum Entropy Method[J].Journal of Electronic Measurement and Instrument,2005,42(8):5-8.

[8] YU X M,TONG L. Evaluation of Uncertainty in Measurement Based on the Maximum Entropy Method[J].Journal of Electronic Measurement and Instrument,2004,17(6):98-103.

[9] 葛陽, 寧劍平. 間斷工作導彈裝備維修備件需求量仿真模型[J].起重運輸機械, 2010,22(12):42-48.

GE Yang, NING Jian-ping. Demand Simulation Model of the Missile Equipment Maintenance Parts in Discontinuous Work[J].Lift Transportation Mechanism, 2010,22(12):42-48.

[10] 金星, 文明, 李俊美. 壽命服從指數分布產品相關失效解析分析[J].裝備指揮技術學院學報, 2002, 13(4):37-39.

JIN Xing, WEN Ming, LI Jun-mei. Correlative Failure Analysis on Exponential Distribution Manufacture[J].Journal of Academy of Equipment, 2002, 13(4):37-39.

[11] 程文鑫, 陳立強, 龔沈光. 基于蒙特卡洛法的艦船裝備戰備完好性仿真[J].兵工學報, 2006, 27(6):1090-1094.

CHENG Wen-xin, CHEN Li-qiang, GONG Shen-guang. Ship Equipment Operational Readiness Simulation Based on Monte Carlo Method[J].Acta Armamentarii, 2006, 27(6):1090-1094.

[12] 曹晉華, 程侃. 可靠性數學引論[M].北京:高等教育出版社, 2006.

CAO Jin-hua, CHENG Kan. Reliability Mathematics Generality Induction[M].Beijing: Higher Education Press, 2006.