二元共載系統(tǒng)退化相關(guān)的兩階段退化模型及可靠性分析

張帆，程伯晗，王鵬，*，董磊

1.中國(guó)民航大學(xué) 民航航空器適航審定技術(shù)重點(diǎn)實(shí)驗(yàn)室，天津 300300

2.中國(guó)民航大學(xué) 科技創(chuàng)新研究院，天津 300300

3.中國(guó)民航大學(xué) 安全科學(xué)與工程學(xué)院，天津 300300

隨著航空系統(tǒng)中裝備的可靠性不斷提高，國(guó)產(chǎn)機(jī)載系統(tǒng)呈現(xiàn)高可靠、長(zhǎng)壽命等特征，且缺乏失效數(shù)據(jù)，因此基于退化的建模方法被認(rèn)為是評(píng)估系統(tǒng)可靠性有效的方法［1-5］。同時(shí)為保證飛行安全，航空系統(tǒng)多具有余度設(shè)計(jì)，系統(tǒng)中存在多組件共同支持系統(tǒng)運(yùn)行。傳統(tǒng)對(duì)于多組件系統(tǒng)的研究基于組件的獨(dú)立性，然而由于組件的結(jié)構(gòu)相關(guān)性、運(yùn)行條件一致性等原因，其退化性能可能會(huì)存在隨機(jī)依賴關(guān)系，這種隨機(jī)依賴關(guān)系主要存在于兩方面［6］：①退化進(jìn)程中的相互影響關(guān)系。組件當(dāng)前的退化量可能會(huì)影響其他組件的退化進(jìn)程；② 組件失效的影響。最典型的為共載系統(tǒng)中組件負(fù)載分擔(dān)關(guān)系。目前已存在大量文獻(xiàn)分別對(duì)這兩種隨機(jī)依賴關(guān)系進(jìn)行研究。

在對(duì)多元相關(guān)退化的建模中，Copula 函數(shù)是一類(lèi)運(yùn)用廣泛的方法。文獻(xiàn)［7］以Wiener 過(guò)程與Gamma 過(guò)程描述退化失效過(guò)程，選擇合適的Copula 函數(shù)對(duì)多元退化失效相關(guān)性進(jìn)行退化建模；文獻(xiàn)［8］由漂移布朗運(yùn)動(dòng)得到目標(biāo)退化模型，利用基于Vine-Copula 的多變量損傷耦合建模方法，建立了失效行為建模的相關(guān)描述模型；文獻(xiàn)［9］改進(jìn)了Vine-Copula 函數(shù)來(lái)描述性能參數(shù)之間的相關(guān)性；文獻(xiàn)［10］在Wiener 過(guò)程與Gamma 過(guò)程的基礎(chǔ)上，將退化過(guò)程分為線性與非線性，加入了外部因素，后用Copula 函數(shù)表述性能參數(shù)之間的相關(guān)性。然而如何選擇合適的Copula 函數(shù)以及如何驗(yàn)證其準(zhǔn)確性仍是待解決的難題。

將退化增量分為自身退化增量部分和相互作用影響增量部分為另一種常用的方法［11-13］，文獻(xiàn)［12-13］將相互作用影響增量部分的相關(guān)性運(yùn)用了一個(gè)固定冪律公式表述。文獻(xiàn)［14］利用因子分析法來(lái)描述相互作用影響增量部分的隨機(jī)依賴關(guān)系。此外，還有部分學(xué)者通過(guò)退化速率對(duì)退化相互作用的影響進(jìn)行建模［15-18］，由相互作用系數(shù)表述影響程度，可以得出退化模型的數(shù)學(xué)表達(dá)式。研究［19-20］將這種依賴關(guān)系歸結(jié)為共因退化進(jìn)行分析。但這類(lèi)模型沒(méi)有建立關(guān)于部件退化過(guò)程相互影響的解析模型，也未得到相應(yīng)的退化過(guò)程顯式模型，這限制了對(duì)該方法的深入研究和使用。

在負(fù)載分擔(dān)模型方面，大多注重于負(fù)載分配機(jī)制問(wèn)題，在負(fù)載分擔(dān)系統(tǒng)可靠性評(píng)估和壽命預(yù)測(cè)方面則采用傳統(tǒng)的壽命統(tǒng)計(jì)方法［21-22］。針對(duì)缺失壽命數(shù)據(jù)的情況，將退化模型引入負(fù)載分擔(dān)系統(tǒng)是一種可行的方法。文獻(xiàn)［23］研究負(fù)載分擔(dān)系統(tǒng)組件故障后前后狀態(tài)之間的依賴關(guān)系，采用了Wiener 退化過(guò)程，并加入了鏈接函數(shù)來(lái)體現(xiàn)應(yīng)力改變對(duì)漂移參數(shù)的影響。文獻(xiàn)［24］通過(guò)設(shè)置應(yīng)力系數(shù)，將其加入Wiener 退化模型分析了共載系統(tǒng)的可靠性和壽命預(yù)測(cè)。文獻(xiàn)［25］在組件退化的基礎(chǔ)上，考慮了隨機(jī)沖擊的影響，分析了共載系統(tǒng)的可靠性。文獻(xiàn)［26］研究了共載系統(tǒng)中退化速率與沖擊載荷的關(guān)系，說(shuō)明退化速率的改變使系統(tǒng)可靠性顯著降低。值得注意的是，這類(lèi)基于退化的共載系統(tǒng)研究，僅針對(duì)單一性能參數(shù)或組件進(jìn)行退化建模。建模思路與多階段退化模型類(lèi)似，如文獻(xiàn)［24］中的圖1，存活組件呈現(xiàn)多階段特性，退化拐點(diǎn)為組件失效時(shí)間。

圖1 雙組件共載系統(tǒng)退化過(guò)程Fig.1 Two-component load-sharing system degradation process

