符弘嵐，張皓，高揚(yáng)

1. 中國(guó)科學(xué)院 空間應(yīng)用工程與技術(shù)中心，北京 100094 2. 中國(guó)科學(xué)院大學(xué) 計(jì)算機(jī)與控制學(xué)院，北京 100049

立方星是一種以1U (10 cm×10 cm×10 cm)為基本單元的可擴(kuò)展小型衛(wèi)星，與傳統(tǒng)衛(wèi)星相比其開(kāi)發(fā)周期短且研制成本低[1-3]。立方星體積小重量輕，配備的傳感器、執(zhí)行機(jī)構(gòu)等硬件精度差，需要設(shè)計(jì)專(zhuān)門(mén)的控制算法抵消低精度帶來(lái)的定姿、控制誤差[4]，這導(dǎo)致立方星在功能上比不上大型衛(wèi)星，因此在開(kāi)發(fā)應(yīng)用初期立方星主要用于技術(shù)驗(yàn)證和作為教育平臺(tái)使用[5]。為了更好地發(fā)揮立方星的優(yōu)勢(shì)，人們使用大量立方星以較低成本形成分布式空間傳感器網(wǎng)絡(luò)，實(shí)現(xiàn)傳統(tǒng)單個(gè)衛(wèi)星無(wú)法實(shí)現(xiàn)的功能[5]。構(gòu)建立方星星座以執(zhí)行商業(yè)、科學(xué)、技術(shù)驗(yàn)證任務(wù)已經(jīng)率先在近地軌道實(shí)現(xiàn)。例如，Planet Labs建立了由120多顆Dove衛(wèi)星組成的星座[6]。噴氣推進(jìn)實(shí)驗(yàn)室(JPL)設(shè)計(jì)在低地球軌道(LEO)上部署雨立方星座(Raincube Constellation)形成空間雷達(dá)網(wǎng)絡(luò)[7]。除了執(zhí)行近地空間任務(wù)，立方星星座在輔助深空探測(cè)任務(wù)領(lǐng)域也具有巨大潛能。如Chen 等提出在地月L2點(diǎn)部署4顆立方星為月球背面的用戶(hù)提供定位服務(wù)[8]。

立方星通常作為搭載載荷與其他主要載荷一同發(fā)射。到達(dá)預(yù)定軌道后，由部署系統(tǒng)利用彈簧機(jī)構(gòu)彈射出去[9-10]，發(fā)射和彈射過(guò)程的振動(dòng)會(huì)影響立方星的功能而降低任務(wù)成功率，統(tǒng)計(jì)數(shù)據(jù)顯示立方星失敗任務(wù)中約20%發(fā)生在發(fā)射和部署階段[11]。另外可以用軟袋將立方星打包發(fā)射到天基平臺(tái)，如國(guó)際空間站(ISS)上[12]，等待部署命令。ISS作為立方星部署的中轉(zhuǎn)站，可以降低立方星在發(fā)射過(guò)程中遭受的振動(dòng)級(jí)別，提高任務(wù)成功率。

分布式立方星網(wǎng)絡(luò)立方星數(shù)目多，相較于傳統(tǒng)衛(wèi)星星座呈現(xiàn)出新的特點(diǎn)，立方星的發(fā)射成功率、部署成功率、故障發(fā)生頻率等不確定性都與傳統(tǒng)衛(wèi)星不同。截至2018年5月31日，總的立方星任務(wù)失敗率仍有25%，且立方星在早期運(yùn)行階段易發(fā)生故障的狀況仍未得到明顯改善[11]。傳統(tǒng)的衛(wèi)星星座(如導(dǎo)航衛(wèi)星星座、通信衛(wèi)星星座、遙感衛(wèi)星星座)一般會(huì)有備份星，在部分衛(wèi)星失效后啟動(dòng)，代替故障衛(wèi)星工作[13-16]。對(duì)于分布式立方星網(wǎng)絡(luò)，相較于直接從地面發(fā)射新的立方星，在軌道上儲(chǔ)存?zhèn)溆昧⒎叫强梢蕴岣呦到y(tǒng)對(duì)故障的響應(yīng)性[14]。傳統(tǒng)備份星一般直接被部署于空間軌道，但是立方星更易受到惡劣空間環(huán)境的影響，直接暴露于軌道上會(huì)降低其壽命而影響備份效果，因此考慮使用在軌倉(cāng)庫(kù)儲(chǔ)存?zhèn)浞萘⒎叫?。在軌倉(cāng)庫(kù)與運(yùn)載火箭、服務(wù)航天器、備用立方星組成在軌維護(hù)架構(gòu)。然而由于立方星故障隨機(jī)發(fā)生，用于維護(hù)立方星網(wǎng)絡(luò)的備用立方星的需求不確定，在軌倉(cāng)庫(kù)應(yīng)提前存儲(chǔ)多少備用立方星難以確定。此外備用立方星還應(yīng)適時(shí)從地面補(bǔ)充。因此，一個(gè)重要的問(wèn)題是確定備用立方星的發(fā)射頻率和發(fā)射數(shù)量，從而在保證系統(tǒng)性能的同時(shí)，降低維護(hù)架構(gòu)的運(yùn)行成本。由于故障的隨機(jī)性，這個(gè)問(wèn)題具有極大的挑戰(zhàn)性。

