1. 項(xiàng)目概述當(dāng)配電網(wǎng)遇上“不確定”的分布式電源搞電力系統(tǒng)研究或者做配網(wǎng)規(guī)劃的朋友最近幾年肯定沒少被“分布式電源”和“不確定性”這兩個(gè)詞刷屏。風(fēng)光這些新能源是好綠色環(huán)保但它們的出力完全看天吃飯今天艷陽高照光伏滿發(fā)明天陰雨連綿功率驟降這種間歇性和波動(dòng)性給傳統(tǒng)的配電網(wǎng)運(yùn)行帶來了巨大挑戰(zhàn)。想象一下你精心設(shè)計(jì)好的電網(wǎng)運(yùn)行方式可能因?yàn)橐黄骑h過整個(gè)系統(tǒng)的潮流、電壓就全亂套了嚴(yán)重時(shí)甚至導(dǎo)致設(shè)備過載、供電中斷。所以我們面臨的核心問題就是如何在分布式電源出力“說不準(zhǔn)”的情況下依然能保證配電網(wǎng)安全、可靠、經(jīng)濟(jì)地運(yùn)行這就是“考慮分布式電源不確定性的配電網(wǎng)魯棒動(dòng)態(tài)重構(gòu)”這個(gè)課題要啃的硬骨頭。它不是一個(gè)靜態(tài)的、一次性的優(yōu)化而是動(dòng)態(tài)的、適應(yīng)性的決策過程。重構(gòu)簡單說就是通過調(diào)整配電網(wǎng)中分段開關(guān)和聯(lián)絡(luò)開關(guān)的狀態(tài)改變網(wǎng)絡(luò)的拓?fù)浣Y(jié)構(gòu)就像給交通網(wǎng)絡(luò)動(dòng)態(tài)調(diào)整紅綠燈和單行線以達(dá)到降低網(wǎng)損、平衡負(fù)荷、消除過載、提高供電質(zhì)量等目的。而“魯棒”Robust這個(gè)詞是精髓所在它追求的不是在某種理想場景下的最優(yōu)而是在各種可能出現(xiàn)的糟糕情況比如分布式電源出力最低或最高下系統(tǒng)方案都“扛得住”性能不會(huì)變得太差保證最基本的安全運(yùn)行底線。這個(gè)項(xiàng)目用Matlab來實(shí)現(xiàn)再合適不過。Matlab強(qiáng)大的矩陣運(yùn)算能力、豐富的優(yōu)化工具箱像YALMIP、CPLEX/Gurobi接口以及便捷的數(shù)據(jù)可視化功能讓它成為解決這類復(fù)雜優(yōu)化問題的利器。復(fù)現(xiàn)一篇EI級(jí)別的論文不僅能讓我們深入理解魯棒優(yōu)化理論在電力系統(tǒng)的應(yīng)用更能掌握一套從建模、求解到分析的全流程實(shí)戰(zhàn)技能。接下來我就把自己在復(fù)現(xiàn)過程中的思路、關(guān)鍵步驟、踩過的坑以及一些心得毫無保留地分享給大家。2. 核心思路與魯棒優(yōu)化模型拆解傳統(tǒng)的配電網(wǎng)重構(gòu)優(yōu)化通常假設(shè)所有參數(shù)如負(fù)荷、電源出力都是確定已知的建立一個(gè)確定性的數(shù)學(xué)模型比如混合整數(shù)非線性規(guī)劃模型來求解最優(yōu)開關(guān)組合。但一旦引入分布式電源的不確定性這招就不靈了。因?yàn)槟銦o法預(yù)測下一刻風(fēng)光的具體出力值。應(yīng)對(duì)不確定性主流思路有兩大類隨機(jī)優(yōu)化和魯棒優(yōu)化。隨機(jī)優(yōu)化需要知道不確定量的精確概率分布比如風(fēng)速、光照強(qiáng)度的歷史統(tǒng)計(jì)分布然后去優(yōu)化期望成本或風(fēng)險(xiǎn)成本。但很多時(shí)候獲取精確的概率分布非常困難或者即使有計(jì)算“期望”也意味著接受了小概率極端事件可能帶來的巨大風(fēng)險(xiǎn)。魯棒優(yōu)化則采用了另一種哲學(xué)它不假設(shè)概率分布而是定義一個(gè)不確定集這個(gè)集合包含了所有可能的不確定性實(shí)現(xiàn)比如光伏出力在預(yù)測值的±30%范圍內(nèi)波動(dòng)。優(yōu)化目標(biāo)是在這個(gè)不確定集內(nèi)的最壞情況下我的方案仍然可行且性能最好或者說最壞情況下的損失最小。這就像為你的系統(tǒng)穿上了一層“盔甲”專防各種意外襲擊。2.1 不確定性的數(shù)學(xué)描述盒式不確定集在配電網(wǎng)動(dòng)態(tài)重構(gòu)中我們主要關(guān)心分布式電源出力的不確定性。最常用也最直觀的模型是盒式不確定集。假設(shè)第i個(gè)分布式電源在t時(shí)刻的預(yù)測出力為P_{DG,i,t}^{forecast}其實(shí)際出力P_{DG,i,t}在一個(gè)區(qū)間內(nèi)波動(dòng)[ P_{DG,i,t} \in [P_{DG,i,t}^{forecast} - \hat{P}{DG,i,t}, \quad P{DG,i,t}^{forecast} \hat{P}_{DG,i,t}] ]其中(\hat{P}_{DG,i,t}) 是最大預(yù)測偏差。所有分布式電源在所有時(shí)段的偏差組合起來就構(gòu)成了一個(gè)多維的“盒子”。但如果我們?cè)试S所有電源同時(shí)、同方向地達(dá)到偏差極限這個(gè)不確定集會(huì)過于保守求出的方案可能成本極高卻不實(shí)用。因此通常會(huì)引入一個(gè)預(yù)算參數(shù)Γ它限制了所有偏差同時(shí)達(dá)到極端的程度使不確定集更符合現(xiàn)實(shí)。例如Γ2表示在所有的不確定變量中最多只有2個(gè)能同時(shí)取到其偏差邊界值。這被稱為預(yù)算不確定集是平衡保守性與經(jīng)濟(jì)性的關(guān)鍵。2.2 魯棒動(dòng)態(tài)重構(gòu)模型框架我們的目標(biāo)是在一個(gè)調(diào)度周期內(nèi)比如24小時(shí)以1小時(shí)為間隔動(dòng)態(tài)地決定每個(gè)時(shí)段配電網(wǎng)的開關(guān)狀態(tài)即網(wǎng)絡(luò)結(jié)構(gòu)以最小化總運(yùn)行成本主要是網(wǎng)損成本同時(shí)必須滿足在各種可能的風(fēng)光出力場景下系統(tǒng)的物理約束都不被破壞。這形成了一個(gè)兩階段魯棒優(yōu)化問題第一階段“這里-現(xiàn)在”決策開關(guān)操作決策。這些是整數(shù)變量0或1代表開關(guān)分合必須在不確定性揭示之前做出并且要考慮開關(guān)動(dòng)作次數(shù)限制設(shè)備壽命和操作成本。第二階段“等待-觀望”決策潮流平衡、電壓、電流等連續(xù)變量。這些可以在不確定性實(shí)現(xiàn)后即知道了實(shí)際風(fēng)光出力后再調(diào)整但必須應(yīng)對(duì)最壞情況下的不確定性。模型的核心約束包括輻射狀網(wǎng)絡(luò)約束配電網(wǎng)必須保持開環(huán)、輻射狀運(yùn)行即無環(huán)、全連通。這通常通過虛擬流、生成樹等方法來建模。潮流約束采用DistFlow配電潮流方程或其線性化版本如LinDistFlow來描述功率平衡和電壓降落關(guān)系。魯棒優(yōu)化中這些約束要對(duì)所有不確定集內(nèi)的場景成立。運(yùn)行安全約束節(jié)點(diǎn)電壓必須在允許范圍內(nèi)如0.95~1.05 p.u.支路電流不能超過熱穩(wěn)定極限。開關(guān)操作約束相鄰時(shí)段開關(guān)狀態(tài)變化次數(shù)有限制模擬實(shí)際操作頻率。目標(biāo)函數(shù)最小化最壞情況下的總網(wǎng)損成本。網(wǎng)損與電流平方成正比是一個(gè)非線性項(xiàng)通常需要線性化處理比如分段線性化或使用二階錐松弛技術(shù)將其轉(zhuǎn)化為可求解的形式。注意這里有一個(gè)重要的權(quán)衡。使用精確的DistFlow模型是非凸非線性的求解極其困難。而LinDistFlow模型忽略了支路損耗對(duì)電壓的影響是一個(gè)線性模型大大降低了求解難度在電壓水平較高的配網(wǎng)中近似效果很好。在魯棒優(yōu)化框架下為了能求解我們通常首選LinDistFlow或其改進(jìn)版本作為基礎(chǔ)潮流模型。2.3 為什么用魯棒優(yōu)化而非隨機(jī)優(yōu)化在我復(fù)現(xiàn)和對(duì)比的過程中深刻體會(huì)到兩者思維的不同。隨機(jī)優(yōu)化像是在“賭概率”期望值好看但萬一碰上百年一遇的極端天氣系統(tǒng)可能直接崩潰。而魯棒優(yōu)化是“做最壞的打算爭取最好的結(jié)果”它給出的方案天生具有可靠性溢價(jià)。對(duì)于電網(wǎng)這種關(guān)系到國計(jì)民生的關(guān)鍵基礎(chǔ)設(shè)施很多時(shí)候可靠性比經(jīng)濟(jì)性那一點(diǎn)點(diǎn)的優(yōu)化更重要。魯棒優(yōu)化方案可能平均成本略高但它為你規(guī)避了災(zāi)難性風(fēng)險(xiǎn)這對(duì)于高比例新能源接入的配電網(wǎng)來說是規(guī)劃者和運(yùn)行人員更看重的屬性。3. 基于Matlab的求解算法列與約束生成CCG兩階段魯棒優(yōu)化模型是一個(gè)min-max-min結(jié)構(gòu)的三層優(yōu)化問題外層最小化成本中層最大化不確定性以制造最壞場景內(nèi)層最小化運(yùn)行調(diào)整直接求解是NP難問題。目前最有效的求解算法之一是列與約束生成算法。CCG算法將原問題分解為主問題和子問題通過迭代求解來逼近原問題的最優(yōu)解。理解這個(gè)算法的流程對(duì)于用Matlab實(shí)現(xiàn)至關(guān)重要。3.1 算法流程詳解初始化設(shè)定一個(gè)初始的最壞場景比如所有分布式電源出力均取預(yù)測值設(shè)定迭代次數(shù)k0上界UB∞下界LB-∞。求解主問題Master Problem, MP輸入歷史迭代中子問題發(fā)現(xiàn)的一系列“最壞場景”集合。模型在主問題中針對(duì)每一個(gè)已知的最壞場景都有一套完整的第二階段潮流、電壓等變量和約束。但第一階段變量開關(guān)狀態(tài)是唯一的、共同的。目標(biāo)尋找一組開關(guān)狀態(tài)使得在所有這些已知的最壞場景下系統(tǒng)的最大運(yùn)行成本網(wǎng)損最小。同時(shí)滿足所有場景下的約束。輸出得到當(dāng)前最優(yōu)的開關(guān)狀態(tài)決策以及該決策下針對(duì)已知場景的最小最大成本這個(gè)值作為新的下界LB。因?yàn)橹鲉栴}只考慮了部分場景所以它的解是原問題的一個(gè)樂觀估計(jì)下界。求解子問題Subproblem, SP輸入主問題求出的固定的開關(guān)狀態(tài)。模型在給定的網(wǎng)絡(luò)結(jié)構(gòu)下子問題去尋找一個(gè)最壞的不確定性場景即分布式電源的出力組合使得該場景下的系統(tǒng)運(yùn)行成本網(wǎng)損最大化同時(shí)還要檢查在這個(gè)最壞場景下是否存在可行的潮流解即第二階段變量能否滿足所有約束。這通常是一個(gè)雙線性規(guī)劃或混合整數(shù)線性規(guī)劃問題取決于不確定集和模型線性化程度。輸出得到一個(gè)新的“最壞場景”及其對(duì)應(yīng)的最大成本。這個(gè)成本值作為新的上界UB的候選因?yàn)檫@是某個(gè)可行方案下的實(shí)際最大成本。同時(shí)如果這個(gè)最大成本超過了當(dāng)前LB說明這個(gè)新場景對(duì)主問題構(gòu)成了威脅需要將其“場景”信息添加到主問題的場景集合中。收斂判斷計(jì)算間隙( gap (UB - LB) / |UB| )。如果 ( gap \leq \epsilon )ε是一個(gè)很小的正數(shù)如1e-3則算法收斂輸出當(dāng)前主問題的解作為魯棒最優(yōu)解。否則將子問題產(chǎn)生的新最壞場景添加到主問題的場景列表中令kk1返回步驟2。這個(gè)算法的精妙之處在于它通過不斷將子問題發(fā)現(xiàn)的“威脅性”場景加入主問題迫使主問題開關(guān)決策去適應(yīng)這些越來越壞的情況從而逐步逼近真正的魯棒解。3.2 Matlab實(shí)現(xiàn)中的關(guān)鍵技巧在Matlab中實(shí)現(xiàn)CCG核心是高效地構(gòu)建和求解主問題與子問題的優(yōu)化模型。我強(qiáng)烈推薦使用YALMIP建模語言 Gurobi或CPLEX求解器。YALMIP它讓你可以用近乎數(shù)學(xué)公式的方式描述優(yōu)化問題自動(dòng)處理變量定義、約束添加和目標(biāo)函數(shù)設(shè)置大大降低了建模復(fù)雜度。它就像一個(gè)翻譯官把你的數(shù)學(xué)模型“翻譯”成求解器能聽懂的語言。Gurobi/CPLEX商業(yè)級(jí)的高性能數(shù)學(xué)規(guī)劃求解器尤其擅長處理混合整數(shù)線性規(guī)劃問題速度和穩(wěn)定性遠(yuǎn)超Matlab自帶的intlinprog。實(shí)現(xiàn)步驟簡述數(shù)據(jù)準(zhǔn)備讀取配電網(wǎng)參數(shù)節(jié)點(diǎn)、支路、阻抗、負(fù)荷、分布式電源預(yù)測出力及波動(dòng)范圍、開關(guān)操作成本、電價(jià)等。定義變量二進(jìn)制變量開關(guān)狀態(tài)每個(gè)時(shí)段、每個(gè)開關(guān)。連續(xù)變量節(jié)點(diǎn)電壓幅值平方、支路功率流、分布式電源實(shí)際出力不確定變量等。構(gòu)建主問題函數(shù)函數(shù)輸入是“已知最壞場景集合”輸出是優(yōu)化模型對(duì)象。在函數(shù)內(nèi)部使用YALMIP為集合中的每個(gè)場景創(chuàng)建一套第二階段的變量和約束但共享第一階段的開關(guān)變量。構(gòu)建子問題函數(shù)函數(shù)輸入是“固定的開關(guān)狀態(tài)”輸出是優(yōu)化模型對(duì)象。這里需要利用對(duì)偶理論或Karush-Kuhn-Tucker條件將內(nèi)層的min-max問題轉(zhuǎn)化為一個(gè)單層的最大化問題這是求解的關(guān)鍵也是難點(diǎn)。對(duì)于線性模型子問題通常可以轉(zhuǎn)化為一個(gè)混合整數(shù)線性規(guī)劃。迭代循環(huán)編寫一個(gè)while循環(huán)依次調(diào)用主問題函數(shù)和子問題函數(shù)更新上下界和場景集合直到滿足收斂條件。結(jié)果解析與可視化提取最優(yōu)開關(guān)序列、各時(shí)段網(wǎng)損、電壓分布并繪制成圖。實(shí)操心得在編寫子問題時(shí)處理不確定集約束特別是預(yù)算約束Γ要格外小心。一種有效的方法是將分布式電源的實(shí)際出力表示為P_forecast ζ * P_hat其中ζ是一個(gè)在[-1,1]之間變化的連續(xù)變量然后通過約束sum(abs(ζ)) Γ來施加預(yù)算限制。YALMIP可以直接處理這種絕對(duì)值約束。4. 完整復(fù)現(xiàn)流程與代碼模塊解析下面我將以一個(gè)標(biāo)準(zhǔn)的33節(jié)點(diǎn)配電網(wǎng)系統(tǒng)為例拆解完整的Matlab復(fù)現(xiàn)流程。假設(shè)系統(tǒng)接入了3個(gè)光伏電源調(diào)度周期為24小時(shí)。4.1 數(shù)據(jù)準(zhǔn)備與預(yù)處理模塊首先需要一個(gè)清晰的數(shù)據(jù)結(jié)構(gòu)。我習(xí)慣創(chuàng)建一個(gè)network_data.m腳本定義所有基礎(chǔ)參數(shù)。% network_data.m % 1. 網(wǎng)絡(luò)拓?fù)?bus_data importdata(33bus_data.txt); % 節(jié)點(diǎn)數(shù)據(jù)編號(hào)類型負(fù)荷Pd, Qd branch_data importdata(33branch_data.txt); % 支路數(shù)據(jù)首端節(jié)點(diǎn)末端節(jié)點(diǎn)電阻R電抗X最大電流Imax % 2. 開關(guān)信息 % 假設(shè)有5個(gè)常閉分段開關(guān)在支路上和1個(gè)常開聯(lián)絡(luò)開關(guān) switch_locations [2, 3, 4, 5, 6; ... % 分段開關(guān)所在支路編號(hào) 33]; % 聯(lián)絡(luò)開關(guān)所在支路編號(hào)虛擬支路連接節(jié)點(diǎn)33和1 initial_status [1;1;1;1;1;0]; % 初始狀態(tài)1閉合0斷開 % 3. 分布式電源數(shù)據(jù) DG_bus [6, 18, 33]; % DG接入的節(jié)點(diǎn) DG_forecast importdata(DG_forecast_24h.csv); % 24小時(shí)預(yù)測出力3列 DG_uncertainty 0.3; % 最大波動(dòng)為預(yù)測值的±30% DG_hat DG_forecast * DG_uncertainty; % 計(jì)算最大偏差 % 4. 時(shí)間與成本參數(shù) T 24; % 時(shí)段數(shù) rho 0.5; % 電價(jià)元/kWh switch_cost 10; % 單次開關(guān)操作成本元/次 max_switch_ops 3; % 單個(gè)開關(guān)在調(diào)度周期內(nèi)最大操作次數(shù) % 5. 魯棒參數(shù) Gamma 2; % 不確定預(yù)算可根據(jù)保守程度調(diào)整 convergence_tol 1e-3; % 收斂容忍度4.2 主問題MP構(gòu)建模塊主問題是一個(gè)大規(guī)?；旌险麛?shù)線性規(guī)劃。我們需要?jiǎng)討B(tài)地添加場景。function [master_model, x_opt, LB] build_master_problem(scenarios, network_data) % scenarios: 結(jié)構(gòu)體數(shù)組每個(gè)元素包含一個(gè)最壞場景的DG出力數(shù)據(jù) % network_data: 包含所有網(wǎng)絡(luò)參數(shù)的結(jié)構(gòu)體 T network_data.T; num_switches length(network_data.initial_status); % 定義第一階段變量開關(guān)狀態(tài)二進(jìn)制 x binvar(num_switches, T, full); % x(i,t) 開關(guān)i在時(shí)段t的狀態(tài) % 初始化目標(biāo)函數(shù)和約束 obj 0; constraints []; % 為每個(gè)已知場景創(chuàng)建第二階段的變量和約束 for s 1:length(scenarios) scenario scenarios(s); % 為該場景定義第二階段的連續(xù)變量如電壓V功率流Pij, Qij V{s} sdpvar(network_data.num_buses, T, full); Pij{s} sdpvar(network_data.num_branches, T, full); Qij{s} sdpvar(network_data.num_branches, T, full); % 為該場景添加所有約束輻射狀約束、線性化潮流約束、電壓電流安全約束 % 注意潮流約束中的DG出力使用該場景的固定值 scenario.DG_output constraints [constraints, ... get_radial_constraints(x, network_data), ... get_lin_distflow_constraints(Pij{s}, Qij{s}, V{s}, x, network_data, scenario.DG_output), ... get_safety_constraints(V{s}, Pij{s}, Qij{s}, network_data)]; end % 目標(biāo)函數(shù)最小化所有場景中最大的網(wǎng)損成本 開關(guān)操作成本 % 網(wǎng)損近似為支路電阻上的損耗之和: sum(R * I^2) ≈ sum(R * (Pij^2 Qij^2) / V^2) % 線性化后可近似為 sum(2 * R * Pij) 等這里用簡化形式 max_scenario_cost max(network_data.rho * sum(sum(network_data.R * (Pij{s}.^2 Qij{s}.^2) ./ V{s}(branch_from,:))) ...); % 實(shí)際操作中網(wǎng)損線性化需要更精細(xì)的處理例如使用分段線性化或錐松弛 switch_ops_cost network_data.switch_cost * sum(sum(abs(diff(x, 1, 2)), 2)); % 操作次數(shù)成本 obj max_scenario_cost switch_ops_cost; % 構(gòu)建主問題模型 master_model optimizer(constraints, obj, sdpsettings(solver, gurobi), {scenarios}, {x, obj}); % 注意這里用optimizer封裝便于迭代中更新scenarios后快速重構(gòu)問題。 % 首次求解或獲取當(dāng)前解 if ~isempty(scenarios) diagnostics optimize(constraints, obj); x_opt value(x); LB value(obj); else x_opt []; LB -inf; end end4.3 子問題SP構(gòu)建與求解模塊子問題是算法的核心也是性能瓶頸。我們需要固定主問題給出的x然后尋找使網(wǎng)損最大化的DG出力場景。function [worst_scenario, UB_candidate, feasible] solve_subproblem(x_fixed, network_data) % x_fixed: 主問題給出的當(dāng)前開關(guān)狀態(tài)方案 % 返回最壞場景數(shù)據(jù)該場景下的成本以及該場景下潮流是否可行 T network_data.T; num_DG length(network_data.DG_bus); % 定義不確定變量 ζ ∈ [-1, 1] zeta sdpvar(num_DG, T, full); % 實(shí)際DG出力 預(yù)測值 ζ * 最大偏差 P_DG_actual network_data.DG_forecast zeta .* network_data.DG_hat; % 定義第二階段變量在給定x_fixed和P_DG_actual后 V sdpvar(network_data.num_buses, T, full); Pij sdpvar(network_data.num_branches, T, full); Qij sdpvar(network_data.num_branches, T, full); % 構(gòu)建約束包括潮流約束、安全約束以及不確定集約束 constraints []; % 1. 不確定集約束盒式約束 預(yù)算約束 constraints [constraints, -1 zeta 1]; constraints [constraints, sum(sum(abs(zeta))) network_data.Gamma]; % 2. 配電網(wǎng)運(yùn)行約束此時(shí)x_fixed是已知參數(shù) constraints [constraints, get_lin_distflow_constraints(Pij, Qij, V, x_fixed, network_data, P_DG_actual)]; constraints [constraints, get_safety_constraints(V, Pij, Qij, network_data)]; % 目標(biāo)函數(shù)最大化網(wǎng)損即運(yùn)行成本 % 同樣需要線性化處理 loss network_data.rho * sum(sum( ... )); % 計(jì)算總網(wǎng)損成本 obj -loss; % 由于YALMIP默認(rèn)最小化我們最大化loss等價(jià)于最小化-loss % 求解子問題 ops sdpsettings(solver, gurobi, verbose, 0); diagnostics optimize(constraints, obj, ops); if diagnostics.problem 0 feasible true; worst_scenario.DG_output value(P_DG_actual); UB_candidate -value(obj); % 恢復(fù)為正的成本值 else feasible false; worst_scenario []; UB_candidate inf; end end4.4 主循環(huán)與收斂判斷模塊這是整個(gè)程序的驅(qū)動(dòng)引擎。% main_robust_reconfiguration.m clear; clc; network network_data(); % 加載數(shù)據(jù) % 初始化 scenarios []; % 最壞場景集合 UB inf; % 上界 LB -inf; % 下界 iter 0; max_iter 20; x_optimal []; gap_history []; fprintf(開始魯棒動(dòng)態(tài)重構(gòu)求解...\n); while iter max_iter iter iter 1; fprintf(\n--- 迭代 %d ---\n, iter); % 步驟1求解主問題 [~, x_current, LB] build_master_problem(scenarios, network); fprintf(下界 LB %.4f\n, LB); % 步驟2求解子問題固定當(dāng)前開關(guān)方案 [new_scenario, UB_current, feasible] solve_subproblem(x_current, network); if ~feasible fprintf(子問題不可行當(dāng)前開關(guān)方案無法應(yīng)對(duì)某些不確定性。\n); % 處理不可行情況可能需要調(diào)整Gamma或約束 break; end fprintf(子問題找到場景候選上界 UB_candidate %.4f\n, UB_current); % 更新上界取歷史所有子問題結(jié)果的最小值因?yàn)樯辖缡亲詈谩白顗那闆r”成本 UB min(UB, UB_current); % 步驟3收斂判斷 if abs(UB - LB) convergence_tol * abs(UB) fprintf(收斂最優(yōu)解找到。\n); x_optimal x_current; break; end % 步驟4添加新場景到主問題 scenarios [scenarios; new_scenario]; fprintf(添加新場景到主問題場景總數(shù)%d\n, length(scenarios)); gap (UB - LB) / abs(UB); gap_history [gap_history; gap]; fprintf(當(dāng)前間隙%.4f%%\n, gap*100); end if isempty(x_optimal) fprintf(未在最大迭代次數(shù)內(nèi)收斂。最終方案使用最后一次迭代結(jié)果。\n); x_optimal x_current; end % 輸出與可視化結(jié)果 plot_results(x_optimal, scenarios, network, gap_history);4.5 結(jié)果可視化與分析模塊可視化是理解結(jié)果的關(guān)鍵。function plot_results(x_opt, scenarios, network, gap_history) figure(Position, [100, 100, 1200, 800]); % 子圖1開關(guān)狀態(tài)時(shí)序圖 subplot(2,2,1); imagesc(1:network.T, 1:size(x_opt,1), x_opt); colorbar; colormap(jet); xlabel(時(shí)段 (h)); ylabel(開關(guān)編號(hào)); title(最優(yōu)開關(guān)狀態(tài)動(dòng)態(tài)變化); set(gca, YTick, 1:size(x_opt,1)); % 子圖2各時(shí)段系統(tǒng)總網(wǎng)損在最壞場景下 subplot(2,2,2); total_loss zeros(network.T, 1); % 計(jì)算最終采用方案在最壞場景下的各時(shí)段網(wǎng)損 % ... (計(jì)算代碼) plot(1:network.T, total_loss, b-o, LineWidth, 1.5); xlabel(時(shí)段 (h)); ylabel(網(wǎng)損 (kW)); title(最壞場景下系統(tǒng)總網(wǎng)損); grid on; % 子圖3關(guān)鍵節(jié)點(diǎn)電壓幅值最壞場景下最低電壓曲線 subplot(2,2,3); % ... (計(jì)算并繪制電壓最低的節(jié)點(diǎn)在各時(shí)段的電壓曲線) plot(1:network.T, V_min, r-s, LineWidth, 1.5); hold on; yline(0.95, k--, 電壓下限); yline(1.05, k--, 電壓上限); xlabel(時(shí)段 (h)); ylabel(電壓 (p.u.)); title(最壞場景下關(guān)鍵節(jié)點(diǎn)電壓); legend(最低電壓, Location, best); grid on; % 子圖4算法收斂過程間隙變化 subplot(2,2,4); plot(1:length(gap_history), gap_history*100, m-^, LineWidth, 1.5); xlabel(迭代次數(shù)); ylabel(間隙 (%)); title(CCG算法收斂過程); grid on; fprintf(優(yōu)化完成。\n); fprintf(最終方案開關(guān)操作次數(shù)%d\n, sum(sum(abs(diff(x_opt,1,2))))); fprintf(最壞情況下總成本%.2f 元\n, UB); end5. 復(fù)現(xiàn)過程中的典型問題與調(diào)試心得復(fù)現(xiàn)這類復(fù)雜算法不可能一帆風(fēng)順。下面是我踩過的一些坑和總結(jié)的排查技巧。5.1 問題1模型不可行Infeasible這是最常見的問題。主問題或子問題報(bào)告“無可行解”。可能原因1約束過緊。比如電壓允許范圍設(shè)得太窄如[0.99, 1.01]或者支路容量設(shè)得太小在分布式電源波動(dòng)下確實(shí)無法找到滿足所有約束的運(yùn)行點(diǎn)。排查先放松所有安全約束比如電壓范圍放到[0.9, 1.1]看問題是否變得可行。如果可行再逐步收緊約束找到臨界點(diǎn)。技巧在調(diào)試初期可以嘗試用確定性模型即固定DG出力為預(yù)測值先跑通確保基礎(chǔ)模型和網(wǎng)絡(luò)數(shù)據(jù)無誤?？赡茉?輻射狀約束建模錯(cuò)誤。這是配電網(wǎng)重構(gòu)的難點(diǎn)。確保你的建模方法如虛擬流法、生成樹法正確實(shí)現(xiàn)了“連通且無環(huán)”。排查手動(dòng)驗(yàn)證一個(gè)簡單的開關(guān)組合計(jì)算其對(duì)應(yīng)的節(jié)點(diǎn)-支路關(guān)聯(lián)矩陣檢查是否滿足輻射狀條件?？梢詫懸粋€(gè)小函數(shù)來自動(dòng)檢查任意開關(guān)狀態(tài)下的網(wǎng)絡(luò)拓?fù)洹？赡茉?線性化潮流模型誤差過大。LinDistFlow在重負(fù)載或高阻抗線路上誤差顯著可能導(dǎo)致基于其的約束在實(shí)際非線性潮流下無法滿足。排查將優(yōu)化得到的開關(guān)序列和DG場景代入更精確的潮流計(jì)算工具如Matpower的runpf函數(shù)進(jìn)行校驗(yàn)。如果電壓越限嚴(yán)重需要考慮使用更精確的線性化方法如Second-Order Cone Relaxation或者將非線性潮流約束以近似方式加入魯棒框架。5.2 問題2算法不收斂或收斂慢CCG迭代幾十次后上下界間隙仍然很大?？赡茉?不確定集過于保守。預(yù)算參數(shù)Γ設(shè)置過大導(dǎo)致不確定集包含了許多物理上極不可能發(fā)生的極端場景迫使算法去尋找應(yīng)對(duì)這些“奇葩”場景的方案使得上下界難以接近。調(diào)整根據(jù)歷史數(shù)據(jù)或經(jīng)驗(yàn)合理設(shè)置Γ?？梢詮妮^小的值如0.5*DG數(shù)量開始測試觀察結(jié)果保守程度再逐步增加。可能原因2子問題求解不精確。子問題是最大化問題可能存在多個(gè)局部最優(yōu)解。如果求解器沒有找到全局最優(yōu)的“最壞場景”那么返回給主問題的威脅就不夠“致命”導(dǎo)致收斂緩慢。技巧對(duì)于子問題嘗試使用全局優(yōu)化求解器如果問題規(guī)模允許或者多次從不同的初始點(diǎn)求解取目標(biāo)函數(shù)最大的那個(gè)解作為最壞場景。檢查子問題模型的對(duì)偶間隙是否很小?？赡茉?數(shù)值問題。模型中存在數(shù)值量級(jí)差異巨大的系數(shù)如電阻是0.001級(jí)別功率是1000級(jí)別可能導(dǎo)致求解器數(shù)值不穩(wěn)定。處理對(duì)數(shù)據(jù)進(jìn)行標(biāo)幺化處理將所有物理量轉(zhuǎn)換到相近的數(shù)量級(jí)。這是電力系統(tǒng)計(jì)算中的標(biāo)準(zhǔn)做法能極大提高數(shù)值穩(wěn)定性。5.3 問題3求解時(shí)間過長對(duì)于大規(guī)模配網(wǎng)上百節(jié)點(diǎn)或長調(diào)度周期問題規(guī)模激增求解主問題混合整數(shù)規(guī)劃耗時(shí)很長。優(yōu)化策略1有效不等式。向主問題添加一些能加快整數(shù)解收斂的額外約束例如基于網(wǎng)絡(luò)流或割平面的不等式。優(yōu)化策略2并行計(jì)算。主問題和子問題的求解在每次迭代中是獨(dú)立的可以考慮并行求解多個(gè)場景的子問題如果存在多個(gè)極端場景候選。優(yōu)化策略3啟發(fā)式或簡化。對(duì)于實(shí)時(shí)性要求高的場景可以考慮使用魯棒優(yōu)化的近似算法或者將動(dòng)態(tài)重構(gòu)問題分解為多個(gè)靜態(tài)重構(gòu)問題來降低復(fù)雜度。Matlab特定技巧使用YALMIP的optimizer對(duì)象和solvemp函數(shù)針對(duì)多參數(shù)規(guī)劃有時(shí)能更高效地處理這類問題結(jié)構(gòu)。確保你的Gurobi/CPLEX許可證支持多線程并在sdpsettings中設(shè)置合適的線程數(shù)。5.4 一份實(shí)用的調(diào)試檢查清單[ ]數(shù)據(jù)校驗(yàn)網(wǎng)絡(luò)拓?fù)涫欠襁B通基礎(chǔ)潮流計(jì)算結(jié)果是否合理[ ]模型驗(yàn)證固定所有開關(guān)和DG你的LinDistFlow約束是否能還原出一個(gè)合理的潮流解[ ]確定性測試設(shè)置不確定集為0即無波動(dòng)你的魯棒模型是否退化為標(biāo)準(zhǔn)的確定性動(dòng)態(tài)重構(gòu)模型結(jié)果是否合理[ ]單次迭代手動(dòng)執(zhí)行一次CCG迭代檢查主問題和子問題的輸入輸出是否正確傳遞。[ ]可視化中間結(jié)果在每次迭代后繪制當(dāng)前開關(guān)方案下的潮流分布直觀感受方案是否合理。[ ]求解器日志打開求解器詳細(xì)輸出sdpsettings(verbose, 1)觀察每次迭代的求解狀態(tài)、對(duì)偶間隙等信息。6. 擴(kuò)展思考與項(xiàng)目價(jià)值成功復(fù)現(xiàn)這個(gè)模型只是一個(gè)起點(diǎn)。在實(shí)際研究和應(yīng)用中還有更多可以深入和擴(kuò)展的方向考慮網(wǎng)絡(luò)損耗的精確建模本項(xiàng)目為了求解效率使用了高度線性化的模型。一個(gè)重要的改進(jìn)方向是采用二階錐規(guī)劃來更精確地描述支路損耗和電壓關(guān)系實(shí)現(xiàn)魯棒凸優(yōu)化能在保證求解性的同時(shí)提高精度。多類型不確定性聯(lián)合考慮除了分布式電源出力負(fù)荷本身也具有不確定性?？梢詫⒇?fù)荷預(yù)測誤差也建模到不確定集中形成“源-荷”雙側(cè)不確定性模型使魯棒性更強(qiáng)。與儲(chǔ)能系統(tǒng)協(xié)同優(yōu)化在模型中引入儲(chǔ)能系統(tǒng)其充放電行為可以作為應(yīng)對(duì)不確定性的靈活資源。這時(shí)的決策變量還包括儲(chǔ)能的充放電計(jì)劃問題會(huì)變得更復(fù)雜但實(shí)用性大大增強(qiáng)。數(shù)據(jù)驅(qū)動(dòng)的魯棒優(yōu)化傳統(tǒng)的盒式不確定集可能過于粗糙?？梢岳脷v史數(shù)據(jù)通過機(jī)器學(xué)習(xí)方法如聚類、支持向量機(jī)構(gòu)建更精細(xì)、更貼合實(shí)際概率分布的“數(shù)據(jù)驅(qū)動(dòng)不確定集”在不過度保守的前提下提升魯棒性。從離線規(guī)劃到在線滾動(dòng)將24小時(shí)動(dòng)態(tài)重構(gòu)模型改為在線滾動(dòng)優(yōu)化模式每15分鐘或1小時(shí)根據(jù)最新的超短期預(yù)測更新一次未來數(shù)小時(shí)的決策更貼合實(shí)際調(diào)度運(yùn)行需求。復(fù)現(xiàn)這個(gè)項(xiàng)目最大的收獲不僅僅是學(xué)會(huì)了一段Matlab代碼或者一個(gè)算法。更重要的是它訓(xùn)練了一種系統(tǒng)性的思維框架如何將物理世界的復(fù)雜不確定性用嚴(yán)謹(jǐn)?shù)臄?shù)學(xué)語言進(jìn)行描述和約束如何在“最優(yōu)”與“可靠”之間尋找工程上的平衡點(diǎn)如何將一個(gè)復(fù)雜的min-max-min問題拆解成可迭代求解的模塊。這種能力對(duì)于從事電力系統(tǒng)優(yōu)化、能源管理乃至任何涉及決策與不確定性的領(lǐng)域都是極為寶貴的。