綜上所述，目前基于退化模型的共載系統(tǒng)的研究中，大多集中于負(fù)載變換后不同階段退化的相關(guān)性以及內(nèi)、外部應(yīng)力對(duì)退化的影響方面；而基于退化模型進(jìn)行多組件系統(tǒng)可靠性分析時(shí)，多元相關(guān)性退化模型得到廣泛引用，但組件之間有除串聯(lián)外的多種關(guān)聯(lián)方式以及其他因素（負(fù)載分擔(dān)）對(duì)組件退化進(jìn)程的影響，使得多元相關(guān)退化模型在進(jìn)行如共載系統(tǒng)的可靠性分析時(shí)會(huì)出現(xiàn)偏差。因此針對(duì)共載系統(tǒng)中存在二元退化相關(guān)的可靠性分析問(wèn)題，提出了一種兩階段退化模型，并給出了系統(tǒng)可靠性和壽命評(píng)估方法，能夠?qū)崿F(xiàn)兩種依賴關(guān)系影響下系統(tǒng)的可靠性分析與剩余壽命預(yù)測(cè)。

1 問(wèn)題描述

在多組件系統(tǒng)退化可靠性分析中，對(duì)組件間退化關(guān)聯(lián)性的研究是十分重要的。存在一種并聯(lián)系統(tǒng)，系統(tǒng)內(nèi)組件按照一定的規(guī)律各自承擔(dān)整個(gè)系統(tǒng)的一部分負(fù)載。其中一個(gè)組件失效時(shí)，存活組件會(huì)按照一定的規(guī)律重新分配負(fù)載，當(dāng)失效組件達(dá)到某一數(shù)量時(shí)，系統(tǒng)會(huì)表現(xiàn)出故障狀態(tài)，這種系統(tǒng)稱為共載系統(tǒng)。分析可知共載系統(tǒng)中組件會(huì)因負(fù)載分配策略而改變退化速率，從而呈現(xiàn)出階段性退化的特點(diǎn)。

圖1 為某雙組件共載系統(tǒng)退化過(guò)程。第一階段時(shí)，二者共同承擔(dān)任務(wù)負(fù)載。當(dāng)其中退化較快的組件到達(dá)失效閾值發(fā)生故障后退出系統(tǒng)，存活組件因承擔(dān)所有負(fù)載導(dǎo)致工作負(fù)荷加重，退化速率改變，出現(xiàn)變點(diǎn)。第二階段由存活組件承擔(dān)負(fù)載直至其到達(dá)失效閾值，系統(tǒng)失效。

因此在多組件系統(tǒng)退化可靠性分析中，若在構(gòu)建相關(guān)性退化模型時(shí)忽略了階段性變化，將使可靠性預(yù)測(cè)結(jié)果偏離實(shí)際。在考慮退化關(guān)聯(lián)性的情況下對(duì)共載系統(tǒng)進(jìn)行退化建模和可靠性分析，為組件退化過(guò)程和共載系統(tǒng)設(shè)立以下基本假設(shè)。

1.1 退化過(guò)程基本假設(shè)

1）組件僅存在退化失效，組件退化過(guò)程服從Wiener 過(guò)程，退化過(guò)程為增長(zhǎng)趨勢(shì)，初始退化量為0。

2）組件的退化水平超過(guò)其失效閾值時(shí)，認(rèn)為該組件失效，且不再影響其他組件的退化進(jìn)程。

3）組件只在工作應(yīng)力下開(kāi)始退化，工作應(yīng)力恒定的情況下，自身退化速率不變。

4）組件兩階段退化變點(diǎn)處的退化進(jìn)程是連續(xù)的，不存在躍遷。

5）退化進(jìn)程是可觀測(cè)的，即可知組件失效的先后順序。

1.2 二元共載系統(tǒng)假設(shè)

1）共載系統(tǒng)中包含有兩個(gè)組件，故障閾值相同，采用并聯(lián)的連接方式。當(dāng)兩個(gè)組件都失效時(shí)，系統(tǒng)失效。

2）系統(tǒng)載荷采用均分的方式分配給各個(gè)部件。一旦有組件失效，系統(tǒng)的載荷將由存活組件承擔(dān)。

3）系統(tǒng)組件出現(xiàn)階段性特征的條件是存在組件失效，使負(fù)載重新分配，在不考慮組件同時(shí)失效的情況下，可認(rèn)為階段變點(diǎn)前系統(tǒng)是可靠的。

4）若存在組件替換，則組件失效即可替換且只替換一次，忽略更換時(shí)間。替換組件具有相同的故障閾值，初始退化量為0。

2 雙組件共載系統(tǒng)兩階段退化過(guò)程建模

在分析雙組件共載系統(tǒng)的退化過(guò)程時(shí)，需同時(shí)考慮組件之間的相互影響對(duì)于退化進(jìn)程的影響和共載狀態(tài)下組件失效導(dǎo)致剩余組件退化速率的變化，基于此建立退化速率相關(guān)性的雙組件共載系統(tǒng)兩階段退化模型。

2.1 基于Wiener 過(guò)程的退化模型

2.1.1 一元Wiener 退化模型

根據(jù)Wiener 過(guò)程的定義，首先建立單個(gè)組件的退化模型，設(shè)組件的退化量為Xi(t)，i=1，2。定義組件i 的退化過(guò)程為非線性Wiener 退化過(guò)程［27］：

式中：ui(t)為組件的退化軌跡數(shù)學(xué)表達(dá)式；Bi(t)為布朗運(yùn)動(dòng)，表示退化過(guò)程中的不確定性；σi為擴(kuò)散參數(shù)；函數(shù)φi(t；θ)為組件的退化速率，θ 表示影響退化過(guò)程相關(guān)參數(shù)的集合。可以發(fā)現(xiàn)當(dāng)φi(t；θ)為一常數(shù)時(shí)，該模型為線性Wiener 退化過(guò)程。

2.1.2 二元Wiener 退化模型

設(shè)兩個(gè)組件的退化量由相關(guān)的二元Wiener過(guò)程X(t)=[ X1(t)，X2(t)]T描述，結(jié)合式（1），退化模型為

2.2 退化速率相關(guān)性模型建立

設(shè)組件退化速率包括：

1）βi，i=1，2，表示組件在工作應(yīng)力下的退化率參數(shù)，當(dāng)工作應(yīng)力改變時(shí)，βi會(huì)隨之改變。

2）g(Xj(t)，μij)i≠j，i，j=1，2，其 中g(shù)(?) 表示系統(tǒng)中兩個(gè)組件退化速率相互作用的數(shù)學(xué)表達(dá)式，μij表示組件j 對(duì)組件i 退化的影響程度，設(shè)為影響相關(guān)參數(shù)。μij的大小表示著影響程度的大小，μij為0 時(shí)表示無(wú)影響，μij為負(fù)數(shù)時(shí)表示退化影響存在抑制作用。