Gu等提出基于航空部件壽命分布的部件采購(gòu)策略以處理部件需求的不確定性，并且假定部件壽命服從正態(tài)分布[17]。Du Jonchay和Ho描述了一種由運(yùn)載火箭、機(jī)器人服務(wù)設(shè)備和軌道倉(cāng)庫(kù)組成的機(jī)器人在軌服務(wù)基礎(chǔ)設(shè)施，這篇文章假定故障發(fā)生服從泊松過(guò)程[14]。然而航天器故障行為并不一定遵守正態(tài)分布或者泊松過(guò)程[18]，真實(shí)的模型應(yīng)該來(lái)自真實(shí)的數(shù)據(jù)。部分學(xué)者對(duì)小衛(wèi)星數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析，結(jié)果表明Weibull分布適合用于建模衛(wèi)星以及衛(wèi)星子系統(tǒng)的故障行為[19]。

本工作旨在為維護(hù)架構(gòu)優(yōu)化提供更實(shí)際的解決方案，考慮隨機(jī)故障行為，優(yōu)化備用立方星的補(bǔ)給時(shí)刻和補(bǔ)給量，以降低維護(hù)架構(gòu)的運(yùn)行成本。提出了分布式立方星網(wǎng)絡(luò)維護(hù)架構(gòu)成本模型。接著收集了立方星故障數(shù)據(jù)，設(shè)計(jì)了最大化擬合優(yōu)度參數(shù)估計(jì)方法，比較了不同分布類(lèi)型對(duì)故障數(shù)據(jù)的擬合效果，最終選擇Weibull分布來(lái)建模立方星的壽命分布[19-21]。最后結(jié)合成本模型和壽命分布模型采用蒙特卡羅仿真和遺傳算法分析了成本模型參數(shù)對(duì)最優(yōu)補(bǔ)給策略的影響。

1 維護(hù)架構(gòu)與成本模型

1.1 維護(hù)架構(gòu)

在合適時(shí)刻發(fā)射備用立方星替換故障立方星可以對(duì)立方星網(wǎng)絡(luò)進(jìn)行維護(hù)。分布式立方星網(wǎng)絡(luò)維護(hù)架構(gòu)如圖1所示，由4部分組成：運(yùn)載火箭、在軌立方星倉(cāng)庫(kù)、服務(wù)航天器和備用立方星(包括地面?zhèn)浞莺驮谲墏浞?。該架構(gòu)進(jìn)行維護(hù)的過(guò)程如下：一定數(shù)量的備用立方星適時(shí)搭載運(yùn)載火箭到達(dá)在軌倉(cāng)庫(kù)，當(dāng)立方星網(wǎng)絡(luò)中有立方星發(fā)生故障時(shí)，由?？坑谠谲墏}(cāng)庫(kù)的服務(wù)航天器攜帶備用立方星/或者立方星采用自身的軌道機(jī)動(dòng)系統(tǒng)轉(zhuǎn)移至故障立方星的軌道，替換掉故障立方星，然后服務(wù)航天器返回在軌倉(cāng)庫(kù)，等待下次任務(wù)。

圖1 立方星網(wǎng)絡(luò)維護(hù)架構(gòu)Fig.1 Maintenance architecture of CubeSat networks

