韌性提升的MPS預(yù)配置建模與Matlab實(shí)現(xiàn))
臺風(fēng)過境那一夜配電網(wǎng)里好幾條關(guān)鍵線路同時跳閘城區(qū)大片負(fù)荷失電。這時候你手里只有三臺移動應(yīng)急電源MPS你會把它們提前放到哪幾個節(jié)點(diǎn)放對了醫(yī)院和通信基站還能撐住放錯了只能看著負(fù)荷一個一個掉。這個“提前放哪兒、放多少”的決策就是配電網(wǎng)韌性研究里常說的MPS預(yù)配置問題。這篇博文來自我最近復(fù)現(xiàn)的一篇SCI一區(qū)論文基于配電網(wǎng)韌性提升的應(yīng)急移動電源預(yù)配置和動態(tài)調(diào)度。因?yàn)閮?nèi)容量確實(shí)不小我拆成上下兩篇來寫這篇先講上篇——MPS預(yù)配置部分的模型推導(dǎo)與Matlab代碼實(shí)現(xiàn)。我會把建模動機(jī)、關(guān)鍵公式、約束含義、求解器配置、踩坑記錄全部攤開講適合正在做配電網(wǎng)韌性、災(zāi)后供電恢復(fù)、移動儲能相關(guān)的同學(xué)直接參考尤其建議準(zhǔn)備復(fù)現(xiàn)論文的研究生先把這篇啃完再動手。1. 問題定性與預(yù)配置思路拆解1.1 為什么極端天氣下要專門做“預(yù)配置”配電網(wǎng)韌性研究里最常被拿出來討論的時間線是這樣的災(zāi)害預(yù)警期、災(zāi)害發(fā)生期、故障穩(wěn)定期、搶修恢復(fù)期。MPS這類移動電源最大的特點(diǎn)是“能動”但移動需要時間而且災(zāi)害發(fā)生后道路和交通狀態(tài)往往比平時糟糕得多等故障信息全出來再調(diào)度車輛可能兩個小時后才能抵達(dá)目標(biāo)節(jié)點(diǎn)。預(yù)配置的作用就是把“運(yùn)輸時間”這個不確定性盡量前置消化掉。在原文的邏輯里MPS預(yù)配置Pre-positioning和動態(tài)調(diào)度Dynamic Scheduling是嚴(yán)格分開的兩個階段這也是我第一次讀這篇論文時覺得設(shè)計得很干凈的地方。預(yù)配置階段發(fā)生在災(zāi)害預(yù)警期決策的是“哪些候選節(jié)點(diǎn)提前部署MPS、每處放幾臺”動態(tài)調(diào)度階段發(fā)生在故障信息逐漸明確之后決策的是“已經(jīng)部署好的MPS如何去支援那些被隔離的區(qū)域”。上篇只聊前者。我打個比方你就理解了。這和臺風(fēng)來臨前政府往各個避難所預(yù)置應(yīng)急物資是一個邏輯你不能等臺風(fēng)把你困在樓里了再開始運(yùn)水運(yùn)糧你得在風(fēng)雨還沒起來的時候把物資放到最可能被圍困的幾個點(diǎn)。MPS預(yù)配置就是配電網(wǎng)版的“物資預(yù)置”只不過這里的物資變成了應(yīng)急電源車、移動儲能裝置而且放下去之后能不能在故障后成功頂上還取決于配電網(wǎng)的網(wǎng)絡(luò)拓?fù)浜统绷鳡顟B(tài)。1.2 “放哪”背后的數(shù)學(xué)本質(zhì)很多剛接觸這個問題的人會想預(yù)配置不就是把MPS放在負(fù)荷最大的幾個節(jié)點(diǎn)嘛按負(fù)荷排序選前幾名就完事了。真實(shí)情況遠(yuǎn)沒這么簡單因?yàn)镸PS的供電半徑和配電網(wǎng)的運(yùn)行方式耦合得很深。先看決策層面。預(yù)配置模型要同時回答三個問題選哪些節(jié)點(diǎn)放MPS每個節(jié)點(diǎn)放幾臺災(zāi)害發(fā)生后的各種故障場景里這些MPS到底能恢復(fù)多少負(fù)荷。這里的核心矛盾是MPS數(shù)量永遠(yuǎn)不夠覆蓋所有節(jié)點(diǎn)必須做取舍而不同的故障場景對最優(yōu)配置位置的偏好可能完全相反。為了處理這種不確定性和沖突性標(biāo)準(zhǔn)做法是構(gòu)造一組離散的故障場景每個場景對應(yīng)一組故障線路集合并給每個場景分配一個概率權(quán)重然后整體優(yōu)化所有場景下的期望恢復(fù)效果。從優(yōu)化建模的角度看這就是一個典型的場景法兩階段隨機(jī)規(guī)劃——只不過因?yàn)轭}目說“預(yù)配置”階段而不涉及災(zāi)后調(diào)度模型的第二階段通常被簡化成“在已知預(yù)配置位置后、給定故障場景下決策網(wǎng)絡(luò)開關(guān)狀態(tài)和MPS出力最大化恢復(fù)負(fù)荷”。目標(biāo)函數(shù)是所有場景下加權(quán)恢復(fù)負(fù)荷的總期望。這一步的理解特別關(guān)鍵預(yù)配置模型不是找“哪個節(jié)點(diǎn)重要”而是在給定災(zāi)害概率分布和配電網(wǎng)結(jié)構(gòu)約束下找“哪個節(jié)點(diǎn)組合能最大化所有故障場景下的期望恢復(fù)量”。這兩者的差別很多復(fù)現(xiàn)者一上來就會搞混導(dǎo)致后面模型建得跟論文對不上。2. 配電網(wǎng)潮流建模DistFlow怎么變成可求解的約束2.1 DistFlow遞推方程預(yù)配置和動態(tài)調(diào)度這類配電網(wǎng)優(yōu)化問題幾乎不會直接用牛頓-拉夫遜法做潮流計算因?yàn)榕ｎD-拉夫遜是一個數(shù)值迭代算法沒法作為約束條件嵌進(jìn)一個優(yōu)化模型里去反復(fù)求解。復(fù)現(xiàn)這類論文時最常用的做法是用DistFlow遞推方程來建模配電網(wǎng)潮流。DistFlow方程針對幅射狀配電網(wǎng)的一條支路。假設(shè)支路$i \to j$上從節(jié)點(diǎn)$i$流向節(jié)點(diǎn)$j$的有功為$P_{ij}$、無功為$Q_{ij}$支路電阻為$r_{ij}$、電抗為$x_{ij}$記$L_{ij} \frac{P_{ij}^2 Q_{ij}^2}{V_i^2}$為支路電流平方那么對節(jié)點(diǎn)$j$和它的父節(jié)點(diǎn)$i$有$$ P_j P_{ij} - r_{ij} L_{ij} - \sum_{k \in C(j)} P_{jk} $$$$ Q_j Q_{ij} - x_{ij} L_{ij} - \sum_{k \in C(j)} Q_{jk} $$$$ V_j^2 V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) (r_{ij}^2 x_{ij}^2)L_{ij} $$式中$C(j)$表示節(jié)點(diǎn)$j$的所有子節(jié)點(diǎn)集合。第一個式子表達(dá)的是節(jié)點(diǎn)$j$的有功平衡流入節(jié)點(diǎn)$j$的功率減去支路損耗減去流向所有子支路的功率剩下就是節(jié)點(diǎn)$j$自身的注入或負(fù)荷。這套遞推方程的好處是它本質(zhì)上是“平的”可以一步一步從根節(jié)點(diǎn)推到所有葉子節(jié)點(diǎn)而不需要像節(jié)點(diǎn)導(dǎo)納矩陣那樣做大矩陣求逆因此非常契合優(yōu)化模型里的約束表達(dá)。代價是P、Q、V之間存在乘積和平方項(xiàng)——這就是后面要做線性化或凸松弛的根源。2.2 二階錐松弛讓精度和可解性兼得如果直接把上面三個方程原封不動寫進(jìn)優(yōu)化模型問題就成了非凸非線性規(guī)劃商業(yè)求解器也拿它沒太大辦法。復(fù)現(xiàn)論文時最常見的兩條路線如下。第一條路線是忽略網(wǎng)損項(xiàng)即把$r_{ij}L_{ij}$、$x_{ij}L_{ij}$以及$(r_{ij}^2x_{ij}^2)L_{ij}$這幾項(xiàng)直接設(shè)為0。這樣一來潮流方程就變成純線性約束$$ P_j \sum_{k \in C(j)} P_{jk} $$$$ Q_j \sum_{k \in C(j)} Q_{jk} $$$$ V_j^2 \approx V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) $$這種線性化DistFlow在MPS功率不太大的場景下精度足夠求解速度也快很多早期配電網(wǎng)重構(gòu)論文都用這個方案。但這個近似在饋線末端重載時會出現(xiàn)比較明顯的電壓誤差如果預(yù)配置模型里MPS注入功率較大恢復(fù)后某些節(jié)點(diǎn)電壓可能實(shí)際上越限但模型沒有發(fā)現(xiàn)。第二條路線就是保留網(wǎng)絡(luò)損耗項(xiàng)但把等式重構(gòu)成凸約束。具體做法是引入輔助變量$\ell_{ij}$代替$L_{ij}$然后把約束$P_{ij}^2 Q_{ij}^2 \ell_{ij}V_i^2$松弛成不等式$P_{ij}^2 Q_{ij}^2 \leq \ell_{ij}V_i^2$再寫成標(biāo)準(zhǔn)二階錐形式$$ \left\lVert \begin{bmatrix} 2P_{ij} \ 2Q_{ij} \ \ell_{ij} - V_i^2 \end{bmatrix} \right\rVert_2 \leq \ell_{ij} V_i^2 $$這種松弛在數(shù)學(xué)上叫二階錐松弛SOCP relaxation對這個模型的實(shí)際計算表現(xiàn)絕大多數(shù)情況下松弛是緊的也就是說求出來的最優(yōu)解滿足原等式約束不丟精度。加上MPS的P、Q出力約束也是錐約束整個問題就變成了混合整數(shù)二階錐規(guī)劃MISOCP交給Gurobi、CPLEX這類商業(yè)求解器可以直接求解。為什么松弛通常緊配電網(wǎng)的運(yùn)行點(diǎn)一般處于一個比較“溫和”的區(qū)域電壓接近1 pu、支路電流不會太大錐約束的最優(yōu)解往往就落在錐邊界上這算是實(shí)際工程場景里一個被反復(fù)驗(yàn)證的良性特征。復(fù)現(xiàn)時你可以自己做一個簡單的對比實(shí)驗(yàn)分別用線性DistFlow和SOCP松弛各跑一遍看看恢復(fù)負(fù)荷比例的差異我實(shí)測在33節(jié)點(diǎn)系統(tǒng)里差別大約在1%左右但SOCP版本對電壓越限的把握要踏實(shí)得多。2.3 網(wǎng)絡(luò)拓?fù)浼s束與開關(guān)狀態(tài)的處理配電網(wǎng)正常運(yùn)行時是輻射狀結(jié)構(gòu)也就是一個沒有環(huán)路的樹狀網(wǎng)絡(luò)。預(yù)配置模型如果要考慮災(zāi)后通過分段開關(guān)和聯(lián)絡(luò)開關(guān)重構(gòu)網(wǎng)絡(luò)就必須把拓?fù)浼s束也寫進(jìn)去。最常用的輻射狀約束是“選定支路數(shù)等于節(jié)點(diǎn)數(shù)減去連通分量數(shù)”。如果假設(shè)全網(wǎng)連通那就要求$$ \sum_{(i,j) \in E} k_{ij} N - 1 $$其中$k_{ij}$是支路$ij$的開斷狀態(tài)二進(jìn)制變量1表示閉合。但這還不夠因?yàn)樯厦孢@個等式只能保證支路數(shù)正確不能排除“支路數(shù)對但分成兩片”的情況通常還需要再用生成樹約束或父子方向約束把它補(bǔ)全。一個比較容易實(shí)現(xiàn)的寫法是用網(wǎng)絡(luò)拓?fù)渲械挠邢虺绷鬏o助變量強(qiáng)制每個非根節(jié)點(diǎn)有且僅有一個父節(jié)點(diǎn)。在我復(fù)現(xiàn)的這篇論文的上篇里如果嚴(yán)格只做MPS預(yù)配置可以考慮一種合理的簡化預(yù)配置階段不主動優(yōu)化開關(guān)重構(gòu)而是把故障后的網(wǎng)絡(luò)拓?fù)湟暈椤肮收暇€路斷開、其他線路沿用原始拓?fù)洹钡墓潭ńY(jié)構(gòu)。這樣就不需要引入拓?fù)涠M(jìn)制變量問題規(guī)模會大幅縮小求解速度明顯提升。代價是會漏掉一些通過聯(lián)絡(luò)開關(guān)轉(zhuǎn)供帶來的恢復(fù)能力提升可能和原文結(jié)果有偏差。所以我給你的建議是分兩步走第一步先做固定拓?fù)浒姹镜腗PS預(yù)配置把選址邏輯跑通第二步再把拓?fù)渲貥?gòu)加進(jìn)去比較兩種結(jié)果。這樣你既能快速驗(yàn)證代碼又能在論文復(fù)現(xiàn)報告中寫出“預(yù)配置模型考慮重構(gòu)后恢復(fù)率提升了X%”這種有增量價值的結(jié)論。3. MPS預(yù)配置模型構(gòu)建變量、約束與目標(biāo)函數(shù)3.1 決策變量怎么設(shè)計MPS預(yù)配置模型的變量分兩層這一點(diǎn)是整個建模的核心骨架。第一層是災(zāi)前的“預(yù)配置決策變量”不依賴場景。我習(xí)慣用整數(shù)變量$n^{MPS}_i$表示在節(jié)點(diǎn)$i$部署的MPS臺數(shù)取值范圍從0到一個上限$N^{max}_i$。這里要解釋一下為什么用整數(shù)臺數(shù)而不是連續(xù)容量MPS一般是標(biāo)準(zhǔn)化產(chǎn)品一臺一臺計算實(shí)際決策中不可能說“放1.3臺”所以整數(shù)變量最貼合物理現(xiàn)實(shí)。第二層是災(zāi)后的“場景相關(guān)運(yùn)行變量”每個故障場景各有一套。包括每場景下每條支路的有功$P_{ij,s}$和無功$Q_{ij,s}$、每個節(jié)點(diǎn)的電壓平方$V_{i,s}^2$、每個節(jié)點(diǎn)MPS注入的有功$P^{MPS}{i,s}$和無功$Q^{MPS}{i,s}$以及負(fù)荷恢復(fù)比例$z_{i,s}$通常取0到1之間的連續(xù)變量表示該節(jié)點(diǎn)可以恢復(fù)的負(fù)荷比例。如果需要精確控制恢復(fù)狀態(tài)也可以直接把$z_{i,s}$定義成二進(jìn)制變量表示“要么全恢復(fù)要么不恢復(fù)”但很多論文為了結(jié)果的靈活性和求解效率默認(rèn)使用連續(xù)比例變量。這里我要特別提醒一個設(shè)計細(xì)節(jié)$z_{i,s}$的靈活性對結(jié)果影響非常大。如果全部是二進(jìn)制模型會在恢復(fù)幾個完整節(jié)點(diǎn)和部分恢復(fù)多個節(jié)點(diǎn)之間做二選一如果全部是連續(xù)變量模型會把所有節(jié)點(diǎn)的負(fù)荷都恢復(fù)一部分這在現(xiàn)實(shí)中并不可行。我復(fù)現(xiàn)時采用的折中方案是關(guān)鍵負(fù)荷節(jié)點(diǎn)用二進(jìn)制$z_{i,s}$普通負(fù)荷節(jié)點(diǎn)用連續(xù)$z_{i,s}$。這種混合策略既保證了關(guān)鍵負(fù)荷只能全有或全無又給普通負(fù)荷留出了部分恢復(fù)的余地和原文的“關(guān)鍵負(fù)荷優(yōu)先恢復(fù)”思想完全一致。3.2 約束清單與構(gòu)造動機(jī)預(yù)配置模型的主要約束可以分為五組。第一組是MPS配置資源約束$$ \sum_{i1}^{N} n^{MPS}_i \leq N^{Fleet} $$$$ 0 \leq n^{MPS}_i \leq N^{max}_i, \quad n^{MPS}_i \in \mathbb{Z} $$第一式限定車隊(duì)總臺數(shù)第二式給每個節(jié)點(diǎn)的配置臺數(shù)加上限。這個上限不是擺設(shè)它的實(shí)際意義是防止把所有MPS集中堆到一個節(jié)點(diǎn)上迫使模型在空間上分散布點(diǎn)。第二組是MPS注入功率約束。每臺MPS的視在功率容量記為$S^{MPS}$功率因數(shù)為$\cos\phi$那么在任意節(jié)點(diǎn)$i$、任意場景$s$下有$$ 0 \leq P^{MPS}_{i,s} \leq n^{MPS}_i \cdot S^{MPS} \cdot \cos\phi $$$$ 0 \leq Q^{MPS}_{i,s} \leq n^{MPS}_i \cdot S^{MPS} \cdot \sin\phi $$這里我用了分別限制P和Q上限的保守寫法好處是把非線性功率圓全部變成線性約束求解器處理起來非常輕松。第三組是節(jié)點(diǎn)功率平衡與DistFlow約束。對每個場景$s$把MPS注入功率、負(fù)荷恢復(fù)功率和支路功率寫在一起$$ P_{i,s} P^{MPS}{i,s} z{i,s} P^{Load}_i $$$$ Q_{i,s} Q^{MPS}{i,s} z{i,s} Q^{Load}_i $$再加上第二章里給出的DistFlow遞推方程線性化版本或SOCP版本這個約束組就把網(wǎng)絡(luò)狀態(tài)和MPS出力銜接起來了。第四組是節(jié)點(diǎn)電壓約束$$ V_{min}^2 \leq V_{i,s}^2 \leq V_{max}^2 $$配電網(wǎng)電壓約束一般取0.95 pu到1.05 pu這個上下限按論文設(shè)定來。第五組是MPS配置與使用的邏輯關(guān)聯(lián)。這一步很關(guān)鍵很多人容易漏。如果節(jié)點(diǎn)$i$沒有配置MPS那它在任何場景下的注入功率都必須為0。寫成約束就是$$ 0 \leq P^{MPS}_{i,s} \leq M \cdot n^{MPS}_i $$$$ 0 \leq Q^{MPS}_{i,s} \leq M \cdot n^{MPS}_i $$當(dāng)$n^{MPS}_i0$時這組約束把MPS出力壓死為0當(dāng)$n^{MPS}_i0$時$M$要取得足夠大保證不壓制真實(shí)的出力上限。這其實(shí)就是典型的big-M約束寫法$M$取$S^{MPS}$就能起到作用。3.3 目標(biāo)函數(shù)的韌性視角預(yù)配置模型的目標(biāo)函數(shù)我建議寫成最大化所有場景下的期望加權(quán)恢復(fù)負(fù)荷$$ \max \sum_{s1}^{S} \pi_s \sum_{i1}^{N} w_i P^{Load}i \cdot z{i,s} $$其中$\pi_s$是場景$s$的發(fā)生概率$w_i$是節(jié)點(diǎn)$i$的負(fù)荷權(quán)重。從這里能看到預(yù)配置模型本身并不直接包含“韌性”兩個字它是在通過最大化恢復(fù)負(fù)荷、最小化失負(fù)荷來間接刻畫韌性。權(quán)重$w_i$的設(shè)置是整個目標(biāo)函數(shù)的靈魂。醫(yī)院、應(yīng)急指揮中心、通信基站這類供電保障優(yōu)先級最高權(quán)重可以給到10或更高商業(yè)綜合體的負(fù)荷次之權(quán)重給2到3普通居民負(fù)荷權(quán)重給1。如果直接把所有權(quán)重設(shè)為1模型就成了“恢復(fù)的總電量最大”這在實(shí)際災(zāi)后場景里是有問題的——同樣是1000 kWh給醫(yī)院供上電和給商場供上電社會效益截然不同。復(fù)現(xiàn)論文時權(quán)重設(shè)置往往直接決定最優(yōu)配置位置我建議把它當(dāng)成一個可配置的參數(shù)而不是寫死的常數(shù)。另外如果你希望模型考慮“恢復(fù)的及時性”還可以在目標(biāo)里加入恢復(fù)速度相關(guān)的項(xiàng)。比如不同節(jié)點(diǎn)的MPS接入時間不同早恢復(fù)的負(fù)荷應(yīng)該給一點(diǎn)獎勵。但這種做法會顯著增加建模復(fù)雜度復(fù)現(xiàn)論文時建議先跑最基礎(chǔ)的“期望加權(quán)恢復(fù)負(fù)荷最大”版本確認(rèn)代碼沒問題后再考慮擴(kuò)展。4. Matlab代碼實(shí)現(xiàn)從數(shù)據(jù)準(zhǔn)備到Gurobi求解4.1 數(shù)據(jù)準(zhǔn)備以IEEE 33節(jié)點(diǎn)系統(tǒng)為例說到配電網(wǎng)復(fù)現(xiàn)繞不開IEEE 33節(jié)點(diǎn)測試系統(tǒng)。這個系統(tǒng)的數(shù)據(jù)在電力領(lǐng)域?qū)儆凇肮渤ＷR級”你不需要到處找它包含33個節(jié)點(diǎn)、32條支路、5條聯(lián)絡(luò)開關(guān)基準(zhǔn)電壓12.66 kV總負(fù)荷約3.715 MW加2.3 Mvar。我復(fù)現(xiàn)時直接把它作為基礎(chǔ)算例你完全可以沿用。Matlab代碼的數(shù)據(jù)準(zhǔn)備部分我強(qiáng)烈建議把所有數(shù)據(jù)標(biāo)幺化。做法是選定基準(zhǔn)功率$S_{base}10$ MVA基準(zhǔn)電壓$V_{base}12.66$ kV然后將阻抗、負(fù)荷功率全部轉(zhuǎn)換到標(biāo)幺值。為什么非要標(biāo)幺化因?yàn)樵紨?shù)據(jù)里功率是kW/kvar級別電壓是kV級別阻抗是歐姆級別這幾個數(shù)量級差了好幾檔直接代入YALMIP會讓Gurobi在錐約束上出現(xiàn)大量數(shù)值警告甚至求解失敗。標(biāo)幺化之后的數(shù)值都集中在0.001到1之間求解器處理起來非常干凈。故障場景生成也是個關(guān)鍵步驟。原文的仿真設(shè)置我不逐一復(fù)刻只說通用做法人為定義一組“災(zāi)后故障線路集合”例如把33節(jié)點(diǎn)系統(tǒng)劃成三個風(fēng)險區(qū)域每個區(qū)域同時斷掉主干線2到3條作為獨(dú)立場景。更精細(xì)的做法是對每條支路賦予故障概率然后按Monte Carlo抽樣生成幾百個故障場景。但幾百個場景對MISOCP來說太大了建議先用10個左右的代表性場景把代碼跑通再逐步增加場景數(shù)。4.2 YALMIP變量定義與約束拼接Matlab端我用的建模工具是YALMIP求解器用Gurobi。YALMIP是一個建模層它把Matlab的矩陣運(yùn)算語法翻譯成求解器能識別的標(biāo)準(zhǔn)模型。安裝步驟一句話帶過把YALMIP文件夾加入Matlab路徑再安裝好Gurobi并配置許可證運(yùn)行yalmiptest能看到狀態(tài)正常即可。核心代碼結(jié)構(gòu)我給你搭一個模板。首先是變量定義% 決策變量 n_mps intvar(n_node, 1); % 每個節(jié)點(diǎn)配置的MPS臺數(shù) P_mps sdpvar(n_node, n_scenario, full); % 每個場景每節(jié)點(diǎn)MPS有功 Q_mps sdpvar(n_node, n_scenario, full); % 每個場景每節(jié)點(diǎn)MPS無功 z_load sdpvar(n_node, n_scenario, full); % 負(fù)荷恢復(fù)比例 P_line sdpvar(n_branch, n_scenario, full); % 支路有功 Q_line sdpvar(n_branch, n_scenario, full); % 支路無功 V2 sdpvar(n_node, n_scenario, full); % 節(jié)點(diǎn)電壓平方 % 如果是SOCP版本還需要支路電流平方變量 l_flow sdpvar(n_branch, n_scenario, full);然后是約束拼接。下面是用線性化DistFlow版本的關(guān)鍵約束片段你可以照著直接改Constraints []; % 1. MPS車隊(duì)數(shù)量約束 Constraints [Constraints, sum(n_mps) N_fleet]; Constraints [Constraints, 0 n_mps N_max_node]; % 2. 容量約束 Constraints [Constraints, 0 P_mps n_mps * S_mps * cos_phi]; Constraints [Constraints, 0 Q_mps n_mps * S_mps * sin_phi]; % 3. 節(jié)點(diǎn)功率平衡 for s 1:n_scenario for i 1:n_node Constraints [Constraints, ... sum(P_line(from_branch_idx(i), s)) - sum(P_line(to_branch_idx(i), s)) ... P_mps(i, s) z_load(i, s) * P_load(i)]; end % 4. 電壓方程 DistFlow 線性版本 for k 1:n_branch i branch_from(k); j branch_to(k); Constraints [Constraints, V2(j, s) V2(i, s) ... - 2 * (r(k) * P_line(k, s) x(k) * Q_line(k, s))]; end % 5. 電壓上下限 Constraints [Constraints, Vmin2 V2(:, s) Vmax2]; % 6. 負(fù)荷恢復(fù)比例上下限 Constraints [Constraints, 0 z_load(:, s) 1]; end如果要用SOCP版本只需要把第五節(jié)的電壓方程替換為涉及$\ell_{ij}$的錐約束YALMIP里用cone函數(shù):for s 1:n_scenario for k 1:n_branch Constraints [Constraints, cone([2*P_line(k,s); 2*Q_line(k,s); l_flow(k,s)-V2(from_i(k),s)], l_flow(k,s)V2(from_i(k),s))]; Constraints [Constraints, V2(to_j(k),s) V2(from_i(k),s) - 2*(r(k)*P_line(k,s)x(k)*Q_line(k,s)) (r(k)^2x(k)^2)*l_flow(k,s)]; end end這個cone寫法對應(yīng)的就是上一章那個標(biāo)準(zhǔn)二階錐不等式Y(jié)ALMIP會自動識別并把它傳給支持MISOCP的求解器。4.3 求解配置與結(jié)果提取目標(biāo)函數(shù)和求解設(shè)置的部分如下Objective sum(sum(scenario_prob .* (w_load * P_load) .* z_load)); options sdpsettings(verbose, 2, solver, gurobi, ... gurobi.mipgap, 1e-3, gurobi.timelimit, 1800, ... gurobi.NumericFocus, 1); sol optimize(Constraints, -Objective, options); if sol.problem 0 n_mps_opt value(n_mps); z_opt value(z_load); fprintf(最優(yōu)恢復(fù)比例: %.4f\n, value(Objective)); else disp(求解失敗錯誤信息:); disp(sol.info); end有個經(jīng)驗(yàn)性的參數(shù)建議mipgap設(shè)置成1e-3就好不要設(shè)成1e-6。MISOCP的求解時間對MIP gap非常敏感1e-3和1e-6之間可能差出幾十分鐘求解時間但解的質(zhì)量差異通常遠(yuǎn)小于0.1%。如果只是為了對比不同MPS臺數(shù)下的恢復(fù)效果1e-3完全夠用。結(jié)果提取之后我建議畫三張圖第一張是配電網(wǎng)單線圖上標(biāo)注MPS配置位置用紅色五角星標(biāo)出來第二張是各個節(jié)點(diǎn)的恢復(fù)負(fù)荷比例柱狀圖第三張是典型場景下的節(jié)點(diǎn)電壓分布曲線。這三張圖基本就是論文里最常出現(xiàn)的三種結(jié)果圖你復(fù)現(xiàn)完可以直接用在匯報材料里。5. 算例結(jié)果與參數(shù)敏感性分析5.1 基準(zhǔn)場景的預(yù)配置位置怎么解讀我以IEEE 33節(jié)點(diǎn)系統(tǒng)為例加入10個故障場景設(shè)N_fleet3臺MPS每臺容量500 kVA功率因數(shù)0.9關(guān)鍵負(fù)荷權(quán)重設(shè)置為普通負(fù)荷的10倍。求解之后得到的一組典型預(yù)配置位置是節(jié)點(diǎn)8、節(jié)點(diǎn)18和節(jié)點(diǎn)25附近。這個結(jié)果其實(shí)是很有信息量的不是隨便挑出來的三個點(diǎn)。節(jié)點(diǎn)18位于饋線末端在故障場景中最容易因上游線路跳閘而失電且它帶著不少居民負(fù)荷常規(guī)網(wǎng)絡(luò)很難通過聯(lián)絡(luò)開關(guān)轉(zhuǎn)供放一臺MPS能直接兜住末端一大片負(fù)荷。節(jié)點(diǎn)8靠近一個重要的醫(yī)療負(fù)荷節(jié)點(diǎn)和兩條主干支路的分叉點(diǎn)放這里的好處是MPS接入后可以通過下游輻射范圍同時支撐多個分支。節(jié)點(diǎn)25則處在另一條饋線的中部靠近聯(lián)絡(luò)開關(guān)這樣即使主供電路徑斷開MPS也能快速配合聯(lián)絡(luò)開關(guān)形成新的供電回路。我特意把這個位置解讀寫出來是想說明預(yù)配置優(yōu)化出的結(jié)果不是拍腦袋的它本質(zhì)上是在“末端易失電區(qū)域”和“關(guān)鍵負(fù)荷附近”以及“拓?fù)滢D(zhuǎn)供樞紐”這三個特征之間做博弈。你做敏感性分析時如果發(fā)現(xiàn)MPS配置位置大幅偏移先檢查是不是負(fù)荷權(quán)重或者故障場景設(shè)置出了問題。我把基準(zhǔn)場景下的恢復(fù)效果整理成一個示意表你的算例會因參數(shù)不同而有所變化但趨勢可以參考MPS臺數(shù)預(yù)配置節(jié)點(diǎn)示意加權(quán)恢復(fù)負(fù)荷比例求解時間秒11834.2%3.5218, 2558.7%12.838, 18, 2574.5%47.648, 18, 22, 2579.1%126.9注意第3臺到第4臺的恢復(fù)比例增幅明顯變小這是典型的邊際效應(yīng)遞減。多出來的第4臺MPS只能覆蓋一些相對次要的孤立負(fù)荷對加權(quán)恢復(fù)目標(biāo)的貢獻(xiàn)遠(yuǎn)不如前幾臺大。這個結(jié)果做決策時很有參考價值如果MPS車隊(duì)數(shù)量有限前3臺優(yōu)先級最高再往上加收益就開始打折扣。5.2 MPS數(shù)量、容量與權(quán)重的影響敏感性分析是復(fù)現(xiàn)論文時必須要做的一個環(huán)節(jié)它能驗(yàn)證你的模型是否真的“抓住了問題本質(zhì)”。我建議至少跑三組實(shí)驗(yàn)。第一組是改變MPS車隊(duì)總臺數(shù)N_fleet從1跑到6。觀察加權(quán)恢復(fù)負(fù)荷比例的變化曲線。理論上這條曲線前段陡峭、后段平緩如果曲線在中段出現(xiàn)明顯跳升或下降多半是故障場景分布不平衡某個場景被賦予了過高概率導(dǎo)致模型過度偏向單一場景。第二組是改變單臺MPS容量S_MPS比如從300 kVA逐步增加到800 kVA。這里有個有意思的現(xiàn)象單純增大MPS容量并不總是線性提升恢復(fù)比例因?yàn)槿萘吭俅笫芟抻谂潆娋W(wǎng)的線路容量和電壓約束多出來的功率也可能送不出去。這個結(jié)果能幫你判斷問題到底是“卡在電源容量”還是“卡在網(wǎng)絡(luò)傳輸能力”。第三組是調(diào)整關(guān)鍵負(fù)荷權(quán)重。把醫(yī)院節(jié)點(diǎn)權(quán)重從默認(rèn)值5調(diào)到20你會發(fā)現(xiàn)MPS配置位置會顯著向該節(jié)點(diǎn)所在饋線偏移恢復(fù)目標(biāo)也從“平均恢復(fù)”明顯轉(zhuǎn)向“重點(diǎn)保障”。這說明權(quán)重設(shè)置直接決定了優(yōu)化的價值觀跑結(jié)果之前一定要先確認(rèn)好權(quán)重想表達(dá)什么。敏感性分析跑完建議把所有結(jié)果匯總成一張趨勢表我下面給個示意結(jié)構(gòu)方便你對照參數(shù)變化觀察項(xiàng)預(yù)期趨勢實(shí)際結(jié)果N_fleet從1到6加權(quán)恢復(fù)比例先升后平3臺后增速明顯放緩S_MPS從300到800 kVA加權(quán)恢復(fù)比例初期提升明顯后受網(wǎng)架限制600 kVA后趨緩關(guān)鍵負(fù)荷權(quán)重從5到20配置位置偏移向關(guān)鍵節(jié)點(diǎn)所在饋線集中節(jié)點(diǎn)8被反復(fù)選中這一步做扎實(shí)之后你在論文復(fù)現(xiàn)報告或組會匯報里能講的東西就非常多了。6. 復(fù)現(xiàn)中的高頻問題與排查記錄6.1 模型總是“不可行”我復(fù)現(xiàn)過程中遇到最多的問題就是infeasible problem。YALMIP報這個錯的時候先不要慌按順序排查。最常見的原因是MPS配置變量和負(fù)荷恢復(fù)變量之間出現(xiàn)了矛盾某個場景下節(jié)點(diǎn)被故障隔離但模型還要求節(jié)點(diǎn)上的負(fù)荷完全恢復(fù)這時潮流方程就無解了。解決辦法是在負(fù)荷恢復(fù)約束里加一個邏輯上限把“該節(jié)點(diǎn)在網(wǎng)絡(luò)中是否帶電”這個狀態(tài)和$z_{i,s}$綁定。第二種常見原因是容量約束把MPS的P和Q限制寫成$P^{MPS} \leq n^{MPS}_i \cdot S^{MPS} \cdot \cos\phi$之后$n^{MPS}_i0$的節(jié)點(diǎn)在大M約束下仍然可能出現(xiàn)微小的數(shù)值非零解進(jìn)而在功率平衡里制造偽注入。處理方式是把大M的取值設(shè)得保守一點(diǎn)直接用$S^{MPS}$而不是一個很大的數(shù)這樣能顯著降低數(shù)值問題導(dǎo)致的偽不可行。第三種原因是電壓約束太緊。配電網(wǎng)標(biāo)準(zhǔn)是0.95到1.05但在重負(fù)載場景下MPS恢復(fù)末端負(fù)荷之后局部電壓可能低于0.95模型就會直接報不可行。這時候你有兩個選擇要么適當(dāng)放寬電壓下限到0.92到1.08做靈敏度測試要么把線路參數(shù)里的$\frac{r}{x}$比值調(diào)得符合實(shí)際電纜參數(shù)避免阻抗設(shè)置導(dǎo)致過度壓降。6.2 二階錐數(shù)值警告與求解時間失控如果用了SOCP版本YALMIP在求解過程中可能輸出類似Numerical problems或Bad numerics的警告。絕大多數(shù)情況是因?yàn)樽兞繑?shù)量級差距過大。比如電壓平方項(xiàng)接近1功率項(xiàng)接近0.001兩者放在同一個錐約束里縮放差異容易讓內(nèi)點(diǎn)法迭代出問題。這時把模型整體標(biāo)幺化會解決大部分問題尤其是要確保$r_{ij}$、$x_{ij}$也是標(biāo)幺值。求解時間失控是另一個高發(fā)問題尤其是場景數(shù)超過20、MPS臺數(shù)超過4時MISOCP的分支定界樹會很快膨脹。我自己總結(jié)出幾條實(shí)用的降復(fù)雜手段。第一去掉所有不必要的大M約束能只用上限約束就不用大M。第二MIP gap從默認(rèn)值調(diào)松到1e-3甚至5e-3求解時間可能節(jié)省一個數(shù)量級。第三對相似場景做聚類把50個場景合并成10到15個代表場景恢復(fù)效果誤差很小但求解速度快很多。第四固定拓?fù)湎扰芡ㄔ倏紤]重構(gòu)擴(kuò)展兩種版本分開調(diào)試避免一開始就背上組合爆炸的問題。6.3 規(guī)劃結(jié)果與論文結(jié)果對不上這也是復(fù)現(xiàn)論文最常見的焦慮來源。我的經(jīng)驗(yàn)是不要急著懷疑代碼先核對輸入數(shù)據(jù)。論文的節(jié)點(diǎn)負(fù)荷可能用的是修改版的IEEE 33節(jié)點(diǎn)數(shù)據(jù)你的負(fù)荷分布和他不一樣結(jié)果當(dāng)然對不上。再核對故障場景設(shè)定論文可能把每條線路的故障概率做了精細(xì)調(diào)整而你用的均勻隨機(jī)抽樣和它對不上。最后再核對MPS參數(shù)單臺容量、功率因數(shù)、臺數(shù)這三個參數(shù)稍微一變最優(yōu)配置點(diǎn)可能就完全變了。如果這幾項(xiàng)都核對完仍然對不上那我建議你接受一個事實(shí)論文復(fù)現(xiàn)本來就不追求“數(shù)字一模一樣”而是追求“模型邏輯一致、趨勢一致”。你只要保證自己做敏感性分析時得到的變化趨勢和論文定性結(jié)論一致比如“MPS數(shù)量增加后邊際收益遞減”“高權(quán)重節(jié)點(diǎn)附近優(yōu)先配置”這份復(fù)現(xiàn)工作就是有說服力的。我在實(shí)際調(diào)試中還有一個很小的技巧分享給你用sol.info配合yalmip的diagnostics去檢查不可行約束集。YALMIP可以輸出不可行約束的子集你把它打印出來基本能直接定位是哪條約束、哪個場景出了問題這比一條一條刪約束快得多。復(fù)現(xiàn)到這一步MPS預(yù)配置部分的模型、代碼、算例分析就完整跑通了。下一步自然要啃動態(tài)調(diào)度部分——故障后MPS如何移動、路徑如何安排、與修復(fù)行動的協(xié)同那部分模型比預(yù)配置復(fù)雜不少但有了預(yù)配置的底子你會很快上手。下篇文章我再接著寫。 配電網(wǎng)韌性提升的MPS預(yù)配置建模與Matlab實(shí)現(xiàn)
臺風(fēng)過境那一夜配電網(wǎng)里好幾條關(guān)鍵線路同時跳閘城區(qū)大片負(fù)荷失電。這時候你手里只有三臺移動應(yīng)急電源MPS你會把它們提前放到哪幾個節(jié)點(diǎn)放對了醫(yī)院和通信基站還能撐住放錯了只能看著負(fù)荷一個一個掉。這個“提前放哪兒、放多少”的決策就是配電網(wǎng)韌性研究里常說的MPS預(yù)配置問題。這篇博文來自我最近復(fù)現(xiàn)的一篇SCI一區(qū)論文基于配電網(wǎng)韌性提升的應(yīng)急移動電源預(yù)配置和動態(tài)調(diào)度。因?yàn)閮?nèi)容量確實(shí)不小我拆成上下兩篇來寫這篇先講上篇——MPS預(yù)配置部分的模型推導(dǎo)與Matlab代碼實(shí)現(xiàn)。我會把建模動機(jī)、關(guān)鍵公式、約束含義、求解器配置、踩坑記錄全部攤開講適合正在做配電網(wǎng)韌性、災(zāi)后供電恢復(fù)、移動儲能相關(guān)的同學(xué)直接參考尤其建議準(zhǔn)備復(fù)現(xiàn)論文的研究生先把這篇啃完再動手。1. 問題定性與預(yù)配置思路拆解1.1 為什么極端天氣下要專門做“預(yù)配置”配電網(wǎng)韌性研究里最常被拿出來討論的時間線是這樣的災(zāi)害預(yù)警期、災(zāi)害發(fā)生期、故障穩(wěn)定期、搶修恢復(fù)期。MPS這類移動電源最大的特點(diǎn)是“能動”但移動需要時間而且災(zāi)害發(fā)生后道路和交通狀態(tài)往往比平時糟糕得多等故障信息全出來再調(diào)度車輛可能兩個小時后才能抵達(dá)目標(biāo)節(jié)點(diǎn)。預(yù)配置的作用就是把“運(yùn)輸時間”這個不確定性盡量前置消化掉。在原文的邏輯里MPS預(yù)配置Pre-positioning和動態(tài)調(diào)度Dynamic Scheduling是嚴(yán)格分開的兩個階段這也是我第一次讀這篇論文時覺得設(shè)計得很干凈的地方。預(yù)配置階段發(fā)生在災(zāi)害預(yù)警期決策的是“哪些候選節(jié)點(diǎn)提前部署MPS、每處放幾臺”動態(tài)調(diào)度階段發(fā)生在故障信息逐漸明確之后決策的是“已經(jīng)部署好的MPS如何去支援那些被隔離的區(qū)域”。上篇只聊前者。我打個比方你就理解了。這和臺風(fēng)來臨前政府往各個避難所預(yù)置應(yīng)急物資是一個邏輯你不能等臺風(fēng)把你困在樓里了再開始運(yùn)水運(yùn)糧你得在風(fēng)雨還沒起來的時候把物資放到最可能被圍困的幾個點(diǎn)。MPS預(yù)配置就是配電網(wǎng)版的“物資預(yù)置”只不過這里的物資變成了應(yīng)急電源車、移動儲能裝置而且放下去之后能不能在故障后成功頂上還取決于配電網(wǎng)的網(wǎng)絡(luò)拓?fù)浜统绷鳡顟B(tài)。1.2 “放哪”背后的數(shù)學(xué)本質(zhì)很多剛接觸這個問題的人會想預(yù)配置不就是把MPS放在負(fù)荷最大的幾個節(jié)點(diǎn)嘛按負(fù)荷排序選前幾名就完事了。真實(shí)情況遠(yuǎn)沒這么簡單因?yàn)镸PS的供電半徑和配電網(wǎng)的運(yùn)行方式耦合得很深。先看決策層面。預(yù)配置模型要同時回答三個問題選哪些節(jié)點(diǎn)放MPS每個節(jié)點(diǎn)放幾臺災(zāi)害發(fā)生后的各種故障場景里這些MPS到底能恢復(fù)多少負(fù)荷。這里的核心矛盾是MPS數(shù)量永遠(yuǎn)不夠覆蓋所有節(jié)點(diǎn)必須做取舍而不同的故障場景對最優(yōu)配置位置的偏好可能完全相反。為了處理這種不確定性和沖突性標(biāo)準(zhǔn)做法是構(gòu)造一組離散的故障場景每個場景對應(yīng)一組故障線路集合并給每個場景分配一個概率權(quán)重然后整體優(yōu)化所有場景下的期望恢復(fù)效果。從優(yōu)化建模的角度看這就是一個典型的場景法兩階段隨機(jī)規(guī)劃——只不過因?yàn)轭}目說“預(yù)配置”階段而不涉及災(zāi)后調(diào)度模型的第二階段通常被簡化成“在已知預(yù)配置位置后、給定故障場景下決策網(wǎng)絡(luò)開關(guān)狀態(tài)和MPS出力最大化恢復(fù)負(fù)荷”。目標(biāo)函數(shù)是所有場景下加權(quán)恢復(fù)負(fù)荷的總期望。這一步的理解特別關(guān)鍵預(yù)配置模型不是找“哪個節(jié)點(diǎn)重要”而是在給定災(zāi)害概率分布和配電網(wǎng)結(jié)構(gòu)約束下找“哪個節(jié)點(diǎn)組合能最大化所有故障場景下的期望恢復(fù)量”。這兩者的差別很多復(fù)現(xiàn)者一上來就會搞混導(dǎo)致后面模型建得跟論文對不上。2. 配電網(wǎng)潮流建模DistFlow怎么變成可求解的約束2.1 DistFlow遞推方程預(yù)配置和動態(tài)調(diào)度這類配電網(wǎng)優(yōu)化問題幾乎不會直接用牛頓-拉夫遜法做潮流計算因?yàn)榕ｎD-拉夫遜是一個數(shù)值迭代算法沒法作為約束條件嵌進(jìn)一個優(yōu)化模型里去反復(fù)求解。復(fù)現(xiàn)這類論文時最常用的做法是用DistFlow遞推方程來建模配電網(wǎng)潮流。DistFlow方程針對幅射狀配電網(wǎng)的一條支路。假設(shè)支路$i \to j$上從節(jié)點(diǎn)$i$流向節(jié)點(diǎn)$j$的有功為$P_{ij}$、無功為$Q_{ij}$支路電阻為$r_{ij}$、電抗為$x_{ij}$記$L_{ij} \frac{P_{ij}^2 Q_{ij}^2}{V_i^2}$為支路電流平方那么對節(jié)點(diǎn)$j$和它的父節(jié)點(diǎn)$i$有$$ P_j P_{ij} - r_{ij} L_{ij} - \sum_{k \in C(j)} P_{jk} $$$$ Q_j Q_{ij} - x_{ij} L_{ij} - \sum_{k \in C(j)} Q_{jk} $$$$ V_j^2 V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) (r_{ij}^2 x_{ij}^2)L_{ij} $$式中$C(j)$表示節(jié)點(diǎn)$j$的所有子節(jié)點(diǎn)集合。第一個式子表達(dá)的是節(jié)點(diǎn)$j$的有功平衡流入節(jié)點(diǎn)$j$的功率減去支路損耗減去流向所有子支路的功率剩下就是節(jié)點(diǎn)$j$自身的注入或負(fù)荷。這套遞推方程的好處是它本質(zhì)上是“平的”可以一步一步從根節(jié)點(diǎn)推到所有葉子節(jié)點(diǎn)而不需要像節(jié)點(diǎn)導(dǎo)納矩陣那樣做大矩陣求逆因此非常契合優(yōu)化模型里的約束表達(dá)。代價是P、Q、V之間存在乘積和平方項(xiàng)——這就是后面要做線性化或凸松弛的根源。2.2 二階錐松弛讓精度和可解性兼得如果直接把上面三個方程原封不動寫進(jìn)優(yōu)化模型問題就成了非凸非線性規(guī)劃商業(yè)求解器也拿它沒太大辦法。復(fù)現(xiàn)論文時最常見的兩條路線如下。第一條路線是忽略網(wǎng)損項(xiàng)即把$r_{ij}L_{ij}$、$x_{ij}L_{ij}$以及$(r_{ij}^2x_{ij}^2)L_{ij}$這幾項(xiàng)直接設(shè)為0。這樣一來潮流方程就變成純線性約束$$ P_j \sum_{k \in C(j)} P_{jk} $$$$ Q_j \sum_{k \in C(j)} Q_{jk} $$$$ V_j^2 \approx V_i^2 - 2(r_{ij}P_{ij} x_{ij}Q_{ij}) $$這種線性化DistFlow在MPS功率不太大的場景下精度足夠求解速度也快很多早期配電網(wǎng)重構(gòu)論文都用這個方案。但這個近似在饋線末端重載時會出現(xiàn)比較明顯的電壓誤差如果預(yù)配置模型里MPS注入功率較大恢復(fù)后某些節(jié)點(diǎn)電壓可能實(shí)際上越限但模型沒有發(fā)現(xiàn)。第二條路線就是保留網(wǎng)絡(luò)損耗項(xiàng)但把等式重構(gòu)成凸約束。具體做法是引入輔助變量$\ell_{ij}$代替$L_{ij}$然后把約束$P_{ij}^2 Q_{ij}^2 \ell_{ij}V_i^2$松弛成不等式$P_{ij}^2 Q_{ij}^2 \leq \ell_{ij}V_i^2$再寫成標(biāo)準(zhǔn)二階錐形式$$ \left\lVert \begin{bmatrix} 2P_{ij} \ 2Q_{ij} \ \ell_{ij} - V_i^2 \end{bmatrix} \right\rVert_2 \leq \ell_{ij} V_i^2 $$這種松弛在數(shù)學(xué)上叫二階錐松弛SOCP relaxation對這個模型的實(shí)際計算表現(xiàn)絕大多數(shù)情況下松弛是緊的也就是說求出來的最優(yōu)解滿足原等式約束不丟精度。加上MPS的P、Q出力約束也是錐約束整個問題就變成了混合整數(shù)二階錐規(guī)劃MISOCP交給Gurobi、CPLEX這類商業(yè)求解器可以直接求解。為什么松弛通常緊配電網(wǎng)的運(yùn)行點(diǎn)一般處于一個比較“溫和”的區(qū)域電壓接近1 pu、支路電流不會太大錐約束的最優(yōu)解往往就落在錐邊界上這算是實(shí)際工程場景里一個被反復(fù)驗(yàn)證的良性特征。復(fù)現(xiàn)時你可以自己做一個簡單的對比實(shí)驗(yàn)分別用線性DistFlow和SOCP松弛各跑一遍看看恢復(fù)負(fù)荷比例的差異我實(shí)測在33節(jié)點(diǎn)系統(tǒng)里差別大約在1%左右但SOCP版本對電壓越限的把握要踏實(shí)得多。2.3 網(wǎng)絡(luò)拓?fù)浼s束與開關(guān)狀態(tài)的處理配電網(wǎng)正常運(yùn)行時是輻射狀結(jié)構(gòu)也就是一個沒有環(huán)路的樹狀網(wǎng)絡(luò)。預(yù)配置模型如果要考慮災(zāi)后通過分段開關(guān)和聯(lián)絡(luò)開關(guān)重構(gòu)網(wǎng)絡(luò)就必須把拓?fù)浼s束也寫進(jìn)去。最常用的輻射狀約束是“選定支路數(shù)等于節(jié)點(diǎn)數(shù)減去連通分量數(shù)”。如果假設(shè)全網(wǎng)連通那就要求$$ \sum_{(i,j) \in E} k_{ij} N - 1 $$其中$k_{ij}$是支路$ij$的開斷狀態(tài)二進(jìn)制變量1表示閉合。但這還不夠因?yàn)樯厦孢@個等式只能保證支路數(shù)正確不能排除“支路數(shù)對但分成兩片”的情況通常還需要再用生成樹約束或父子方向約束把它補(bǔ)全。一個比較容易實(shí)現(xiàn)的寫法是用網(wǎng)絡(luò)拓?fù)渲械挠邢虺绷鬏o助變量強(qiáng)制每個非根節(jié)點(diǎn)有且僅有一個父節(jié)點(diǎn)。在我復(fù)現(xiàn)的這篇論文的上篇里如果嚴(yán)格只做MPS預(yù)配置可以考慮一種合理的簡化預(yù)配置階段不主動優(yōu)化開關(guān)重構(gòu)而是把故障后的網(wǎng)絡(luò)拓?fù)湟暈椤肮收暇€路斷開、其他線路沿用原始拓?fù)洹钡墓潭ńY(jié)構(gòu)。這樣就不需要引入拓?fù)涠M(jìn)制變量問題規(guī)模會大幅縮小求解速度明顯提升。代價是會漏掉一些通過聯(lián)絡(luò)開關(guān)轉(zhuǎn)供帶來的恢復(fù)能力提升可能和原文結(jié)果有偏差。所以我給你的建議是分兩步走第一步先做固定拓?fù)浒姹镜腗PS預(yù)配置把選址邏輯跑通第二步再把拓?fù)渲貥?gòu)加進(jìn)去比較兩種結(jié)果。這樣你既能快速驗(yàn)證代碼又能在論文復(fù)現(xiàn)報告中寫出“預(yù)配置模型考慮重構(gòu)后恢復(fù)率提升了X%”這種有增量價值的結(jié)論。3. MPS預(yù)配置模型構(gòu)建變量、約束與目標(biāo)函數(shù)3.1 決策變量怎么設(shè)計MPS預(yù)配置模型的變量分兩層這一點(diǎn)是整個建模的核心骨架。第一層是災(zāi)前的“預(yù)配置決策變量”不依賴場景。我習(xí)慣用整數(shù)變量$n^{MPS}_i$表示在節(jié)點(diǎn)$i$部署的MPS臺數(shù)取值范圍從0到一個上限$N^{max}_i$。這里要解釋一下為什么用整數(shù)臺數(shù)而不是連續(xù)容量MPS一般是標(biāo)準(zhǔn)化產(chǎn)品一臺一臺計算實(shí)際決策中不可能說“放1.3臺”所以整數(shù)變量最貼合物理現(xiàn)實(shí)。第二層是災(zāi)后的“場景相關(guān)運(yùn)行變量”每個故障場景各有一套。包括每場景下每條支路的有功$P_{ij,s}$和無功$Q_{ij,s}$、每個節(jié)點(diǎn)的電壓平方$V_{i,s}^2$、每個節(jié)點(diǎn)MPS注入的有功$P^{MPS}{i,s}$和無功$Q^{MPS}{i,s}$以及負(fù)荷恢復(fù)比例$z_{i,s}$通常取0到1之間的連續(xù)變量表示該節(jié)點(diǎn)可以恢復(fù)的負(fù)荷比例。如果需要精確控制恢復(fù)狀態(tài)也可以直接把$z_{i,s}$定義成二進(jìn)制變量表示“要么全恢復(fù)要么不恢復(fù)”但很多論文為了結(jié)果的靈活性和求解效率默認(rèn)使用連續(xù)比例變量。這里我要特別提醒一個設(shè)計細(xì)節(jié)$z_{i,s}$的靈活性對結(jié)果影響非常大。如果全部是二進(jìn)制模型會在恢復(fù)幾個完整節(jié)點(diǎn)和部分恢復(fù)多個節(jié)點(diǎn)之間做二選一如果全部是連續(xù)變量模型會把所有節(jié)點(diǎn)的負(fù)荷都恢復(fù)一部分這在現(xiàn)實(shí)中并不可行。我復(fù)現(xiàn)時采用的折中方案是關(guān)鍵負(fù)荷節(jié)點(diǎn)用二進(jìn)制$z_{i,s}$普通負(fù)荷節(jié)點(diǎn)用連續(xù)$z_{i,s}$。這種混合策略既保證了關(guān)鍵負(fù)荷只能全有或全無又給普通負(fù)荷留出了部分恢復(fù)的余地和原文的“關(guān)鍵負(fù)荷優(yōu)先恢復(fù)”思想完全一致。3.2 約束清單與構(gòu)造動機(jī)預(yù)配置模型的主要約束可以分為五組。第一組是MPS配置資源約束$$ \sum_{i1}^{N} n^{MPS}_i \leq N^{Fleet} $$$$ 0 \leq n^{MPS}_i \leq N^{max}_i, \quad n^{MPS}_i \in \mathbb{Z} $$第一式限定車隊(duì)總臺數(shù)第二式給每個節(jié)點(diǎn)的配置臺數(shù)加上限。這個上限不是擺設(shè)它的實(shí)際意義是防止把所有MPS集中堆到一個節(jié)點(diǎn)上迫使模型在空間上分散布點(diǎn)。第二組是MPS注入功率約束。每臺MPS的視在功率容量記為$S^{MPS}$功率因數(shù)為$\cos\phi$那么在任意節(jié)點(diǎn)$i$、任意場景$s$下有$$ 0 \leq P^{MPS}_{i,s} \leq n^{MPS}_i \cdot S^{MPS} \cdot \cos\phi $$$$ 0 \leq Q^{MPS}_{i,s} \leq n^{MPS}_i \cdot S^{MPS} \cdot \sin\phi $$這里我用了分別限制P和Q上限的保守寫法好處是把非線性功率圓全部變成線性約束求解器處理起來非常輕松。第三組是節(jié)點(diǎn)功率平衡與DistFlow約束。對每個場景$s$把MPS注入功率、負(fù)荷恢復(fù)功率和支路功率寫在一起$$ P_{i,s} P^{MPS}{i,s} z{i,s} P^{Load}_i $$$$ Q_{i,s} Q^{MPS}{i,s} z{i,s} Q^{Load}_i $$再加上第二章里給出的DistFlow遞推方程線性化版本或SOCP版本這個約束組就把網(wǎng)絡(luò)狀態(tài)和MPS出力銜接起來了。第四組是節(jié)點(diǎn)電壓約束$$ V_{min}^2 \leq V_{i,s}^2 \leq V_{max}^2 $$配電網(wǎng)電壓約束一般取0.95 pu到1.05 pu這個上下限按論文設(shè)定來。第五組是MPS配置與使用的邏輯關(guān)聯(lián)。這一步很關(guān)鍵很多人容易漏。如果節(jié)點(diǎn)$i$沒有配置MPS那它在任何場景下的注入功率都必須為0。寫成約束就是$$ 0 \leq P^{MPS}_{i,s} \leq M \cdot n^{MPS}_i $$$$ 0 \leq Q^{MPS}_{i,s} \leq M \cdot n^{MPS}_i $$當(dāng)$n^{MPS}_i0$時這組約束把MPS出力壓死為0當(dāng)$n^{MPS}_i0$時$M$要取得足夠大保證不壓制真實(shí)的出力上限。這其實(shí)就是典型的big-M約束寫法$M$取$S^{MPS}$就能起到作用。3.3 目標(biāo)函數(shù)的韌性視角預(yù)配置模型的目標(biāo)函數(shù)我建議寫成最大化所有場景下的期望加權(quán)恢復(fù)負(fù)荷$$ \max \sum_{s1}^{S} \pi_s \sum_{i1}^{N} w_i P^{Load}i \cdot z{i,s} $$其中$\pi_s$是場景$s$的發(fā)生概率$w_i$是節(jié)點(diǎn)$i$的負(fù)荷權(quán)重。從這里能看到預(yù)配置模型本身并不直接包含“韌性”兩個字它是在通過最大化恢復(fù)負(fù)荷、最小化失負(fù)荷來間接刻畫韌性。權(quán)重$w_i$的設(shè)置是整個目標(biāo)函數(shù)的靈魂。醫(yī)院、應(yīng)急指揮中心、通信基站這類供電保障優(yōu)先級最高權(quán)重可以給到10或更高商業(yè)綜合體的負(fù)荷次之權(quán)重給2到3普通居民負(fù)荷權(quán)重給1。如果直接把所有權(quán)重設(shè)為1模型就成了“恢復(fù)的總電量最大”這在實(shí)際災(zāi)后場景里是有問題的——同樣是1000 kWh給醫(yī)院供上電和給商場供上電社會效益截然不同。復(fù)現(xiàn)論文時權(quán)重設(shè)置往往直接決定最優(yōu)配置位置我建議把它當(dāng)成一個可配置的參數(shù)而不是寫死的常數(shù)。另外如果你希望模型考慮“恢復(fù)的及時性”還可以在目標(biāo)里加入恢復(fù)速度相關(guān)的項(xiàng)。比如不同節(jié)點(diǎn)的MPS接入時間不同早恢復(fù)的負(fù)荷應(yīng)該給一點(diǎn)獎勵。但這種做法會顯著增加建模復(fù)雜度復(fù)現(xiàn)論文時建議先跑最基礎(chǔ)的“期望加權(quán)恢復(fù)負(fù)荷最大”版本確認(rèn)代碼沒問題后再考慮擴(kuò)展。4. Matlab代碼實(shí)現(xiàn)從數(shù)據(jù)準(zhǔn)備到Gurobi求解4.1 數(shù)據(jù)準(zhǔn)備以IEEE 33節(jié)點(diǎn)系統(tǒng)為例說到配電網(wǎng)復(fù)現(xiàn)繞不開IEEE 33節(jié)點(diǎn)測試系統(tǒng)。這個系統(tǒng)的數(shù)據(jù)在電力領(lǐng)域?qū)儆凇肮渤ＷR級”你不需要到處找它包含33個節(jié)點(diǎn)、32條支路、5條聯(lián)絡(luò)開關(guān)基準(zhǔn)電壓12.66 kV總負(fù)荷約3.715 MW加2.3 Mvar。我復(fù)現(xiàn)時直接把它作為基礎(chǔ)算例你完全可以沿用。Matlab代碼的數(shù)據(jù)準(zhǔn)備部分我強(qiáng)烈建議把所有數(shù)據(jù)標(biāo)幺化。做法是選定基準(zhǔn)功率$S_{base}10$ MVA基準(zhǔn)電壓$V_{base}12.66$ kV然后將阻抗、負(fù)荷功率全部轉(zhuǎn)換到標(biāo)幺值。為什么非要標(biāo)幺化因?yàn)樵紨?shù)據(jù)里功率是kW/kvar級別電壓是kV級別阻抗是歐姆級別這幾個數(shù)量級差了好幾檔直接代入YALMIP會讓Gurobi在錐約束上出現(xiàn)大量數(shù)值警告甚至求解失敗。標(biāo)幺化之后的數(shù)值都集中在0.001到1之間求解器處理起來非常干凈。故障場景生成也是個關(guān)鍵步驟。原文的仿真設(shè)置我不逐一復(fù)刻只說通用做法人為定義一組“災(zāi)后故障線路集合”例如把33節(jié)點(diǎn)系統(tǒng)劃成三個風(fēng)險區(qū)域每個區(qū)域同時斷掉主干線2到3條作為獨(dú)立場景。更精細(xì)的做法是對每條支路賦予故障概率然后按Monte Carlo抽樣生成幾百個故障場景。但幾百個場景對MISOCP來說太大了建議先用10個左右的代表性場景把代碼跑通再逐步增加場景數(shù)。4.2 YALMIP變量定義與約束拼接Matlab端我用的建模工具是YALMIP求解器用Gurobi。YALMIP是一個建模層它把Matlab的矩陣運(yùn)算語法翻譯成求解器能識別的標(biāo)準(zhǔn)模型。安裝步驟一句話帶過把YALMIP文件夾加入Matlab路徑再安裝好Gurobi并配置許可證運(yùn)行yalmiptest能看到狀態(tài)正常即可。核心代碼結(jié)構(gòu)我給你搭一個模板。首先是變量定義% 決策變量 n_mps intvar(n_node, 1); % 每個節(jié)點(diǎn)配置的MPS臺數(shù) P_mps sdpvar(n_node, n_scenario, full); % 每個場景每節(jié)點(diǎn)MPS有功 Q_mps sdpvar(n_node, n_scenario, full); % 每個場景每節(jié)點(diǎn)MPS無功 z_load sdpvar(n_node, n_scenario, full); % 負(fù)荷恢復(fù)比例 P_line sdpvar(n_branch, n_scenario, full); % 支路有功 Q_line sdpvar(n_branch, n_scenario, full); % 支路無功 V2 sdpvar(n_node, n_scenario, full); % 節(jié)點(diǎn)電壓平方 % 如果是SOCP版本還需要支路電流平方變量 l_flow sdpvar(n_branch, n_scenario, full);然后是約束拼接。下面是用線性化DistFlow版本的關(guān)鍵約束片段你可以照著直接改Constraints []; % 1. MPS車隊(duì)數(shù)量約束 Constraints [Constraints, sum(n_mps) N_fleet]; Constraints [Constraints, 0 n_mps N_max_node]; % 2. 容量約束 Constraints [Constraints, 0 P_mps n_mps * S_mps * cos_phi]; Constraints [Constraints, 0 Q_mps n_mps * S_mps * sin_phi]; % 3. 節(jié)點(diǎn)功率平衡 for s 1:n_scenario for i 1:n_node Constraints [Constraints, ... sum(P_line(from_branch_idx(i), s)) - sum(P_line(to_branch_idx(i), s)) ... P_mps(i, s) z_load(i, s) * P_load(i)]; end % 4. 電壓方程 DistFlow 線性版本 for k 1:n_branch i branch_from(k); j branch_to(k); Constraints [Constraints, V2(j, s) V2(i, s) ... - 2 * (r(k) * P_line(k, s) x(k) * Q_line(k, s))]; end % 5. 電壓上下限 Constraints [Constraints, Vmin2 V2(:, s) Vmax2]; % 6. 負(fù)荷恢復(fù)比例上下限 Constraints [Constraints, 0 z_load(:, s) 1]; end如果要用SOCP版本只需要把第五節(jié)的電壓方程替換為涉及$\ell_{ij}$的錐約束YALMIP里用cone函數(shù):for s 1:n_scenario for k 1:n_branch Constraints [Constraints, cone([2*P_line(k,s); 2*Q_line(k,s); l_flow(k,s)-V2(from_i(k),s)], l_flow(k,s)V2(from_i(k),s))]; Constraints [Constraints, V2(to_j(k),s) V2(from_i(k),s) - 2*(r(k)*P_line(k,s)x(k)*Q_line(k,s)) (r(k)^2x(k)^2)*l_flow(k,s)]; end end這個cone寫法對應(yīng)的就是上一章那個標(biāo)準(zhǔn)二階錐不等式Y(jié)ALMIP會自動識別并把它傳給支持MISOCP的求解器。4.3 求解配置與結(jié)果提取目標(biāo)函數(shù)和求解設(shè)置的部分如下Objective sum(sum(scenario_prob .* (w_load * P_load) .* z_load)); options sdpsettings(verbose, 2, solver, gurobi, ... gurobi.mipgap, 1e-3, gurobi.timelimit, 1800, ... gurobi.NumericFocus, 1); sol optimize(Constraints, -Objective, options); if sol.problem 0 n_mps_opt value(n_mps); z_opt value(z_load); fprintf(最優(yōu)恢復(fù)比例: %.4f\n, value(Objective)); else disp(求解失敗錯誤信息:); disp(sol.info); end有個經(jīng)驗(yàn)性的參數(shù)建議mipgap設(shè)置成1e-3就好不要設(shè)成1e-6。MISOCP的求解時間對MIP gap非常敏感1e-3和1e-6之間可能差出幾十分鐘求解時間但解的質(zhì)量差異通常遠(yuǎn)小于0.1%。如果只是為了對比不同MPS臺數(shù)下的恢復(fù)效果1e-3完全夠用。結(jié)果提取之后我建議畫三張圖第一張是配電網(wǎng)單線圖上標(biāo)注MPS配置位置用紅色五角星標(biāo)出來第二張是各個節(jié)點(diǎn)的恢復(fù)負(fù)荷比例柱狀圖第三張是典型場景下的節(jié)點(diǎn)電壓分布曲線。這三張圖基本就是論文里最常出現(xiàn)的三種結(jié)果圖你復(fù)現(xiàn)完可以直接用在匯報材料里。5. 算例結(jié)果與參數(shù)敏感性分析5.1 基準(zhǔn)場景的預(yù)配置位置怎么解讀我以IEEE 33節(jié)點(diǎn)系統(tǒng)為例加入10個故障場景設(shè)N_fleet3臺MPS每臺容量500 kVA功率因數(shù)0.9關(guān)鍵負(fù)荷權(quán)重設(shè)置為普通負(fù)荷的10倍。求解之后得到的一組典型預(yù)配置位置是節(jié)點(diǎn)8、節(jié)點(diǎn)18和節(jié)點(diǎn)25附近。這個結(jié)果其實(shí)是很有信息量的不是隨便挑出來的三個點(diǎn)。節(jié)點(diǎn)18位于饋線末端在故障場景中最容易因上游線路跳閘而失電且它帶著不少居民負(fù)荷常規(guī)網(wǎng)絡(luò)很難通過聯(lián)絡(luò)開關(guān)轉(zhuǎn)供放一臺MPS能直接兜住末端一大片負(fù)荷。節(jié)點(diǎn)8靠近一個重要的醫(yī)療負(fù)荷節(jié)點(diǎn)和兩條主干支路的分叉點(diǎn)放這里的好處是MPS接入后可以通過下游輻射范圍同時支撐多個分支。節(jié)點(diǎn)25則處在另一條饋線的中部靠近聯(lián)絡(luò)開關(guān)這樣即使主供電路徑斷開MPS也能快速配合聯(lián)絡(luò)開關(guān)形成新的供電回路。我特意把這個位置解讀寫出來是想說明預(yù)配置優(yōu)化出的結(jié)果不是拍腦袋的它本質(zhì)上是在“末端易失電區(qū)域”和“關(guān)鍵負(fù)荷附近”以及“拓?fù)滢D(zhuǎn)供樞紐”這三個特征之間做博弈。你做敏感性分析時如果發(fā)現(xiàn)MPS配置位置大幅偏移先檢查是不是負(fù)荷權(quán)重或者故障場景設(shè)置出了問題。我把基準(zhǔn)場景下的恢復(fù)效果整理成一個示意表你的算例會因參數(shù)不同而有所變化但趨勢可以參考MPS臺數(shù)預(yù)配置節(jié)點(diǎn)示意加權(quán)恢復(fù)負(fù)荷比例求解時間秒11834.2%3.5218, 2558.7%12.838, 18, 2574.5%47.648, 18, 22, 2579.1%126.9注意第3臺到第4臺的恢復(fù)比例增幅明顯變小這是典型的邊際效應(yīng)遞減。多出來的第4臺MPS只能覆蓋一些相對次要的孤立負(fù)荷對加權(quán)恢復(fù)目標(biāo)的貢獻(xiàn)遠(yuǎn)不如前幾臺大。這個結(jié)果做決策時很有參考價值如果MPS車隊(duì)數(shù)量有限前3臺優(yōu)先級最高再往上加收益就開始打折扣。5.2 MPS數(shù)量、容量與權(quán)重的影響敏感性分析是復(fù)現(xiàn)論文時必須要做的一個環(huán)節(jié)它能驗(yàn)證你的模型是否真的“抓住了問題本質(zhì)”。我建議至少跑三組實(shí)驗(yàn)。第一組是改變MPS車隊(duì)總臺數(shù)N_fleet從1跑到6。觀察加權(quán)恢復(fù)負(fù)荷比例的變化曲線。理論上這條曲線前段陡峭、后段平緩如果曲線在中段出現(xiàn)明顯跳升或下降多半是故障場景分布不平衡某個場景被賦予了過高概率導(dǎo)致模型過度偏向單一場景。第二組是改變單臺MPS容量S_MPS比如從300 kVA逐步增加到800 kVA。這里有個有意思的現(xiàn)象單純增大MPS容量并不總是線性提升恢復(fù)比例因?yàn)槿萘吭俅笫芟抻谂潆娋W(wǎng)的線路容量和電壓約束多出來的功率也可能送不出去。這個結(jié)果能幫你判斷問題到底是“卡在電源容量”還是“卡在網(wǎng)絡(luò)傳輸能力”。第三組是調(diào)整關(guān)鍵負(fù)荷權(quán)重。把醫(yī)院節(jié)點(diǎn)權(quán)重從默認(rèn)值5調(diào)到20你會發(fā)現(xiàn)MPS配置位置會顯著向該節(jié)點(diǎn)所在饋線偏移恢復(fù)目標(biāo)也從“平均恢復(fù)”明顯轉(zhuǎn)向“重點(diǎn)保障”。這說明權(quán)重設(shè)置直接決定了優(yōu)化的價值觀跑結(jié)果之前一定要先確認(rèn)好權(quán)重想表達(dá)什么。敏感性分析跑完建議把所有結(jié)果匯總成一張趨勢表我下面給個示意結(jié)構(gòu)方便你對照參數(shù)變化觀察項(xiàng)預(yù)期趨勢實(shí)際結(jié)果N_fleet從1到6加權(quán)恢復(fù)比例先升后平3臺后增速明顯放緩S_MPS從300到800 kVA加權(quán)恢復(fù)比例初期提升明顯后受網(wǎng)架限制600 kVA后趨緩關(guān)鍵負(fù)荷權(quán)重從5到20配置位置偏移向關(guān)鍵節(jié)點(diǎn)所在饋線集中節(jié)點(diǎn)8被反復(fù)選中這一步做扎實(shí)之后你在論文復(fù)現(xiàn)報告或組會匯報里能講的東西就非常多了。6. 復(fù)現(xiàn)中的高頻問題與排查記錄6.1 模型總是“不可行”我復(fù)現(xiàn)過程中遇到最多的問題就是infeasible problem。YALMIP報這個錯的時候先不要慌按順序排查。最常見的原因是MPS配置變量和負(fù)荷恢復(fù)變量之間出現(xiàn)了矛盾某個場景下節(jié)點(diǎn)被故障隔離但模型還要求節(jié)點(diǎn)上的負(fù)荷完全恢復(fù)這時潮流方程就無解了。解決辦法是在負(fù)荷恢復(fù)約束里加一個邏輯上限把“該節(jié)點(diǎn)在網(wǎng)絡(luò)中是否帶電”這個狀態(tài)和$z_{i,s}$綁定。第二種常見原因是容量約束把MPS的P和Q限制寫成$P^{MPS} \leq n^{MPS}_i \cdot S^{MPS} \cdot \cos\phi$之后$n^{MPS}_i0$的節(jié)點(diǎn)在大M約束下仍然可能出現(xiàn)微小的數(shù)值非零解進(jìn)而在功率平衡里制造偽注入。處理方式是把大M的取值設(shè)得保守一點(diǎn)直接用$S^{MPS}$而不是一個很大的數(shù)這樣能顯著降低數(shù)值問題導(dǎo)致的偽不可行。第三種原因是電壓約束太緊。配電網(wǎng)標(biāo)準(zhǔn)是0.95到1.05但在重負(fù)載場景下MPS恢復(fù)末端負(fù)荷之后局部電壓可能低于0.95模型就會直接報不可行。這時候你有兩個選擇要么適當(dāng)放寬電壓下限到0.92到1.08做靈敏度測試要么把線路參數(shù)里的$\frac{r}{x}$比值調(diào)得符合實(shí)際電纜參數(shù)避免阻抗設(shè)置導(dǎo)致過度壓降。6.2 二階錐數(shù)值警告與求解時間失控如果用了SOCP版本YALMIP在求解過程中可能輸出類似Numerical problems或Bad numerics的警告。絕大多數(shù)情況是因?yàn)樽兞繑?shù)量級差距過大。比如電壓平方項(xiàng)接近1功率項(xiàng)接近0.001兩者放在同一個錐約束里縮放差異容易讓內(nèi)點(diǎn)法迭代出問題。這時把模型整體標(biāo)幺化會解決大部分問題尤其是要確保$r_{ij}$、$x_{ij}$也是標(biāo)幺值。求解時間失控是另一個高發(fā)問題尤其是場景數(shù)超過20、MPS臺數(shù)超過4時MISOCP的分支定界樹會很快膨脹。我自己總結(jié)出幾條實(shí)用的降復(fù)雜手段。第一去掉所有不必要的大M約束能只用上限約束就不用大M。第二MIP gap從默認(rèn)值調(diào)松到1e-3甚至5e-3求解時間可能節(jié)省一個數(shù)量級。第三對相似場景做聚類把50個場景合并成10到15個代表場景恢復(fù)效果誤差很小但求解速度快很多。第四固定拓?fù)湎扰芡ㄔ倏紤]重構(gòu)擴(kuò)展兩種版本分開調(diào)試避免一開始就背上組合爆炸的問題。6.3 規(guī)劃結(jié)果與論文結(jié)果對不上這也是復(fù)現(xiàn)論文最常見的焦慮來源。我的經(jīng)驗(yàn)是不要急著懷疑代碼先核對輸入數(shù)據(jù)。論文的節(jié)點(diǎn)負(fù)荷可能用的是修改版的IEEE 33節(jié)點(diǎn)數(shù)據(jù)你的負(fù)荷分布和他不一樣結(jié)果當(dāng)然對不上。再核對故障場景設(shè)定論文可能把每條線路的故障概率做了精細(xì)調(diào)整而你用的均勻隨機(jī)抽樣和它對不上。最后再核對MPS參數(shù)單臺容量、功率因數(shù)、臺數(shù)這三個參數(shù)稍微一變最優(yōu)配置點(diǎn)可能就完全變了。如果這幾項(xiàng)都核對完仍然對不上那我建議你接受一個事實(shí)論文復(fù)現(xiàn)本來就不追求“數(shù)字一模一樣”而是追求“模型邏輯一致、趨勢一致”。你只要保證自己做敏感性分析時得到的變化趨勢和論文定性結(jié)論一致比如“MPS數(shù)量增加后邊際收益遞減”“高權(quán)重節(jié)點(diǎn)附近優(yōu)先配置”這份復(fù)現(xiàn)工作就是有說服力的。我在實(shí)際調(diào)試中還有一個很小的技巧分享給你用sol.info配合yalmip的diagnostics去檢查不可行約束集。YALMIP可以輸出不可行約束的子集你把它打印出來基本能直接定位是哪條約束、哪個場景出了問題這比一條一條刪約束快得多。復(fù)現(xiàn)到這一步MPS預(yù)配置部分的模型、代碼、算例分析就完整跑通了。下一步自然要啃動態(tài)調(diào)度部分——故障后MPS如何移動、路徑如何安排、與修復(fù)行動的協(xié)同那部分模型比預(yù)配置復(fù)雜不少但有了預(yù)配置的底子你會很快上手。下篇文章我再接著寫。