因此根據(jù)式（1），退化速率φi(t；θ)可表示為φi(t；βi，g(Xj(t)，μij)i≠j)，此時(shí)單 組件的 退化模型中退化路徑模型可改寫(xiě)為

分析可知得出uip(t)的解析形式是推導(dǎo)出退化進(jìn)程Xi(t)退化模型的關(guān)鍵，但Wiener 退化模型中含有布朗運(yùn)動(dòng)，存在隨機(jī)性，所以很難推導(dǎo)出精確的Xi(t)解析表達(dá)形式。結(jié)合式（1）由Wiener 退化過(guò)程的特性可知退化路徑代表著退化的趨勢(shì)，可將隨機(jī)布朗運(yùn)動(dòng)視作退化過(guò)程中的噪聲［27］，在分析推導(dǎo)退化路徑函數(shù)時(shí)去掉該類(lèi)隨機(jī)項(xiàng)的影響，因此可分別近似得到組件1和組件2的退化速率表達(dá)式（4）。

在已知各函數(shù)形式的情況下，可通過(guò)一階微分方程組得出退化率相互作用模型的解析形式：

式中：φ1、φ2為組件1和組件2 的退化速率函數(shù)。

根據(jù)式（4）可以看出，組件之間的相關(guān)性是通過(guò)uip(t)來(lái)表述的，當(dāng)給定一個(gè)時(shí)刻時(shí)，uip(t)是確定值。同時(shí)根據(jù)式（4）可得到退化速率中不包含uip(t)的退化模型的解析解形式，因此可認(rèn)為得到解析解后的退化模型是相互獨(dú)立的［18］，在二元Wiener 退化模型的基礎(chǔ)上結(jié)合式（4），可得到雙組件退化速率相關(guān)性退化模型：

式中：up(t)=[u1p(t)，u2p(t)]T為退化路徑向量；∑p=diag(σ12，σ22)為兩個(gè)組件的退化量方差矩陣。

2.3 雙組件共載系統(tǒng)分階段退化模型建立

1）考慮一般情況

假設(shè)組件i 在τi時(shí)刻因負(fù)載變化導(dǎo)致退化速率發(fā)生變化，組件退化過(guò)程即分成兩階段，若考慮隨機(jī)因素與負(fù)載之間的相互作用，則模型中的擴(kuò)散參數(shù)也隨之改變。因此首先給出單組件兩階段的退化模型，結(jié)合式（1）可得Xi(t)的表達(dá)式為

式中：ui2為第二階段退化路徑表達(dá)式；σi2為第二階段的擴(kuò)散參數(shù)；τi為變點(diǎn)時(shí)刻。

對(duì)于雙組件共載系統(tǒng)，由分析可知當(dāng)雙組件均處于正常工作狀態(tài)時(shí)，共同分擔(dān)負(fù)載，組件間存在相互影響關(guān)系，影響著彼此的退化進(jìn)程。當(dāng)其中一個(gè)組件失效時(shí)，存活組件負(fù)載變化，其退化會(huì)出現(xiàn)兩階段的特性，同時(shí)可知，階段變點(diǎn)時(shí)刻即為組件失效時(shí)刻，此時(shí)變點(diǎn)時(shí)刻τ1為組件1的壽命。設(shè)組件1 最先失效，其失效時(shí)間為t1，由上述分析有t1=τ1，結(jié)合式（5）和式（6）可得到雙組件共載兩階段退化模型X(t)的表達(dá)式為

式中：u22(t-t1)為組件2 第二階段的退化軌跡；σ22為組件2 第二階段的擴(kuò)散參數(shù)。

2）考慮更為復(fù)雜的情況

圖2 為考慮組件替換時(shí)的退化過(guò)程，以式（7）建立的模型為基礎(chǔ)進(jìn)行擴(kuò)展分析。組件1 失效后立即替換，此時(shí)雙組件退化進(jìn)程仍保持相互影響的狀態(tài)，共同分擔(dān)負(fù)載。故退化分為3 個(gè)階段：

第一階段在組件1 到達(dá)壽命τ1之前，雙組件退化模型符合式（5）。

第二階段更換失效組件后，設(shè)更換件的退化模型為

式中：X1s(t)為替換組件1 的退化量；u1s(t-t1)為替換組件1 的退化路徑；σ1s為換組件1 的擴(kuò)散參數(shù)；B1s(t-t1)為替換組件1 的布朗運(yùn)動(dòng)。此時(shí)系統(tǒng)仍處于負(fù)載分擔(dān)階段且相互影響，由式（4）知組件間的影響因素與組件的退化量相關(guān)，替換后組件退化量減小使組件2 的退化進(jìn)程減緩，因此組件2 會(huì)在替換時(shí)刻出現(xiàn)階段特性，結(jié)合式（3），該階段的退化路徑模型為

式中：β1s表示替換組件1 工作應(yīng)力下的退化率參數(shù)；u2s(t-t1)表示替換組件后，組件2 的退化路徑；X2s(x)表示替換組件后，組件2 的退化量；φ1s為替換組件1 的退化速率函數(shù)；φ2s為替換組件后，組件2 的退化速率函數(shù)。

第三階段替換組件1 后設(shè)組件2 先失效，替換組件1承擔(dān)全部負(fù)載，退化會(huì)出現(xiàn)兩階段的特性，階段變點(diǎn)為組件2失效時(shí)刻，失效時(shí)間為t2，由上述分析有t2=τ2，結(jié)合式（5）、式（8）和式（9）得到考慮組件替換后，雙組件共載兩階段退化模型X(t)表達(dá)式為

式中：us(t-t1)=[u1s(t-t1)，u2s(t-t1)]T為退化路徑向量；∑s=diag(σ1s2，σ22)為兩個(gè)組件的退化量方差矩陣；u1s2為替換組件1 第二階段退化路徑表達(dá)式；σ1s2為替換組件1 第二階段的擴(kuò)散參數(shù)；τ2為組件2 的壽命。

3 基于退化速率相關(guān)性的共載系統(tǒng)可靠性分析

