做電力系統(tǒng)風(fēng)險(xiǎn)評(píng)估的同行應(yīng)該都遇到過(guò)這個(gè)尷尬蒙特卡洛抽樣跑了幾萬(wàn)次期望損失算出來(lái)漂漂亮亮可一旦想回答“哪幾組初始故障組合能把系統(tǒng)搞癱”統(tǒng)計(jì)結(jié)果往往是“查無(wú)此案”。原因很簡(jiǎn)單大面積停電屬于典型的低概率高損失事件在均勻隨機(jī)抽樣里可能幾千次都輪不到一次。我當(dāng)時(shí)接手一個(gè)級(jí)聯(lián)故障風(fēng)險(xiǎn)分析的需求也是先撞上這面墻評(píng)估結(jié)果始終“偏樂(lè)觀”后來(lái)把思路從“被動(dòng)隨機(jī)抽樣”翻轉(zhuǎn)為“主動(dòng)搜索最壞情景”用隨機(jī)化學(xué)算法去直接狩獵高風(fēng)險(xiǎn)故障組合評(píng)估效率和數(shù)據(jù)質(zhì)量都上了一個(gè)臺(tái)階。這篇是完整的Matlab實(shí)現(xiàn)復(fù)盤既是算法講解也把我在建模、仿真、調(diào)參路上踩過(guò)的坑一并交代清楚。想拿直流潮流做連鎖故障分析又不想成天跑幾十萬(wàn)次抽樣的朋友這篇可以直接當(dāng)參考。1. 為什么蒙特卡洛抽樣會(huì)低估級(jí)聯(lián)故障風(fēng)險(xiǎn)組合爆炸是第一道坎1.1 從N-1到N-k搜索空間漲得比想象中快先看一個(gè)具體數(shù)字。一個(gè)39節(jié)點(diǎn)的標(biāo)準(zhǔn)測(cè)試系統(tǒng)支路大約46條N-1故障只有46種。N-2按排列組合算是C(46,2)1035個(gè)雙回線組合。到N-3就變成15180組N-4直接到163185組。還沒(méi)算上故障順序、保護(hù)動(dòng)作時(shí)限、隱性故障這些因素光是“初始斷開哪幾條線路”這個(gè)維度就已經(jīng)到了枚舉吃不消的規(guī)模。蒙特卡洛的辦法是把這些組合當(dāng)成均勻樣本去抽。抽得多了自然能覆蓋一部分組合但覆蓋≠有效覆蓋。大部分隨機(jī)抽到的故障可能連鎖反應(yīng)一兩步就結(jié)束了真正能引發(fā)大范圍停電的組合在組合空間中占比很低。我在測(cè)試中跑過(guò)3萬(wàn)次隨機(jī)抽樣最嚴(yán)重的一次系統(tǒng)失負(fù)荷也只到了總負(fù)荷的20%左右因?yàn)槲膊渴录某霈F(xiàn)概率實(shí)在不高。不是說(shuō)蒙特卡洛一無(wú)是處它在算期望值、生成統(tǒng)計(jì)分布方面仍然可靠只是用來(lái)找“最壞情景”的時(shí)候搜索效率很不夠看。1.2 級(jí)聯(lián)過(guò)程不是“斷開一條線”那么簡(jiǎn)單的級(jí)聯(lián)故障評(píng)估還有一個(gè)難點(diǎn)在于潮流會(huì)重分布。一條線路跳閘后原本由它承載的功率會(huì)轉(zhuǎn)移到剩余支路上其他線路的負(fù)載率可能瞬間上升保護(hù)裝置動(dòng)作切除過(guò)載線路然后功率再次轉(zhuǎn)移形成連鎖反應(yīng)。如果過(guò)程中出現(xiàn)孤島某些區(qū)域發(fā)電機(jī)出力與負(fù)荷功率嚴(yán)重不平衡就只能切負(fù)荷保電網(wǎng)。這一串事件在每一次潮流計(jì)算之后都要重新判斷所以評(píng)估一次故障組合的代價(jià)遠(yuǎn)不是“解一次潮流”那么簡(jiǎn)單而是要把整個(gè)連鎖過(guò)程走完。這也是為什么風(fēng)險(xiǎn)評(píng)估不能只看初始故障本身。同樣一條線路斷開在不同的系統(tǒng)負(fù)荷水平下后續(xù)連鎖路徑可能差別很大。我見過(guò)某組N-2故障在低負(fù)荷時(shí)只損了少量負(fù)荷負(fù)荷抬高10%之后直接發(fā)展成系統(tǒng)性停電。所以模型里的負(fù)荷參數(shù)必須能調(diào)節(jié)才有可比性。1.3 風(fēng)險(xiǎn)指標(biāo)怎么定期望值之外還要盯住尾部風(fēng)險(xiǎn)量化最常用的是失負(fù)荷概率PLC和期望失電量EENS這兩個(gè)指標(biāo)對(duì)規(guī)劃決策很有用。但做風(fēng)險(xiǎn)研究時(shí)我習(xí)慣額外關(guān)注條件風(fēng)險(xiǎn)值CVaR和最大失負(fù)荷量Max Loss。因?yàn)殡娏ο到y(tǒng)安全評(píng)估的特殊之處在于小概率大后果事件恰恰是運(yùn)行人員最需要提前防范的。用期望值評(píng)價(jià)容易被海量“無(wú)事發(fā)生”的樣本稀釋掩蓋掉真正的隱患。所以在本文的方案里目標(biāo)函數(shù)不是單用失負(fù)荷概率而是把“故障后果的嚴(yán)重度”作為核心采用失負(fù)荷量、受影響節(jié)點(diǎn)數(shù)、連鎖切除的線路數(shù)加權(quán)組合來(lái)定義。這樣隨機(jī)化學(xué)算法在搜索時(shí)有明確的目標(biāo)可以優(yōu)化而不是去擬合一條毫無(wú)區(qū)分度的平均曲線。2. 隨機(jī)化學(xué)算法原理把候選解當(dāng)分子把風(fēng)險(xiǎn)值當(dāng)勢(shì)能2.1 用化學(xué)反應(yīng)的語(yǔ)言描述優(yōu)化過(guò)程隨機(jī)化學(xué)算法Stochastic Chemical Reaction Optimization是化學(xué)反應(yīng)優(yōu)化框架下的隨機(jī)變體思路很直觀把優(yōu)化問(wèn)題里的每個(gè)候選解想象成一個(gè)在反應(yīng)容器里運(yùn)動(dòng)的分子。分子有結(jié)構(gòu)、有勢(shì)能、有動(dòng)能這其中的映射關(guān)系如下。分子結(jié)構(gòu)一組故障組合編碼。對(duì)于級(jí)聯(lián)故障問(wèn)題可以是一組初始斷開線路的編號(hào)。勢(shì)能PE候選解對(duì)應(yīng)的目標(biāo)函數(shù)值。勢(shì)能越低說(shuō)明這個(gè)故障組合造成的后果越輕微勢(shì)能越高說(shuō)明后果越嚴(yán)重。動(dòng)能KE分子的活躍程度?？梢岳斫鉃檫@個(gè)解還“有多少余力”被改造動(dòng)能大代表它容易跳出當(dāng)前局部區(qū)域。最低勢(shì)能minPE和最低結(jié)構(gòu)minStruct記錄該分子歷史上見過(guò)的最優(yōu)解相當(dāng)于帶記憶功能的個(gè)體最優(yōu)。化學(xué)反應(yīng)優(yōu)化的核心邏輯是讓分子通過(guò)碰撞、分解、合成等反應(yīng)算子不斷改變結(jié)構(gòu)同時(shí)在勢(shì)能和動(dòng)能之間做能量交換最終讓整個(gè)系統(tǒng)趨向低勢(shì)能狀態(tài)。放到我們的場(chǎng)景里就是不斷搜索“能導(dǎo)致更大停電范圍”的初始故障組合。2.2 四種反應(yīng)算子各司其職隨機(jī)化學(xué)算法一般包含四類反應(yīng)恰好對(duì)應(yīng)四種搜索行為。第一類單分子無(wú)效碰撞。單個(gè)分子撞到容器壁或內(nèi)部微小擾動(dòng)結(jié)構(gòu)發(fā)生小范圍突變比如翻轉(zhuǎn)編碼中的一個(gè)位把某條線路從故障集合里換成另一條。如果新結(jié)構(gòu)勢(shì)能更低就保留勢(shì)能變高則可能需要消耗動(dòng)能才允許接受這是典型的局部搜索行為。第二類分解反應(yīng)。一個(gè)分子分裂成兩個(gè)新分子新分子結(jié)構(gòu)與原分子完全不同通常伴隨較大幅度的隨機(jī)重組。這一算子是全局探測(cè)的主力用來(lái)跳出局部最優(yōu)代價(jià)是需要系統(tǒng)提供額外能量也就是要消耗緩沖區(qū)里的能量。第三類分子間無(wú)效碰撞。兩個(gè)分子碰撞交換部分信息類似遺傳算法里的交叉保留各自片段再組合成兩個(gè)新分子。這個(gè)算子在中等范圍上進(jìn)行搜索讓比較優(yōu)秀的故障片段在種群中擴(kuò)散。第四類合成反應(yīng)。兩個(gè)分子合二為一新分子繼承雙方的部分特征。合成往往在兩個(gè)分子的動(dòng)能都很低、都趨向收斂時(shí)觸發(fā)是加速局部收斂的手段。2.3 隨機(jī)性體現(xiàn)在哪些地方隨機(jī)化學(xué)算法和標(biāo)準(zhǔn)化學(xué)優(yōu)化算法的區(qū)別主要體現(xiàn)在三個(gè)地方一是反應(yīng)通道的選擇帶隨機(jī)性每次迭代按概率決定發(fā)生哪類反應(yīng)二是結(jié)構(gòu)擾動(dòng)的幅度由隨機(jī)因子控制不像遺傳算法那樣固定比例變異三是能量交換不完全按確定性規(guī)則執(zhí)行新分子的動(dòng)能常常從區(qū)間內(nèi)隨機(jī)取值這樣種群能保持一個(gè)比較自然的多樣性水平。我自己的體會(huì)是這種隨機(jī)性對(duì)于級(jí)聯(lián)故障搜索特別合適。因?yàn)楣收辖M合的目標(biāo)函數(shù)非常不平滑往往改變一條線路后果就從“很小”跳到“極大”存在大量斷層基于梯度的算法根本無(wú)從下手而基于化學(xué)反應(yīng)的隨機(jī)搜索能靠分解反應(yīng)跨過(guò)這些斷層。2.4 算法主流程偽代碼隨機(jī)化學(xué)算法的主循環(huán)相當(dāng)簡(jiǎn)潔偽代碼如下。初始化種群N個(gè)分子隨機(jī)生成N組故障組合 初始化中央緩沖區(qū)能量 buffer for iter 1:maxIter 按概率選擇單分子反應(yīng)或雙分子反應(yīng) if 單分子反應(yīng): if 分子能量較高且隨機(jī)條件滿足: 分解 else: 單分子無(wú)效碰撞(鄰域擾動(dòng)) else: if 兩分子動(dòng)能都低于閾值: 合成 else: 分子間無(wú)效碰撞(信息交換) 評(píng)估新分子目標(biāo)函數(shù)(級(jí)聯(lián)故障模擬) 更新分子勢(shì)能、動(dòng)能、最優(yōu)結(jié)構(gòu)和種群最優(yōu) 按規(guī)則向緩沖區(qū)存取能量 返回全局最優(yōu)故障組合寫完這段邏輯再看圖就非常直接了每次迭代都在“擾動(dòng)故障組合→算一次級(jí)聯(lián)模擬→看損失大不大→根據(jù)能量調(diào)整是否接受”這個(gè)圈子里轉(zhuǎn)。目標(biāo)函數(shù)評(píng)估頻率高所以級(jí)聯(lián)模擬器必須足夠輕量這也是后文選擇直流潮流的核心原因。3. 級(jí)聯(lián)故障模擬器直流潮流與過(guò)載切除的循環(huán)博弈3.1 為什么用直流潮流而不是交流潮流級(jí)聯(lián)故障評(píng)估要在每個(gè)候選解上跑多輪潮流一次連鎖過(guò)程往往要算好幾輪而隨機(jī)化學(xué)算法一代就要評(píng)估幾十次整個(gè)優(yōu)化算下來(lái)潮流的求解次數(shù)是幾千甚至上萬(wàn)次。這個(gè)計(jì)算量用交流潮流硬扛不太現(xiàn)實(shí)尤其是還要考慮牛頓-拉夫遜法收斂性在系統(tǒng)解列時(shí)可能出問(wèn)題所以我在這類風(fēng)險(xiǎn)搜索研究中統(tǒng)一用直流潮流近似。直流潮流把無(wú)功和電壓?jiǎn)栴}全部忽略只保留有功功率和相角的關(guān)系網(wǎng)絡(luò)方程為PBθB是由線路電抗倒數(shù)構(gòu)成的導(dǎo)納矩陣。求解只需要一次線性消元速度快數(shù)值穩(wěn)定性好在分析“哪條線路過(guò)載、哪個(gè)孤島缺電”這類問(wèn)題時(shí)精度足夠。它當(dāng)然算不出電壓崩潰但級(jí)聯(lián)故障研究中更關(guān)心拓?fù)溥B通性和有功潮流轉(zhuǎn)移這個(gè)精度換速度的買賣做得過(guò)。3.2 級(jí)聯(lián)模擬主循環(huán)的完整步驟級(jí)聯(lián)模擬器是評(píng)估函數(shù)的核心它的職責(zé)是給定一組初始斷開的線路完整播放一次“斷開→潮流轉(zhuǎn)移→過(guò)載→再斷開”的連鎖過(guò)程。實(shí)現(xiàn)邏輯其實(shí)不復(fù)雜重要的是循環(huán)出口條件設(shè)計(jì)合理。第一步把初始故障線路的開關(guān)狀態(tài)置0。第二步對(duì)剩余網(wǎng)格形成B矩陣并求解直流潮流。第三步檢查所有仍在運(yùn)行線路的功率與其容量之比超過(guò)過(guò)載閾值通常取1.0即負(fù)載率達(dá)到100%以上的線路進(jìn)入候選切除列表。第四步如果沒(méi)有過(guò)載線路連鎖停止否則按過(guò)載率從高到低排序切除最嚴(yán)重的一條或一批回到第二步重新求潮流。第五步重復(fù)直到?jīng)]有新過(guò)載、或者潮流矩陣因孤島而奇異、或者切除線路數(shù)達(dá)到上限。這套結(jié)構(gòu)里最需要注意的是防止振蕩。某些網(wǎng)架結(jié)構(gòu)下切了一條線后另一條線過(guò)載切掉后又導(dǎo)致之前那條線過(guò)載兩邊反復(fù)橫跳。解決辦法是維護(hù)一個(gè)“已斷開”集合任何線路一旦切除就不再恢復(fù)從邏輯上杜絕反復(fù)投切的循環(huán)。3.3 失負(fù)荷量的折算方法級(jí)聯(lián)結(jié)束后要統(tǒng)計(jì)后果嚴(yán)重程度。對(duì)于直流潮流模型我的做法是在裂解成孤島之后對(duì)每個(gè)連通區(qū)域分別計(jì)算發(fā)電機(jī)總出力和負(fù)荷總量。如果發(fā)電機(jī)出力小于負(fù)荷說(shuō)明這個(gè)島缺電缺額部分作為切負(fù)荷量計(jì)入損失如果發(fā)電機(jī)出力大于負(fù)荷通常還會(huì)設(shè)置一個(gè)小的調(diào)節(jié)余量默認(rèn)認(rèn)為可以調(diào)低出力不計(jì)損失。于是總失負(fù)荷量就是所有孤島切負(fù)荷量的總和。另外我會(huì)把“系統(tǒng)解列后的受影響母線數(shù)”和“連鎖過(guò)程中動(dòng)作切除的線路總數(shù)”也記錄成結(jié)構(gòu)化字段用于最終目標(biāo)函數(shù)的加權(quán)。孤島內(nèi)功率失衡的計(jì)算依賴圖連通分量判斷Matlab里可以用graph和conncomp函數(shù)實(shí)現(xiàn)代碼量不大。3.4 模擬器的接口設(shè)計(jì)為了讓后面的優(yōu)化算法調(diào)用方便我在Matlab里把模擬器封裝成了一個(gè)函數(shù)function res cascade_simulator(init_status, bus, branch) % init_status: 邏輯向量1表示該線路初始斷開 % bus: 節(jié)點(diǎn)表包含負(fù)荷有功、發(fā)電有功 % branch: 支路表包含首末端節(jié)點(diǎn)、電抗、容量 res.loss_mw 0; res.lost_buses 0; res.tripped_lines 0; ... end所有候選解都會(huì)調(diào)用這個(gè)函數(shù)返回的結(jié)果就是目標(biāo)函數(shù)的輸入。這里有一個(gè)工程上的好習(xí)慣級(jí)聯(lián)模擬器內(nèi)部不要寫任何隨機(jī)數(shù)必須保持確定性。同一個(gè)故障組合任何時(shí)候調(diào)用結(jié)果都應(yīng)該一致否則隨機(jī)化學(xué)算法會(huì)在優(yōu)化時(shí)被噪聲干擾很難穩(wěn)定收斂。4. Matlab代碼骨架從分子編碼到四種反應(yīng)算子的實(shí)現(xiàn)4.1 故障組合的編碼方式級(jí)聯(lián)故障的初始故障集合我用一組整數(shù)編碼表示假設(shè)系統(tǒng)有L條線路初始故障集合就是從1到L中選出k個(gè)互不相同的整數(shù)。種群中的每個(gè)分子結(jié)構(gòu)就是一個(gè)1×k的整型向量比如[7, 23, 41]表示第7、23、41條線路同時(shí)斷開。這個(gè)編碼比二進(jìn)制向量更緊湊在實(shí)現(xiàn)分解和合成算子時(shí)也更方便。不過(guò)要小心隨機(jī)改動(dòng)結(jié)構(gòu)時(shí)可能導(dǎo)致重復(fù)元素所以每次擾動(dòng)之后都要做一次unique去重再補(bǔ)齊長(zhǎng)度。4.2 分子的數(shù)據(jù)結(jié)構(gòu)設(shè)計(jì)在Matlab里我習(xí)慣直接用struct數(shù)組表示分子這樣代碼清晰調(diào)試也方便。mol struct(... struct, initial_structure, ... PE, 0, ... KE, 0, ... minPE, inf, ... minStruct, [], ... numHit, 0);PE初始值要設(shè)為極大值因?yàn)閯?shì)能越低越好第一輪評(píng)估之后會(huì)更新。minStruct用來(lái)記錄這個(gè)分子歷史最優(yōu)結(jié)構(gòu)。KE的初始值通常設(shè)為一個(gè)與問(wèn)題規(guī)模相關(guān)的正數(shù)太小則分解難以觸發(fā)太大則搜索過(guò)于發(fā)散。4.3 反應(yīng)算子的實(shí)現(xiàn)要點(diǎn)單分子無(wú)效碰撞的核心是鄰域擾動(dòng)。對(duì)于整數(shù)編碼的故障集合我每次以0.2概率隨機(jī)替換一條線路以0.3概率交換兩條線路位置其余情況做微量偏移。新結(jié)構(gòu)評(píng)估后如果勢(shì)能更小直接接受如果勢(shì)能變大而且能量差ΔPE大于當(dāng)前動(dòng)能則拒絕否則接受同時(shí)將動(dòng)能減去ΔPE。分解算子的實(shí)現(xiàn)是產(chǎn)生兩個(gè)新分子第一個(gè)保留原結(jié)構(gòu)的一半另一半重新隨機(jī)生成第二個(gè)反過(guò)來(lái)取原結(jié)構(gòu)的另一半和新的隨機(jī)部分。分解后新分子的動(dòng)能分配按照能量守恒思路來(lái)一個(gè)偏高一個(gè)偏低用randi做隨機(jī)劃分。分子間無(wú)效碰撞可以用類似交叉的方式實(shí)現(xiàn)。兩個(gè)分子結(jié)構(gòu)分別是a和b生成兩個(gè)新分子a_new和b_newa_new取a的前半段加b的后半段b_new取b的前半段加a的后半段同樣要做去重并與原長(zhǎng)度對(duì)齊。如果a的勢(shì)能和b的勢(shì)能都屬于種群上游這種碰撞更容易保留優(yōu)秀片段。合成反應(yīng)相對(duì)簡(jiǎn)單把兩個(gè)分子合并成一個(gè)新結(jié)構(gòu)隨機(jī)從兩個(gè)父結(jié)構(gòu)中抽取元素組成。觸發(fā)條件是兩個(gè)分子的動(dòng)能都小于一個(gè)閾值說(shuō)明它們已經(jīng)沒(méi)有多少探索空間了合并有助于集中搜索力量。4.4 算法主循環(huán)的結(jié)構(gòu)主循環(huán)骨架如下我把它寫成單文件腳本方便實(shí)驗(yàn)時(shí)反復(fù)調(diào)整參數(shù)。for iter 1:maxIter for idx 1:popsize if rand() 0.6 % 單分子反應(yīng) if mol(idx).KE buffer decompose_threshold ... rand() decompose_rate mol decompose(mol(idx), ...); else mol(idx) ineffective_collision(mol(idx), ...); end else % 雙分子反應(yīng) j random_select_other(idx); if mol(idx).KE synth_threshold mol(j).KE synth_threshold mol(idx) synthesize(mol(idx), mol(j), ...); else [mol(idx), mol(j)] inter_inelastic_collision(mol(idx), mol(j)); end end [mol(idx), global_best] evaluate_and_update(mol(idx), global_best); end % 更新緩沖區(qū) if mod(iter, 10) 0 buffer buffer 5.0 * randn(); end end這里evaluate_and_update負(fù)責(zé)調(diào)用級(jí)聯(lián)模擬器計(jì)算勢(shì)能并按規(guī)則更新分子的PE、KE和minStruct。我在實(shí)際項(xiàng)目里把評(píng)估函數(shù)單獨(dú)拆開使用MATLAB函數(shù)句柄傳入這樣萬(wàn)一要換成交流潮流或者引入隱性故障模型主循環(huán)一行都不用改。4.5 參數(shù)經(jīng)驗(yàn)值參考隨機(jī)化學(xué)算法的參數(shù)不像遺傳算法那么多但幾個(gè)關(guān)鍵量還是要重點(diǎn)調(diào)的。我給出自己在39節(jié)點(diǎn)系統(tǒng)上測(cè)試后比較穩(wěn)定的區(qū)間供剛上手的同學(xué)參考。參數(shù)含義經(jīng)驗(yàn)值區(qū)間popsize種群規(guī)模20~40maxIter迭代代數(shù)100~200decompose_rate分解觸發(fā)系數(shù)0.2~0.4synth_threshold合成觸發(fā)動(dòng)能閾值2~8initial_KE分子初始動(dòng)能5~15buffer中央緩沖區(qū)初始能量10~30需要注意的是這些參數(shù)是和目標(biāo)函數(shù)數(shù)值尺度耦合的。如果目標(biāo)函數(shù)是失負(fù)荷量兆瓦級(jí)勢(shì)能數(shù)量級(jí)就很大如果目標(biāo)函數(shù)是歸一化指標(biāo)勢(shì)能就在0到1之間。合成和分解的閾值必須跟著勢(shì)能尺度調(diào)整否則會(huì)長(zhǎng)時(shí)間不觸發(fā)分解算法退化成局部搜索。5. 39節(jié)點(diǎn)系統(tǒng)實(shí)測(cè)SCRO找出了蒙特卡洛漏掉的“高影響場(chǎng)景”5.1 一次同規(guī)模對(duì)比實(shí)驗(yàn)的過(guò)程與結(jié)果為了驗(yàn)證隨機(jī)化學(xué)算法的實(shí)用性我在39節(jié)點(diǎn)標(biāo)準(zhǔn)測(cè)試系統(tǒng)上做了一組對(duì)比隨機(jī)抽樣3萬(wàn)次隨機(jī)化學(xué)算法跑150代、種群規(guī)模20兩種方法都用同一個(gè)級(jí)聯(lián)模擬器指標(biāo)都看最大失負(fù)荷量、平均失負(fù)荷量和失負(fù)荷超過(guò)20%總負(fù)荷的場(chǎng)景數(shù)量。結(jié)果很有參考價(jià)值。蒙特卡洛3萬(wàn)次抽樣中最嚴(yán)重的故障組合導(dǎo)致了約19.8%的總負(fù)荷損失失負(fù)荷超過(guò)15%的場(chǎng)景數(shù)量只有個(gè)位數(shù)大部分樣本的失負(fù)荷都在5%以下。而隨機(jī)化學(xué)算法在150代中找到的最壞組合能把接近27%的負(fù)荷甩掉而且搜索結(jié)束時(shí)還給出了好幾組后果相當(dāng)嚴(yán)重的候選故障這些在蒙特卡洛抽樣里基本都被漏掉了。這說(shuō)明不是蒙特卡洛“不夠準(zhǔn)”而是它的抽樣策略對(duì)罕見高影響事件天然不敏感剛好切中了級(jí)聯(lián)故障風(fēng)險(xiǎn)的痛點(diǎn)越需要防范的事故越難被隨機(jī)看見。這也是啟發(fā)式搜索在風(fēng)險(xiǎn)評(píng)估里真正的價(jià)值所在。5.2 踩坑記錄一潮流矩陣奇異導(dǎo)致連鎖中斷第一次把模擬器接進(jìn)算法時(shí)我遇到的第一個(gè)bug是某些故障組合切掉線路后系統(tǒng)被切成兩個(gè)互不連通的孤島孤立子系統(tǒng)對(duì)應(yīng)的直流潮流B矩陣變成奇異矩陣直接用\求解會(huì)報(bào)NaN。我當(dāng)時(shí)的第一反應(yīng)是把奇異矩陣單獨(dú)處理成“連鎖結(jié)束”統(tǒng)計(jì)這個(gè)瞬間的失負(fù)荷量。后來(lái)發(fā)現(xiàn)這個(gè)處理太粗暴了。因?yàn)槿绻麅蓚€(gè)孤島各自還能獨(dú)立平衡潮流連鎖過(guò)程未必停止甚至某些孤島內(nèi)部可能還有拓?fù)湔{(diào)整空間。正確的做法是用圖論的連通性分析先找到所有連通分量對(duì)每個(gè)分量單獨(dú)判斷功率平衡再?zèng)Q定是否需要切負(fù)荷。5.3 踩坑記錄二過(guò)載切除順序影響模擬結(jié)果另一個(gè)讓我印象很深的坑是過(guò)載切除策略。最開始我按“一次切除所有過(guò)載線路”去寫結(jié)果在環(huán)網(wǎng)結(jié)構(gòu)較多的系統(tǒng)里一下子切掉多條線路會(huì)導(dǎo)致潮流矩陣劇烈變化連鎖過(guò)程一步就跳過(guò)去了模擬結(jié)果要么失負(fù)荷量特別大要么出現(xiàn)一些明顯不合理的“假停電影像”。后來(lái)我改成“每次只切除過(guò)載率最高的那一條”每切一條就重新求解一次潮流模擬結(jié)果變得合理多了也更貼近實(shí)際保護(hù)動(dòng)作的序列化特征。代價(jià)是級(jí)聯(lián)模擬器的計(jì)算量增加了幾倍但換來(lái)的是結(jié)果的解釋性強(qiáng)很多畢竟風(fēng)險(xiǎn)評(píng)估里任何一個(gè)異常偏高的數(shù)值都要能說(shuō)清它來(lái)源。5.4 踩坑記錄三目標(biāo)函數(shù)權(quán)重分配不合理時(shí)算法“鉆空子”我最早把目標(biāo)函數(shù)設(shè)為“失負(fù)荷量50×切除線路數(shù)”意圖是讓算法兼顧停電損失和連鎖范圍。結(jié)果算法確實(shí)收斂得很快但最優(yōu)解全是能夠制造大量線路跳閘、卻沒(méi)什么真實(shí)失負(fù)荷的故障組合。說(shuō)白了搜索器發(fā)現(xiàn)了我的權(quán)重漏洞在拼命刷線路斷開次數(shù)。這也是啟發(fā)式優(yōu)化的通用問(wèn)題目標(biāo)函數(shù)怎么寫算法就朝哪個(gè)方向走。級(jí)聯(lián)故障研究里后果評(píng)價(jià)體系必須貼近真實(shí)的保電價(jià)號(hào)。我最后的處理是把切除線路數(shù)的權(quán)重降得很低只作為輔助項(xiàng)主目標(biāo)完全落在失負(fù)荷量和受影響的負(fù)荷節(jié)點(diǎn)數(shù)上。5.5 進(jìn)一步驗(yàn)證最優(yōu)解的可信度啟發(fā)式算法找到的最優(yōu)解不能直接當(dāng)成真實(shí)系統(tǒng)的最壞故障去上報(bào)。我的標(biāo)準(zhǔn)流程是在算法結(jié)束后把排名前10的故障組合全部單獨(dú)提交給原始交流潮流程序重新跑一遍確認(rèn)它們?cè)诟?xì)模型下仍然會(huì)引發(fā)大規(guī)模停電。這一步很有必要直流潮流下的高風(fēng)險(xiǎn)組合交流模型下未必同樣嚴(yán)重但反過(guò)來(lái)凡是直流模型下穩(wěn)定復(fù)現(xiàn)的高后果場(chǎng)景絕大多數(shù)在實(shí)際系統(tǒng)中都是危險(xiǎn)情況值得深挖。6. 參數(shù)調(diào)優(yōu)與工程化建議從跑通Demo到真正常態(tài)化使用6.1 先用小網(wǎng)絡(luò)校準(zhǔn)參數(shù)再上大系統(tǒng)我剛拿到算法時(shí)直接在39節(jié)點(diǎn)系統(tǒng)上調(diào)參每次實(shí)驗(yàn)要跑一個(gè)多小時(shí)效率極低。后來(lái)學(xué)著先用一個(gè)14節(jié)點(diǎn)的簡(jiǎn)化測(cè)試系統(tǒng)調(diào)試模擬器和基本參數(shù)確認(rèn)主循環(huán)、評(píng)估函數(shù)、四種反應(yīng)算子全部工作正常后再切換到39節(jié)點(diǎn)系統(tǒng)做完整實(shí)驗(yàn)時(shí)間成本一下降了下來(lái)。小網(wǎng)絡(luò)的價(jià)值在于能快速驗(yàn)證趨勢(shì)分解率太高會(huì)讓搜索變得像隨機(jī)抽樣合成率太高會(huì)讓種群過(guò)早抱團(tuán)。找到一組在小網(wǎng)絡(luò)上的合理參數(shù)遷移到大規(guī)模系統(tǒng)時(shí)往往只在小范圍內(nèi)微調(diào)即可。如果有條件可以再做一個(gè)參數(shù)網(wǎng)格搜索結(jié)合目標(biāo)函數(shù)的收斂曲線判斷參數(shù)是否合適。6.2 計(jì)算加速的三個(gè)簡(jiǎn)單手段隨機(jī)化學(xué)算法計(jì)算瓶頸幾乎全部在級(jí)聯(lián)模擬器上所以最劃算的優(yōu)化一定是讓模擬器變輕。我這里有三條實(shí)踐經(jīng)驗(yàn)第一把潮流矩陣的稀疏結(jié)構(gòu)預(yù)計(jì)算好直流潮流的B矩陣因子化后后續(xù)重算只需要更新被切除線路對(duì)應(yīng)的元素第二級(jí)聯(lián)模擬時(shí)及時(shí)移除孤島節(jié)點(diǎn)縮小矩陣維數(shù)第三使用MATLAB的parfor并行評(píng)估種群因?yàn)槊總€(gè)候選解之間互不依賴并行收益非常直接。我還在模擬器里加了一層緩存如果某條線路組合之前已經(jīng)評(píng)估過(guò)直接從緩存表讀取結(jié)果不再重算。隨機(jī)化學(xué)算法在后期經(jīng)常會(huì)在最優(yōu)解附近反復(fù)試探這一層緩存省下的計(jì)算量相當(dāng)可觀。6.3 風(fēng)險(xiǎn)搜索的下游應(yīng)用從評(píng)估到加固決策搜索到的高風(fēng)險(xiǎn)故障組合本身不是項(xiàng)目終點(diǎn)。我在實(shí)踐中會(huì)把這些組合轉(zhuǎn)成交互性分析統(tǒng)計(jì)被反復(fù)觸發(fā)的高頻故障線路把這些線路對(duì)應(yīng)的斷面做成風(fēng)險(xiǎn)熱區(qū)圖再用隨機(jī)化學(xué)算法跑“加固決策”的反向問(wèn)題比如只允許加固三條線路哪種加固方案能讓最壞情景下的失負(fù)荷量降到最低。這個(gè)反向問(wèn)題同樣適合用隨機(jī)化學(xué)算法求解因?yàn)榧庸探M合也是典型的組合優(yōu)化問(wèn)題整數(shù)編碼和四種反應(yīng)算子基本可以照搬只需要把目標(biāo)函數(shù)改成“加固后最壞故障組合的剩余風(fēng)險(xiǎn)”。如果你手頭已經(jīng)有級(jí)聯(lián)故障模擬器加一個(gè)透?jìng)鲄?shù)就能把加固邏輯嵌進(jìn)去實(shí)現(xiàn)成本很低。6.4 從39節(jié)點(diǎn)走向更大系統(tǒng)的擴(kuò)展思路電網(wǎng)規(guī)模變大以后故障組合搜索空間會(huì)漲得更夸張。我的建議是不要一味增加迭代次數(shù)而是引入“故障深度限制”和“關(guān)鍵斷面篩選”兩層前置約束。先用直流靈敏度分析挑出那些功率傳輸占比高、斷開后影響面大的關(guān)鍵線路集合讓算法只在關(guān)鍵斷面內(nèi)搜索初始故障把N-k中的k從全系統(tǒng)縮小到關(guān)鍵線路集合范圍內(nèi)。這樣計(jì)算量大幅下降結(jié)果卻往往更有工程決策參考價(jià)值。從39節(jié)點(diǎn)擴(kuò)展到幾百節(jié)點(diǎn)系統(tǒng)時(shí)我還把隨機(jī)化學(xué)算法的單分子變異概率隨迭代代數(shù)做了指數(shù)衰減讓算法在早期大范圍探索、晚期集中精修。這個(gè)策略對(duì)大系統(tǒng)效果很突出收斂速度明顯加快而且沒(méi)有明顯犧牲最優(yōu)解質(zhì)量。最后再分享一點(diǎn)個(gè)人體會(huì)隨機(jī)化學(xué)算法在整個(gè)項(xiàng)目里并不是取代蒙特卡洛而是補(bǔ)上蒙特卡洛最不擅長(zhǎng)的“尾部搜索”環(huán)節(jié)。比較穩(wěn)妥的評(píng)估流程是蒙特卡洛負(fù)責(zé)統(tǒng)計(jì)期望風(fēng)險(xiǎn)隨機(jī)化學(xué)算法負(fù)責(zé)挖掘最壞情景兩者結(jié)論合在一起做風(fēng)險(xiǎn)管理決策。我在實(shí)際使用中發(fā)現(xiàn)這類元啟發(fā)式算法最怕的不是參數(shù)難調(diào)而是目標(biāo)函數(shù)定義得含糊一旦后果評(píng)價(jià)標(biāo)準(zhǔn)自己都說(shuō)不清算法一定會(huì)鉆空子。所以新手入門先把級(jí)聯(lián)模擬器寫扎實(shí)把失負(fù)荷統(tǒng)計(jì)的口徑定清楚再談優(yōu)化算法。模擬器靠譜了算法即便參數(shù)次優(yōu)出來(lái)的結(jié)果也不會(huì)離譜模擬器本身有問(wèn)題再花哨的算法也只是把錯(cuò)誤答案搜得更穩(wěn)定而已。希望這套實(shí)現(xiàn)記錄能幫你少踩一些我踩過(guò)的坑。