峰配置方案與經(jīng)濟性分析的Matlab復(fù)現(xiàn)指南)
剛收到“參與調(diào)峰的儲能系統(tǒng)配置方案及經(jīng)濟性分析”這個題目時我原以為就是把儲能容量、功率套進一個優(yōu)化函數(shù)在Matlab里調(diào)出結(jié)果就完事。真做下來才發(fā)現(xiàn)調(diào)峰場景牽扯到目標(biāo)函數(shù)怎么建、約束條件怎么表達、充放電時機怎么建模每個環(huán)節(jié)都有一堆細節(jié)。這篇文章就把我復(fù)現(xiàn)EI論文時踩過的坑、摸過的路完整整理出來適合正在做論文復(fù)現(xiàn)、畢業(yè)設(shè)計或者儲能項目前期方案論證的同行參考。我會把建模邏輯、經(jīng)濟性核算、Matlab實現(xiàn)路徑、結(jié)果敏感性分析和常見報錯一并講清楚盡量做到看完就能動手。1. 調(diào)峰場景下儲能配置的建模邏輯從負(fù)荷曲線到?jīng)Q策變量1.1 調(diào)峰需求從哪來儲能又扮演什么角色電力系統(tǒng)的負(fù)荷一天之內(nèi)波動很大夜間低谷、白天高峰。新能源大規(guī)模接入之后情況變得更復(fù)雜光伏在午間大發(fā)把凈負(fù)荷負(fù)荷減去新能源出力壓得很低到了傍晚光伏出力歸零凈負(fù)荷又快速爬升形成所謂的“鴨子曲線”。這樣一來系統(tǒng)對調(diào)峰資源的需求不再只是簡單的峰谷差而是要求機組具備快速爬坡、深度下調(diào)的能力。傳統(tǒng)火電機組調(diào)峰深度有限頻繁啟停和深調(diào)不僅增加煤耗還帶來機組壽命損耗。儲能系統(tǒng)的優(yōu)勢正好體現(xiàn)在這里充電時相當(dāng)于增加負(fù)荷放電時相當(dāng)于減少負(fù)荷配合PCS的快速響應(yīng)能力可以在幾分鐘內(nèi)完成充放電狀態(tài)切換。在配置層面儲能參與調(diào)峰主要回答三個問題裝多大功率、裝多大容量、充放電策略怎么安排。功率決定了單位時間能“挪”多少電量容量決定了能持續(xù)“挪”多久兩者互相耦合。1.2 為什么配置方案要當(dāng)成優(yōu)化問題來解很多人一上來就套“兩小時儲能”的固定比例比如負(fù)荷峰值的10%、持續(xù)2小時這是工程估算思路不能用于論文復(fù)現(xiàn)或精細方案論證。原因是功率和容量之間存在最優(yōu)配比關(guān)系功率大了容量不足高峰時段放一會兒就沒電容量大了功率不足高峰時段派不上用場。同時儲能參與調(diào)峰還面臨“在哪里充電、在哪里放電”的時序決策這是一個典型的規(guī)劃-運行聯(lián)合優(yōu)化問題。EI論文里常見做法是把問題寫成一個混合整數(shù)規(guī)劃MILP上層決策儲能的額定功率和額定容量運行層決策每個時段的充放電功率和0-1狀態(tài)變量。兩層可以解耦迭代也可以用KKT條件合并成單層更多人直接用單層MILP表達。判別標(biāo)準(zhǔn)很簡單如果目標(biāo)函數(shù)里同時出現(xiàn)投資成本項和運行收益項且運行收益依賴逐時段的充放電變量就需要放在同一個模型里聯(lián)立求解而不是先定規(guī)模再算運行收益那樣會漏掉規(guī)模與運行策略之間的耦合關(guān)系。1.3 典型日與全時間尺度的取舍復(fù)現(xiàn)這類論文時最先要決定的事情就是時間尺度。直接用全年8760小時建模變量維度大求解時間長而且很多論文數(shù)據(jù)只給典型日。實用做法是選取典型日來代表不同季節(jié)的運行特性比如夏季典型日、冬季典型日、過渡季典型日每個典型日按24時段或96時段離散。目標(biāo)函數(shù)中把典型日天數(shù)作為權(quán)重乘進去折算全年收益。典型日的選取不是隨手挑一天平均曲線我建議從全年負(fù)荷數(shù)據(jù)里按峰谷差、總電量、早晚高峰形態(tài)做聚類或者直接選取峰谷差最大的保守日。保守日算出來的儲能規(guī)模偏大適合可研階段聚類典型日算出來的規(guī)模更貼近實際運營適合方案比選。如果原文沒寫清楚典型日怎么來的可以在復(fù)現(xiàn)說明里注明“采用聚類得到的三個典型日各代表一個季節(jié)”這樣審稿人或?qū)煵粫|(zhì)疑數(shù)據(jù)來源。2. 經(jīng)濟性目標(biāo)函數(shù)與約束條件凈年值法下的成本-收益核算2.1 成本項投資成本怎么折算成每年的費用儲能系統(tǒng)配置的經(jīng)濟性核心是凈年值NAV或等年值而不是靜態(tài)投資回收。不同方案壽命不同、規(guī)模不同只有折算成年值才能放到同一個尺度上比較??偼顿Y成本通常分成兩部分C_inv c_p × P_ess c_e × E_essc_p是單位功率成本包含PCS、變壓器、場地等常用單位是元/kWc_e是單位容量成本包含電池本體、BMS、溫控等常用單位是元/kWh。投資發(fā)生在建設(shè)期運行收益發(fā)生在運營期所以要把投資換算成等年值換算系數(shù)是資金回收系數(shù)CRF r(1r)^n / ((1r)^n - 1)r是折現(xiàn)率n是儲能系統(tǒng)壽命。年化投資成本 C_inv × CRF。運行維護成本一般按投資額的一定比例估算常見取值是每年2%到3%也可以用單位電量運維成本乘年充放電量。最后成本項表達式為C_total C_inv × CRF C_om我復(fù)現(xiàn)時習(xí)慣把折現(xiàn)率設(shè)為8%、壽命15年這兩個參數(shù)對結(jié)果影響很大后面敏感性分析會專門講。2.2 收益項峰谷套利和調(diào)峰價值的量化口徑收益項是經(jīng)濟性分析里最容易出問題的地方。參與調(diào)峰場景下儲能收益主要有兩筆口徑必須分清楚。第一筆是峰谷套利收益。儲能低價時段充電、高價時段放電每時段收益等于放電量乘放電電價減去充電量乘充電電價。表達式為R_arb Σ [ p_dis(t) × P_dis(t) - p_ch(t) × P_ch(t) ] × Δt第二筆是調(diào)峰價值。儲能削峰填谷后系統(tǒng)原本需要調(diào)用深度調(diào)峰火電機組現(xiàn)在這部分電量由儲能承擔(dān)替代的深調(diào)成本就是儲能的調(diào)峰價值。如果按調(diào)峰補償機制算可以寫成R_reg c_reg × E_total_disc_reg是單位調(diào)峰電量補償價格。這里要特別提醒峰谷套利和調(diào)峰補償不能對同一度電重復(fù)計算。如果儲能在高峰時段放電已經(jīng)賺了峰谷價差又把同一段放電量拿一份調(diào)峰補償收益就虛高了。我的處理原則是把儲能放電量分成兩類——價差套利部分按分時電價結(jié)算額外深度調(diào)峰支持部分按調(diào)峰補償結(jié)算兩者通過時段劃分和系統(tǒng)調(diào)峰缺口來界定。2.3 約束條件從SOC遞推到充放電互斥約束條件直接決定模型是否物理可行。復(fù)現(xiàn)時至少要有以下五組約束。儲能SOC荷電狀態(tài)遞推約束是核心SOC(t1) SOC(t) η_ch × P_ch(t) × Δt / E_ess - P_dis(t) × Δt / (η_dis × E_ess)η_ch是充電效率η_dis是放電效率。注意放電時效率要放在分母意思是電池釋放出的電能要大于送入電網(wǎng)的電能這個位置放反了會直接把儲能效率算到120%經(jīng)濟性結(jié)果完全失真。功率上下限約束0 ≤ P_ch(t) ≤ P_ess 0 ≤ P_dis(t) ≤ P_essSOC上下限約束SOC_min ≤ SOC(t) ≤ SOC_max充放電互斥約束用0-1變量u(t)表示P_ch(t) ≤ M × u(t) P_dis(t) ≤ M × (1-u(t))M是足夠大的常數(shù)但不能太大否則數(shù)值計算會病態(tài)。最后還要加一個周期約束讓典型日結(jié)束時的SOC回到初始值否則儲能會“免費”消耗初始電量收益被高估SOC(T1) SOC_0如果是調(diào)峰場景還可能需要加入系統(tǒng)凈負(fù)荷峰谷差約束或火電最小出力約束具體形式取決于原文采用的機制模型。3. Matlab實現(xiàn)路徑從參數(shù)初始化到優(yōu)化求解的完整流程3.1 數(shù)據(jù)準(zhǔn)備負(fù)荷曲線歸一化與時段劃分我復(fù)現(xiàn)時第一步不是寫代碼而是先把負(fù)荷數(shù)據(jù)整理成標(biāo)準(zhǔn)格式。假設(shè)一個典型日24時段負(fù)荷單位為MW需要轉(zhuǎn)換成與儲能功率統(tǒng)一量綱。如果原文給的是標(biāo)幺值要乘以基準(zhǔn)負(fù)荷。分時電價按時段填入向量峰平谷時段要跟負(fù)荷曲線對應(yīng)上。一個容易被忽略的問題是時段步長。如果按24點離散Δt 1小時能量單位就是MWh如果按96點離散Δt 0.25小時公式里的Δt就不能省略。很多報錯都源于步長沒乘導(dǎo)致SOC遞推錯得離譜。我習(xí)慣把數(shù)據(jù)組織成結(jié)構(gòu)體data.T 24; % 時段數(shù) data.dt 1; % 步長小時 data.P_load [..]; % 典型日負(fù)荷MW data.price_ch [..]; % 充電電價元/MWh data.price_dis [..]; % 放電電價元/MWh data.eta_ch 0.95; data.eta_dis 0.95; data.SOC_min 0.1; data.SOC_max 0.9; data.SOC_0 0.5;3.2 變量定義與約束建模YALMIP還是intlinprogMatlab里處理MILP有兩條路一是直接用內(nèi)置的intlinprog二是裝YALMIP工具箱讓YALMIP把模型翻譯成求解器能識別的形式。我強烈建議用YALMIP因為論文復(fù)現(xiàn)要頻繁修改約束條件YALMIP維護變量的維度關(guān)系更方便不容易出現(xiàn)索引錯位。變量聲明如下P_ch sdpvar(data.T, 1); P_dis sdpvar(data.T, 1); u binvar(data.T, 1); E_ess sdpvar(1, 1); P_ess sdpvar(1, 1);SOC變量可以顯式聲明也可以直接用遞推公式約束表達SOC sdpvar(data.T1, 1); SOC(1) data.SOC_0;約束寫成數(shù)組形式Y(jié)ALMIP會自動處理每個時段的約束Constraints []; for t 1:data.T Constraints [Constraints, SOC(t1) SOC(t) ... data.eta_ch*P_ch(t)*data.dt/E_ess - ... P_dis(t)*data.dt/(data.eta_dis*E_ess)]; Constraints [Constraints, P_ch(t) data.P_max_scale * u(t)]; Constraints [Constraints, P_dis(t) data.P_max_scale * (1-u(t))]; end注意互斥約束里M的取值我用的是P_ess的上界比如取負(fù)荷峰值的2倍這個上界在所有可行解范圍內(nèi)都成立又不會大到引起數(shù)值問題。目標(biāo)函數(shù)按凈年值定義C_inv c_p * P_ess c_e * E_ess; C_annual C_inv * CRF c_om * C_inv; R_arb sum(price_dis .* P_dis - price_ch .* P_ch) * year_days * data.dt; R_reg c_reg * sum(P_dis) * year_days * data.dt; Objective C_annual - R_arb - R_reg; optimize(Constraints, Objective);這里的year_days是典型日在全年出現(xiàn)的天數(shù)三個典型日各算各的再累加。3.3 求解器選型與參數(shù)設(shè)置如果只裝Matlab可以用內(nèi)置的intlinprog。YALMIP底層默認(rèn)調(diào)用的求解器如果不是專門針對MILP優(yōu)化過的大模型會跑得很慢。我復(fù)現(xiàn)時用Gurobi或者Cplex速度和穩(wěn)定性都明顯好于內(nèi)置求解器。如果沒有外部求解器許可證先用intlinprog也能跑通小算例但迭代次數(shù)要調(diào)大。幾個關(guān)鍵求解參數(shù)供參考o(jì)ptions sdpsettings(solver, gurobi, ... gurobi.MIPGap, 0.01, ... gurobi.TimeLimit, 300, ... verbose, 2);MIPGap設(shè)到1%就足夠論文精度再往下壓會顯著增加求解時間。如果模型規(guī)模大建議先把某個典型日的時段數(shù)從96降到24做調(diào)試跑通再加密。3.4 結(jié)果驗證與輸出一套模型跑完不能只看目標(biāo)函數(shù)值就收工還要做三層驗證。第一層檢查可行性SOC序列是否在上下限內(nèi)充放電功率是否同時為正。第二層檢查能量平衡充電電量乘以效率減去放電電量應(yīng)等于24小時前后SOC差折合的能量。第三層檢查經(jīng)濟性指標(biāo)單位容量年收益、靜態(tài)回收期、內(nèi)部收益率這些指標(biāo)應(yīng)處于行業(yè)合理區(qū)間。我最常用的驗證手段是畫圖觀察SOC曲線。如果SOC曲線出現(xiàn)鋸齒狀頻繁震蕩大概率是電價時段劃分和儲能響應(yīng)邏輯不匹配如果SOC長期頂在上限說明容量配置過大儲能大部分時間用不上。這些直覺判斷比單純看報表有效得多。4. 典型算例演算配置結(jié)果隨邊界條件的變化規(guī)律4.1 算例參數(shù)與基礎(chǔ)電價環(huán)境為了把方法落地我構(gòu)造一個典型算例參數(shù)取當(dāng)前工程常見量級。負(fù)荷峰值100MW峰谷差率40%分時電價峰、平、谷分別為1.0、0.6、0.3元/kWh。儲能單位功率成本1500元/kW單位容量成本800元/kWh運維費率2%效率95%SOC范圍0.1到0.9壽命15年折現(xiàn)率8%。這里要說明參數(shù)不是從某一篇特定論文抄來的而是綜合多個文獻和行業(yè)可研報告取的中間值。不同論文原文參數(shù)差別很大復(fù)現(xiàn)時第一步就是把原文參數(shù)表完整落到代碼里再按我這里的流程跑。4.2 最優(yōu)配置結(jié)果與運行特性按上述參數(shù)求解典型的優(yōu)化結(jié)果會是20MW/80MWh左右對應(yīng)4小時持續(xù)放電。年放電量約2600萬kWh峰谷套利收益和調(diào)峰補償合計年收益約900萬元年化成本約500萬元凈年收益約400萬元。這個量級下靜態(tài)回收期約7到8年符合當(dāng)前儲能項目可研常用判斷區(qū)間。運行曲線有幾個特征值得注意。儲能會在谷段滿充在峰段滿放平段基本不出力。SOC曲線在谷段從0.1爬到0.9峰段從0.9掉回0.1整個循環(huán)干干凈凈。如果優(yōu)化結(jié)果里出現(xiàn)SOC只用了0.3到0.7的情況說明電價差還沒有大到值得用滿儲能容量配置規(guī)模相對電價環(huán)境偏大了。4.3 敏感性分析單位成本、價差與調(diào)峰補償?shù)挠绊懡?jīng)濟性分析不能只給一個方案結(jié)果必須有敏感性分析。我通常掃三個變量儲能單位容量成本、峰谷價差、調(diào)峰補償價格。單位容量成本從1000元/kWh降到600元/kWh時最優(yōu)配置容量大約增加30%到50%。這是因為成本下降直接拉低邊際成本原本不經(jīng)濟的多配置部分變得可接受。峰谷價差從0.5元/kWh擴大到0.9元/kWh時儲能運行收益提升最優(yōu)功率和容量都會上升且儲能更傾向于“滿充滿放”。調(diào)峰補償價格從0.1元/kWh提高到0.3元/kWh時儲能收益曲線變得更陡配置結(jié)果也會變大。有一個所有算例都成立的經(jīng)濟規(guī)律最優(yōu)配置落點一定在邊際收益等于邊際成本的位置。儲能容量每增加1MWh帶來額外的年收益遞減因為峰谷電量空間有限、充放電次數(shù)有上限。一旦邊際收益小于年化邊際成本繼續(xù)加容量就是負(fù)優(yōu)化。論文里的“經(jīng)濟拐點”指的就是這個位置畫成曲線就是凈年值隨容量先升后降的倒U型。5. 復(fù)現(xiàn)中的常見坑位與實戰(zhàn)排查思路5.1 SOC遞推與維度對不齊這是我在復(fù)現(xiàn)時踩過的第一個坑。Matlab索引從1開始t1對應(yīng)第一個時段SOC變量如果長度是T1第一行是初始值循環(huán)里SOC(t1)和SOC(t)的對應(yīng)關(guān)系一旦寫錯整個SOC曲線就會發(fā)生平移或跳變。最直接的排查方法打印t1、t2時段的SOC與前一步計算人工手算驗證。還有一個隱蔽問題SOC遞推式子里E_ess做分母如果E_ess是決策變量且初始值為0求解器在計算約束時會出現(xiàn)除以0的問題。解決方法是給E_ess加一個極小下界比如0.001或者在約束里做等價變形SOC(t1) × E_ess SOC(t) × E_ess ...這在實際求解時更穩(wěn)。5.2 充放電互斥條件遺漏導(dǎo)致“既充又放”這是一個特別典型的模型錯誤。如果不加0-1變量約束只寫功率上下限優(yōu)化器會利用“低價充電的同時高價放電”這種虛假的套利空間把收益計算到極高結(jié)果完全失真。檢查辦法是看結(jié)果里max(min(P_ch, P_dis))是否接近0只要大于求解容差就說明約束沒起作用。我遇到過一次YALMIP模型里M取值過大導(dǎo)致求解器在數(shù)值上繞過了互斥約束。M取P_ess上界的1.1倍就夠不要順手寫10000。另外如果用了商業(yè)化求解器的高精度模式數(shù)值容差變小也會緩解這個問題。5.3 收益重復(fù)計算與口徑混淆復(fù)現(xiàn)論文時收益項口徑不一致是審稿人最愛挑的問題。有些原文同時寫了削峰填谷收益和調(diào)峰輔助服務(wù)收益表面上兩類實際卻存在放電量的重疊部分。同一個放電時段不能既算峰谷價差收益又算調(diào)峰補償收益。我的處理方式是在代碼注釋里明確每個時段的收益來源先按放電量拆分分類匯總后跟原文收益表逐項核對。具體拆分邏輯我這樣寫峰段放電量屬于套利電量按當(dāng)?shù)胤骞葍r差計收益在系統(tǒng)調(diào)峰缺口時段比如火電最小出力已經(jīng)壓到下限仍然無法平衡的時段放電的部分才計入調(diào)峰補償。代碼實現(xiàn)上就是增加一個調(diào)峰缺口標(biāo)記向量gap(t)收益函數(shù)改為R_arb sum(price_dis .* P_dis .* (1-gap) - price_ch .* P_ch) * year_days; R_reg c_reg * sum(P_dis .* gap) * year_days;這個gap向量來自系統(tǒng)的調(diào)峰需求計算雖然增加了建模復(fù)雜度但能讓經(jīng)濟性分析站得住腳。5.4 從論文復(fù)現(xiàn)走向?qū)嶋H項目時的三個方向擴展如果做的是真實項目而不是單純復(fù)現(xiàn)我建議在基礎(chǔ)模型之上做三方面擴展。第一把單典型日改成多季節(jié)加權(quán)場景每個場景按不同天數(shù)權(quán)重計入目標(biāo)函數(shù)這能反映儲能全年的真實利用率。第二計及儲能壽命衰減每充放一個循環(huán)容量輕微退坡長期運行成本會更貼近實際。第三如果是新能源配儲把棄電率約束或新能源消納收益納入目標(biāo)函數(shù)否則模型會傾向于“只做套利不做消納”與實際政策目標(biāo)不一致。這三個擴展在YALMIP框架下實現(xiàn)成本并不高無非是給SDP變量加下標(biāo)約束批量生成。我后來做工程項目時就是在復(fù)現(xiàn)模型基礎(chǔ)上加了多場景聚類和壽命衰減模塊整個代碼結(jié)構(gòu)沒有推倒重來。復(fù)現(xiàn)這類EI論文最大的體會是先把經(jīng)濟機制想清楚再寫代碼。目標(biāo)函數(shù)每一項對應(yīng)什么物理過程、約束每一條對應(yīng)什么運行邊界全部落到紙上Matlab代碼只是一個翻譯過程。先用小規(guī)模算例手算驗證再逐步加復(fù)雜度中間結(jié)果每步都要打印出來看合理性。按照這個流程即便原文沒有公開代碼也能通過復(fù)現(xiàn)理解每一處細節(jié)的來龍去脈拿到一份可信的儲能配置方案與經(jīng)濟性評估結(jié)論。