雙組件系統(tǒng)在分析壽命分布和可靠度時(shí)，由于組件之間的相互影響以及組件失效帶來(lái)的負(fù)載變化，導(dǎo)致剩余組件退化速率發(fā)生變化，系統(tǒng)的壽命以及可靠度無(wú)法由簡(jiǎn)單的獨(dú)立關(guān)系得出。因此建立基于退化速率相關(guān)性的共載系統(tǒng)可靠性模型，并建立剩余壽命模型以用于剩余壽命預(yù)測(cè)。

3.1 可靠性模型

3.1.1 單組件可靠性模型

對(duì)于式（1），當(dāng)φi(t；θ)為一常數(shù)βi時(shí)，式（1）為標(biāo)準(zhǔn)的線性Wiener 過(guò)程，設(shè)故障閾值為D，其壽命分布由逆高斯分布［27］給出：

其概率密度函數(shù)為

然而在一般情況下，式（1）是一個(gè)非線性Wiener 過(guò)程，在這種情況下，時(shí)間-尺度變換是一種線性變換方法，比如et、tb等，通過(guò)線性化簡(jiǎn)便壽命分布的推導(dǎo)并減小后續(xù)參數(shù)估計(jì)的困難。但并非所有函數(shù)都能做線性化處理，使得很難推導(dǎo)出精確的非線性Wiener 過(guò)程首達(dá)時(shí)間分布函數(shù)的一般形式。對(duì)于非線性Wiener 過(guò)程文獻(xiàn)［28］基于Kolmogorov 方程推導(dǎo)出了Xi(t)首次達(dá)到失效閾值D 時(shí)間的近似概率密度函數(shù)：

可靠性模型為Ri(t)=1-Fi(t)。

3.1.2 單組件兩階段可靠性模型

對(duì)于服從Wiener 過(guò)程的單組件兩階段退化模型，一般可靠性模型也采用分段的描述，以線性Wiener 過(guò)程為例，在僅考慮在第二階段失效的情況下結(jié)合式（6），令階段變點(diǎn)為τi，Xi(τi)=xi，組件的壽命分布［29-30］為

然而式（14）中xi作為第二階段退化進(jìn)程的初始退化量，是隨機(jī)的且服從正態(tài)分布，忽略這一特性可能會(huì)限制兩階段系統(tǒng)可靠性評(píng)估的準(zhǔn)確性［31］。為了克服這一局限性，構(gòu)建新的兩階段系統(tǒng)可靠性模型，首先建立僅考慮第二階段失效的情況下的壽命分布：

3.1.3 雙組件共載系統(tǒng)可靠性模型

1）考慮一般情況：設(shè)雙組件共載系統(tǒng)的系統(tǒng)壽命累積分布函數(shù)為F(t)，組件1 與組件2 的壽命分別為τ1、τ2，則系統(tǒng)壽命分布表示為：F(t)=P{τ1≤t，τ2≤t}。在第一階段可靠的條件下，設(shè)τ1<τ2。在t1=τ1時(shí)，變點(diǎn)與失效組件的失效點(diǎn)一致，存活組件的變點(diǎn)是由于組件失效產(chǎn)生的，故變點(diǎn)是隨組件的失效時(shí)間而改變，因此可認(rèn)為變點(diǎn)的概率密度函數(shù)fτ1(t1)與失效組件的壽命概率密度函數(shù)一致，結(jié)合式（15）得系統(tǒng)的壽命分布：

式中：x2為存活組件X2在變點(diǎn)τ1處的退化量值；(x2，t1)為變點(diǎn)時(shí)間和存活組件退化量的聯(lián)合概率分布。

分析P{τ1≤t，τ2≤t|τ1=t1，X2(τ1)=x2}可發(fā)現(xiàn)，在給定τ1的情況下，有X2(τ1)=x2，該條件符合初始值為x2的退化過(guò)程，可描述為

對(duì)比式（7）可知該條件概率所描繪的為存活組件第二階段的退化過(guò)程，同時(shí)這也證明了在并聯(lián)系統(tǒng)中，影響系統(tǒng)壽命的是系統(tǒng)中壽命更長(zhǎng)的組件，即可得出第二階段存活組件的壽命累積分布函數(shù)，P{τ1≤t，τ2≤t|τ1=t1，X2(τ1)=x2}可改寫(xiě) 為P{X22(t-t1)≥D-x2}=F2(t-t1)，F(xiàn)2(t-t1)為組件2的壽命分布，結(jié)合式（14）即可得出其表達(dá)式。

因組件之間不是獨(dú)立的，在同時(shí)考慮x2與t1的隨機(jī)性時(shí)，兩個(gè)條件概率分布(x2)和(t1)之間存在相關(guān)性，在統(tǒng)計(jì)學(xué)上，需考慮為二者的聯(lián)合概率分布(x2，t1)。通過(guò)2.2節(jié)的描述可知，在引入退化速率相互作用模型，得出的新的退化模型解析式后，可認(rèn)為兩個(gè)組件是獨(dú)立的，其相關(guān)系數(shù)為0，因此可得(x2，t1)=(x2) fτ1(t1)。其中，fτ1(t1)為組件1 的壽命 概率密度函數(shù)，結(jié)合式（12）或式（13）可得(x2)為X2退化量在τ1處的退化量分布，由Wiener 過(guò)程的數(shù)學(xué)特性，該變量服從正態(tài)分布：

綜上所述，在式（16）的基礎(chǔ)上，共載系統(tǒng)的系統(tǒng)壽命分布為

系統(tǒng)的壽命密度函數(shù)為

系統(tǒng)可靠度R(t)為

式中：f2(t-t1)為第二階段退化過(guò)程的壽命分布密度函數(shù)。

2）考慮組件替換的情況：

與無(wú)組件替換的分析方法類(lèi)似，分為三階段討論：①通過(guò)圖2 分析可知，替換組件的起始退化時(shí)間與失效組件1 的失效點(diǎn)一致，服從組件1 的壽命分布；②因退化量改變，組件2 出現(xiàn)退化變點(diǎn)，該變點(diǎn)由于組件替換產(chǎn)生的，故變點(diǎn)是隨組件1 的失效時(shí)間而改變。在組件2 失效前，替換組件與組件2 相互影響；③替換組件的變點(diǎn)是由于組件2 失效產(chǎn)生的，故變點(diǎn)是隨組件2 的失效時(shí)間而改變。設(shè)替換組件壽命為τs，結(jié)合式（16）得系統(tǒng)的壽命分布：

