基于蒙特卡羅仿真的全自動運行系統(tǒng)彈性評估方法研究

郭子薇，燕 飛

（北京交通大學電子信息工程學院，北京 100044）

隨著中國城軌的高速發(fā)展，全自動運行系統(tǒng)(Fully Automatic Operation，F(xiàn)AO)以其安全、高效的特點引發(fā)了新的發(fā)展需求。當前，北京、上海等城市日均客流已超千萬，運營負荷巨大，安全運營面臨嚴峻挑戰(zhàn)，采用全自動運行技術是有效解決途徑和必然選擇[1]。

彈性的概念經(jīng)Holling[2]引入生態(tài)領域后，經(jīng)過幾十年的發(fā)展，彈性研究已擴展到組織管理、經(jīng)濟、基礎設施建設以及電網(wǎng)等諸多領域，以評估目標系統(tǒng)對威脅其正常運行的中斷作出的反應能力。為進一步提升全自動運行系統(tǒng)抵抗擾動以及從非正常態(tài)中快速恢復的能力，研究人員必須能夠定量評估在不同的設計或維護管理方案下系統(tǒng)的彈性水平，以便做出最佳決策。因此，提出圍繞全自動運行系統(tǒng)的彈性評估方法是對全自動運行系統(tǒng)的彈性工程實踐的基礎環(huán)節(jié)，同時對彈性在工程領域的應用也起到了重要作用。目前，彈性的度量方法并沒有標準化，彈性度量仍是一個具有挑戰(zhàn)的問題，從彈性的定義出發(fā)，目前已有的系統(tǒng)彈性度量通常圍繞系統(tǒng)性能降級程度與快速恢復性展開，一般可分為確定型度量和概率型度量兩類[3]。

確定型度量主要是對系統(tǒng)在某次確定的擾動情形下從性能降級到性能恢復這個過程的彈性定量評估，由于系統(tǒng)所經(jīng)受的擾動、產(chǎn)生的性能降級和恢復時間都具有隨機性，基于單次擾動下的系統(tǒng)彈性能力的確定型度量方式本身具備了隨機性，因此較多的研究人員開始研究概率型系統(tǒng)彈性能力度量方式。Jin等[3]對道路運輸網(wǎng)絡的各條線路受擾動情況以及降級水平都以概率的形式給出，以故障后與故障前的網(wǎng)絡最大數(shù)據(jù)流的比值作為彈性度量值，經(jīng)過多次的隨機抽樣，統(tǒng)計出了系統(tǒng)的平均彈性水平。

以北京全自動運行線路燕房線為例，采用蒙特卡羅仿真方法模擬系統(tǒng)受擾動以及性能下降的隨機性，以故障前后系統(tǒng)數(shù)據(jù)流權值的損失比量化系統(tǒng)性能下降水平，經(jīng)過多次仿真統(tǒng)計出系統(tǒng)的平均彈性水平，為全自動運行系統(tǒng)的彈性評估提供一種更為可靠、能規(guī)避擾動隨機性問題的概率型彈性評估方法。

1 彈性計算方法

彈性通常被定義為系統(tǒng)抵御擾動并快速恢復的能力，一般采用“彈性三角”曲線來表征系統(tǒng)性能的變化過程，同時將由于遭遇擾動導致曲線下降而缺失的面積比作為衡量系統(tǒng)彈性的指標。但對于全自動運行系統(tǒng)這樣一個龐大而復雜的系統(tǒng)來說，很難將系統(tǒng)各設備故障前后的性能變化用曲線進行準確描述，因此采用Zobel[4]提出的基于幾何數(shù)學模型的彈性度量方法，在已知故障導致的性能下降程度，以及故障恢復時間后，該方法假設故障以一個恒定的速度進行恢復，降低了彈性的量化難度。歸一化剩余性能值以系統(tǒng)故障后和故障前的系統(tǒng)數(shù)據(jù)流權值的比值來衡量。

為了便于對全自動運行系統(tǒng)的彈性展開定量計算，做出以下假設。

1) 擾動：每次只考慮系統(tǒng)遭受一次擾動，每次擾動只影響系統(tǒng)的一個組成部分，中斷發(fā)生在第i個設備的概率為qi；