使用在軌倉(cāng)庫(kù)臨時(shí)存儲(chǔ)備用立方星有兩點(diǎn)優(yōu)勢(shì)：① 當(dāng)網(wǎng)絡(luò)損失立方星后，在軌倉(cāng)庫(kù)可以迅速釋放備用立方星，對(duì)故障的響應(yīng)性高；② 在軌倉(cāng)庫(kù)可以提供相對(duì)穩(wěn)定的儲(chǔ)存環(huán)境，降低惡劣的空間環(huán)境對(duì)備用星的影響，提高立方星的壽命。

1.2 成本模型

假定立方星網(wǎng)絡(luò)維護(hù)架構(gòu)的運(yùn)行成本包括3部 分：固定成本、儲(chǔ)存成本和短缺成本[17]。固定成本是指設(shè)計(jì)和制造立方星的成本；儲(chǔ)存成本是用于為在軌倉(cāng)庫(kù)中的備用立方星提供穩(wěn)定的儲(chǔ)存環(huán)境所產(chǎn)生的成本；短缺成本是指空間傳感器網(wǎng)絡(luò)由于未能及時(shí)補(bǔ)充立方星而造成的經(jīng)濟(jì)損失。為了簡(jiǎn)化分析，本文不考慮地面?zhèn)浞荨⑦\(yùn)載火箭發(fā)射成功率和備用立方星軌道轉(zhuǎn)移成功率、備用立方星從在軌倉(cāng)庫(kù)轉(zhuǎn)移到目標(biāo)位置所消耗的時(shí)間。

成本模型用到的記號(hào)如下：s為單位時(shí)間短缺成本；h為單位時(shí)間儲(chǔ)存成本；c為立方星的固定成本；T為規(guī)劃周期；t為備用立方星到達(dá)在軌倉(cāng)庫(kù)時(shí)刻；Q為到達(dá)在軌倉(cāng)庫(kù)的備用立方星數(shù)量；n為故障立方星數(shù)量；j為t時(shí)刻前故障立方星總數(shù)量；tj為第j個(gè)故障發(fā)生時(shí)刻；tk為第k個(gè)故障發(fā)生時(shí)刻。

假定新的立方星在t時(shí)刻到達(dá)在軌倉(cāng)庫(kù)，tj≤t≤tj+1。同時(shí)假定在一個(gè)規(guī)劃周期內(nèi)共有n顆立方星發(fā)生故障。成本模型可以根據(jù)Q與j的大小分成5種情況(Q

圖2 成本計(jì)算示意圖(j+1

假設(shè)Q顆備用立方星在t時(shí)刻到達(dá)在軌倉(cāng)庫(kù)。其中j顆立即用于替換故障立方星。剩余的Q-j顆備用立方星從t時(shí)刻起被儲(chǔ)存在在軌倉(cāng)庫(kù)中。每出現(xiàn)一個(gè)新的故障事件就消耗1顆備用立方星，直到第Q個(gè)故障發(fā)生并消耗掉最后一顆備用立方星，相應(yīng)的儲(chǔ)存成本為：h(tk-t),k=j+1,j+2, …,Q；短缺成本為：s(t-tk),k=1, 2, …,j。在這之后發(fā)生故障的n-Q顆立方星因?yàn)槿鄙賯溆昧⒎叫菬o(wú)法被替換，相應(yīng)的短缺成本：s(T-tk),k=Q+1,Q+2, …,n。此種情況的總成本計(jì)算公式為

(1)

另外4種情況總成本分析方法與第4種情況類(lèi)似，不再贅述。總的成本模型可統(tǒng)一寫(xiě)為

M=

(2)

2 立方星可靠性模型

可靠性分析方法包括非參數(shù)分析方法和參數(shù)分析方法。非參數(shù)分析不需要假定壽命服從某一分布，適合于不同類(lèi)型的故障數(shù)據(jù)，也可以用于評(píng)估參數(shù)分析的結(jié)果。參數(shù)分析首先假定壽命服從某一分布，然后估計(jì)該分布的參數(shù)。

由于立方星網(wǎng)絡(luò)維護(hù)架構(gòu)成本分析需要考慮立方星的壽命分布，因此本文選擇參數(shù)分析方法，并用非參數(shù)分析方法評(píng)估參數(shù)估計(jì)結(jié)果。

2.1 立方星壽命數(shù)據(jù)

收集了立方星壽命數(shù)據(jù)，包括立方星的發(fā)射時(shí)刻、故障發(fā)生時(shí)刻或項(xiàng)目結(jié)束時(shí)刻等[22]。如果立方星一直正常工作到任務(wù)完成，即沒(méi)有故障發(fā)生，這樣的數(shù)據(jù)稱(chēng)為右刪失數(shù)據(jù)。本文獲得的數(shù)據(jù)包含右刪失數(shù)據(jù)。

由于一些立方星沒(méi)有公開(kāi)任務(wù)運(yùn)行狀態(tài)，最終共收集到111個(gè)有效數(shù)據(jù)，包括84個(gè)故障數(shù)據(jù)，27個(gè)刪失數(shù)據(jù)，數(shù)據(jù)時(shí)間跨度為2000—2018年，足以進(jìn)行統(tǒng)計(jì)分析并用于擬合參數(shù)模型。用于可靠性分析的數(shù)據(jù)示例如表1所示,其中刪失標(biāo)記0表示故障，1表示刪失。

表1 立方星壽命數(shù)據(jù)示例Table 1 Lifetime data of CubeSats

2.2 非參數(shù)分析

Kaplan-Meier(KM)估計(jì)方法可以從包含刪失數(shù)據(jù)的故障數(shù)據(jù)中推導(dǎo)出經(jīng)驗(yàn)可靠性函數(shù)?？煽啃院瘮?shù)又稱(chēng)為生存函數(shù)，對(duì)于立方星而言，其可靠性函數(shù)表示立方星故障發(fā)生在t時(shí)刻之后的概率，計(jì)算公式為

