搞這個方向的人應(yīng)該都有過類似的體驗看到一篇EI論文摘要里寫著“提出了一種改進粒子群算法求解建筑集成光儲系統(tǒng)規(guī)劃運行綜合優(yōu)化模型”理論上思路很清晰可真到自己動手用Python復(fù)現(xiàn)的時候光是把這個模型從論文公式變成能跑通的代碼就夠折騰好幾天。我這次復(fù)現(xiàn)的就是這類模型——基于改進粒子群算法的建筑集成光儲系統(tǒng)規(guī)劃運行綜合優(yōu)化目標(biāo)是在同一個框架里同時回答兩個問題光伏和儲能到底該配多大容量規(guī)劃層以及配好之后儲能每天怎么充放電、建筑跟電網(wǎng)怎么交互運行層。這篇文章會把完整的復(fù)現(xiàn)過程拆開講先分析模型到底在優(yōu)化什么再把目標(biāo)函數(shù)和約束條件翻譯成代碼邏輯然后說清楚標(biāo)準(zhǔn)粒子群算法PSO為什么不夠用、改進在哪里最后給出可運行的Python實現(xiàn)要點和幾個復(fù)現(xiàn)時最容易踩的坑。適合正在做光儲優(yōu)化、微電網(wǎng)調(diào)度方向論文復(fù)現(xiàn)的讀者也適合工程師想快速驗證一個新算法方案在真實場景下是否可行。1. 這個題到底在優(yōu)化什么規(guī)劃決策與運行策略的雙層耦合1.1 規(guī)劃層在決定什么建筑集成光儲系統(tǒng)通俗講就是在建筑屋頂裝光伏板、在配電房或者地下室裝儲能電池讓建筑既可以自己發(fā)電又可以存電還能跟電網(wǎng)買賣電。規(guī)劃層的決策變量比較直白光伏的裝機容量kW、儲能電池的額定容量kWh、儲能變流器的額定功率kW。有時候還包含逆變器選型、變壓器擴容與否這類0-1變量。這些變量一旦定下來后續(xù)十幾年的收益上限基本就被鎖死了。光伏裝少了白天發(fā)的電不夠用儲能也沒多少電量可存儲能容量配大了初始投資高得嚇人一年到頭可能就閑置在那里。1.2 運行層在決定什么運行層解決的是“怎么用”的問題。在一個典型日里光伏出力曲線是一個正弦波樣子的序列負荷曲線是建筑里空調(diào)、照明、辦公設(shè)備疊加出來的另一個序列儲能則需要在每個時段決定充電還是放電、充多少放多少。同時建筑還要決定跟電網(wǎng)買電還是賣電以及各自多少功率。所以運行層的決策變量是全天每個時段通常按1小時或15分鐘劃分的儲能充放電功率、購電功率、售電功率。如果按4個典型日代表全年四季那么運行變量的維度就是4典型日× 24小時× 若干個動作變量直接撐起了優(yōu)化問題里最大的變量部分。1.3 為什么“規(guī)劃運行”必須一起算這倆階段聽著像是可以分開做的先粗算一下容量需求再在固定容量下做運行調(diào)度。但這樣做的結(jié)果往往是偏保守甚至互相矛盾。舉個例子如果我一開始把儲能容量設(shè)得過大運行調(diào)度時會發(fā)現(xiàn)最優(yōu)策略是每天最多用掉一小半電量因為電價差不足以覆蓋深度充放電帶來的效率損耗。反過來如果運行策略只追求最大化自消納率那就會傾向于把儲能耗盡長期滿充滿放會加速電池衰減運行維護成本上升又反過來影響規(guī)劃段的經(jīng)濟性評估。這就是典型的強耦合問題必須在一個統(tǒng)一的優(yōu)化模型里同時找容量和運行策略的“最優(yōu)組合”。2. 目標(biāo)函數(shù)與約束條件怎么轉(zhuǎn)成代碼能算的數(shù)學(xué)模型2.1 目標(biāo)函數(shù)經(jīng)濟性指標(biāo)怎么取舍復(fù)現(xiàn)這類論文時最常見的目標(biāo)函數(shù)是年綜合費用最小化。它由四部分組成初始投資的年化值、年運行維護費用、年購電費用、減去年售電收入。初始投資要用資金回收系數(shù)CRF折算成年值公式是CRF r × (1r)^L / ((1r)^L - 1)其中r是貼現(xiàn)率L是設(shè)備壽命。光伏和儲能壽命不一樣光伏通常按25年、儲能按10年左右計所以折算的時候得分開算。這一步很多人會漏直接拿總投資除以壽命年限誤差會很大。用公式表達就是min C_total CRF_pv × C_pv × P_pv CRF_bat × C_bat × E_bat C_power × P_pcs C_om Σ_t (price_buy(t) × P_buy(t) × Δt) ? Σ_t (price_sell(t) × P_sell(t) × Δt)運行維護費一般按初始投資的百分比比如1%~2%來估算購電分時電價和售電電價則是兩個不同的價格序列。2.2 等式約束與不等式約束的矩陣化壓縮模型里約束條件不少但歸納起來就幾類我在代碼里全部用numpy矩陣和向量處理避免寫for循環(huán)串行判斷導(dǎo)致速度太慢。第一類是功率平衡約束每個時段都要滿足P_pv(t) P_dis(t) P_buy(t) P_load(t) P_ch(t) P_sell(t)四個變量里有兩個自由度代碼里最好另外兩個變量表達P_buy(t) P_load(t) P_ch(t) P_sell(t) ? P_pv(t) ? P_dis(t)這樣平衡約束就被自動消去了優(yōu)化的搜索空間能縮小一圈。第二類是儲能系統(tǒng)約束包括SOC遞推公式SOC(t1) SOC(t) η_ch × P_ch(t) × Δt / E_bat ? P_dis(t) × Δt / (η_dis × E_bat)以及SOC上下限0.1~0.9、充放電功率上下限不超過P_pcs、同一時段不能同時充放電。同一時段不能同時充放電這種約束簡化處理可以用一個很小的正數(shù)罰項嚴(yán)謹處理則引入0-1變量但那樣會變成混合整數(shù)問題粒子群就得配合修復(fù)策略來解。第三類是電網(wǎng)交互約束購電功率上限和售電功率上限通常設(shè)為同一個值比如200kW。第四類是決策變量邊界約束光伏裝機容量不超過屋頂可安裝面積折算的上限儲能容量有離散檔位的話就做就近取整處理。2.3 典型參數(shù)的取值參考復(fù)現(xiàn)時我用了下面這組參數(shù)基本符合多數(shù)EI論文的量級跑出來的結(jié)果趨勢上是合理的參數(shù)取值說明光伏單位造價4500元/kW含安裝儲能單位容量造價1800元/kWh磷酸鐵鋰PCS單位功率造價1200元/kW雙向變流器貼現(xiàn)率5%年化折算用光伏壽命25年CRF折算儲能壽命10年CRF折算分時購電電價峰1.2/平0.75/谷0.38元/kWh三費率售電電價0.65元/kWh固定采樣間隔1小時24個時段/天典型日數(shù)量4個按季度聚類這些參數(shù)不需要跟原論文完全一致因為原始數(shù)據(jù)拿不到但量級一定要對否則結(jié)果很容易出現(xiàn)“儲能配了0容量”這種退化結(jié)果——那往往不是算法問題而是分時電價差不夠覆蓋儲能投資成本。3. 標(biāo)準(zhǔn)PSO為什么不夠用改進粒子群的三個關(guān)鍵改動3.1 慣性權(quán)重與學(xué)習(xí)因子的動態(tài)調(diào)整標(biāo)準(zhǔn)粒子群的核心公式大家都知道v w × v c1 × r1 × (pbest ? x) c2 × r2 × (gbest ? x) x x v但標(biāo)準(zhǔn)PSO在光儲規(guī)劃這個問題上有兩個先天不足一是搜索后期粒子速度衰減很快容易還沒找到好位置就全體聚集在局部最優(yōu)附近二是運行變量維度動輒上百維標(biāo)準(zhǔn)PSO幾乎不具備全局搜索能力。所以論文里的“改進”通常第一步就是慣性權(quán)重非線性遞減。線性遞減雖然簡單但前期全局搜索時間不夠后期局部收斂力度也不夠。我自己試下來效果最好的是冪函數(shù)遞減w (w_max ? w_min) × ((k_max ? k) / k_max)^alpha w_minw_max取0.9w_min取0.4alpha取1.5。迭代前40%期間w始終維持在0.7以上保證粒子有足夠的速度去探索不同容量組合到后期w快速跌向0.4讓粒子群體集中打磨最優(yōu)區(qū)域。學(xué)習(xí)因子也做了互補變化。c1表示向個體歷史最優(yōu)學(xué)習(xí)的程度c2表示向全局最優(yōu)學(xué)習(xí)的程度。前期c1應(yīng)該大一點讓粒子自己做主后期c2應(yīng)該大起來讓粒子更信任群體經(jīng)驗。實測c1從2.5線性降到0.5、c2從0.5升到2.5的效果比固定雙因子好很多。3.2 約束處理方式罰函數(shù)還是修復(fù)策略改進PSO與標(biāo)準(zhǔn)PSO最大的區(qū)別往往藏在約束處理里。規(guī)劃運行優(yōu)化模型中的等式約束和不等式約束非常多如果不處理就直接用粒子群去搜大概率在前幾十次迭代里幾乎找不到一個可行解。罰函數(shù)法是最通用的方案把每個違反約束的量求平方后乘以一個大罰因子加進目標(biāo)函數(shù)里。這樣做的好處是簡單不用改動粒子群的結(jié)構(gòu)缺點是罰因子太小約束會放松罰因子太大會形成“硬墻”粒子飛不進去收斂曲線會非常難看。我自己復(fù)現(xiàn)時采用的是“動態(tài)罰函數(shù)儲能SOC修復(fù)”的混合策略。具體來說SOC越界不直接罰而是把SOC序列重新縮放到[0.1, 0.9]區(qū)間再做功率平衡修正。這個操作能大幅提高粒子群搜索效率。功率平衡通過2.2節(jié)的變量替換自動滿足。充放電功率越限壓縮到上下限邊界。對仍然無法恢復(fù)的約束再計入罰項。這個混合策略的思路是能用物理邏輯修復(fù)的就先修復(fù)修復(fù)不了的才交給算法去罰這樣粒子群不會把大量算力浪費在尋找可行域入口上。3.3 混合離散-連續(xù)變量的特殊編碼規(guī)劃變量里的儲能容量在實際采購時是有離散檔位的比如50kWh一檔但光伏容量和功率變量是連續(xù)的。如果全部按連續(xù)變量搜索最后取整會產(chǎn)生誤差如果全部按離散變量處理維度和搜索復(fù)雜度會爆炸。我采用的方式是把粒子位置向量拆成三段來編碼前2維是光伏容量和儲能容量連續(xù)值第三維是PCS功率連續(xù)值后面的維度是各個典型日每個時段的充放電功率和購售電功率連續(xù)值。每次迭代結(jié)束以后只對儲能容量做取整到最近檔位的操作其他變量保持連續(xù)。這樣既保持了粒子群的連續(xù)搜索特性又保證了最終給出來的方案是工程上能落地采購的離散容量。4. Python實現(xiàn)的分模塊拆解從數(shù)據(jù)處理到迭代求解4.1 數(shù)據(jù)準(zhǔn)備與典型日場景提取原始數(shù)據(jù)通常是8760小時的建筑負荷序列和光伏出力理論值序列。復(fù)現(xiàn)時我沒有直接拿全年8760小時去做優(yōu)化那樣不僅計算量爆炸而且相鄰時段的高度相關(guān)性會讓粒子群陷入維度災(zāi)難。先用KMeans做典型日聚類。把每天24個負荷值24個光伏出力值拼成48維的特征向量然后用sklearn.cluster.KMeans聚成4類代表春夏秋冬四種典型日。聚類完成后每類的樣本數(shù)量占比就是該典型日在全年優(yōu)化中的權(quán)重。需要注意KMeans對初始聚類中心敏感我固定了random_state42并且跑了多次驗證聚類結(jié)果的穩(wěn)定性。from sklearn.cluster import KMeans import numpy as np # X_daily shape: (365, 48)前24維是負荷后24維是光伏 kmeans KMeans(n_clusters4, random_state42, n_init20).fit(X_daily) labels kmeans.labels_ weights np.bincount(labels, minlength4) / 365.0 typical_days np.array([X_daily[kmeans.labels_ i].mean(axis0) for i in range(4)])聚類完成之后4個典型日的24小時負荷序列和24小時光伏出力序列就是后續(xù)運行優(yōu)化的基礎(chǔ)輸入。權(quán)重weights則用來把4個典型日的運行費用折算成年費用這樣目標(biāo)函數(shù)才能跟“全年”的經(jīng)濟賬對得上。4.2 粒子編碼與種群初始化粒子編碼是整個代碼最核心的數(shù)據(jù)結(jié)構(gòu)。我定義每個粒子為一個一維numpy數(shù)組長度由三部分組成3規(guī)劃變量 4典型日× 24時段× 3儲能充電功率、放電功率、購電功率。實際算下來長度是3 4×24×3 291維每個粒子291維種群規(guī)模設(shè)40個這個維度下粒子群勉強還能跑但已經(jīng)能明顯感到收斂速度變慢。初始化的時候不能純隨機不然大量粒子初始解嚴(yán)重違反功率平衡罰函數(shù)值會大到把所有粒子都拖在不可行域邊緣。我的做法是先根據(jù)當(dāng)天光伏出力和負荷曲線用啟發(fā)式規(guī)則生成一個粗糙的初始運行策略光伏不足且電價高就放電光伏富余且電價低就充電。把這個粗糙解作為其中一個粒子的初始位置其余粒子在這個解周圍加隨機擾動生成。這樣做的好處是種群一開始就有幾個粒子落在可行域附近收斂速度明顯加快。4.3 迭代主循環(huán)與多目標(biāo)處理主循環(huán)就是一個典型的PSO迭代框架但每次迭代里要依次完成計算所有粒子的目標(biāo)函數(shù)與罰項、更新個體最優(yōu)pbest、更新全局最優(yōu)gbest、按動態(tài)參數(shù)的公式更新速度和位置、對位置做邊界處理和儲能SOC修復(fù)。多目標(biāo)處理方面如果原論文是同時優(yōu)化經(jīng)濟性和碳排放我復(fù)現(xiàn)時用的策略是加權(quán)和法把碳排放量乘以一個碳價因子折算成費用加進目標(biāo)函數(shù)。這樣粒子群只需要維護單目標(biāo)的排序關(guān)系實現(xiàn)簡單而且改動目標(biāo)函數(shù)不影響算法主體邏輯。追求帕累托前沿的話需要引入非支配排序代碼量會增加不少不利于快速復(fù)現(xiàn)驗證。核心迭代代碼大致是這樣的結(jié)構(gòu)for k in range(max_iter): w (w_max - w_min) * ((max_iter - k) / max_iter) ** 1.5 w_min c1 2.5 - 2.0 * k / max_iter c2 0.5 2.0 * k / max_iter for i in range(pop_size): # 計算目標(biāo)值含罰函數(shù) cost[i] objective_function(positions[i], typical_days, weights) # 更新個體最優(yōu) if cost[i] pbest_cost[i]: pbest_cost[i] cost[i] pbest_pos[i] positions[i].copy() # 更新全局最優(yōu) if cost[i] gbest_cost: gbest_cost cost[i] gbest_pos positions[i].copy() for i in range(pop_size): r1, r2 np.random.rand(2) velocities[i] w * velocities[i] c1 * r1 * (pbest_pos[i] - positions[i]) \ c2 * r2 * (gbest_pos - positions[i]) positions[i] velocities[i] positions[i] repair_and_bounds(positions[i], typical_days)每一輪迭代結(jié)束后把gbest_cost記錄到列表里后期直接畫收斂曲線。建議每一輪都打印一次當(dāng)前最優(yōu)配置和最優(yōu)費用方便肉眼判斷是否卡在局部最優(yōu)。4.4 結(jié)果可視化與敏感性分析跑完之后要輸出的結(jié)果至少有三樣一是收斂曲線這是證明改進算法有效性的第一手證據(jù)二是最優(yōu)規(guī)劃結(jié)果表包含光伏容量、儲能容量、PCS功率和年化總費用三是典型日的儲能充放電曲線和購售電曲線用來檢查運行策略是否符合直覺比如應(yīng)該在光伏富余時充電、在傍晚電價峰值時放電。我習(xí)慣額外加一個敏感性分析腳本把光伏單位造價從4500元/kWh降到2000元/kWh看最優(yōu)配置里光伏容量會不會上升儲能容量會不會因此改變。這種分析寫論文時非常好用能直觀體現(xiàn)模型的工程指導(dǎo)意義而且實現(xiàn)成本也就一個for循環(huán)的事。5. 復(fù)現(xiàn)EI論文時最容易翻車的幾個點5.1 典型日聚類的坑最容易翻車的環(huán)節(jié)其實是數(shù)據(jù)預(yù)處理而不是算法本身。KMeans聚類之前如果負荷和光伏出力沒有做歸一化光伏出力的數(shù)值范圍0到幾百kW會完全壓過負荷曲線的形狀特征幾百到上千kW的變化聚出來的“典型日”基本只看負荷大小光伏特征全被淹沒了。解決辦法是先把每個特征維度做標(biāo)準(zhǔn)化再聚類聚類完成后再把典型日的數(shù)值還原回原始物理量綱。另外聚幾類也是個問題K取太大典型日代表性好但計算量成倍上漲K取太小儲能充放電的策略會被平均化到不倫不類。我試下來4個典型日是精度和速度平衡得比較好的選擇部分論文用3類也可以接受。5.2 收斂判據(jù)與早熟現(xiàn)象粒子群最典型的早熟現(xiàn)象就是收斂曲線在迭代前30%就變成一條水平線之后的迭代完全沒有任何下降。遇到這種情況第一反應(yīng)不應(yīng)該是加大迭代次數(shù)而是檢查是不是粒子多樣性丟失太快。我在調(diào)試時試過幾個方案在速度更新里加入隨機擾動項、每隔一部分迭代把部分粒子重新隨機初始化、引入模擬退火式的Metropolis準(zhǔn)則接受劣解來跳出局部最優(yōu)。效果最穩(wěn)定的是“精英重啟”策略——每50輪迭代保留gbest粒子隨機初始化其他30%的粒子同時把這些新粒子的速度重置。這個策略不復(fù)雜但對打破早熟很有效。另外一定要固定隨機種子做多次重復(fù)實驗。EI論文的復(fù)現(xiàn)結(jié)果如果每次跑都不一樣很難讓審稿人相信這個算法是穩(wěn)定的。我一般固定np.random.seed(0)跑5次獨立重復(fù)看最優(yōu)解的方差如果方差太大說明算法穩(wěn)定性有問題。5.3 結(jié)果對比要對比什么復(fù)現(xiàn)論文時最讓人頭疼的是“怎么我的結(jié)果跟論文對不上”。坦白講除非你能拿到作者的原始數(shù)據(jù)和超參數(shù)否則數(shù)值完全一致幾乎不可能審稿人預(yù)期也不是這個。有效的對比方式是跑三個版本的算法標(biāo)準(zhǔn)PSO、改進PSO、再加一個常規(guī)的遺傳算法GA做參照畫收斂曲線對比收斂速度和最終目標(biāo)值。如果改進PSO的最終年費用比標(biāo)準(zhǔn)PSO低2%~5%同時收斂速度快50%以上這已經(jīng)是一個很好看的結(jié)果了。我在自己的復(fù)現(xiàn)里改進PSO比標(biāo)準(zhǔn)PSO的年費用降低約3.8%收斂到穩(wěn)定值的迭代次數(shù)從420輪降到260輪左右這說明動態(tài)權(quán)重和約束修復(fù)策略確實起了作用。最后分享一個復(fù)現(xiàn)完這個模型后的感受改進粒子群算法在這類問題里真正的價值不是追求全局最優(yōu)的數(shù)學(xué)保證而是在模型足夠復(fù)雜、商用求解器不一定吃得住的時候給你一個能在可接受時間內(nèi)找到好解的務(wù)實手段。調(diào)試過程中如果發(fā)現(xiàn)某一組參數(shù)下改進PSO還不如標(biāo)準(zhǔn)PSO別急著懷疑算法先回看約束處理是不是逼得太狠了。我自己的習(xí)慣是始終保留一個“基礎(chǔ)參數(shù)對照”腳本每改一個改進機制就跑一遍同一組隨機種子這樣到底哪個改動起了作用一目了然。