式中：x1s為存活的替換組件1 在變點(diǎn)τ2處的退化量值(x1s，t2)為變點(diǎn)時(shí)間和替換組件退化量的聯(lián)合概率分布。

與式（16）的分析方法一致，對(duì)式（23）進(jìn)行分析，P{τ2≤t，τs≤t|τ2=t2，X1s(τ2-τ1)=x1s} 可改寫(xiě)為P{X1s(t-t2)≥D-x1s}=F1s(t-t2)。獨(dú)立性(x1s，t2)=(x1s)(t2)，其中(x1s)為X1s退化量在τ2處的退化量分布，由Wiener 過(guò)程的數(shù)學(xué)特性，該變量服從正態(tài)分布。(t2)為組件2 的壽命概率密度函數(shù)，由上述分析可知其存在階段性，與組件1 失效有關(guān)，故(t2)為組件2 的剩余壽命概率密度函數(shù)，結(jié)合式（22）的分析方法，可得其表達(dá)式。

綜上，組件替換情形作為雙組件共載系統(tǒng)分階段退化模型的拓展，其退化模型建立與可靠性分析均與雙組件共載系統(tǒng)分階段退化模型類(lèi)似，因此后續(xù)工作將以雙組件共載系統(tǒng)分階段退化模型為重點(diǎn)進(jìn)行分析。

3.1.4 模型數(shù)值求解

式（21）是復(fù)雜二重積分，無(wú)法得到積分原函數(shù)，難以推導(dǎo)出可靠性的解析解形式。因此采用蒙特卡洛模擬的方法，具體的流程如圖3所示。

圖3 積分近似計(jì)算過(guò)程Fig.3 Integral approximation process

由圖3 所求f 為各離散點(diǎn)處密度函數(shù)之和，令積分域?yàn)镾，則各離散化面積為S/(n1×n2)，可得該二重積分的數(shù)值解為(S×f)/(n1×n2)，采樣的值越多，計(jì)算結(jié)果越精確。

3.2 模型參數(shù)估計(jì)

由式（3）給出的退化模型中，對(duì)于組件i，待估參數(shù)有｛βi，μij，σi｝。為不失一般性，假設(shè)雙組件的初始退化模型均為線性模型，且組件間相互影響關(guān)系也為線性關(guān)系，則可得退化軌跡表達(dá)式：在可觀 測(cè)的情況下已知組件退化順序，設(shè)組件1 先達(dá)到失效閾值，雙組件的初始退化量為均0。根據(jù)式（7）可得退化模型X(t)為

式中：uL(t)=[u1L(t)，u2L(t)]T；方差矩陣∑L=diag(σ12，σ22)；β22為組件2 第二階 段的退化率參數(shù)。

采用極大似然方法，估計(jì)式（24）退化模型的參數(shù)。一般的退化速率相互影響估計(jì)方法采用定時(shí)截尾退化數(shù)據(jù)，而該方法需根據(jù)組件失效閾值對(duì)退化數(shù)據(jù)進(jìn)行預(yù)處理，以變點(diǎn)為時(shí)間節(jié)點(diǎn)將退化數(shù)據(jù)分成兩部分，用以估計(jì)不同階段模型的參數(shù)。設(shè)觀測(cè)間隔時(shí)間為一定值Δt，有8 個(gè)待估未知參數(shù){ β1，β2，β22，μ12，μ21，σ1，σ2，σ22}。此模型類(lèi)似于兩階段退化模型，兩階段退化模型含有退化拐點(diǎn)，一般認(rèn)為拐點(diǎn)是隨機(jī)變量且缺乏完全數(shù)據(jù)，參數(shù)估計(jì)十分困難。而該模型中的退化變點(diǎn)與組件1 的失效時(shí)間有關(guān)，變點(diǎn)分布與組件1 的壽命分布一致，所以該模型可根據(jù)組件1 的壽命分布，對(duì)變點(diǎn)及其分布參數(shù)進(jìn)行估計(jì)。

由Wiener 退化過(guò)程的數(shù)學(xué)特性可知，Wiener退化過(guò)程具有獨(dú)立增量，且滿足正態(tài)分布。則ΔX 在不同階段下的概率密度函數(shù)［18］為

其中，平均增量Δui為

式中：Δxi為退化增量。

可以計(jì)算得到觀測(cè)數(shù)據(jù)的似然函數(shù)中各項(xiàng)，第一階段組件1和組件2 的似然函數(shù)和第二階段組件2 的似然函數(shù)如式（26）所示。

將同一類(lèi)型多個(gè)組件實(shí)際的退化增量數(shù)據(jù)代入后，對(duì)該類(lèi)型似然函數(shù)求和，即可得到同類(lèi)型的總體似然函數(shù)，再對(duì)待估計(jì)參數(shù)求偏導(dǎo)，從而得到各參數(shù)的點(diǎn)估計(jì)值。為滿足參數(shù)的無(wú)偏性，需對(duì)各參數(shù)進(jìn)行無(wú)偏性檢驗(yàn)。根據(jù)正態(tài)分布的特性，在均值參數(shù)與方差參數(shù)均未知時(shí)，均值參數(shù)是無(wú)偏的，方差參數(shù)需根據(jù)樣本數(shù)進(jìn)行修正。

根據(jù)Cayley-Hamilton 定理，結(jié) 合式（4）和式（24），可推導(dǎo)出第一階段退化速率相互作用模型退化路徑解 析解式（27）。將參數(shù){β1，β2，μ12，μ21，σ1，σ2}代入式（27）中，結(jié)合式（24）即可得到第一階段退化模型。再將估計(jì)的參數(shù){ β22，σ22}代入式（24）中，可得第二階段退化模型。

3.3 剩余壽命預(yù)測(cè)

若系統(tǒng)運(yùn)行到tk時(shí)刻還沒(méi)有失效，則其剩余壽命可定義為