(3)

式中：ti為故障時(shí)刻；ai為ti時(shí)刻故障立方星總數(shù)量；ni為ti時(shí)刻前正常運(yùn)行的立方星數(shù)量，等于立方星總數(shù)量減去ti時(shí)刻前完成在軌任務(wù)的立方星數(shù)量(刪失數(shù)量)和故障立方星數(shù)量。

該可靠性函數(shù)是階梯狀函數(shù)，在觀察到的故障時(shí)刻處不連續(xù)，在刪失時(shí)刻函數(shù)值不變。

2.3 參數(shù)分析

參數(shù)分析首先需要選擇合適的模型，然后采用參數(shù)估計(jì)算法估計(jì)模型參數(shù)。對(duì)數(shù)正態(tài)分布適合用于建模故障率隨時(shí)間先增加后降低的對(duì)象的壽命分布，Weibull分布可以用于建模故障率隨時(shí)間降低、增加、或不變的對(duì)象的壽命分布[21]。采用最大化擬合優(yōu)度估計(jì)上述分布的參數(shù)。

1) 最大化擬合優(yōu)度參數(shù)估計(jì)方法

統(tǒng)計(jì)模型的擬合優(yōu)度表明觀測(cè)數(shù)據(jù)與模型預(yù)測(cè)值的差異，擬合優(yōu)度函數(shù)值越大表明模型擬合效果越好[19]，因此可以對(duì)擬合優(yōu)度函數(shù)進(jìn)行優(yōu)化，確定模型參數(shù)值。

對(duì)實(shí)際觀測(cè)故障數(shù)據(jù)進(jìn)行非參數(shù)分析得到的可靠性記為R(ti)，由參數(shù)模型得到的可靠性記為R′(ti,α)，α為模型參數(shù)，則擬合優(yōu)度R2計(jì)算公式為

(4)

(5)

(6)

(7)

最大化擬合優(yōu)度參數(shù)估計(jì)建模為

α*=arg maxR2(α)

(8)

采用遺傳算法對(duì)式(8)進(jìn)行求解。

2) 對(duì)數(shù)正態(tài)分布

假定立方星壽命x服從對(duì)數(shù)正態(tài)分布，即lnx～N(μ,σ2)，則立方星的可靠性函數(shù)為

(9)

式中：erf為誤差函數(shù)。

3) Weibull分布

Weibull分布廣泛應(yīng)用于可靠性分析中，能靈活描述不同故障行為，其PDF和CDF分別為

f(x;β,η)=

(10)

F(x;β,η)=1-exp(-(x/η)β)

(11)

式中：β為形狀參數(shù)；η為尺度參數(shù)。

假定立方星壽命x服從Weibull分布，則立方星的可靠性函數(shù)為

(12)

3 仿真方法