2) 性能下降：設備的性能下降水平遵循離散分布(對于初始數(shù)據(jù)量為Di的設備i，其性能下降后的數(shù)據(jù)量可能降為Di,1，Di,2，…，Di,k，每個值對應的概率為pi,k)；

3) 恢復：設備的恢復時間服從對數(shù)正態(tài)分布(對于設備i，其恢復時間為

假設1是為了簡化研究的問題，同時也是彈性分析中常用的假設，假設2是對設備遭受擾動后的性能變化進行了狀態(tài)抽象，降低了運算的難度，而假設 3將設備的恢復時間采用對數(shù)正態(tài)分布進行描述，是因為通過查閱相關文獻[5-10]發(fā)現(xiàn)，對數(shù)正態(tài)分布式系統(tǒng)修復時間是最廣泛使用的分布，且參考文獻中有交通事故的持續(xù)時間服從對數(shù)正態(tài)分布的分析，因此，本文采取對數(shù)正態(tài)分布對設備的恢復時間進行描述。

Zobel[4]在其彈性度量方法中考慮到由于在系統(tǒng)設計過程中無法確定系統(tǒng)遭遇擾動后的性能變化過程，因此其假設系統(tǒng)在t0時刻遭遇擾動后，性能下降X，隨后以恒定的恢復速率進行恢復，記恢復時間為T，如圖1所示。T*是T的可能值集合中的嚴格上界，系統(tǒng)遭受的任何擾動行為均可在T*時間內(nèi)恢復?？梢钥闯?，三角形的面積就是系統(tǒng)在某種特定擾動后的性能損失量，因此，用歸一化的性能損失量表達彈性可用下列公式表示：

圖1 基于恒定恢復速度的彈性度量Figure 1 Resilience measurement based on constant recovery speed

全自動運行系統(tǒng)設備組成較為復雜，找到統(tǒng)一的衡量各設備狀態(tài)的指標較困難，同時也難以保證準確，因此用狀態(tài)虛擬值來量化表示設備性能，設備的輸出數(shù)據(jù)流的權值與狀態(tài)虛擬值相同。設備正常工作時，狀態(tài)虛擬值最大，當設備遭受擾動后虛擬值降低，統(tǒng)計設備間的數(shù)據(jù)流權值總和，用系統(tǒng)故障前后總的數(shù)據(jù)流權值的比值來表征系統(tǒng)性能的下降水平。因此考慮數(shù)據(jù)流權值的設備i的彈性可計算為：

考慮設備降級服從離散分布，恢復時間服從對數(shù)正態(tài)分布，因此，設備i的基于數(shù)據(jù)流權值的彈性期望可以計算為：

2 全自動運行系統(tǒng)參數(shù)確定

2.1 FAO系統(tǒng)原理及架構

全自動運行系統(tǒng)較傳統(tǒng)的基于通信的列車控制系統(tǒng)(CBTC)具有更加完善的自動控制功能，以行車為核心，信號、車輛、綜合監(jiān)控、通信、站臺門等多系統(tǒng)深度互聯(lián)，信號系統(tǒng)與綜合監(jiān)控系統(tǒng)深度集成(行車綜合自動化系統(tǒng))，提升了城市軌道交通運行系統(tǒng)的整體自動化水平。

全自動運行系統(tǒng)架構核心包括以下部分：

1) 控制中心核心設備：列車自動監(jiān)控系統(tǒng)(ATS)、中心電力調(diào)度系統(tǒng)(SCADA)以及綜合監(jiān)控系統(tǒng)(ISCS)高度集成，形成行車綜合自動化系統(tǒng)(TIAS)，在運營組織方面，加強緊急情況處理，對信號、通信、綜合監(jiān)控等系統(tǒng)深度集成與整合，保證系統(tǒng)內(nèi)部的協(xié)調(diào)工作。

2) 車站控制設備：以車站TIAS系統(tǒng)為核心，集成控制車站內(nèi)與乘客乘降相關的設備，如乘客信息系統(tǒng)(PIS)等，響應中心命令，負責站內(nèi)引導與應急處理。

