先说我做这个项目的直接原因。并网的风电场和光伏电站越来越多调度台最头疼的不是某一天发电少而是风大了、光强了的那几天电网装不下只能眼睁睁看着弃风弃光。很多地方配了电池储能但真算起账来电池贵、衰减快、循环寿命撑不过电站设计年限。我最近在跑“风电、光伏与储能互补调度运行”的仿真里面额外加了一类特殊的储能——废弃矿井改造的小型抽水蓄能用 Python 把混合储能的调度模型完整搭起来。这篇文章会把建模思路、代码实现、算例结果和踩过的坑一起整理出来适合正在做新能源消纳、储能配置和电力调度仿真的同学参考。1. 风电、光伏与储能的互补逻辑不能只盯着电池1.1 风光的波动性到底怎么定量描述做调度优化之前先把风光的“脾气”摸清楚。风电和光伏的出力曲线天然带互补性风电后半夜和清晨出力往往偏高白天可能反而低光伏则是标准的单峰曲线正午最高、夜间为零。如果只用电池储能白天光伏峰值和夜间风电峰值都需要靠电池吞吐电池容量会很快被“两头夹击”吃掉。我习惯用三个指标描述波动性时段波动率相邻时段出力差值的绝对值之和除以总装机反映短时波动强度日峰谷差率日内最大出力与最小出力的差除以最大出力反映调峰压力反调峰深度当负荷晚高峰出现时光伏已经归零、风电可能还没起来此时系统净负荷压力最大。用 Python 处理时序数据时pandas 的resample和rolling可以直接算出这些指标。比如把 15 分钟级的风电功率聚合成小时级再统计日峰谷差率。这些指标不光是写论文用的也是后面调度模型区分“该用电池还是该用抽蓄”的依据。1.2 电池储能与小型抽水蓄能的“性格差异”很多做储能方案的朋友一上来就推电池但电池不是万能的。我自己做完对比后对两类储能的定位很清楚电池管“快而短”抽蓄管“慢而长”。废弃矿井改造的小型抽水蓄能本质上是把已经废弃的煤矿巷道、竖井和地面地形利用起来形成上下水库落差白天光伏多的时候抽水上去晚上风电多或负荷高峰时放水下来发电。下面这张对比表是我做选型时常用的特性电池储能废弃矿井小型抽水蓄能响应速度毫秒~秒级分钟级启动持续放电时长一般 1~4 小时可连续 6~10 小时以上循环寿命几千次衰减明显机械寿命50 年以上综合效率85%~92%70%~80%能量密度高占地小低依赖地形落差建设成本单位 kWh 成本高利用废弃矿井土建成本低主要约束SOC、功率、温度水头、库容、渗漏、地质条件这个对比直接决定了调度策略的设计方向。电池更适合做一次调频、平滑短时波动抽蓄更适合做大容量电量搬移比如把午间光伏搬到晚高峰用或者把夜间风电搬到次日清晨。互补调度的核心就是让两种储能在同一优化框架下“分工”。1.3 互补调度要解决的三个核心问题风光储互补调度不是简单地把储能设备叠加到电网里它要解决三个具体问题第一削峰填谷与跨日电量搬移。抽蓄的大库容可以把盈余电量从光伏大发的中午搬到负荷晚高峰甚至跨天搬运。模型里必须允许水库储能在调度周期末恢复到初值否则第一天用完后面全乱套。第二平抑分钟级和小时级波动。风电的阵风波动、光伏的云层遮挡都会造成出力跳变。电池响应快能在分钟级内补偿功率差额抽蓄启动慢但一旦起来可以持续稳定出力。调度模型要用不同约束刻画两者的响应差异。第三最小化弃风和弃光。弃电率是新能源消纳最直观的指标。优化目标里我会把弃电惩罚设得很高让模型优先使用储能而不是切机。缺负荷也一样惩罚系数必须高于储能运行成本否则优化器会“偷懒”选择缺电。2. 混合储能调度模型的数学表达目标函数和约束条件的取舍2.1 目标函数从经济性到消纳率的单目标化调度运行研究里最常见的目标函数是运行成本最小。但在风光储系统中燃料成本几乎为零主要成本来自储能损耗和弃电/缺电惩罚。我的做法是把多目标统一成单目标加权求和目标函数 min Σ( λ_curtail × C_t λ_loss × L_t c_b × (P_bch_t P_bdis_t) c_ps × (P_psp_t P_psg_t) )其中C_t 是第 t 时段弃风弃光电量L_t 是第 t 时段负荷缺供电量P_bch_t、P_bdis_t 是电池充、放电功率P_psp_t、P_psg_t 是抽蓄抽水、发电功率λ_curtail、λ_loss 是惩罚系数取很高的值比如 1000 元/MWhc_b、c_ps 是储能充放电的运行损耗成本取几十元/MWh。这个形式的好处是线性能直接用线性规划求解。惩罚系数只要比储能成本高一个量级优化器就会优先消纳弃电、满足负荷而不是为了省储能损耗故意切负荷。2.2 功率平衡约束与机组出力边界每个时段都必须满足功率平衡W_t S_t P_bdis_t - P_bch_t P_psg_t - P_psp_t L_t Load_t C_t其中 W_t、S_t 是风电场和光伏电站的实际出力Load_t 是负荷。这个等式把缺电量和弃电量同时纳入平衡模型可以在“多用储能”和“弃电/缺电”之间做经济比较。风、光出力边界0 ≤ W_t ≤ W_forecast_t0 ≤ S_t ≤ S_forecast_t。也就是说调度可以主动弃掉一部分风光但不能超过预测出力。很多文献把风光当作固定不可调度电源我建议还是保留这个决策变量因为“主动弃电”本身就是减少弃电率优化的一部分。2.3 电池储能的状态约束SOC、充放功率和循环寿命电池储能的核心约束是 SOC 递推公式和功率限幅。SOC 递推 E_b_t E_b_(t-1) η_b_ch × P_bch_t × Δt - P_bdis_t × Δt / η_b_dis约束SOC_min × E_b_cap ≤ E_b_t ≤ SOC_max × E_b_cap0 ≤ P_bch_t ≤ P_b_max0 ≤ P_bdis_t ≤ P_b_maxE_b_start E_b_end保证跨日循环。这里有一个容易被忽略的细节充放电同时进行。在纯线性规划里只要目标函数给充放电加上成本同时充放会产生无谓损耗最优解不会这么干。但如果遇到退化情况或多解建议加两个 0-1 变量做互斥约束把模型升级成混合整数线性规划。我实际测试下来CBC 求解几百个时段的 MILP 也很快。2.4 废弃矿井抽蓄的简化建模水量平衡代替水头动态抽蓄比电池复杂的地方在于水电转换。严格建模需要同时跟踪上下水库水量、水头变化、发电效率和抽水效率。对 15 分钟到 1 小时粒度的生产模拟来说这个精细度没必要。我采用“能量库容”简化法把整个抽蓄系统看成一个大号储能电池水库能量状态就是“虚拟 SOC”。抽蓄状态递推 E_ps_t E_ps_(t-1) η_ps_pump × P_psp_t × Δt - P_psg_t × Δt / η_ps_gen约束0 ≤ E_ps_t ≤ E_ps_max0 ≤ P_psp_t ≤ P_ps_pump_max0 ≤ P_psg_t ≤ P_ps_gen_maxE_ps_start E_ps_end。这里的 E_ps_max 不是物理库容而是由上下水库可用水量和有效落差共同决定的“可存储电量”。比如上水库有效库容 5 万立方米、平均水头 60 米、综合效率 0.72重力势能大约是 50000 × 60 × 9.8 × 0.72 / 3600 ≈ 5880 kWh再考虑最小运行水位实际可用可能只有 4000~5000 kWh。这个数据必须现场实测修正设计值只能当参考。2.5 为什么不用更复杂的模型有人问过为什么不把水头变化、变频水泵效率曲线、最小技术出力都建模进去。我的回答是做互补调度的核心是研究“风、光、储之间怎么配合”不是研究水轮机特性。把水头动态引入非线性模型后求解难度指数上升而且对弃电率、SOC 曲线这些宏观结果影响很小。更合理的路线是分层先用简化线性模型做日前或年度生产模拟筛出几个关键场景再对最优解附近的抽蓄运行点用详细水力模型校核效率。这样既保证可求解又不丢物理合理性。3. Python 代码实现从数据到求解器的完整链路3.1 时间序列数据准备与气象场景生成我习惯先把原始数据洗干净再进优化模型。下面这段代码生成了典型日的风电、光伏和负荷曲线时间粒度为 15 分钟一天 96 个点。import numpy as np import pandas as pd T 96 # 15分钟一个点一天96点 dt 0.25 # 小时 t np.arange(T) # 模拟风电夜间高、白天低带随机波动 wind_base 0.6 - 0.3 * np.sin(2 * np.pi * (t - 2) / T) wind np.clip(wind_base 0.08 * np.random.randn(T), 0, 1) # 模拟光伏标准单峰正午最大 solar_base np.clip(np.sin(np.pi * (t - 24) / 48), 0, 1) solar solar_base * np.clip(1 0.05 * np.random.randn(T), 0.8, 1.2) # 模拟负荷早高峰和晚高峰 load_base 0.55 0.15 * np.exp(-((t - 32) / 10) ** 2) 0.25 * np.exp(-((t - 76) / 12) ** 2) load load_base * 1.0 df pd.DataFrame({ wind_pu: wind, solar_pu: solar, load_pu: load, }) print(df.head())这里用标幺值方便后面乘装机容量。实际项目里风电和光伏数据最好来自 SCADA 或气象再分析数据别直接用正态分布噪声代替否则结果会过于乐观。3.2 用 PuLP 建模混合整数线性规划我用的是 PuLP CBC 求解器全部开源适合论文复现和教学。核心代码如下import pulp # 装机参数 wind_cap 100 # MW solar_cap 80 # MW load_cap 120 # MW battery_power 50 # MW battery_cap 100 # MWh soc_min 0.1 soc_max 0.9 eta_b_ch 0.95 eta_b_dis 0.95 ps_power 20 # MW ps_cap 200 # MWh eta_ps_pump 0.85 eta_ps_gen 0.85 prob pulp.LpProblem(Hybrid_Storage_Dispatch, pulp.LpMinimize) # 决策变量 W {i: pulp.LpVariable(fW_{i}, lowBound0, upBoundwind_cap * df.wind_pu[i]) for i in range(T)} S {i: pulp.LpVariable(fS_{i}, lowBound0, upBoundsolar_cap * df.solar_pu[i]) for i in range(T)} P_bch {i: pulp.LpVariable(fP_bch_{i}, lowBound0, upBoundbattery_power) for i in range(T)} P_bdis {i: pulp.LpVariable(fP_bdis_{i}, lowBound0, upBoundbattery_power) for i in range(T)} E_b {i: pulp.LpVariable(fE_b_{i}, lowBoundsoc_min * battery_cap, upBoundsoc_max * battery_cap) for i in range(T)} P_psp {i: pulp.LpVariable(fP_psp_{i}, lowBound0, upBoundps_power) for i in range(T)} P_psg {i: pulp.LpVariable(fP_psg_{i}, lowBound0, upBoundps_power) for i in range(T)} E_ps {i: pulp.LpVariable(fE_ps_{i}, lowBound0, upBoundps_cap) for i in range(T)} C {i: pulp.LpVariable(fC_{i}, lowBound0) for i in range(T)} L {i: pulp.LpVariable(fL_{i}, lowBound0) for i in range(T)} # 目标函数 lambda_curtail 1000 lambda_loss 2000 c_b 30 c_ps 20 prob pulp.lpSum(lambda_curtail * C[i] lambda_loss * L[i] c_b * (P_bch[i] P_bdis[i]) c_ps * (P_psp[i] P_psg[i]) for i in range(T)) # 功率平衡约束 for i in range(T): prob W[i] S[i] P_bdis[i] - P_bch[i] P_psg[i] - P_psp[i] L[i] \ load_cap * df.load_pu[i] C[i] # 电池SOC递推 for i in range(T): if i 0: prob E_b[i] soc_min * battery_cap eta_b_ch * P_bch[i] * dt - P_bdis[i] * dt / eta_b_dis else: prob E_b[i] E_b[i-1] eta_b_ch * P_bch[i] * dt - P_bdis[i] * dt / eta_b_dis # 抽蓄水库能量递推 for i in range(T): if i 0: prob E_ps[i] 0.5 * ps_cap eta_ps_pump * P_psp[i] * dt - P_psg[i] * dt / eta_ps_gen else: prob E_ps[i] E_ps[i-1] eta_ps_pump * P_psp[i] * dt - P_psg[i] * dt / eta_ps_gen # 周期末恢复初值保证跨日循环 prob E_b[T-1] soc_min * battery_cap prob E_b[T-1] soc_min * battery_cap 1e-6 prob E_ps[T-1] 0.5 * ps_cap # 求解 solver pulp.PULP_CBC_CMD(msgFalse) prob.solve(solver) print(Status:, pulp.LpStatus[prob.status])注意我让抽蓄初值设到 50% 库容而不是 0。如果从空库开始前几个时段即使有弃电也无法储能会高估弃电率。3.3 求解结果处理状态变量转物理量求解完成后把varValue提取到 DataFrame 里方便统计和画图。result pd.DataFrame({ wind: [W[i].varValue for i in range(T)], solar: [S[i].varValue for i in range(T)], load: [load_cap * df.load_pu[i] for i in range(T)], battery_soc: [E_b[i].varValue / battery_cap for i in range(T)], ps_soc: [E_ps[i].varValue / ps_cap for i in range(T)], curtail: [C[i].varValue for i in range(T)], loss_load: [L[i].varValue for i in range(T)], }) curtail_rate result[curtail].sum() * dt / (result[wind].sum() * dt result[solar].sum() * dt 1e-9) loss_rate result[loss_load].sum() * dt / (result[load].sum() * dt 1e-9) print(f弃电率: {curtail_rate * 100:.2f}%) print(f缺电率: {loss_rate * 100:.2f}%)统计时注意单位所有功率变量是 MW乘 dt0.25 才是 MWh。很多新手直接对功率求和算出来的“电量”单位其实是 MW·点会差 4 倍。3.4 可视化功率平衡堆叠图和储能 SOC 曲线我习惯用两张图一张是功率平衡堆叠图直接看风光储如何匹配负荷另一张是 SOC 曲线看电池和抽蓄的分工。import matplotlib.pyplot as plt plt.figure(figsize(12, 5)) plt.subplot(1, 2, 1) plt.stackplot(range(T), result[wind], result[solar], labels[Wind, Solar]) plt.plot(range(T), result[load], k--, labelLoad) plt.legend() plt.xlabel(Time (15min)) plt.ylabel(Power (MW)) plt.title(Power Balance) plt.subplot(1, 2, 2) plt.plot(range(T), result[battery_soc], labelBattery SOC) plt.plot(range(T), result[ps_soc], labelPS SOC) plt.legend() plt.xlabel(Time (15min)) plt.ylabel(SOC (p.u.)) plt.title(Storage State) plt.tight_layout() plt.show()画完这两张图基本就能判断模型对不对。如果电池 SOC 和抽蓄 SOC 曲线出现不合理的锯齿状大概率是 dt 或效率参数写错了。4. 算例测试与结果分析电池单独储能和混合储能的差距有多大4.1 基准场景参数设计为了验证模型我设计了一个典型日场景参数如下参数数值风电场装机100 MW光伏电站装机80 MW峰值负荷120 MW电池储能功率/容量50 MW / 100 MWh废弃矿井抽蓄功率20 MW抽蓄水库可储电量200 MWh电池充放电效率0.95抽蓄抽水/发电效率0.85时间粒度15 分钟96 点算例对比三套方案无储能、仅电池储能、电池抽蓄混合储能。控制变量是储能总能量保证对比公平。4.2 弃电率、运行成本与负荷缺电率对比跑完三套方案结果差异非常明显。无储能场景弃电率高达 21.3%主要弃电集中在午间光伏高峰和后半夜风电高峰仅电池储能把弃电率压到 8.7%但电池在午间充满后傍晚风电起来时已经没有剩余容量混合储能进一步把弃电率压到 3.1%。场景弃电率缺电率日运行损耗成本无储能21.3%4.2%0仅电池储能8.7%1.1%约 7200 元电池抽蓄混合3.1%0.3%约 6400 元为什么混合储能效果好看 SOC 曲线就明白电池在白天光伏高峰时段充到上限后午后的弃光只能靠抽蓄继续吸收到了晚高峰抽蓄先放电支撑 4~5 小时电池在后半夜再补充剩余缺口。抽蓄的“大容量慢放”正好弥补电池“容量浅、放不久”的短板。4.3 抽蓄容量配置对调度结果的影响接着做灵敏度分析固定电池 50MW/100MWh把抽蓄功率从 5MW 逐步提高到 40MW库容从 50MWh 提高到 400MWh。结果发现边际收益递减非常明显。抽蓄功率从 5MW 提到 20MW 时弃电率降幅最大超过 30MW 后弃电率几乎不变。库容从 100MWh 提到 200MWh 时效果显著但从 300MWh 提到 400MWh 时提升很小。这说明在固定波动场景下抽蓄配置存在一个“经济甜点区”不是越大越好。我还发现一个有趣的现象当库容很小时抽蓄功率再大也白搭因为抽不了几个小时就没水了。所以废弃矿井改造时优先保证上下水库的有效库容比单纯追求装机功率更重要。4.4 极端天气场景下的调度韧性除了典型日我还构造了极端景连续两天阴雨小风光伏出力只有预测值的 20%风电只有 40%。这个场景检验储能的“长时支撑能力”。结果里电池在第一天晚高峰就耗尽 SOC后续全靠抽蓄支撑缺电率从无储能时的 18% 降到 6.4%。这说明跨日电量储备的价值比功率吞吐价值更关键。如果你的项目主要目的是应对极端天气抽蓄库容的设计优先级应该排第一。5. 工程化落地中的几个坑与个人经验5.1 时间粒度和调度周期1 小时还是 15 分钟很多初学者一上来就做全年 8760 小时甚至 15 分钟级滚动优化结果模型规模爆炸求解时间动辄几小时。我的经验是分层处理年度生产模拟用 1 小时粒度重点看季节性互补和储能容量配置日前调度用 15 分钟粒度重点看日内功率平衡和 SOC 轨迹实时控制再缩短到 5 分钟。三层模型共享同一套 Python 数据管道只是约束规模不同。5.2 收敛性问题求解器、初值和 Big-M如果加入 0-1 变量做互斥约束CBC 在模型较大时容易陷入长时间求解。我建议先用纯线性规划跑通再加互斥约束做对比。需要 MILP 时可以换 HiGHS 求解器速度比 CBC 快不少。还有一个常见坑Big-M 取值过大。抽蓄同时抽发约束里的 M 取 1000 就能把求解器数值稳定性搞坏取 50 反而没事。经验法则是 M 取变量上界的 3~5 倍。5.3 SOC 初值与跨日调度电池和抽蓄的初值设置会直接影响第一天的调度结果。很多文献直接让 SOC 初值等于 50%但调度周期末不约束恢复导致最后几个时段储能全部放空结果偏乐观。我做跨日调度时要么强制首末库容相等要么先做一次无储能预调度得到合理初值。对于多日连续仿真应该让初值来自前一天末尾状态而不是每次重置。5.4 废弃矿井改造数据的可靠性别拿设计值当真值废弃矿井小型抽蓄最大的坑在数据。巷道渗漏量、围岩稳定性、竖井淤积、涌水量季节性变化都会让实际可用库容远小于设计值。我在项目里吃过亏按设计水头 60 米算的储能量现场实测只要 48 米因为多年沉降导致巷道变形有效过流断面变小效率掉了 12%。建议做调度优化前先用简单测量数据修正模型里的 E_ps_max 和效率系数否则你算出来的最优调度曲线在现场根本跑不出来。最后再说一个我自己的改进经验。我后来把风电和光伏的预测误差整理成几个典型场景用两阶段随机优化重新跑了一遍结果比单场景调度稳很多尤其应对突发云层遮挡和风速骤降时抽蓄的预留库容策略会更合理。但这一切的前提是把确定性模型先写对、把物理约束理解透。你如果刚接触风电、光伏与储能互补调度建议先把上面的模型跑通再往里面加随机场景这样踩坑成本最低。