由1.2節(jié)分析可知，分布式立方星網(wǎng)絡(luò)維護(hù)架構(gòu)運(yùn)行成本與立方星的故障時(shí)間序列、備用立方星的發(fā)射時(shí)刻和補(bǔ)給數(shù)量有關(guān)。因此分析成本模型首先需要生成故障時(shí)間序列。以?xún)赡?730天)為一個(gè)規(guī)劃周期T，假定分布式網(wǎng)絡(luò)由Qm顆相同的立方星構(gòu)成，用Weibull壽命分布模型生成Qm個(gè)隨機(jī)數(shù)，用以代表這Qm顆立方星的壽命，刪除大于730的數(shù)，將剩余的隨機(jī)數(shù)從小到大排序， 即為一個(gè)規(guī)劃周期內(nèi)的故障時(shí)間序列，記作{t}。

將故障時(shí)間序列，備用立方星補(bǔ)給時(shí)刻、補(bǔ)給數(shù)量代入式(2)，可得到維護(hù)架構(gòu)的運(yùn)行成本。但是由于故障時(shí)刻是隨機(jī)的，單個(gè)故障時(shí)間序列對(duì)應(yīng)的總成本隨著補(bǔ)給時(shí)刻延后而波動(dòng)變化，因此采用蒙特卡羅仿真,求N個(gè)故障時(shí)間序列對(duì)應(yīng)成本的平均值對(duì)成本進(jìn)行分析。

為了求解最小總成本對(duì)應(yīng)的最優(yōu)補(bǔ)給策略，給出優(yōu)化模型。采用遺傳算法求解可得到最優(yōu)補(bǔ)給時(shí)刻t*和最優(yōu)補(bǔ)給數(shù)量Q*，即

(13)

式中：{t}l表示第l個(gè)故障時(shí)間序列。

4 結(jié)果分析

4.1 立方星壽命分布模型

采用最大化擬合優(yōu)度方法分別估計(jì)了對(duì)數(shù)正態(tài)分布和Weibull分布的參數(shù)，同時(shí)對(duì)比文獻(xiàn)[19] 中結(jié)合馬爾可夫鏈-蒙特卡羅(MCMC)仿真的貝葉斯方法，對(duì)比了這2種參數(shù)估計(jì)算法對(duì)于立方星故障數(shù)據(jù)的擬合效果。

非參數(shù)分析方法Kaplan-Meier估計(jì)得到立方星的可靠性如圖3中的階梯狀曲線所示，可以看出立方星可靠性在初期迅速降低，隨后緩慢降低。

圖3 不同參數(shù)估計(jì)方法對(duì)比Fig.3 Contrast of different parameter estimation methods

參數(shù)估計(jì)結(jié)果列于表2。擬合優(yōu)度R2越接近于1說(shuō)明擬合優(yōu)度越好。對(duì)于Weibull分布，最大化擬合優(yōu)度參數(shù)估計(jì)方法所得結(jié)果的R2比MCMC方法所得結(jié)果的R2更大，說(shuō)明前者的擬合效果更好。如圖3所示，在早期，最大化擬合優(yōu)度參數(shù)估計(jì)方法獲得的可靠性函數(shù)曲線顯然更靠近KM非參數(shù)分析得到的曲線，即該方法能更好地?cái)M合立方星的早期失效率高的特征，也說(shuō)明了最大化擬合優(yōu)度參數(shù)估計(jì)方法擬合精度更高。Weibull分布擬合結(jié)果的R2大于對(duì)數(shù)正態(tài)分布擬合結(jié)果的R2，說(shuō)明Weibull分布更適合用于立方星壽命分布的擬合。

表2 參數(shù)估計(jì)結(jié)果Table 2 Parameter estimation results

下面分析最大擬合優(yōu)度方法Weibull分布擬合結(jié)果。β=0.476 8，β<1，表明立方星具有早期失效率高的特點(diǎn)。η=291.900 4，表示立方星的特征壽命為291.900 4天。當(dāng)t=η時(shí)，代入Weibull分布的累計(jì)分布函數(shù)(CFD)，得到F(t)=0.632，表示當(dāng)立方星網(wǎng)絡(luò)運(yùn)行到第291.900 4天時(shí)，可能有63.2%的立方星發(fā)生故障。