3) 軌旁設備：包括列車自動防護系統(tǒng)(ATP)、聯(lián)鎖系統(tǒng)(CI)、輔助定位LEU系統(tǒng)等，可根據(jù)運營線路實際需要增加軌旁障礙物檢測設備。

4) 車載設備：以車載列控系統(tǒng)(ATP/ATO)為核心，配備牽引制定、輔助測速等車輛設備，在車頭加裝視頻監(jiān)控設備(CCTV)，增加車輛信息管理系統(tǒng)(TCMS)，作為傳輸行車信息的專用系統(tǒng)。

為了降低計算量，將系統(tǒng)架構進行了適當?shù)暮喕?，如圖 2所示，同時對重要的設備之間的交互關系進行了梳理。其中，箭頭的方向表示數(shù)據(jù)流的信息傳遞方向。數(shù)據(jù)流的權值用輸出該數(shù)據(jù)的設備狀態(tài)虛擬值來表示。

圖2 全自動運行系統(tǒng)架構Figure 2 Fully automatic operation system architecture

2.2 設備參數(shù)

全自動運行系統(tǒng)的各設備的故障模式較多，同時很難用單一的指標準確描述設備的性能，因此為了運算的統(tǒng)一與方便，根據(jù)設備故障后果嚴重程度不同進行性能降級等級劃分。對每一等級設置“狀態(tài)虛擬值”來定量表征設備性能，當設備正常工作時，該虛擬值最大，性能降級越嚴重，其虛擬值越低。

在城市軌道交通列控系統(tǒng)中，用SIL等級表示設備應達到的相應安全水平，SIL等級越高，說明該設備的正常運行對保障系統(tǒng)安全越重要。參考《城市軌道交通全自動運行系統(tǒng)運營需求導則》中對全自動運行系統(tǒng)各設備的SIL等級要求，SIL等級越高的設備，狀態(tài)虛擬值設置越高，以此表示其設備故障對系統(tǒng)的安全影響越大。SIL4級的設備，未故障時狀態(tài)虛擬值為8，下降性能共設3級，相對狀態(tài)虛擬值依次降為4、2、0；SIL2級的設備正常工作時狀態(tài)虛擬值為4，各性能降級等級對應狀態(tài)虛擬值依次降為 2、1、0；對沒有SIL等級要求的設備，均設置為正常運行時狀態(tài)虛擬值為2，性能降級等級對應狀態(tài)虛擬值依次為1、0.5、0。

采用上述設備狀態(tài)虛擬值的設置規(guī)則，為圖2全自動運行架構中所涉及的主要設備確定狀態(tài)虛擬值。

調(diào)查了北京燕房線2017年12月到2020年10月的信號系統(tǒng)和綜合監(jiān)控系統(tǒng)的故障數(shù)據(jù)記錄，對CI、ATO、ATP、ATS、ZC、ISCS、轉轍機以及軌道電路這8種重要設備的故障數(shù)據(jù)進行了統(tǒng)計，對故障記錄模糊的數(shù)據(jù)進行了剔除，其中發(fā)生在計算機聯(lián)鎖設備(CI)的故障共有 67條，發(fā)生在自動駕駛(ATO)子系統(tǒng)的故障共有59條，發(fā)生在自動防護(ATP)子系統(tǒng)的故障共有57條，發(fā)生在綜合監(jiān)控系統(tǒng)(ISCS)的故障記錄共有367條，而自動監(jiān)控系統(tǒng)(ATS)的故障記錄最多，共有420條，轉轍機和軌道電路分別為47條和19條，發(fā)生在區(qū)域控制器(ZC)子系統(tǒng)的故障記錄相對較少，只有12條。由于燕房線開通時間較短，所統(tǒng)計到的故障數(shù)據(jù)有限，在之后的研究中，可以持續(xù)監(jiān)測燕房線的故障數(shù)據(jù)，不斷糾正研究結果。根據(jù)故障記錄描述，統(tǒng)計了每條故障從發(fā)生到故障恢復的時間，對各設備對應的故障恢復時間進行了對數(shù)正態(tài)分布參數(shù)計算，同時統(tǒng)計了各設備的故障頻率，將設備的故障頻率作為設備故障的概率，另外根據(jù)設備不同，性能下降等級對應不同的狀態(tài)虛擬值，以各性能等級出現(xiàn)的頻率確定了各虛擬值對應的發(fā)生概率，具體數(shù)據(jù)見表1。