由于性能參數(shù)退化過(guò)程的不確定性，導(dǎo)致產(chǎn)品的剩余壽命具有隨機(jī)性。從實(shí)際應(yīng)用的角度來(lái)看，需要得到剩余壽命的概率密度函數(shù)。同時(shí)基于雙階段模型，變點(diǎn)前后退化模型不同，因此需根據(jù)觀測(cè)點(diǎn)建立不同的剩余壽命的概率密度函數(shù)，分為變點(diǎn)前觀測(cè)與變點(diǎn)后觀測(cè)兩種情況。此處對(duì)于雙組件系統(tǒng)假設(shè)τ1<τ2。

設(shè)定系統(tǒng)從t1到tn時(shí)刻線性能退化數(shù)據(jù)為

情況1當(dāng)X1(tk)

結(jié)合式（13）和式（18），可分別得出tk后組件1 的變點(diǎn)處的剩余壽命分布和退化量分布：

第二階段的剩余壽命分布不變，結(jié)合式（29）和式（30），可最終得出系統(tǒng)的剩余壽命分布的概率密度函數(shù)：

由式（19）和式（22）得X1(tk)

情況2當(dāng)X1(tk)>D時(shí)

系統(tǒng)剩余壽命即為存活組件2 的剩余壽命。結(jié)合式（7）和式（14），可推出系統(tǒng)的剩余壽命分布為

4 算例分析

4.1 系統(tǒng)概述

4.1.1 案例描述

利用實(shí)際的退化數(shù)據(jù)驗(yàn)證提出雙階段退化模型和可靠性模型。一般的飛機(jī)系統(tǒng)中，多采用多余度設(shè)計(jì)，以保障飛機(jī)飛行的安全性。兩個(gè)航空液壓組件構(gòu)成飛機(jī)內(nèi)部的傳動(dòng)裝置，共同驅(qū)動(dòng)同一組件工作，因此兩個(gè)液壓組件構(gòu)成了一個(gè)共載系統(tǒng)。液壓組件因工作應(yīng)力或振動(dòng)等可能產(chǎn)生疲勞裂紋，當(dāng)其疲勞裂紋到達(dá)某一閾值時(shí)，伺服燃油發(fā)生泄漏，液壓組件故障，因此兩個(gè)組件的疲勞裂紋退化信號(hào)可看作該共載系統(tǒng)的兩個(gè)性能參數(shù)。案例數(shù)值數(shù)據(jù)取自文獻(xiàn)［32］，將每百萬(wàn)循環(huán)沖擊次數(shù)看作時(shí)間。退化路徑如圖4 所示，展示了10 組共載系統(tǒng)中組件的退化過(guò)程。

圖4 疲勞裂紋退化數(shù)據(jù)Fig.4 Fatigue crack degradation data

具有兩個(gè)液壓組件的飛機(jī)傳動(dòng)裝置共載系統(tǒng)故障是由兩個(gè)階段退化過(guò)程觸發(fā)的，從其中一個(gè)組件故障，到兩個(gè)組件故障使系統(tǒng)失效。在第一階段，由于共同應(yīng)力環(huán)境，兩個(gè)性能參數(shù)的退化過(guò)程是相關(guān)的。在第二階段，生存組件承擔(dān)所有傳動(dòng)任務(wù)，因此，組件的退化速率改變，這種現(xiàn)象可以從圖4 中觀察到。當(dāng)兩個(gè)組件的退化程度都達(dá)到它們的閾值時(shí)，飛機(jī)傳動(dòng)共載系統(tǒng)完全失效［33］。因此，式（7）適用于對(duì)該案例進(jìn)行退化建模。

這10 組退化數(shù)據(jù)作為歷史數(shù)據(jù)用以訓(xùn)練模型進(jìn)行極大似然估計(jì)，得到模型參數(shù)的點(diǎn)估計(jì)值；再選取d 組的樣本用于檢驗(yàn)退化擬合程度的實(shí)際數(shù)據(jù)，所選取的樣本具有明顯的階段性特征，符合所研究的背景情況。

4.1.2 系統(tǒng)可靠性分析總體思路

共載系統(tǒng)可靠性分析的總體流程包括：①基于Wiener 退化過(guò)程建立考慮負(fù)載分擔(dān)的退化速率相互作用模型，得出退化解析形式；②根據(jù)退化模型與依賴耦合關(guān)系建立共載系統(tǒng)的可靠性模型，并基于極大似然估計(jì)與蒙特卡洛仿真的方法求解。整體思路如圖5 所示。

圖5 雙組件共載系統(tǒng)可靠性分析方法Fig.5 Reliability analysis method for two-component load-sharing system

4.2 退化過(guò)程分析

4.2.1 建立退化模型

為了說(shuō)明所建立模型的有效性，將進(jìn)一步開(kāi)展對(duì)比研究。設(shè)考慮負(fù)載分擔(dān)的雙組件退化速率相關(guān)性模型為M1；不考慮負(fù)載分擔(dān)，僅考慮雙組件退化速率相關(guān)性的模型為M2［16］。對(duì)于組件退化速率相互作用的連續(xù)函數(shù)模型g(Xj(tj)，μij)i≠j，因?yàn)榻M件間多種隨機(jī)影響g(Xj(tj)，μij)i≠j十分復(fù)雜，故采用泰勒展開(kāi)取第一項(xiàng)，M2 與M1第一階段退化模型一致。在假設(shè)M1 中組件1 先達(dá)到閾值的情況下，設(shè)組件1 失效時(shí)刻為τ1，其失效時(shí)間為t1，由上述分析有t1=τ1，M1 的退化模型為式（24），M2 的退化模型為

式中：β1′、β2′為組件1和組件2 工作應(yīng)力下的退化速率；μ12′、μ21′為組件1和組件2 退化影響相關(guān)參數(shù)；σ1′、σ2′為組件1和組件2 的擴(kuò)散參數(shù)。

4.2.2 參數(shù)估計(jì)