早期立方星主要用于技術(shù)驗(yàn)證或科學(xué)實(shí)驗(yàn)，因此收集的故障數(shù)據(jù)前一部分以驗(yàn)證或?qū)嶒?yàn)立方星為主。隨著技術(shù)發(fā)展，執(zhí)行商業(yè)任務(wù)的業(yè)務(wù)立方星出現(xiàn)[6]，對(duì)可靠性要求更高。有必要將數(shù)據(jù)按發(fā)射時(shí)間等分成早期和近期兩部分，分別進(jìn)行參數(shù)估計(jì)，得到兩類(lèi)可靠性，結(jié)果如表3和圖4所示。由表3可知，近期發(fā)射的立方星的特征壽命約為415 天，明顯高于早期發(fā)射的立方星，但是形狀參數(shù)β差別不大，都呈現(xiàn)出早期失效的特點(diǎn)；由圖4 可知近期發(fā)射的立方星的可靠性也高于早期的立方星。該對(duì)比在一定程度上說(shuō)明了驗(yàn)證實(shí)驗(yàn)立方星與業(yè)務(wù)立方星的差異性，因此在優(yōu)化業(yè)務(wù)立方星網(wǎng)絡(luò)補(bǔ)給策略時(shí)，需要考慮這種差異性。

表3 早期、近期數(shù)據(jù)參數(shù)估計(jì)結(jié)果

圖4 早期、近期立方星對(duì)比Fig.4 Contrast of early and recent CubeSats

4.2 參數(shù)對(duì)總成本的影響

分布式立方星網(wǎng)絡(luò)維護(hù)架構(gòu)成本模型中有3個(gè) 未知參數(shù)：c、h、s，通過(guò)調(diào)節(jié)這些參數(shù)的大小，可以涵蓋大多數(shù)的情形。假定立方星網(wǎng)絡(luò)由100顆 立方星組成，分析這3個(gè)參數(shù)對(duì)總成本的影響。

1) 單一成本情況

只考慮固定成本，即c:h:s=1:0:0時(shí)，成本模型簡(jiǎn)化為：M=cQ，總成本與補(bǔ)給時(shí)刻無(wú)關(guān)，僅與補(bǔ)給數(shù)量成正比。

只考慮儲(chǔ)存成本，即c:h:s=0:1:0時(shí)，t、Q對(duì)總成本的影響如圖5所示，總成本隨著補(bǔ)給時(shí)刻延后而降低，隨著補(bǔ)給數(shù)量增加而增加。值得注意的是，當(dāng)補(bǔ)給數(shù)量小于60，補(bǔ)給時(shí)刻延后至第300天之后時(shí)，儲(chǔ)存成本降為零，原因如下：由Weibull分布擬合結(jié)果可知立方星的特征壽命約292天，表示當(dāng)立方星網(wǎng)絡(luò)運(yùn)行到第292天時(shí)，可能有63.2%的立方星發(fā)生故障；立方星網(wǎng)絡(luò)由100顆立方星組成，當(dāng)運(yùn)行到292天時(shí)，可能有63顆 立方星發(fā)生故障，因此在第300天之后補(bǔ)給立方星數(shù)量小于60時(shí)，這些補(bǔ)給的立方星將全部立即用于替代故障立方星，因此總儲(chǔ)存成本為零。

圖5 儲(chǔ)存成本對(duì)總成本的影響(c∶h∶s=0∶1∶0)Fig.5 Effect of holding costs on total costs (c∶h∶s=0∶1∶0)

只考慮短缺成本，即c∶h∶s=0∶0∶1時(shí)，t、Q對(duì)總成本的影響如圖6所示，總成本隨著補(bǔ)給時(shí)刻延后而增加，隨著補(bǔ)給數(shù)量增加而降低。當(dāng)補(bǔ)給數(shù)量超過(guò)80時(shí)，改變補(bǔ)給數(shù)量，總成本隨補(bǔ)給時(shí)間的變化曲線幾乎重合(如圖6放大圖所示，Q=90與Q=100的曲線幾乎是重合的)，這是由于在規(guī)劃周期內(nèi)，總故障立方星數(shù)目小于80的概率非常大，補(bǔ)給立方星數(shù)目大于80時(shí)，所有故障立方星都能被替換，因此總的累積短缺時(shí)間是相同的，又由于沒(méi)有考慮短缺成本和固定成本，所以總成本相同。

圖6 短缺成本對(duì)總成本的影響(c∶h∶s=0∶0∶1)Fig.6 Effect of shortage costs on total costs (c∶h∶s=0∶0∶1)