表1 設備相關參數(shù)Table 1 Device-related parameters

3 基于蒙特卡羅的系統(tǒng)彈性評估

蒙特卡羅方法是一種基于仿真統(tǒng)計的計算方法，其核心思想是建立系統(tǒng)中概率過程的仿真模型，然后使用多次試驗的方法，計算得到系統(tǒng)特征。針對系統(tǒng)彈性這一指標，由于研究對象的復雜性，以及系統(tǒng)可能遭受的擾動、性能下降和恢復過程的隨機性，蒙特卡羅仿真可以作為一種系統(tǒng)彈性評估的有效方法。

采取蒙特卡羅仿真方法來探究基于全自動運行系統(tǒng)狀態(tài)虛擬值的系統(tǒng)彈性，模擬全自動運行系統(tǒng)設備故障、設備故障程度以及恢復時間點的隨機性，通過105次仿真運算，統(tǒng)計出系統(tǒng)彈性的平均值，具體仿真步驟如下。

1) 將設備的輸出數(shù)據(jù)流權值設置為該設備的狀態(tài)虛擬值，計算基于全自動運行系統(tǒng)架構中各設備間的交互關系的數(shù)據(jù)流權值總和，各設備均正常工作，數(shù)據(jù)流權值取最大值Ds。

2) 根據(jù)各設備可能遭受擾動的概率qi，抽樣確定受擾動的設備j。

3) 根據(jù)設備故障后性能下降和恢復時間所服從的分布，抽樣確定設備j受擾動后系統(tǒng)下降信號值和恢復時間tj。

4) 在K次擾動下計算系統(tǒng)彈性

5) 考慮到擾動、性能下降和恢復時間的隨機性，重復2)到4)步，直到迭代次數(shù)M。

根據(jù)燕房線2017年12月到2020年10月的信號及綜合監(jiān)控系統(tǒng)問題跟蹤記錄的故障數(shù)據(jù)，一共獲取了ISCS、ATS、ATO等8個設備的故障統(tǒng)計數(shù)據(jù)，將這8個設備的故障數(shù)據(jù)進行系統(tǒng)的彈性計算，若能獲得更詳細更全面的設備故障數(shù)據(jù)，可擴展計算的范圍。以圖2全自動運行系統(tǒng)架構為基礎，將各設備的輸出數(shù)據(jù)流權值設置為對應的狀態(tài)虛擬值，經(jīng)計算，系統(tǒng)在正常運行時系統(tǒng)總的數(shù)據(jù)流權值為 182，設備故障的邊界時間取24 h，采用Matlab進行蒙特卡羅仿真，得到 105次仿真的全自動運行系統(tǒng)彈性經(jīng)驗估計值：。該系統(tǒng)彈性的概率密度函數(shù)如圖3所示，可以看到全自動運行系統(tǒng)彈性主要集中在0.98之后的區(qū)域，這表明在CI、ATP等論文中，討論的8個主要信號系統(tǒng)設備中的大多數(shù)擾動系統(tǒng)表現(xiàn)的彈性都非常高，在某些特定的擾動下系統(tǒng)也可能體現(xiàn)出較低的彈性。

圖3 全自動運行系統(tǒng)彈性的概率密度分布Figure 3 Probability density distribution of resilience in fully automatic operation system

