1. 為什么微電網(wǎng)調(diào)度必須用多目標(biāo)優(yōu)化從單目標(biāo)到NSGA-III的進(jìn)化邏輯做微電網(wǎng)的人應(yīng)該都有同感早期做經(jīng)濟(jì)調(diào)度大家普遍跑的是單目標(biāo)優(yōu)化目標(biāo)函數(shù)寫“運(yùn)行成本最小”約束加上功率平衡、機(jī)組出力上下限、儲(chǔ)能SOC范圍粒子群或者遺傳算法一跑出個(gè)調(diào)度計(jì)劃表感覺(jué)就算交差了。但實(shí)際并網(wǎng)運(yùn)行之后你會(huì)發(fā)現(xiàn)單目標(biāo)調(diào)度根本不夠用。微電網(wǎng)里的核心矛盾是經(jīng)濟(jì)性和環(huán)保性往往此消彼長(zhǎng)。你想讓燃?xì)廨啓C(jī)少發(fā)電壓低成本光伏和風(fēng)電占比高了出力波動(dòng)帶來(lái)的功率懲罰就大你想最大限度用可再生能源降低碳排放又要頻繁調(diào)節(jié)儲(chǔ)能和柴油機(jī)設(shè)備磨損和維護(hù)成本跟著漲。更麻煩的是電網(wǎng)側(cè)對(duì)聯(lián)絡(luò)線功率的平穩(wěn)性還有硬指標(biāo)你優(yōu)化出來(lái)的調(diào)度方案可能經(jīng)濟(jì)性最好但聯(lián)絡(luò)線功率波動(dòng)太大并網(wǎng)考核根本過(guò)不了。所以近幾年微電網(wǎng)調(diào)度研究基本都轉(zhuǎn)到了多目標(biāo)優(yōu)化框架。學(xué)術(shù)圈和工程界比較常用的目標(biāo)組合一般是三個(gè)運(yùn)行成本最小化包含燃料成本、購(gòu)售電費(fèi)用、設(shè)備啟停成本、維護(hù)成本。碳排放量最小化火電、燃?xì)廨啓C(jī)、柴油機(jī)的單位發(fā)電碳排放系數(shù)不同需要按出力加權(quán)累計(jì)。聯(lián)絡(luò)線功率波動(dòng)最小化或者配網(wǎng)交互功率平抑這條目標(biāo)在并網(wǎng)型微電網(wǎng)里越來(lái)越重要直接關(guān)系到電網(wǎng)考核和電能質(zhì)量。這三個(gè)目標(biāo)放到同一個(gè)模型里問(wèn)題就來(lái)了不存在一個(gè)解能讓三個(gè)目標(biāo)同時(shí)最優(yōu)。比如你讓成本降到最低聯(lián)絡(luò)線功率波動(dòng)指標(biāo)必然惡化你讓碳排放最低就要盡量用光伏風(fēng)電儲(chǔ)能調(diào)度壓力變大成本也可能上升。這時(shí)候單目標(biāo)算法里的“最優(yōu)解”概念失效我們需要的是一組在目標(biāo)空間里互相不支配的折中解也就是Pareto前沿——調(diào)度人員最終根據(jù)實(shí)際情況從這一組解里挑一個(gè)最順眼的方案執(zhí)行。而NSGA-III帶參考點(diǎn)的非支配排序遺傳算法就是用來(lái)找這組Pareto解的。它在NSGA-II的基礎(chǔ)上改進(jìn)了選擇機(jī)制用一組均勻分布的參考點(diǎn)來(lái)引導(dǎo)種群向整個(gè)Pareto前沿逼近特別適合目標(biāo)數(shù)在3個(gè)以上的優(yōu)化場(chǎng)景——微電網(wǎng)調(diào)度恰好就是3目標(biāo)甚至4目標(biāo)的問(wèn)題。我在Matlab里復(fù)現(xiàn)了這個(gè)完整流程從日前負(fù)荷預(yù)測(cè)數(shù)據(jù)輸入到機(jī)組出力曲線輸出整個(gè)算例跑下來(lái)大約兩個(gè)小時(shí)。這篇就來(lái)拆解一下模型怎么建、約束怎么處理、NSGA-III代碼結(jié)構(gòu)怎么搭、以及算出來(lái)之后怎么從Pareto解集里挑工程可用的調(diào)度方案。2. 微電網(wǎng)調(diào)度模型搭建目標(biāo)函數(shù)、約束條件與決策變量設(shè)計(jì)NSGA-III只是一個(gè)求解器算法本身不關(guān)心你的優(yōu)化問(wèn)題是微電網(wǎng)還是生產(chǎn)線排程。真正決定調(diào)度方案質(zhì)量的是你往算法里塞的模型長(zhǎng)什么樣。這一節(jié)先把調(diào)度模型建模的幾個(gè)關(guān)鍵環(huán)節(jié)理清楚。2.1 典型微電網(wǎng)拓?fù)渑c設(shè)備組成我這里建模的微電網(wǎng)是典型交流并網(wǎng)型拓?fù)浒夥l(fā)電單元PV出力由光照強(qiáng)度決定按日前預(yù)測(cè)值給定可調(diào)空間很小通常當(dāng)作負(fù)的負(fù)荷處理或者允許棄光。風(fēng)力發(fā)電單元WT同理出力按預(yù)測(cè)曲線給定可少量棄風(fēng)。微型燃?xì)廨啓C(jī)MT出力可調(diào)響應(yīng)速度快是微電網(wǎng)里最主要的可控單元。柴油發(fā)電機(jī)DE作為后備可控電源成本比燃?xì)廨啓C(jī)高碳排放系數(shù)也更高。儲(chǔ)能系統(tǒng)BESS雙向功率可控充電放電都行SOC有上下限約束。聯(lián)絡(luò)線Grid與上級(jí)配電網(wǎng)交換功率可以買電也可以賣電一般有功率上限約束。調(diào)度周期按24小時(shí)劃分時(shí)間步長(zhǎng)取1小時(shí)——這是目前微電網(wǎng)日前調(diào)度的主流配置步長(zhǎng)取得太小計(jì)算量太大取得太大又沒(méi)法反映光伏中午的高峰和晚間的負(fù)荷高峰。2.2 決策變量編碼方式建模時(shí)最容易踩坑的就是決策變量怎么設(shè)計(jì)。我這里把決策變量分成兩組連續(xù)變量微型燃?xì)廨啓C(jī)在24個(gè)時(shí)段的出力 P_MT(t)儲(chǔ)能系統(tǒng)在24個(gè)時(shí)段的充放電功率 P_BESS(t)。整數(shù)變量?jī)?chǔ)能充放電狀態(tài)標(biāo)志 SOC_state(t)柴油機(jī)啟停狀態(tài) u_DE(t)。這里有個(gè)經(jīng)驗(yàn)問(wèn)題要提醒一下主體優(yōu)化用NSGA-III跑連續(xù)變量整數(shù)變量在實(shí)際編碼時(shí)建議用約束罰函數(shù)處理而不是把所有變量混在一起都交給進(jìn)化算法搞交叉變異。因?yàn)榛旌险麛?shù)變量的交叉操作在進(jìn)化算法里非常容易產(chǎn)生大量不可行解你后面的約束修正邏輯會(huì)變得非常復(fù)雜。2.3 目標(biāo)函數(shù)公式與量綱歸一化處理三個(gè)目標(biāo)函數(shù)定義如下目標(biāo)1運(yùn)行成本最小化F1 Σ [C_fuel_MT(t) C_fuel_DE(t) C_OM(t) C_grid(t)]C_fuel_MT(t) 是燃?xì)廨啓C(jī)燃料成本用二次函數(shù)擬合a*P2 b*P c系數(shù)按實(shí)際機(jī)型取。C_fuel_DE(t) 是柴油機(jī)燃料成本同樣二次函數(shù)但系數(shù)比燃?xì)廨啓C(jī)大。C_OM(t) 是各設(shè)備的運(yùn)行維護(hù)成本按出力線性比例計(jì)算。C_grid(t) 是聯(lián)絡(luò)線購(gòu)售電費(fèi)用買電按分時(shí)電價(jià)計(jì)賣電按上網(wǎng)電價(jià)計(jì)這兩個(gè)電價(jià)在峰谷時(shí)段不一樣。目標(biāo)2碳排放量最小化F2 Σ [EF_MT * P_MT(t) EF_DE * P_DE(t)]光伏和風(fēng)電的碳排放系數(shù)近似為0儲(chǔ)能是能量搬運(yùn)不產(chǎn)生直接碳排放所以這里只算燃?xì)廨啓C(jī)和柴油機(jī)。目標(biāo)3聯(lián)絡(luò)線功率波動(dòng)最小化F3 Σ |P_grid(t1) - P_grid(t)|這個(gè)函數(shù)直接懲罰相鄰時(shí)段聯(lián)絡(luò)線功率的變化量值越小說(shuō)明微電網(wǎng)對(duì)外交互越平穩(wěn)。但把三個(gè)目標(biāo)直接丟給NSGA-III之前必須先做量綱歸一化。原因是F1的量級(jí)通常在幾千到幾萬(wàn)元F2是幾百千克F3可能是幾百千瓦如果不歸一化進(jìn)化算法在計(jì)算參考點(diǎn)距離的時(shí)候量綱大的目標(biāo)會(huì)主導(dǎo)選擇壓力結(jié)果就是你得到的“Pareto解集”里F1方向的解非常密另外兩個(gè)目標(biāo)方向的解稀疏得沒(méi)法用。我的做法是先用單目標(biāo)優(yōu)化分別求出每個(gè)目標(biāo)的最優(yōu)值然后用這個(gè)最優(yōu)值作為歸一化基準(zhǔn)把三個(gè)目標(biāo)都?jí)旱?~1區(qū)間。具體在代碼里就是在評(píng)價(jià)函數(shù)里加一個(gè)歸一化步驟把每個(gè)目標(biāo)的原始值除以對(duì)應(yīng)基準(zhǔn)值。2.4 約束條件與罰函數(shù)設(shè)計(jì)約束條件分三類等式約束功率平衡約束。P_MT(t) P_DE(t) P_PV(t) P_WT(t) P_BESS_discharge(t) P_grid(t) P_load(t) P_BESS_charge(t)這個(gè)約束本質(zhì)上是一個(gè)等式但在進(jìn)化算法的染色體里每個(gè)時(shí)段的各個(gè)變量是獨(dú)立的等式大概率不會(huì)被嚴(yán)格滿足。我這里的處理方式是在評(píng)價(jià)函數(shù)里做功率平衡修正把不平衡量作為罰項(xiàng)加進(jìn)去。不等式約束各機(jī)組出力上下限P_min ≤ P ≤ P_max儲(chǔ)能充放電功率上下限-P_charge_max ≤ P_BESS ≤ P_discharge_max儲(chǔ)能SOC范圍SOC_min ≤ SOC(t) ≤ SOC_max一般取0.1~0.9聯(lián)絡(luò)線傳輸功率上下限-P_grid_max ≤ P_grid ≤ P_grid_max柴油機(jī)爬坡約束|P_DE(t1) - P_DE(t)| ≤ ramp_rate罰函數(shù)設(shè)計(jì)經(jīng)驗(yàn)罰函數(shù)系數(shù)不能給太大否則所有解都被罰得差不多進(jìn)化算法會(huì)喪失選擇壓力也不能給太小否則不可行解一樣可以混進(jìn)下一代。我調(diào)試后的經(jīng)驗(yàn)是罰項(xiàng)的量級(jí)控制在目標(biāo)函數(shù)值的10%~20%范圍內(nèi)這樣進(jìn)化早期能保留一部分接近可行域的“邊界解”幫助搜索進(jìn)化后期逐步提高罰項(xiàng)比重逼種群收斂到可行區(qū)域。3. NSGA-III算法核心機(jī)制拆解參考點(diǎn)、非支配排序與自適應(yīng)歸一化NSGA-III剛上手的時(shí)候容易被它的名字嚇到覺(jué)得是NSGA-II的簡(jiǎn)單加加改改。但實(shí)際上從NSGA-II到NSGA-III最關(guān)鍵的變化是選擇機(jī)制的重構(gòu)。NSGA-II用擁擠距離來(lái)維持解的多樣性這個(gè)方法在2目標(biāo)問(wèn)題上非常好用但在3個(gè)以上目標(biāo)的時(shí)候效果急劇下降——你想象一下三維空間里一堆點(diǎn)用擁擠距離判斷誰(shuí)周圍更“稀疏”容易把邊界上的點(diǎn)和角落里的點(diǎn)誤判成高密度區(qū)域?qū)е陆夥植疾痪鶆颉SGA-III改用參考點(diǎn)引導(dǎo)的生存選擇這也是它名字里“III”的核心價(jià)值。3.1 參考點(diǎn)怎么生成Das-Dennis方法NSGA-III的一大優(yōu)勢(shì)在于不用手動(dòng)設(shè)置參考點(diǎn)位置而是通過(guò)Das-Dennis方法自動(dòng)生成一組在單位超平面上均勻分布的點(diǎn)。如果目標(biāo)數(shù)是M每個(gè)維度等分?jǐn)?shù)是p參考點(diǎn)數(shù)量由組合數(shù) C(Mp-1, p) 決定。舉例來(lái)說(shuō)M3p4時(shí)參考點(diǎn)數(shù)量 C(6, 4) 15個(gè)M3p10時(shí)參考點(diǎn)數(shù)量 C(12, 10) 66個(gè)M4p4時(shí)參考點(diǎn)數(shù)量 C(7, 4) 35個(gè)參考點(diǎn)數(shù)量直接決定了種群規(guī)模的推薦值。一般建議種群大小N取大于參考點(diǎn)數(shù)量的某個(gè)整數(shù)比如參考點(diǎn)66個(gè)的時(shí)候就取N68或者N70略微多個(gè)幾個(gè)個(gè)體不會(huì)影響結(jié)果。注意一點(diǎn)分段數(shù)p并不是越大越好。p越大參考點(diǎn)越多種群分布覆蓋更細(xì)但計(jì)算量和非支配排序壓力也變大。對(duì)3目標(biāo)的微電網(wǎng)調(diào)度p取10~15比較合適再大就可能出現(xiàn)過(guò)擬合調(diào)度曲線的情況解雖然很多但真正工程可用的也就那幾條區(qū)域。3.2 非支配排序和選擇過(guò)程的完整流程N(yùn)SGA-III每一代的進(jìn)化流程可以分成五步種群合并父代種群P_t規(guī)模N通過(guò)模擬二進(jìn)制交叉和多項(xiàng)式變異生成子代種群Q_t規(guī)模N合并得到R_t規(guī)模2N。非支配排序?qū)_t做快速非支配排序得到若干前沿層F1、F2、F3...逐層填充從F1開始逐層填充新種群S_t直到S_t規(guī)模達(dá)到N。問(wèn)題在于最后一層F_l如果全部放進(jìn)去會(huì)超N只放一部分又不知道怎么選。參考點(diǎn)關(guān)聯(lián)對(duì)S_t中所有的解做目標(biāo)空間歸一化然后逐個(gè)關(guān)聯(lián)到最近的參考點(diǎn)。選擇保留對(duì)F_l中的解統(tǒng)計(jì)每個(gè)參考點(diǎn)周圍已有的關(guān)聯(lián)解數(shù)量?jī)?yōu)先保留那些周圍解少的參考點(diǎn)所關(guān)聯(lián)的個(gè)體這樣能保證最終種群覆蓋整個(gè)Pareto前沿的各個(gè)區(qū)域。非支配排序本身怎么做每個(gè)解計(jì)算兩個(gè)值——被支配數(shù)n_p有多少個(gè)解支配它和支配列表S_p它支配哪些解。把n_p0的解放進(jìn)F1然后遍歷F1中每個(gè)解的S_p列表把被支配數(shù)減1減到0的進(jìn)入F2以此類推。3.3 自適應(yīng)歸一化與ideal point的更新機(jī)制我在實(shí)際運(yùn)行時(shí)遇到一個(gè)很容易被論文一句話帶過(guò)的坑——?dú)w一化不是簡(jiǎn)單把每個(gè)目標(biāo)除以其最大最小值就完事的。NSGA-III的歸一化方法是這樣先找到當(dāng)前種群在每個(gè)目標(biāo)方向上的最小值組成“理想點(diǎn)”z_min。用z_min對(duì)所有目標(biāo)值做平移。在平移后的目標(biāo)空間里構(gòu)造一個(gè)極值點(diǎn)集合每個(gè)目標(biāo)各取一個(gè)極值點(diǎn)這些極值點(diǎn)構(gòu)成的超平面。計(jì)算截距a_i用截距做歸一化。這個(gè)做法的好處是歸一化會(huì)隨種群進(jìn)化動(dòng)態(tài)調(diào)整。早期種群分布散極值點(diǎn)位置不穩(wěn)定截距變動(dòng)大后期種群收斂到Pareto前沿附近極值點(diǎn)穩(wěn)定歸一化就準(zhǔn)了。3.4 算例參數(shù)與收斂性觀察我這里微電網(wǎng)算例的設(shè)置參數(shù)如下參數(shù)值種群規(guī)模N200目標(biāo)維度M3參考點(diǎn)等分?jǐn)?shù)p12交叉概率0.9交叉分布指數(shù)eta_c20變異概率1/決策變量數(shù)變異分布指數(shù)eta_m20最大迭代次數(shù)500我建議不要一上來(lái)就加大迭代次數(shù)先跑200代看一眼種群分布如果Pareto前沿已經(jīng)有基本輪廓就加到300代或者500代繼續(xù)跑。微電網(wǎng)調(diào)度這個(gè)問(wèn)題的計(jì)算開銷大頭在約束評(píng)價(jià)和功率平衡修正上一次種群評(píng)估200個(gè)解每個(gè)解24小時(shí)均衡約束要做24次修正500代就是240萬(wàn)次計(jì)算Matlab里不加并行大概要跑15~25分鐘。如果機(jī)器有并行計(jì)算工具箱可以直接用parfor把種群評(píng)估并行化我實(shí)測(cè)跑400個(gè)種群的算例并行能縮到原來(lái)1/4的時(shí)間。4. Matlab代碼實(shí)現(xiàn)結(jié)構(gòu)從主函數(shù)到約束修正的完整拆解NSGA-III在Matlab里的實(shí)現(xiàn)網(wǎng)上有不少版本但很多代碼結(jié)構(gòu)混亂評(píng)價(jià)函數(shù)、進(jìn)化算子、約束處理全堆在一起想改成自己的模型非常痛苦。我這里按模塊化的思路整理了代碼結(jié)構(gòu)每個(gè)模塊的職責(zé)單一這樣你想換目標(biāo)函數(shù)或者加一個(gè)設(shè)備只需要改對(duì)應(yīng)的小函數(shù)就行。4.1 代碼模塊總覽完整代碼按以下文件組織main.m主程序負(fù)責(zé)讀取負(fù)荷數(shù)據(jù)、初始化參數(shù)、循環(huán)迭代、輸出結(jié)果。nsga3_select.mNSGA-III選擇操作包含非支配排序、參考點(diǎn)關(guān)聯(lián)、環(huán)境選擇。nd_sort.m快速非支配排序函數(shù)。reference_points.m生成Das-Dennis參考點(diǎn)。normalization.m自適應(yīng)歸一化輸入當(dāng)前種群的各目標(biāo)值、理想點(diǎn)、極值點(diǎn)。associate.m每個(gè)解關(guān)聯(lián)到最近的參考點(diǎn)返回關(guān)聯(lián)參考點(diǎn)編號(hào)和最短距離。niching.m小生境選擇保持解在參考點(diǎn)周圍的多樣性。crossover.m和mutation.m模擬二進(jìn)制交叉和多項(xiàng)式變異算子。evaluate.m核心評(píng)價(jià)函數(shù)輸入一個(gè)解計(jì)算三個(gè)目標(biāo)值和約束違反量。power_balance_fix.m功率平衡約束修正。init_pop.m種群初始化。4.2 主循環(huán)流程與關(guān)鍵參數(shù)傳遞main.m的主流程如下%% 數(shù)據(jù)加載 load(load_data.mat); % 負(fù)荷曲線 load(pv_wt_data.mat); % 光伏、風(fēng)電預(yù)測(cè)出力 %% 參數(shù)設(shè)置 global data_params; data_params.num_time 24; data_params.pop_size 200; data_params.num_gen 500; data_params.num_obj 3; data_params.num_var 48; % 24小時(shí)MT出力 24小時(shí)儲(chǔ)能功率 %% 生成參考點(diǎn) ref_points reference_points(3, 12); %% 種群初始化 pop init_pop(data_params.num_var, data_params.pop_size); pop_obj zeros(pop_size, num_obj); for i 1:pop_size [pop_obj(i,:), ~, ~] evaluate(pop(i,:), data_params); end %% 進(jìn)化主循環(huán) for gen 1:num_gen % 生成子代 child_pop crossover(pop, crossover_prob, eta_c); child_pop mutation(child_pop, mutation_prob, eta_m); % 評(píng)價(jià)子代 child_obj zeros(pop_size, num_obj); for i 1:pop_size [child_obj(i,:), ~, ~] evaluate(child_pop(i,:), data_params); end % 合并父代和子代 combined_pop [pop; child_pop]; combined_obj [pop_obj; child_obj]; % NSGA-III選擇 [pop, pop_obj] nsga3_select(combined_pop, combined_obj, ref_points, ...); end %% 輸出Pareto前沿解集 pareto_solution pop(1:pareto_count, :); plot_pareto_front(pop_obj);這里要注意一個(gè)細(xì)節(jié)決策變量的染色體編碼順序不要隨便換。比如我設(shè)的前24位是MT出力后24位是儲(chǔ)能功率交叉變異都是作用在這個(gè)染色體上的。如果你中途改了編碼順序必須同步修改evaluate函數(shù)里的解碼邏輯不然后面排查起來(lái)非常折磨。4.3 evaluate函數(shù)三目標(biāo)與罰函數(shù)的整合evaluate.m是整篇代碼的靈魂每個(gè)解都要過(guò)一遍這里。核心邏輯分成四步第一步從染色體中解碼出各時(shí)段變量P_MT solution(1:24); % 燃?xì)廨啓C(jī)24小時(shí)出力 P_BESS solution(25:48); % 儲(chǔ)能24小時(shí)功率正為放電負(fù)為充電 P_DE zeros(1,24); % 柴油機(jī)出力稍后由功率平衡計(jì)算實(shí)際算例中我不把柴油機(jī)出力和聯(lián)絡(luò)線功率直接作為決策變量而是作為功率平衡修正后計(jì)算出來(lái)的量。這樣做的原因是減少?zèng)Q策空間維度讓進(jìn)化算法集中搜索可控單元而把通量計(jì)算交給確定性規(guī)則處理。第二步逐時(shí)段功率平衡修正對(duì)每個(gè)時(shí)段t先根據(jù)負(fù)載和可再生出力計(jì)算不平衡量unbalance P_load(t) - P_MT(t) - P_PV(t) - P_WT(t) - P_BESS(t);如果unbalance 0說(shuō)明供電不足先讓柴油機(jī)補(bǔ)P_DE(t) min(unbalance, P_DE_max)如果unbalance 0說(shuō)明發(fā)電過(guò)剩優(yōu)先減少柴油機(jī)出力如果還過(guò)剩額外扣罰項(xiàng)對(duì)應(yīng)棄風(fēng)棄光如果還不行多余功率賣給電網(wǎng)聯(lián)絡(luò)線功率計(jì)算。第三步目標(biāo)函數(shù)計(jì)算逐時(shí)段累加成本、碳排放、聯(lián)絡(luò)線功率變化量得到三個(gè)原始目標(biāo)值。第四步約束違反量計(jì)算與罰項(xiàng)輸出逐時(shí)段檢查各出力是否越限、SOC是否越限SOC初值設(shè)0.5按儲(chǔ)能效率逐時(shí)段遞推把所有越限量平方和作為約束違反量。最后一個(gè)輸出參數(shù)cv會(huì)被交給選擇模塊作為不可行解篩選的依據(jù)。4.4 約束修正與SOC遞推的坑點(diǎn)SOC遞推公式要特別注意儲(chǔ)能充放電效率不是簡(jiǎn)單的加減SOC(t1) SOC(t) eta_charge * P_charge(t) * delta_t - P_discharge(t) * delta_t / eta_discharge多數(shù)論文取充電效率0.9、放電效率0.9。這個(gè)對(duì)稱假設(shè)簡(jiǎn)化了問(wèn)題但真實(shí)電池是雙向效率不對(duì)稱的如果做工程仿真建議分開取值。一個(gè)常見錯(cuò)誤是效率放錯(cuò)了位置導(dǎo)致SOC在24小時(shí)周期結(jié)束后無(wú)法回到初值調(diào)度方案不可行。另外一個(gè)坑是SOC初值和終值要一致。如果前一天結(jié)束SOC和當(dāng)天開始SOC差距很大儲(chǔ)能調(diào)度就沒(méi)有可持續(xù)性。我在罰函數(shù)里加了一項(xiàng) SOC(end) - SOC(init) 的差值懲罰保證求解出的調(diào)度方案在連續(xù)多日運(yùn)行中不會(huì)“透支”儲(chǔ)能。4.5 小生境保留參考點(diǎn)關(guān)聯(lián)與nigching選擇代碼解析\texttt{nsga3_select.m}內(nèi)部的核心部分如下function next_pop nsga3_select(pop, obj, ref_points, ...) % 快速非支配排序 [fronts, ~] nd_sort(obj); % 逐層填充 next_pop []; i 1; while size(next_pop,1) size(fronts{i},1) pop_size next_pop [next_pop; pop(fronts{i}, :)]; i i 1; end last_front pop(fronts{i}, :); last_obj obj(fronts{i}, :); % 歸一化 [norm_obj, ideal_point] normalization([obj; last_obj], ...); % 關(guān)聯(lián)參考點(diǎn) [ref_assoc, dist_assoc] associate(norm_obj, ref_points); % 小生境選擇 next_pop niching(next_pop, last_front, ref_assoc, dist_assoc, pop_size); end\texttt{associate}函數(shù)的邏輯是對(duì)每個(gè)解計(jì)算它到每個(gè)參考點(diǎn)的垂直距離參考點(diǎn)方向上的長(zhǎng)度差取最小距離對(duì)應(yīng)的參考點(diǎn)作為該解的歸屬。注意這里是“垂直距離”而不是歐氏距離本質(zhì)上是在判定解的方向和哪個(gè)參考點(diǎn)方向最接近。5. 算例結(jié)果分析Pareto前沿、最優(yōu)折中解與調(diào)度計(jì)劃生成模型搭好、代碼跑通之后最關(guān)心的是跑出來(lái)的結(jié)果長(zhǎng)什么樣。這一節(jié)放一版我實(shí)際跑的算例結(jié)果展示從Pareto前沿到最終調(diào)度計(jì)劃的全流程分析思路。5.1 算例基礎(chǔ)數(shù)據(jù)說(shuō)明負(fù)荷曲線取自某園區(qū)微電網(wǎng)冬季典型日最大負(fù)荷約850kW光伏容量300kW、風(fēng)電容量150kW燃?xì)廨啓C(jī)額定出力400kW柴油機(jī)額定出力200kW儲(chǔ)能容量500kWh最大充放電功率100kW聯(lián)絡(luò)線功率上限300kW。分時(shí)電價(jià)峰時(shí)10:00-15:00、18:00-21:001.1元/kWh平時(shí)0.68元/kWh谷時(shí)0.32元/kWh。上網(wǎng)電價(jià)統(tǒng)一0.45元/kWh。5.2 Pareto前沿可視化結(jié)果跑完500代后得到的Pareto前沿在三維目標(biāo)空間中呈清晰的曲面分布。我通常用三張兩兩組合的二維投影圖來(lái)分析成本-碳排放投影最能說(shuō)明經(jīng)濟(jì)性和環(huán)保性的矛盾曲線左下角是最低碳排放的方案對(duì)應(yīng)大量棄風(fēng)棄光限制燃?xì)廨啓C(jī)運(yùn)行右下角是最低成本方案柴油機(jī)幾乎不出力完全依賴聯(lián)絡(luò)線購(gòu)電和燃?xì)廨啓C(jī)碳排放顯著升高。中間區(qū)域的解則是兩者折中后的“甜點(diǎn)區(qū)”。成本-波動(dòng)性投影則直觀反映電網(wǎng)交互壓力越便宜的方案聯(lián)絡(luò)線波動(dòng)越大因?yàn)殡妰r(jià)激勵(lì)下儲(chǔ)能會(huì)大量在谷時(shí)充電、峰時(shí)放電導(dǎo)致聯(lián)絡(luò)線功率大幅漲落。這在單目標(biāo)優(yōu)化中被完全忽略但在NSGA-III的優(yōu)化下可以得到一組波動(dòng)小于20kW、成本僅比最優(yōu)值高5%左右的解工程上很有參考價(jià)值。5.3 折中解選取方法模糊隸屬度函數(shù)Pareto前沿上幾百個(gè)解最終只能選一個(gè)下發(fā)執(zhí)行。工程上常用的方法是模糊隸屬度函數(shù)對(duì)每個(gè)目標(biāo)F_k計(jì)算解的隸屬度u_k (F_k_max - F_k) / (F_k_max - F_k_min)所有目標(biāo)隸屬度的平均值就是該解的“綜合滿意度”取滿意度最大的解作為折中解。這個(gè)方法簡(jiǎn)單直接適合做論文對(duì)比實(shí)驗(yàn)。工程現(xiàn)場(chǎng)如果調(diào)度員有傾向比如更看重經(jīng)濟(jì)性就手動(dòng)調(diào)權(quán)重比如取成本權(quán)重0.6、碳排放權(quán)重0.3、波動(dòng)權(quán)重0.1再做加權(quán)滿意度排序。5.4 最優(yōu)折中解下的各機(jī)組調(diào)度曲線解讀折中解對(duì)應(yīng)的調(diào)度計(jì)劃有幾個(gè)明顯的模式凌晨低谷時(shí)段0:00-6:00負(fù)荷低光伏風(fēng)電出力低但夠用燃?xì)廨啓C(jī)基本不出力儲(chǔ)能谷時(shí)充電聯(lián)絡(luò)線從電網(wǎng)少量購(gòu)電維持平衡。這段時(shí)間成本最低。上午和下午過(guò)渡時(shí)段7:00-10:00、14:00-17:00光伏出力上升燃?xì)廨啓C(jī)緩緩升負(fù)荷儲(chǔ)能處于淺充淺放狀態(tài)。聯(lián)絡(luò)線功率保持平穩(wěn)。晚間高峰時(shí)段18:00-21:00負(fù)荷最高光伏跌到0儲(chǔ)能放電頂上去燃?xì)廨啓C(jī)滿發(fā)柴油機(jī)作為邊際機(jī)組補(bǔ)缺口。這一時(shí)段由于峰時(shí)電價(jià)高盡量自發(fā)自用減少高價(jià)購(gòu)電。這些模式符合微網(wǎng)運(yùn)行直覺(jué)但NSGA-III的價(jià)值在于它在“模式邊界”處給出了量化平衡比如儲(chǔ)能SOC到底在哪個(gè)時(shí)刻開始放電是全部放空還是留10%備用這個(gè)細(xì)節(jié)單靠人工經(jīng)驗(yàn)很難拍板。5.5 與NSGA-II的對(duì)比結(jié)果我用同樣的模型跑了一遍NSGA-II做對(duì)比結(jié)論和文獻(xiàn)一致目標(biāo)數(shù)3個(gè)時(shí)NSGA-III生成的Pareto解集分布明顯更均勻。NSGA-II的解集中在少數(shù)幾個(gè)擁擠區(qū)域參考點(diǎn)法在覆蓋性上有明顯優(yōu)勢(shì)。另外NSGA-II在F3方向波動(dòng)性上經(jīng)常丟失極端解也就是“最低波動(dòng)”和“最高波動(dòng)”兩端都有空洞。這個(gè)對(duì)論文實(shí)驗(yàn)很致命因?yàn)槟惝嫵鰜?lái)的Pareto前沿投影圖不完整審稿人一眼就能看出來(lái)算法多樣性不足。換NSGA-III這個(gè)問(wèn)題基本就消失了。6. 代碼復(fù)現(xiàn)中的常見錯(cuò)誤與調(diào)試心得這部分寫點(diǎn)真正跑代碼才會(huì)遇到的問(wèn)題網(wǎng)上很多開源版本沒(méi)有注釋我踩過(guò)的坑希望能幫你省下幾天的排查時(shí)間。6.1 參考點(diǎn)生成報(bào)錯(cuò)組合數(shù)越界Das-Dennis生成參考點(diǎn)時(shí)用到排列組合計(jì)算如果組合數(shù)過(guò)大比如M5、p10的時(shí)候組合數(shù)C(14,10)1001個(gè)內(nèi)存占用會(huì)瞬間飆升。我的建議是3目標(biāo)問(wèn)題p取10~15參考點(diǎn)數(shù)量66~105個(gè)表現(xiàn)已經(jīng)很好。如果非要跑5目標(biāo)甚至更多就改用兩層的參考點(diǎn)結(jié)構(gòu)邊界層內(nèi)部層這個(gè)技巧在Deb的論文里有詳細(xì)說(shuō)明代碼里也要做相應(yīng)改動(dòng)。6.2 歸一化時(shí)極值點(diǎn)重合導(dǎo)致的截距NaN這是一個(gè)非常隱蔽的bug當(dāng)種群中存在大量相同解進(jìn)化后期常見算法在找每個(gè)目標(biāo)方向的極值點(diǎn)時(shí)可能選到同一個(gè)解導(dǎo)致后續(xù)超平面方程無(wú)解截距計(jì)算出現(xiàn)NaN或Inf。我處理的方法是在\texttt{normalization.m}里加一個(gè)檢查如果極值點(diǎn)重合就退回到目標(biāo)值區(qū)間歸一化用min-max而不是截距歸一化。實(shí)測(cè)這個(gè)兜底方案對(duì)結(jié)果影響很小但能保證程序不崩。6.3 約束罰項(xiàng)過(guò)大導(dǎo)致Pareto前沿“塌陷”如果你發(fā)現(xiàn)跑完的Pareto前沿只剩下幾個(gè)點(diǎn)且都聚集在目標(biāo)空間小角落大概率是罰函數(shù)系數(shù)給太大了。所有不可行解都被罰得遠(yuǎn)遠(yuǎn)的進(jìn)化算法視野里只有極小一片可行區(qū)域多樣性根本起不來(lái)。調(diào)試技巧是先不加罰項(xiàng)跑200代看解集是否分布在多個(gè)方向再逐步加大罰項(xiàng)系數(shù)觀察前沿形態(tài)是否正?！罢归_”。6.4 功率平衡修正后的變量越界怎么處理有些修正在計(jì)算過(guò)程中會(huì)產(chǎn)生越界變量比如柴油機(jī)補(bǔ)不平衡量時(shí)超過(guò)上限必須在修正函數(shù)里把越界部分按比例削掉而不是硬截?cái)?。硬截?cái)鄷?huì)導(dǎo)致功率不平衡量殘留最終負(fù)荷都無(wú)法滿足。我的做法是先讓儲(chǔ)能和柴油機(jī)分別承擔(dān)一部分如果機(jī)組到達(dá)上限還不行就讓聯(lián)絡(luò)線吸收剩余量如果聯(lián)絡(luò)線也到上限就把該時(shí)段標(biāo)記為失負(fù)荷并在罰函數(shù)里加一個(gè)失負(fù)荷懲罰項(xiàng)。6.5 并行計(jì)算的坑隨機(jī)數(shù)一致性問(wèn)題用\texttt{parfor}加速種群評(píng)估時(shí)每個(gè)worker的隨機(jī)數(shù)狀態(tài)是獨(dú)立的。如果交叉變異操作里要用到隨機(jī)數(shù)建議先把父代選定和變異操作的隨機(jī)數(shù)生成放在單個(gè)進(jìn)程中完成只把\texttt{evaluate}函數(shù)放進(jìn)去并行。否則你會(huì)遇到“跑兩次結(jié)果完全一致”的假象——其實(shí)是并行模式下的隨機(jī)流沒(méi)洗干凈進(jìn)化過(guò)程退化成確定性搜索了。7. 從算法到工程落地微電網(wǎng)調(diào)度的實(shí)用擴(kuò)展方向NSGA-III跑通微電網(wǎng)調(diào)度只是第一步。實(shí)際接入工程系統(tǒng)還有幾個(gè)方向值得投入精力。7.1 從日前調(diào)度走向日內(nèi)滾動(dòng)優(yōu)化日前調(diào)度給出的計(jì)劃建立在24小時(shí)預(yù)測(cè)基礎(chǔ)上但光伏和負(fù)荷預(yù)測(cè)誤差隨著時(shí)間推移會(huì)累積。工程上通行做法是滾動(dòng)時(shí)域優(yōu)化每15分鐘做一次優(yōu)化只執(zhí)行未來(lái)1小時(shí)首段的決策到下一個(gè)周期再重新求解。NSGA-III在這種場(chǎng)景下速度夠不夠是個(gè)實(shí)際問(wèn)題目前硬件條件下500代×200種群規(guī)模的優(yōu)化耗時(shí)約5~10秒完全可以滿足15分鐘滾動(dòng)周期的要求。7.2 從靜態(tài)目標(biāo)到動(dòng)態(tài)碳排放電價(jià)目前微電網(wǎng)調(diào)度考慮碳排放是按固定系數(shù)算的。隨著碳市場(chǎng)發(fā)展碳排放成本應(yīng)該按碳價(jià)動(dòng)態(tài)折算進(jìn)目標(biāo)函數(shù)里。這種場(chǎng)景下第三個(gè)目標(biāo)可以改成“碳交易成本最小”和運(yùn)行成本合并成一個(gè)帶碳價(jià)的單目標(biāo)或繼續(xù)保留為多目標(biāo)都行模型只需要改目標(biāo)函數(shù)中的幾個(gè)系數(shù)代碼框架基本不動(dòng)。7.3 考慮多場(chǎng)景魯棒性日前預(yù)報(bào)不確定性帶來(lái)的風(fēng)險(xiǎn)可以通過(guò)多場(chǎng)景法處理生成10個(gè)典型光伏/負(fù)荷場(chǎng)景每個(gè)場(chǎng)景都計(jì)算三個(gè)目標(biāo)值然后用期望值或最壞情況值做目標(biāo)函數(shù)。決策變量仍然是確定性的但評(píng)價(jià)函數(shù)里要做場(chǎng)景循環(huán)。這樣得到的調(diào)度方案對(duì)預(yù)測(cè)誤差有更強(qiáng)的魯棒性代價(jià)是從每代評(píng)估次數(shù)從200×1變成200×10計(jì)算時(shí)間大概增加8~9倍。如果要做這個(gè)方向的研究?jī)?yōu)先優(yōu)化評(píng)價(jià)函數(shù)里的功率平衡修正環(huán)節(jié)把向量化寫到位否則計(jì)算量非常感人。我自己目前做的方向是在NSGA-III基礎(chǔ)上加入場(chǎng)景縮減模塊SCENRED先把原始預(yù)測(cè)誤差的50個(gè)場(chǎng)景用距離矩陣聚合成10個(gè)代表性場(chǎng)景再做多場(chǎng)景魯棒優(yōu)化。跑出來(lái)的調(diào)度曲線比單場(chǎng)景方案平滑很多儲(chǔ)能動(dòng)作次數(shù)明顯減少設(shè)備壽命損耗這一項(xiàng)隱性成本也能估算出來(lái)?；氐阶畛醯膯?wèn)題——為什么微電網(wǎng)調(diào)度研究繞不開多目標(biāo)優(yōu)化因?yàn)檫@個(gè)系統(tǒng)本質(zhì)上就是一個(gè)“省著花、少排碳、穩(wěn)并網(wǎng)”的多目標(biāo)博弈問(wèn)題。NSGA-III的價(jià)值在于把你的偏好選項(xiàng)完整地陳列出來(lái)讓決策者知道“省錢的代價(jià)是什么減排的收獲是什么平穩(wěn)的指標(biāo)要花多少錢去買”。Matlab代碼不是終點(diǎn)把這些代碼接上實(shí)時(shí)數(shù)據(jù)、接上預(yù)測(cè)模塊、接上調(diào)度員的操作界面那才是微電網(wǎng)能量管理系統(tǒng)真正發(fā)揮作用的地方。