2) 全部成本情況

同時(shí)考慮固定成本、儲(chǔ)存成本、短缺成本，即取c∶h∶s=1∶1∶1時(shí)，t、Q對(duì)總成本的影響如圖7所示。補(bǔ)給數(shù)量較少(大約是Q≤60)與補(bǔ)給數(shù)量較多時(shí)對(duì)總成本的影響存在較大差異：當(dāng)補(bǔ)給數(shù)量較少時(shí)，總成本隨著補(bǔ)給數(shù)量增加而減少，隨著補(bǔ)給時(shí)刻延后而增加；當(dāng)補(bǔ)給數(shù)量較多時(shí)，總成本隨著補(bǔ)給數(shù)量增加而增加，隨著補(bǔ)給時(shí)刻延后先降低后增加。這一差異反應(yīng)了短缺成本和儲(chǔ)存成本對(duì)于總成本不同影響效果，這2種成本分別對(duì)應(yīng)備件短缺和備件冗余2種情況。

圖7 3種成本參數(shù)對(duì)總成本的影響(c∶h∶s=1∶1∶1)Fig.7 Effect of three cost parameters on total costs (c∶h∶s=1∶1∶1)

由上述分析可得出下述結(jié)論：只考慮儲(chǔ)存成本，補(bǔ)給時(shí)刻取在立方星的特征壽命292天之后，補(bǔ)給數(shù)量小于網(wǎng)絡(luò)中立方星總量的63.2%時(shí)，總儲(chǔ)存成本為0；只考慮短缺成本，當(dāng)補(bǔ)給數(shù)量Q>80時(shí)，總成本不再隨著補(bǔ)給數(shù)量增加而變化。當(dāng)補(bǔ)給數(shù)量較少時(shí)，短缺成本對(duì)總成本的影響占主導(dǎo)；當(dāng)補(bǔ)給數(shù)量較多時(shí)，儲(chǔ)存成本對(duì)總成本的影響占主導(dǎo)。

4.3 最優(yōu)補(bǔ)給時(shí)刻與最優(yōu)補(bǔ)給數(shù)量

分析不同成本比例對(duì)于最優(yōu)補(bǔ)給策略的影響，這里的不同比例指短缺成本與儲(chǔ)存成本的不同比例，由于立方星成本正在逐步降低，因此暫不考慮固定成本，最優(yōu)補(bǔ)給策略指最優(yōu)補(bǔ)給時(shí)刻t*、最優(yōu)補(bǔ)給數(shù)量Q*；并且對(duì)比了立方星網(wǎng)絡(luò)不同規(guī)模的仿真結(jié)果，考慮網(wǎng)絡(luò)中立方星數(shù)量分別為50、100、150這3種情況。

如圖8所示，隨著s/h比例增大，最優(yōu)補(bǔ)給時(shí)刻均提前，但最優(yōu)補(bǔ)給時(shí)刻與立方星網(wǎng)絡(luò)規(guī)模無(wú)關(guān)，都沒(méi)有超過(guò)150天，均在規(guī)劃早期取得，這與立方星的故障行為有關(guān)，由4.1節(jié)分析可知立方星有早期失效特點(diǎn)，即立方星在軌運(yùn)行早期發(fā)生故障的概率較高，因此在軌運(yùn)行早期發(fā)生故障的立方星數(shù)目比在軌運(yùn)行中期和后期的故障數(shù)目要多，所以應(yīng)該在早期及時(shí)補(bǔ)給備用立方星。

圖8 s/h對(duì)t*的影響(不考慮固定成本)Fig.8 Influence of s/h on t* (without considering fixed costs)

如圖9所示，考慮不同比例的s/h，立方星網(wǎng)絡(luò)最優(yōu)補(bǔ)給數(shù)量大小與規(guī)模正相關(guān)且變化趨勢(shì)相同，隨著短缺成本所占比重增加，最優(yōu)補(bǔ)給數(shù)量也隨之增加，這是由于增大補(bǔ)給數(shù)量能夠緩解備件短缺情況。

圖9 s/h對(duì)Q*的影響(不考慮固定成本)Fig.9 Influence of s/h on Q* (without considering fixed costs)

4.4 不同可靠性模型對(duì)應(yīng)的成本比較