圖4給出了在CI、ATP等8類設備遭受擾動時全自動運行系統(tǒng)彈性的累積概率分布情況，在圖中可以看到CI遭受擾動后系統(tǒng)的彈性跨度最大，系統(tǒng)的彈性可能會低于0.9，其次ATP和ZC遭受擾動后系統(tǒng)的彈性也可能出現(xiàn)較低值，這是由于 CI、ATP、ZC設備本身滿足SIL4級，SIL4級設備的安全性對系統(tǒng)整體的性能影響很大，因此將SIL4級設備的各級性能狀態(tài)值性設置較高；同時這幾個設備在整個系統(tǒng)的拓撲結構中處于關鍵位置，與多個其他設備有交互關系，因此這3個設備故障后對系統(tǒng)整體的性能影響相對更大，也導致其彈性值會出現(xiàn)低值，因此在燕房線的日常運行中要重點關注CI、ATP以及ZC設備的故障情況，要從減少故障發(fā)生頻率以及提升故障處理效率兩方面進行改進。同時，從圖中可以看到發(fā)生在轉轍機和軌道電路的彈性概率累積分布跨度相對更窄，彈性值波動較小，且系統(tǒng)彈性值都在0.97以上，這說明轉轍機和軌道電路故障時，系統(tǒng)的彈性受其影響較小，這是由于與轉轍機和軌道電路有交互關系的設備只有CI，其發(fā)生故障后受其影響的鄰接設備較少，因此其故障對系統(tǒng)整體的影響相對較低，所以當轉轍機和軌道電路發(fā)生故障時，系統(tǒng)的彈性值相對波動較小，彈性值能保持在較高水平。另外，與其他設備相比，發(fā)生在軌道電路以及轉轍機的故障處理過程相對及時，不用等待至夜間列車停運進行系統(tǒng)調(diào)試處理，而本文采用的彈性量化方法與故障的恢復時間密切相關，因此使得因轉轍機和軌道電路故障得到的系統(tǒng)彈性能保持在較高值。

圖4 不同設備遭受擾動下的全自動運行系統(tǒng)彈性Figure 4 Resilience of automatic operation system with different equipment disturbed

4 結論

通過對全自動運行線路燕房線過去 3年間的信號及綜合監(jiān)控系統(tǒng)故障跟蹤記錄表進行數(shù)據(jù)篩選，計算出CI、ATP、ATO、ISCS等8個設備的故障頻率以及故障恢復時間分布函數(shù)，采用基于 Zobel的彈性量化方式以及蒙特卡羅仿真方法進行了系統(tǒng)彈性計算。

1) 通過105次蒙特卡羅仿真計算，求得了燕房線全自動運行系統(tǒng)的彈性估計值為0.989 3。

2) 分別計算了CI、ATP等8個設備遭受隨機擾動時系統(tǒng)彈性值的概率累積分布，其中CI、ATP以及ZC在遭受擾動后系統(tǒng)的彈性值跨度較大，可能會出現(xiàn)較低值，因此在日常的列車運行過程中，應該加強對CI、ATP以及ZC設備的故障管理，盡量縮短其故障處理時間。

3) 轉轍機以及軌道電路隨機故障得到的系統(tǒng)彈性概率累積分布相對比較穩(wěn)定，基本維持在較高水平，這與受轉轍機和軌道電路故障影響的鄰接設備較少有關，導致其故障后系統(tǒng)性能值下降相對較低，因此在計算系統(tǒng)彈性時，彈性值也不會出現(xiàn)較大的下滑。另外，由于處理轉轍機與軌道電路故障無需等待列車停運再進行系統(tǒng)調(diào)試處理，節(jié)約了故障處理時間，這使得計算得到的彈性值也相對較高。

通過采取基于蒙特卡羅仿真方法對全自動運行系統(tǒng)彈性進行計算，將各設備受擾動以及性能下降程度都以概率的形式予以考慮，使彈性計算結果更能貼近系統(tǒng)真實的彈性水平，更能有效地反映系統(tǒng)整體的彈性情況，在實際的運營過程中，彈性計算量可作為綜合反映系統(tǒng)受擾動后重新恢復正常的能力。目前，全自動運行系統(tǒng)的技術已趨于成熟，越來越多新修建的地鐵也開始采用全自動運行系統(tǒng)，但圍繞全自動運行系統(tǒng)的運營管理技術卻還處于起步階段，列車出現(xiàn)的緊急情況，調(diào)度人員只能在遠程進行指導和處理，如何在列車出現(xiàn)緊急情況時，在沒有司機的幫助下還能保證全自動運行系統(tǒng)的高效率，運營管理人員還需要更多的實際經(jīng)驗的摸索，而彈性作為衡量系統(tǒng)從不穩(wěn)定狀態(tài)恢復到正常狀態(tài)的能力，可幫助運營管理人員在實際的運維操作中定量評估其應對措施效率的高低。

致謝：感謝城市軌道交通北京實驗室項目資助。