首先對(duì)這組數(shù)據(jù)做預(yù)處理，令其初始退化量從0 開(kāi)始。假設(shè)當(dāng)退化水平達(dá)到0.3，即各組件的失效閾值D=0.3時(shí)，視為組件故障。首先對(duì)M1進(jìn)行參數(shù)估計(jì)，利用3.2 節(jié)給出的似然估計(jì)方法以及式（26），對(duì)圖4 中的退化數(shù)據(jù)進(jìn)行處理，根據(jù)組件1 的失效時(shí)刻將退化數(shù)據(jù)分為兩個(gè)階段，分別對(duì)式（24）中兩個(gè)階段退化模型的參數(shù)進(jìn)行估計(jì)。M2 則利用一般的定時(shí)截尾退化數(shù)據(jù)［16］，結(jié)合式（26）得到式（34）的參數(shù)估計(jì)值。表1 給出了M1和M2 的參數(shù)點(diǎn)估計(jì)值，結(jié)合式（27）即可得到M1和M2 退化模型中退化路徑的解析形式，再結(jié)合式（24）和式（34），可得到完整的雙組件Wiener退化過(guò)程。

表1 模型未知參數(shù)估計(jì)及模型比較Table 1 Model unknown parameter estimation and model comparison

4.2.3 結(jié)果分析

為了對(duì)比M1和M2 中退化路徑對(duì)實(shí)際退化路線的擬合程度，以赤池信息準(zhǔn)則（Akaike Information Criterion，AIC）和貝葉斯信息準(zhǔn)則（Bayesian Information Criterion，BIC）來(lái)衡量模型擬合數(shù)據(jù)的優(yōu)良性，AIC和BIC 值越小，模型擬合程度越高。M1和M2 的AIC和BIC 值如表1 所示。然后將兩種模型下預(yù)期的退化過(guò)程與實(shí)際的退化過(guò)程進(jìn)行比較，圖6 展示了模型退化路徑的退化過(guò)程曲線對(duì)比情況，圖7 展示了含有隨機(jī)性的完整的Wiener過(guò)程下組件退化對(duì)比的情況。

圖6 退化路徑模型的退化過(guò)程曲線對(duì)比圖Fig.6 Comparison chart of degradation process curve for degradation path model

圖7 Wiener 退化過(guò)程對(duì)比圖Fig.7 Comparison chart of Wiener degradation process

經(jīng)過(guò)M1和M2 參數(shù)的估計(jì)和模型對(duì)比，結(jié)合表1和圖6、圖7 分析可知：

1）M1的AIC和BIC 值更小，說(shuō)明在共載的情況下，M1 更加符合真實(shí)的退化過(guò)程。

2）對(duì)比M1和M2 中βi和μij數(shù)值可發(fā)現(xiàn)，組件1 因應(yīng)力產(chǎn)生的退化速率較組件2 更大，且受到組件2 的影響更劇烈，所以組件1 的退化趨勢(shì)更加迅速，符合圖6和圖7 中所展示的退化趨勢(shì)。

3）結(jié)合圖6 中M1和M2 退化路徑可得，在組件1 到達(dá)失效閾值之前，M1和M2 的退化路徑趨勢(shì)都比較接近于真實(shí)的退化過(guò)程趨勢(shì)。組件1 失效后，觀察M1 與M2 中組件2 的退化進(jìn)程，M1 因考慮了負(fù)載的變化與組件1 失效所帶來(lái)的影響，退化趨勢(shì)改變，其后續(xù)退化過(guò)程能很好的吻合實(shí)際退化過(guò)程，而M2 延續(xù)了之前的退化趨勢(shì)，隨著時(shí)間的推移，與實(shí)際退化進(jìn)程誤差會(huì)越來(lái)越大。

4）圖7 在圖6 的基礎(chǔ)上加入了標(biāo)準(zhǔn)布朗運(yùn)動(dòng)，描繪的是含有隨機(jī)性的Wiener 退化過(guò)程。對(duì)比圖6和圖7 可發(fā)現(xiàn)，在加入隨機(jī)性后，M1 的退化趨勢(shì)仍比M2 更接近實(shí)際退化情況，同時(shí)在Wiener 退化過(guò)程的基礎(chǔ)上可以直接得到系統(tǒng)的壽命分布，這也體現(xiàn)以Wiener 過(guò)程為例的隨機(jī)過(guò)程在退化模型中運(yùn)用的優(yōu)越性。

4.3 可靠性分析

4.3.1 系統(tǒng)可靠性分析

由式（24）和式（27）可知，M1 第一階段為非線性退化過(guò)程，第二階段為線性退化過(guò)程，因此結(jié)合式（15）、式（18）和式（21），代入表1 中的參數(shù)估計(jì)值，即可得到M1 的可靠性模型解析形式。

為了驗(yàn)證式（21）和式（22）的有效性和分析的準(zhǔn)確性，將M1 的解析解和模擬解進(jìn)行了比較。采用蒙特卡洛仿真的方法模擬組件退化過(guò)程，結(jié)合組件間相互影響關(guān)系和負(fù)載分擔(dān)，得到大量系統(tǒng)失效時(shí)間從而得出壽命分布與概率密度函數(shù)仿真解。該案例的退化仿真過(guò)程流程圖如圖8所示。

圖8 蒙特卡洛退化仿真流程Fig.8 Monte Carlo degradation simulation process

通過(guò)退化過(guò)程模擬程序和組件失效閾值生成足夠多故障次數(shù)，將失效時(shí)間頻率壽命直方圖與式（22）給出的解析壽命密度函數(shù)進(jìn)行比較，其中模擬過(guò)程運(yùn)行了5 000次，如圖9 所示。

圖9 基于兩種解的系統(tǒng)壽命概率密度對(duì)比圖Fig.9 System lifetime probability density comparison plot based on two solutions

由結(jié)果圖仿真數(shù)據(jù)和解析分布對(duì)比可以看出，所提的可靠性模型可較為精確地得到系統(tǒng)的壽命分布。同時(shí)基于仿真數(shù)據(jù)可得出系統(tǒng)可靠度的仿真解，將可靠度的解析解和仿真解進(jìn)行比較，如圖10 所示，可見(jiàn)二者結(jié)果吻合較好。

圖10 基于兩種解的系統(tǒng)可靠度對(duì)比圖Fig.10 System reliability comparison chart based on two solutions

4.3.2 結(jié)果分析

根據(jù)式（34）可知M2 雙組件退化過(guò)程均為非線性，同時(shí)結(jié)合2.2 節(jié)的獨(dú)立性分析，可得出M2的可靠度模型為

式中：F1(t)、F2(t)為 服從式（13）的累積分布函數(shù)。