將由故障數(shù)據(jù)和最大化擬合優(yōu)度方法得到的可靠性模型和文獻(xiàn)[19]中的與立方星質(zhì)量相接近的Nano衛(wèi)星、Micro衛(wèi)星的可靠性模型(如表4所示)分別代入到維護(hù)架構(gòu)補(bǔ)給策略?xún)?yōu)化模型中，得到3種補(bǔ)給策略，依次記為立方星模型Nano模型、Micro模型；再將3種補(bǔ)給策略和基于數(shù)據(jù)的可靠性模型代入維護(hù)架構(gòu)，計(jì)算總成本。結(jié)果如表5所示。

表4 文獻(xiàn)[19]模型參數(shù)Table 4 Model parameters of Ref.[19]

表5 不同補(bǔ)給策略比較Table 5 Contrast of different supply strategies

Nano模型與立方星模型的結(jié)果接近，總成本增加了12.23%；而Micro模型相對(duì)于立方星模型結(jié)果總成本增加了181.58%。這說(shuō)明引入適當(dāng)?shù)膶?shí)際故障數(shù)據(jù)對(duì)維護(hù)架構(gòu)進(jìn)行優(yōu)化的必要性。

4.5 早期、近期立方星補(bǔ)給策略對(duì)比

由于早期和近期發(fā)射的立方星均具有早期失效的特點(diǎn)，都應(yīng)該在規(guī)劃早期補(bǔ)給備用立方星，因此只對(duì)補(bǔ)給數(shù)量進(jìn)行對(duì)比分析。用早期和近期的立方星可靠性模型，分別求解由100顆立方星組成的網(wǎng)絡(luò)的最優(yōu)補(bǔ)給數(shù)量。如圖10所示，對(duì)于不同成本比例，近期立方星的補(bǔ)給數(shù)量均小于早期立方星的補(bǔ)給數(shù)量。由4.1節(jié)分析可知，近期發(fā)射的立方星的可靠性更高，因此在規(guī)劃周期內(nèi)立方星網(wǎng)絡(luò)中發(fā)生故障的立方星數(shù)量相對(duì)更少，所以需要的補(bǔ)給量也少。這一結(jié)論說(shuō)明在優(yōu)化立方星網(wǎng)絡(luò)的維護(hù)策略時(shí)，應(yīng)該考慮技術(shù)發(fā)展導(dǎo)致的可靠性模型的變化。

圖10 早期、近期立方星Q*對(duì)比Fig.10 Q* comparison between early and recent CubeSats

5 結(jié) 論

1) 立方星壽命分布的參數(shù)模型可采用Weibull分布或?qū)?shù)正態(tài)分布，進(jìn)而采用最大化擬合優(yōu)度參數(shù)估計(jì)方法進(jìn)行參數(shù)估計(jì)，擬合數(shù)據(jù)表明Weibull模型精度比對(duì)數(shù)正態(tài)分布模型稍高。

2) 立方星具有早期失效率高的特點(diǎn)，由于這一特點(diǎn)，備用立方星最優(yōu)補(bǔ)給時(shí)刻建議為規(guī)劃周期早期。目前立方星總的特征壽命大約為292天，近期發(fā)射的立方星的可靠性高于早期發(fā)射的立方星。

3) 基于本文收集的立方星壽命數(shù)據(jù)，若在軌服務(wù)成本以短缺成本為主，那么建議備用立方星補(bǔ)給時(shí)間不晚于第150天，補(bǔ)給數(shù)量與立方星網(wǎng)絡(luò)規(guī)模有關(guān)。若以?xún)?chǔ)存成本為主，建議備用立方星補(bǔ)給數(shù)量約等于立方星網(wǎng)絡(luò)中立方星總數(shù)量的60%，補(bǔ)給時(shí)刻不晚于第300天。補(bǔ)給時(shí)刻隨短缺成本/儲(chǔ)存成本增加而提前，與立方星網(wǎng)絡(luò)規(guī)模無(wú)明顯相關(guān)關(guān)系；補(bǔ)給數(shù)量隨短缺成本/儲(chǔ)存成本增加而增加。

4) 近期發(fā)射立方星的最優(yōu)補(bǔ)給數(shù)量比早期立方星的少，對(duì)立方星網(wǎng)絡(luò)的補(bǔ)給策略進(jìn)行優(yōu)化時(shí)，不僅需要選擇合適的可靠性模型，還需要考慮由于技術(shù)發(fā)展而導(dǎo)致的立方星可靠性模型的變化。