王普長,孫 輝,陳晉市,張淼淼,何春暉
(1.江蘇徐工工程機(jī)械研究院有限公司,江蘇 徐州 221000; 2.吉林大學(xué) 機(jī)械與航空航天工程學(xué)院,吉林 長春 130025)
近年來,隨著工程機(jī)械行業(yè)的發(fā)展,我國逐漸推進(jìn)國產(chǎn)零件從高產(chǎn)到高質(zhì)的轉(zhuǎn)變;而國產(chǎn)工程機(jī)械可靠性差,適用性不強(qiáng)的問題亟需解決。因此,我國大力支持可靠性領(lǐng)域的研究,推動可靠性理論研究的發(fā)展。載荷譜是整機(jī)零件在典型工況下所受載荷的時(shí)間歷程,經(jīng)數(shù)學(xué)方法處理后,得到其典型載荷的概率頻次特征[1-2]。編制載荷譜可以通過極少的試驗(yàn)樣本獲得大量的樣本數(shù)據(jù)[3],節(jié)約了時(shí)間和人力成本,為后續(xù)的臺架加載提供原始樣本,是可靠性研究的重要環(huán)節(jié)[4-5]。
載荷譜最初是針對航空航天領(lǐng)域提出的,經(jīng)過數(shù)十年的發(fā)展,其研究范圍已超出航空航天領(lǐng)域[6],早已應(yīng)用至汽車、工程機(jī)械、機(jī)床等領(lǐng)域[7]。載荷譜外推分為雨流域外推與時(shí)域外推,時(shí)域外推中主要方法是極值外推法;極值外推法的主要思想就是將對零件損傷產(chǎn)生主要貢獻(xiàn)的載荷進(jìn)行預(yù)測。一般認(rèn)為,超越量服從廣義帕累托(GPD)分布。極值外推法的目的,就是尋找與GPD分布擬合最接近的閾值。因此,準(zhǔn)確的參數(shù)估計(jì)方法,將數(shù)據(jù)對廣義帕累托分布進(jìn)行擬合,是極值外推重要環(huán)節(jié)之一。
針對擬合結(jié)果,國內(nèi)外許多學(xué)者展開了大量研究:沙雪云等[8]分別將貝葉斯與極大似然估計(jì)法進(jìn)行擬合,并對比擬合結(jié)果,得到擬合更好的方法;李景奎等[9]提出一種改進(jìn)極大似然估計(jì)法,對3參數(shù)威布爾分布進(jìn)行擬合,并將該理論應(yīng)用到民航維修,為民航企業(yè)合理制定維修計(jì)劃提供了理論依據(jù);DANESHI N等[10]同樣對極大似然估計(jì)進(jìn)行改進(jìn),并利用到實(shí)際應(yīng)用中,證明改進(jìn)模型的可行性。
在可靠性方面的研究中,國內(nèi)外許多學(xué)者對試驗(yàn)數(shù)據(jù)的處理進(jìn)行探索:王蓉華等[11]提出逆L-矩估計(jì),應(yīng)用到疲勞壽命分析中,并與其他常用參數(shù)估計(jì)法進(jìn)行對比,驗(yàn)證其方法的準(zhǔn)確性;MUSSIE T等[12]提出一種新的壽命分布模型,利用矩估計(jì)與極大似然估計(jì)法進(jìn)行擬合,最后用試驗(yàn)數(shù)據(jù)進(jìn)行驗(yàn)證;崔策等[13]采用3參數(shù)Weibull分布對大型風(fēng)機(jī)液力變槳系統(tǒng)失效率估算結(jié)果,表明3參數(shù)Weibull分布適應(yīng)性更高;YIN F L等[14]基于極大似然估計(jì)法估計(jì)3參數(shù)Weibull分布,提出一種新的參數(shù)估計(jì)模型,實(shí)現(xiàn)更為精確的參數(shù)估計(jì),并應(yīng)用到軸承壽命分析中。
基于以上研究,常認(rèn)為超閾值樣本的分布服從GPD分布。本研究詳述了3種最常用的用來擬合廣義帕累托分布的參數(shù)估計(jì)方法,最大似然估計(jì)法、矩估計(jì)法及概率權(quán)矩估計(jì)法,其過程如圖1所示。總結(jié)了不同參數(shù)估計(jì)方法的過程及注意要點(diǎn),對比不同參數(shù)估計(jì)方法的擬合曲線及相關(guān)系數(shù),并進(jìn)行對比找出最優(yōu)方法。
圖1 參數(shù)估計(jì)過程Fig.1 Parameter estimation process
時(shí)域極值外推法的思想關(guān)鍵在于研究超越量的函數(shù)分布式,并由函數(shù)式隨機(jī)生成新的樣本點(diǎn)。
根據(jù)李昕雪等[15]研究表明,通常對于閾值較大的樣本,超越量服從分布GPD。
GPD累計(jì)函數(shù)表達(dá)式為:
(1)
概率密度表達(dá)式為:
(2)
式中,z=xi-μ(i=1,2,3,…,n)為樣本超越量;μ為閾值;σ為尺度參數(shù);ξ為形狀參數(shù)。
設(shè)x1,x2,…,xn為一組符合GPD分布的超閾值樣本,根據(jù)GPD分布函數(shù)得到相應(yīng)的對數(shù)似然函數(shù)為[16]:
ξ≠0
(3)
式(3)分別對σ和ξ求偏導(dǎo)得:
(4)
令偏導(dǎo)數(shù)等于0,并聯(lián)立方程組,解得:
(5)
通過式(5)可以解得參數(shù)σ和ξ。
MoM是由PEARSON提出的一種參數(shù)估計(jì)方法,根據(jù)辛欽大數(shù)定律可知,簡單隨機(jī)樣本的原點(diǎn)矩依概率收斂于相應(yīng)的總體原點(diǎn)矩,用樣本矩代替總體矩,進(jìn)而計(jì)算出未知參數(shù)的估計(jì)值。
設(shè)x1,x2,…,xn為一組符合GPD分布的超閾值樣本,由式(2)得:
(6)
(7)
式中,E(X)和E(X2)分別表示超閾值樣本的期望和二階矩,σ和ξ分別表示尺度參數(shù)和形狀參數(shù)。
根據(jù)式(6)及式(7)可得:
(8)
式中,D(X)表示超閾值樣本的方差。
用樣本的各階矩替換總體的各階矩,可得到σ和ξ的矩估計(jì),即:
(9)
式中,E(X)和D(X)分別為樣本的均值和方差。
矩法的一般原則是讓所研究的總體分布的各階矩與對應(yīng)的樣本矩相等,但是樣本的二階及高階樣本矩抽樣性質(zhì)不好,因此引入一種新的矩估計(jì)方法,即概率權(quán)矩估計(jì),設(shè)x1,x2,…,xn為一組符合GPD分布的超閾值樣本,其公式如式(10):
ωn(θ)=E[XFn(X;θ)],n∈N0
(10)
將GPD函數(shù)代入到式(10)中,得到GPD分布的n階概率權(quán)矩估計(jì)公式,如式(11):
(11)
為估計(jì)GPD分布中參數(shù)σ和ξ的值,分別令n為0和1,得到相應(yīng)的零階和一階概率權(quán)矩,并用樣本概率權(quán)矩代替總體概率權(quán)矩得到參數(shù)估計(jì):
(12)
式中,σ和ξ分別表示尺度參數(shù)和形狀參數(shù)。
柱塞泵是影響挖掘機(jī)挖掘性能的重要部件[17]。本實(shí)驗(yàn)以某公司20 t級反鏟液壓挖掘機(jī)為研究對象,測試場地為原生土場地且原生場地土壤不低于III級。
為方便駕駛?cè)藛T操作及試驗(yàn)人員采集,本研究利用無線信號傳輸裝置遠(yuǎn)程采集試驗(yàn)數(shù)據(jù),原理如圖2所示。
圖2 壓力傳感器監(jiān)測原理圖Fig.2 Schematic diagram of pressure sensor monitoring
本試驗(yàn)主要測量挖掘機(jī)挖掘期間,主泵出口壓力載荷的時(shí)間歷程。為了方便后續(xù)數(shù)據(jù)處理及工況分割,還同時(shí)監(jiān)測其他參數(shù),如在油缸處安裝壓力傳感器和線位移傳感器,在發(fā)動機(jī)處安裝轉(zhuǎn)速傳感器,監(jiān)測油箱油位等。本試驗(yàn)采用松散土砂作為挖掘物料展開試驗(yàn),挖掘現(xiàn)場如圖3所示。
圖3 物料挖掘現(xiàn)場Fig.3 Material excavation site
挖掘機(jī)工作環(huán)境惡劣,在進(jìn)行數(shù)據(jù)采集時(shí)不可避免的會存在奇異點(diǎn),因此本研究利用小波分析對采集的載荷譜進(jìn)行數(shù)據(jù)降噪及去除奇異點(diǎn)處理,處理后的載荷時(shí)間歷程如圖4所示[18]。
圖4 挖掘機(jī)主泵出口壓力Fig.4 Excavator main pump outlet pressure
為了提高計(jì)算效率,將預(yù)處理后的數(shù)據(jù)按升序進(jìn)行排列,選取該組數(shù)據(jù)中間位置的數(shù)值Z1和98%分位點(diǎn)處的數(shù)值Zm構(gòu)成n個(gè)備選閾值[19]。利用上述3種方法,計(jì)算出各備選閾值所對應(yīng)的GPD分布擬合參數(shù),并以每種參數(shù)估計(jì)方法中各備選閾值所對應(yīng)的GPD分布與超閾值樣本之間的最小均方誤差所對應(yīng)的值作為最優(yōu)閾值[20]。再通過對比各最優(yōu)閾值的GPD分布和超閾值分布的CDF圖,找到其中最優(yōu)的參數(shù)估計(jì)方法。
根據(jù)測試挖掘機(jī)主泵載荷譜數(shù)據(jù),計(jì)算出Z1=203,Z1420=344.9,分別用上述3種參數(shù)估計(jì)方法計(jì)算出1420個(gè)閾值所對應(yīng)的最小均方根誤差,得到各備選閾值和均方根誤差之間的關(guān)系如圖5所示。
根據(jù)圖5所示各備選閾值u和均方根誤差關(guān)系,可得最小均方根誤差所對應(yīng)的閾值,并選取載荷譜中最優(yōu)閾值以上的數(shù)據(jù),通過3種參數(shù)估計(jì)方法進(jìn)行GPD分布擬合,得到3種估計(jì)方法所對應(yīng)的最優(yōu)閾值及GPD分布的參數(shù)估計(jì)結(jié)果,如表1所示。
圖5 備選閾值和均方根誤差關(guān)系Fig.5 Alternative threshold and root mean square error relationship
表1 GPD分布的參數(shù)估計(jì)結(jié)果Tab.1 Parameter estimation results of GPD distribution
為更加直觀地反映每種GPD分布參數(shù)估計(jì)方法的優(yōu)劣,分別繪制ML、MoM及PWM 3種參數(shù)估計(jì)方法擬合的GPD分布函數(shù)的累積分布頻數(shù)圖,如圖6所示。
由表1可知,ML參數(shù)估計(jì)方法所選的最優(yōu)閾值較大,導(dǎo)致較多有效載荷數(shù)據(jù)丟失,MoM參數(shù)估計(jì)方法計(jì)算的最優(yōu)閾值和PWM參數(shù)估計(jì)方法的最優(yōu)閾值相同,但PWM參數(shù)估計(jì)方法的最小均方根誤差略小于MoM參數(shù)估計(jì)的結(jié)果。累積分布函數(shù)P,根據(jù)圖6可知,ML、MoM及PWM 3種參數(shù)估計(jì)方法所得到的樣本值和擬合值之間的相關(guān)系數(shù)分別為0.9966,0.9971及0.9982,均大于0.99,因此3種參數(shù)估計(jì)方法均可以較好地?cái)M合挖掘機(jī)主泵壓力超閾值樣本的GPD分布。相比于PWM參數(shù)估計(jì)方法,當(dāng)超閾值樣本數(shù)據(jù)小于1 MPa時(shí),ML參數(shù)估計(jì)方法得到的CDF擬合圖與超閾值樣本載荷之間略大,當(dāng)超閾值樣本數(shù)據(jù)介于1~1.8 MPa之間時(shí),MoM參數(shù)估計(jì)得到的CDF擬合圖與超閾值樣本載荷的CDF圖之間差異略大。從圖7可知,理論值u1、觀測值u2,當(dāng)超閾值樣本數(shù)據(jù)大于8 MPa,ML參數(shù)估計(jì)的Q-Q圖數(shù)據(jù)點(diǎn)偏離擬合直線,而MoM和PWM參數(shù)估計(jì)方法數(shù)據(jù)點(diǎn)均勻地分布在直線兩側(cè),超閾值樣本數(shù)據(jù)介于2~2.2 MPa之間時(shí),PWM參數(shù)估計(jì)的Q-Q圖偏離程度小于MoM參數(shù)估計(jì)。因此根據(jù)圖6和表1可得,PWM可以作為挖掘機(jī)主泵超閾值樣本GPD參數(shù)估計(jì)的最優(yōu)選擇。
圖6 最優(yōu)閾值擬合CDF圖Fig.6 Optimal threshold fitting CDF graph
圖7 最優(yōu)閾值擬合Q-Q圖Fig.7 Optimal threshold fitting Q-Q diagram
本研究以最小均方根誤差為檢驗(yàn)準(zhǔn)則,首先分別利用3種參數(shù)估計(jì)方法分別計(jì)算挖掘機(jī)主泵載荷譜POT模型的最優(yōu)閾值和超閾值樣本數(shù)據(jù)的GPD擬合參數(shù),其次分別繪制超閾值樣本數(shù)據(jù)及GPD分布的CDF圖,得出以下結(jié)論:
(1) 極大似然估計(jì)、矩估計(jì)及概率加權(quán)估計(jì)所得到的參數(shù)均可以較好地?cái)M合超閾值樣本數(shù)據(jù);
(2) 通過極大似然估計(jì)得到的最優(yōu)閾值相對于矩估計(jì)和概率加權(quán)估計(jì)略大,會造成較多有效載荷數(shù)據(jù)的丟失;
(3) 概率加權(quán)估計(jì)得到的最優(yōu)閾值可以保留較多的有效載荷信息,且通過概率加權(quán)估計(jì)方法得到的GPD分布與超閾值載荷樣本之間的擬合更好。