在4.3.1 節(jié)分析的基礎(chǔ)上，對(duì)比了M1 與M2下系統(tǒng)、組件1和組件2 的可靠度曲線，如圖11所示，其中M1 組件2 的可靠度曲線由式（15）求得。

圖11 M1 與M2 可靠度對(duì)比圖Fig.11 M1 and M2 reliability comparison chart

經(jīng)過(guò)M1和M2 可靠性曲線的對(duì)比分析可知：

1）對(duì)比M2 下組件2和系統(tǒng)的可靠性曲線，二者曲線幾乎重合，因此M2 的可靠性模型無(wú)法體現(xiàn)組件1 失效對(duì)系統(tǒng)可靠性的影響。

2）對(duì)比M1和M2 系統(tǒng)可靠度曲線，可以看出M1 可靠性更低，結(jié)合圖6和圖7 可以得出M1比M2 先失效，因此該差距是合理的。

3）對(duì)比M1 中組件2和系統(tǒng)可靠性曲線可發(fā)現(xiàn)并不重合，這是因?yàn)镸1 可靠度模型式（21）中描述了組件1 失效對(duì)組件2 退化的影響。同時(shí)發(fā)現(xiàn)組件2 的曲線更加陡峭。分析可知，可靠性的陡峭程度反映的為該處壽命密度函數(shù)，越陡峭說(shuō)明其方差越小，壽命集中于期望附近。因?yàn)橄到y(tǒng)可靠度有組件1 失效的影響，隨機(jī)性更強(qiáng)，方差更大。

4.4 剩余壽命分析

首先根據(jù)圖4 的退化數(shù)據(jù)，在除去初始退化時(shí)刻后，列出8 個(gè)觀測(cè)點(diǎn)，根據(jù)組件1 的失效時(shí)刻，剩余壽命建模包括：

階段1：當(dāng)觀測(cè)點(diǎn)在組件1 失效之前時(shí)，結(jié)合式（31）可得每個(gè)觀測(cè)點(diǎn)處的剩余壽命的概率密度函數(shù)。

階段2：當(dāng)觀測(cè)點(diǎn)在組件1 失效之后時(shí)，結(jié)合式（33）可得每個(gè)觀測(cè)點(diǎn)處的剩余壽命的概率密度函數(shù)。

根據(jù)建立的模型，分析了不同觀測(cè)點(diǎn)處剩余壽命的概率密度函數(shù)，并將觀測(cè)點(diǎn)處所預(yù)測(cè)的剩余壽命期望均值與剩余壽命真實(shí)值對(duì)比，并給出了剩余壽命期望均值95%區(qū)間估計(jì)，如圖12和圖13 所示。

圖12 不同預(yù)測(cè)時(shí)刻系統(tǒng)剩余壽命的概率密度函數(shù)Fig.12 Probability density function of remaining lifetime of the system at different prediction moments

圖13 不同預(yù)測(cè)時(shí)刻系統(tǒng)剩余壽命期望均值與區(qū)間估計(jì)Fig.13 Mean and interval estimation of system remaining life expectations at different prediction moments

從圖12 中可以看出，相比而言前幾個(gè)時(shí)刻的預(yù)測(cè)結(jié)果與真實(shí)值相差較大。從實(shí)際出發(fā)可分析到，在系統(tǒng)早期的使用中失效的概率很小，因此早期預(yù)測(cè)結(jié)果不準(zhǔn)確是可以接受的。從整體上看，剩余壽命預(yù)測(cè)均值近似擬合其真實(shí)值，且逐步向真實(shí)值逼近，同時(shí)結(jié)合圖13 可看出95%區(qū)間估計(jì)包含了真實(shí)值。為了具體化實(shí)際值與預(yù)測(cè)值之間的差異性，采用相對(duì)誤差表述，如表2所示。說(shuō)明了該模型可適用性，為后續(xù)維修決策等提供了便利。

表2 各觀測(cè)點(diǎn)剩余壽命相對(duì)誤差值Table 2 Relative error value of the remaining life of each observation point

5 結(jié)論

建立了基于退化相關(guān)的雙組件共載退化模型，在退化進(jìn)程相關(guān)和負(fù)載分擔(dān)的情況下預(yù)測(cè)組件的退化趨勢(shì)，并進(jìn)行系統(tǒng)的可靠性分析。①采用Wiener 退化過(guò)程，模擬由于未觀察到的環(huán)境因素或環(huán)境因素對(duì)退化過(guò)程的未知影響而導(dǎo)致的無(wú)法解釋的隨機(jī)性，同時(shí)具有良好的數(shù)學(xué)運(yùn)算性質(zhì)。②模型具有實(shí)際的物理意義，退化速率相互影響模型能夠描述兩個(gè)組件退化進(jìn)程之間的相關(guān)關(guān)系，且能得出退化模型的解析形式。③基于階段性退化過(guò)程，構(gòu)建了雙組件共載系統(tǒng)的可靠性模型，同時(shí)建立了不同階段下系統(tǒng)的剩余壽命模型。該模型能夠描述組件失效對(duì)于系統(tǒng)可靠度的影響。

基于二元退化相關(guān)的共載系統(tǒng)可靠性分析方法綜合考慮了退化相關(guān)性和負(fù)載分擔(dān)對(duì)于系統(tǒng)的影響，給出了基于退化模型的共載系統(tǒng)可靠性解析式。將蒙特卡羅模擬法應(yīng)用于可靠性分析中，解決了計(jì)算可靠性解析解中復(fù)雜的多重積分運(yùn)算情況，為后續(xù)實(shí)例驗(yàn)證提供便利。

所提方法為雙組件共載系統(tǒng)可靠性分析提供了可行思路，提高系統(tǒng)可靠性預(yù)測(cè)的精度。但更常見(jiàn)的共載系統(tǒng)是基于k-out-of-n 的多組件系統(tǒng)；同時(shí)，針對(duì)不同類(lèi)型的組件，應(yīng)力變化對(duì)于退化速率的影響是不同的。因此針對(duì)多組件共載系統(tǒng)可靠性及應(yīng)力與退化速率關(guān)系的進(jìn)一步分析，提出更通用的退化模型以及更準(zhǔn)確的退化速率表達(dá)形式，是后續(xù)研究的重點(diǎn)。