干过电力系统调度优化的人应该都有体会：风电、光伏和储能放在一起做联合优化调度，难度不是简单叠加。风光的随机性、储能的多时间尺度特性、不同储能形式的性能差异，任何一个环节没处理好，优化结果就会明显失真。我最近用Python把这套“风电+光伏+电池储能+废弃矿井小型抽水蓄能”的互补调度模型完整跑了一遍，从数据生成、数学建模、求解到可视化都做了落地，今天把过程中的思路、关键代码和踩过的坑一并分享出来。
这篇文章适合两类人：一是正在做新能源电力系统优化、储能调度方向研究的学生，论文里需要一套能跑通、能复现的基准模型；二是实际搞微电网或园区综合能源运行的工程师，想用Python把风光储联合调度从概念验证推向工程原型。我会把目标函数、约束条件、代码实现和调试经验都讲透，特别是废弃矿井小型抽水蓄能这个比较少见的储能形式，它在建模上和电池有本质区别，我会单独拿出来讲。
1. 先理清思路：三种电源的“脾气”和互补逻辑
1.1 风电和光伏的出力特征决定了调度为什么难
风电出力的核心特点是随机性和反调峰性。风速本身服从威布尔分布，折算到功率之后，波动幅度会被放大——风速增加一倍，理论上功率可能增加八倍。尤其在我国北方地区，夜间风速往往高于白天，风电出力高峰恰好落在负荷低谷时段，这就是典型的反调峰现象。光伏出力则完全不同，它严格依赖太阳辐射，白天出现高峰，午间甚至可能超过负荷需求，但夜间出力为零，而且受云层遮挡影响，分钟级波动非常剧烈。
这两种电源放到一起，单纯从“量”上看确实有互补性——白天靠光伏，夜间靠风电，但“互补”不等于“稳定”。调度的第一个难点就是：风光出力曲线和负荷曲线在时间尺度上错位，在波动尺度上叠加。负荷高峰通常在上午和晚间，而光伏高峰在正午，风电高峰在深夜，三个曲线凑在一起，如果不做主动的储能调节，要么大量弃风弃光，要么系统需要频繁从电网买电来平衡短时功率缺口。
1.2 储能不是“一个东西”：电池和抽水蓄能扮演的角色完全不同
很多人把储能当成一个统一的概念建模，实际上电池储能和小型抽水蓄能的特性差异非常大，这个差异直接决定了调度模型的约束写法。
电池储能是典型的功率型储能：响应速度快，秒级甚至毫秒级就能完成功率调节，往返效率高（目前锂电系统交流侧往返效率能做到90%以上），但它的能量容量和功率是耦合的——想多存电就得加电芯，成本线性上涨。所以电池适合处理小时级以内的高频波动，相当于系统里的“短期缓冲器”。
抽水蓄能则是典型的能量型储能：它靠上下水库之间的水位差存储势能，能量容量由库容决定，做大容量的边际成本很低，但响应速度慢（从机组启动到满功率通常需要几分钟到十几分钟），往返效率也偏低（综合效率一般70%到80%）。它适合处理跨小时的日内能量搬移，比如把午间过剩的光伏电能转移到晚高峰使用。废弃矿井小型抽蓄就是这个思路的延伸——利用废弃矿井的竖井和巷道高差，建设数十米到数百米水头的分布式抽蓄电站，单机功率一般在0.5MW到10MW量级。
1.3 为什么我选了Python而不是其他工具
调度优化在工业界常用GAMS、MATLAB加YALMIP这些工具，但我在这个项目里坚持用Python，核心原因有三个。第一，Python的数据生态和优化生态都极其完整：pandas处理时间序列、numpy做矩阵运算、scipy直接提供线性规划和混合整数线性规划求解器，matplotlib出图，整个链路不需要换工具。第二，可复现性强，论文审稿或团队协作时，工程的上一段代码在任何机器上能跑出一致结果，这是很多图形化工具做不到的。第三，后续如果要接入机器学习功率预测模型（比如LSTM、Transformer），Python天然无缝衔接，预测输出可以直接作为调度模型的输入，不用做跨语言的接口转换。
当然，Python自带的求解器在大规模问题上性能确实不如商业求解器Gurobi、CPLEX，但在教学演示、方案验证、中小规模微电网调度场景下已经足够。这篇文章的所有代码，只要安装numpy、scipy、matplotlib三个库就能跑通。
2. 调度模型怎么搭：从目标函数到约束条件的完整推导
2.1 数据准备：预测曲线从哪里来，怎么生成测试场景
调度模型的前提是已知未来一段时间（本文以24小时为例）的风电、光伏可用出力和负荷曲线。如果有数值天气预报数据，可以用功率预测模型生成这些序列；没有实测数据时，用统计分布模拟典型场景是学术界和工程上通用的做法。
风电可用出力我用威布尔分布随机生成风速序列，再通过风速-功率特性曲线转成功率。公式是：切入风速以下出力为0，额定风速以上封顶，中间段按三次方关系计算。光伏出力用晴空模型乘以随机云遮系数模拟，负荷曲线则在基础负荷上叠加早晚两个高峰。这样生成的场景虽然不是某一天的真实数据，但保留了风光出力的典型统计特征，足以验证调度模型的正确性。
时间分辨率方面，本文取1小时、24个时段，这对概念验证是合适的。实际工程中建议加密到15分钟甚至5分钟，因为风光分钟级波动对电池充放电策略影响很大，但分辨率越高，变量数量越大，求解时间越长，需要根据场景规模权衡。
2.2 目标函数：运行成本最小化，弃电惩罚必须加
调度模型的目标函数我设计为：购电成本加上弃风弃光惩罚的最小化。公式可以写成：
min Σ [ price(t) × P_import(t) ] + λ × Σ [ (P_wt_fc(t) − P_wt(t)) + (P_pv_fc(t) − P_pv(t)) ]
其中P_import(t)是t时段从电网购入的功率，price(t)是分时电价，P_wt(t)和P_pv(t)是实际消纳的风电和光伏功率，P_wt_fc(t)和P_pv_fc(t)是预测可用出力，λ是弃电惩罚系数。
这个惩罚项非常关键。如果不加它，优化器为了减少购电成本，会在某些时段故意不消纳风电，用网购电代替——这在经济上可能成立，但完全违背了新能源调度的初衷。λ的取值一般要高于最低电价，我在这组参数里取200元/MWh，目的是让优化器只有在“实在消纳不了”的情况下才弃电。
这里有一个实现细节值得注意：由于弃电量等于预测值减去实际消纳值，目标函数展开后会出现常数项λ×Σ(P_wt_fc+P_pv_fc)。常数项不影响优化结果，我在代码里直接把它剥离了，所以在目标函数系数中，P_wt(t)和P_pv(t)对应的系数是负的λ。第一次写的时候如果不理解这一点，很容易把自己绕晕。
2.3 约束条件：功率平衡、储能递推、互斥逻辑
约束条件是调度模型的骨架，我把它分成四组。第一组是系统功率平衡约束：任意时刻，所有电源出力（风电+光伏+电池放电+抽蓄发电+购电）必须等于负荷加上所有充电功率（电池充电+抽蓄泵水）。
第二组是储能能量状态递推。电池的模型是：E_bat(t+1) = E_bat(t) + η_c × P_bat_c(t) × Δt − P_bat_d(t) / η_d × Δt，其中η_c和η_d分别是充放电效率。抽蓄的模型形式上完全一样，只是把充放电换成泵水和发电工况，效率换成抽蓄的泵水效率和发电效率。
第三组是储能功率约束，这里有个常见的坑：电池不能同时充电和放电。这个“互斥逻辑”是0-1整数约束，必须引入二进制变量才能严格表达。如果只用线性规划而忽略互斥，优化器很可能给出“一边充电一边放电”的荒谬结果——因为充电和放电在功率平衡方程里符号相反，同时进行意味着无意义的能量循环，成本上看似没增加，实际是物理上不可能的。
第四组是上下限约束：电池SOC维持在10%到90%之间，抽蓄等效能量状态维持在10%到90%之间，购电功率不超过联络线容量上限。
3. Python实现：数据生成、MILP建模与结果可视化全流程
3.1 工程结构：模块化拆分，别把代码都堆在一个文件里
项目我拆成四个部分：数据生成模块、参数配置模块、优化求解模块、可视化模块。这样做的好处是改参数时不用翻一长串代码，跑实验时也可以单独替换数据生成逻辑。下面这段代码是数据生成和参数配置的核心部分，直接附上。
import numpy as np from scipy.optimize import milp, LinearConstraint, Bounds import matplotlib.pyplot as plt # ---------- 1. 生成测试场景数据 ---------- rng = np.random.default_rng(42) T = 24 # 调度周期24小时 dt = 1 # 时段时长1小时 # 风电功率: 威布尔风速 -> 功率曲线 v = rng.weibull(2.2, T) * 8 v_in, v_r, v_out = 3.0, 12.0, 25.0 P_wt_fc = np.zeros(T) mask1 = (v >= v_in) & (v < v_r) mask2 = (v >= v_r) & (v <= v_out) P_wt_fc[mask1] = 150 * (v[mask1]**3 - v_in**3) / (v_r**3 - v_in**3) P_wt_fc[mask2] = 150 # 光伏功率: 晴空模型 x 云遮系数 h = np.arange(T) clear = np.clip(np.sin(np.pi * (h - 6) / 12), 0, 1) cloud = np.clip(rng.normal(0.8, 0.15, T), 0, 1.2) P_pv_fc = 100 * clear * cloud # 负荷: 基础负荷 + 早晚高峰 + 随机扰动 load_base = 200 load_peak1 = 30 * np.exp(-((h - 8) / 3)**2) load_peak2 = 40 * np.exp(-((h - 19) / 4)**2) P_load = load_base + load_peak1 + load_peak2 + rng.normal(0, 5, T) # 分时电价(元/MWh): 峰平谷三段 price = np.where((h >= 10) & (h <= 15), 900, np.where((h >= 18) & (h <= 22), 800, 400))3.2 核心模型：scipy的MILP怎么搭矩阵
这里我直接给出混合整数线性规划的完整建模代码。变量排列是关键，我按每个时段11个变量组织：风电消纳、光伏消纳、电池充电、电池放电、抽蓄泵水、抽蓄发电、购电功率、电池SOC、抽蓄能量状态、电池互斥标志、抽蓄互斥标志。这样在构建约束矩阵时，通过索引偏移就能精确定位变量。
# ---------- 2. 系统参数 ---------- lambda_curtail = 200.0 # 弃电惩罚系数 P_import_max = 50.0 # 联络线上限 P_bat_rated = 20.0 # 电池额定功率MW E_bat_max = 100.0 # 电池容量MWh E_bat_init = 50.0 # 初始SOC对应能量 eta_bat_c, eta_bat_d = 0.95, 0.95 P_ph_rated = 30.0 # 抽蓄额定功率MW E_ph_max = 300.0 # 抽蓄等效能量容量MWh E_ph_init = 150.0 eta_ph_pump, eta_ph_gen = 0.85, 0.90 # ---------- 3. 构建优化模型 ---------- NVAR = T * 11 c = np.zeros(NVAR) for t in range(T): base = t * 11 # 剥离常数项后, 多消纳风电/光伏相当于减少罚金, 系数为负 c[base] = -lambda_curtail # 风电 c[base + 1] = -lambda_curtail # 光伏 c[base + 6] = price[t] # 购电 rows = [] lb = [] ub = [] for t in range(T): base = t * 11 idx_wt, idx_pv = base, base + 1 idx_bc, idx_bd = base + 2, base + 3 idx_pp, idx_pg = base + 4, base + 5 idx_imp = base + 6 idx_eb, idx_ep = base + 7, base + 8 idx_yb, idx_yp = base + 9, base + 10 # 功率平衡: 源端 = 荷端 row = np.zeros(NVAR) row[[idx_wt, idx_pv, idx_bd, idx_pg, idx_imp]] = 1 row[[idx_bc, idx_pp]] = -1 rows.append(row) lb.append(P_load[t]) ub.append(P_load[t]) # 电池互斥: 充电时不能放电 row = np.zeros(NVAR) row[idx_bc] = 1 row[idx_yb] = -P_bat_rated rows.append(row) lb.append(-np.inf) ub.append(0) row = np.zeros(NVAR) row[idx_bd] = 1 row[idx_yb] = P_bat_rated rows.append(row) lb.append(-np.inf) ub.append(P_bat_rated) # 抽蓄互斥: 泵水工况与发电工况互斥 row = np.zeros(NVAR) row[idx_pp] = 1 row[idx_yp] = -P_ph_rated rows.append(row) lb.append(-np.inf) ub.append(0) row = np.zeros(NVAR) row[idx_pg] = 1 row[idx_yp] = P_ph_rated rows.append(row) lb.append(-np.inf) ub.append(P_ph_rated) # 电池SOC递推: E(t) - E_start - eta_c*P_c + P_d/eta_d = 0 row = np.zeros(NVAR) row[idx_eb] = 1 row[idx_bc] = -eta_bat_c * dt row[idx_bd] = dt / eta_bat_d if t == 0: rhs = E_bat_init else: row[(t - 1) * 11 + 7] = -1 rhs = 0 rows.append(row) lb.append(rhs) ub.append(rhs) # 抽蓄能量递推 row = np.zeros(NVAR) row[idx_ep] = 1 row[idx_pp] = -eta_ph_pump * dt row[idx_pg] = dt / eta_ph_gen if t == 0: rhs = E_ph_init else: row[(t - 1) * 11 + 8] = -1 rhs = 0 rows.append(row) lb.append(rhs) ub.append(rhs) # 变量边界 bounds = [] for t in range(T): base = t * 11 bounds += [(0, P_wt_fc[t]), (0, P_pv_fc[t]), (0, P_bat_rated), (0, P_bat_rated), (0, P_ph_rated), (0, P_ph_rated), (0, P_import_max), (0.1 * E_bat_max, 0.9 * E_bat_max), (0.1 * E_ph_max, 0.9 * E_ph_max), (0, 1), (0, 1)] # 固定初始能量状态 bounds[0 * 11 + 7] = (E_bat_init, E_bat_init) bounds[0 * 11 + 8] = (E_ph_init, E_ph_init) # 整数变量位置: 每个时段的第10、11个变量 integrality = np.zeros(NVAR) integrality[9::11] = 1 integrality[10::11] = 1 # 求解 constraints = LinearConstraint(np.array(rows), np.array(lb), np.array(ub)) res = milp(c=c, integrality=integrality, bounds=Bounds(np.array([b[0] for b in bounds]), np.array([b[1] for b in bounds])), constraints=constraints, options={'time_limit': 10}) assert res.success, f"求解失败: {res.message}" x_opt = res.x.reshape(T, 11) # 结果指标 cost_total = sum(price[t] * x_opt[t, 6] for t in range(T)) curtail_wt = sum(P_wt_fc[t] - x_opt[t, 0] for t in range(T)) curtail_pv = sum(P_pv_fc[t] - x_opt[t, 1] for t in range(T)) print(f"总购电成本: {cost_total:.1f} 元") print(f"弃风量: {curtail_wt:.2f} MWh, 弃光量: {curtail_pv:.2f} MWh")这一段代码有几个地方是新手最容易写错的。第一个是SOC递推公式里系数的符号，充电效率乘在充电功率上，放电效率除在放电功率上，一旦写反，结果里储能会像“永动机”一样越充越多。第二个是互斥约束的方向，p_bd + P_rated × y_bat ≤ P_rated 的含义是：当y_bat=1（充电状态）时，放电功率上限为0；当y_bat=0（非充电）时，放电功率上限为额定值。写反了会导致约束完全失效。
3.3 结果可视化：光看数字不够，调度曲线才能暴露问题
优化结果最终要落到曲线上看，不然很难发现异常。我习惯用三个子图展示：上面是源荷功率曲线，中间是储能功率与能量状态，下面是电网交互与弃电情况。
fig, axes = plt.subplots(3, 1, figsize=(11, 10), sharex=True) ax = axes[0] ax.plot(h, P_wt_fc, label='风电可用出力') ax.plot(h, P_pv_fc, label='光伏可用出力') ax.plot(h, x_opt[:, 0] + x_opt[:, 1], '--', label='实际消纳') ax.plot(h, P_load, color='gray', label='负荷') ax.set_ylabel('功率(MW)') ax.legend() ax = axes[1] ax.step(h, x_opt[:, 2] - x_opt[:, 3], where='mid', label='电池充放(+充/-放)', color='tab:blue') ax.step(h, x_opt[:, 4] - x_opt[:, 5], where='mid', label='抽蓄泵水/发电', color='tab:green') ax.plot(h, x_opt[:, 7], color='tab:blue', linestyle='--', label='电池SOC') ax.plot(h, x_opt[:, 8], color='tab:green', linestyle='--', label='抽蓄能量状态') ax.set_ylabel('功率(MW)/能量(MWh)') ax.legend() ax = axes[2] ax.bar(h, x_opt[:, 6], label='购电功率', color='tab:orange') ax.plot(h, P_wt_fc - x_opt[:, 0], color='tab:red', label='弃风') ax.plot(h, P_pv_fc - x_opt[:, 1], color='tab:purple', label='弃光') ax.set_xlabel('时刻(h)') ax.set_ylabel('功率(MW)') ax.legend() fig.tight_layout() plt.show()正常的调度结果长什么样？简单说：低谷时段电价低，电池和抽蓄都倾向于充电（购电做能量搬移）；高峰时段优先用储能放电，而不是高价购买电网电力；弃风弃光只在风光大发且储能满仓、电网无法消纳的极端时段出现。如果曲线完全违背这些直觉，那大概率是某个约束写错了。
4. 废弃矿井抽水蓄能：特性、建模与协同要点
4.1 为什么废弃矿井能重新变成储能电站
废弃矿井做抽水蓄能，本质上利用的是“高度差”和“地下空间”。煤矿关停后，竖井的深度一般在200米到800米之间，巷道网络形成巨大的可储水空间。工作原理很简单：上水库放在井口附近的地表蓄水池，下水库利用井底巷道空间，需要储能时把上水库的水放到井下，需要发电时把井下蓄水抽回井上，水流经水轮机带动发电机。
相比传统抽水蓄能电站动辄上千万立方米库容，废弃矿井小型抽蓄的优势在于分布式和灵活性。单机功率小，选址限制少，而且把一个废弃矿山改造成储能设施，本身就实现了资源再利用，减少地面沉降和水位变化带来的环境隐患。更重要的是，这种储能形式能量容量由地下空间决定，扩容成本低——矿井巷道挖好了就是现成的库容，不需要像电池一样每增加一度电容量都要采购电芯。
4.2 建模差异：效率、规模和运行约束都和电池不一样
从调度模型的视角，电池和抽蓄虽然都可以抽象成“能量状态+功率输入输出”，但参数和约束的差异很大。第一是效率结构：电池充放电效率可以分别设置，而抽蓄的泵水效率和发电效率也分开设置，但综合往返效率明显更低——我做参数设置时电池往返效率约90%，抽蓄综合往返效率约76%，这意味着抽蓄不适合频繁的短时充放循环，否则效率损耗会把经济性吃掉。
第二是能量容量规模：本文中抽蓄等效能量容量设为300MWh，是电池的3倍，这是刻意体现两者定位的差异——电池负责短时功率调节，抽蓄负责长时能量搬移。第三是运行约束的颗粒度：真实的抽蓄机组有最小技术出力、启停次数限制、工况转换时间等约束，在这些实验代码中没有全加，但做工程落地时必须考虑，否则调度指令下发到现场根本执行不了。
4.3 与电池储能的协同：一个管高频，一个管搬运
两种储能协同调度的价值在结果曲线上体现得非常直观。看优化结果时你会发现：电池的SOC曲线波动频繁，充电和放电的切换次数多，明显在跟踪分钟级到小时级的功率波动，充当“高频缓冲器”的角色；而抽蓄的能量状态呈现平滑的长周期摆动——通常是夜间泵水蓄能、午间也吸收多余光伏、早晚高峰发电，一天只切换一两个完整循环。
这种协同模式背后有明确的物理和经济逻辑：电池响应快但容量有限，让它承担基荷级的长时间能量搬移，会过早耗尽容量且循环次数激增，影响寿命；抽蓄效率低但容量大、每度电的存储成本便宜，让它频繁启停去应对短时波动，不仅响应跟不上，效率损失也不划算。MILP模型天然就能捕捉到这个分工，因为优化器会权衡两种储能在效率、容量和功率上的差异，自动把任务分配给“更擅长”的一方。这就是为什么模型里必须把两种储能分开建模，而不是合并成一个“总储能”变量。
5. 常见报错与结果异常排查：我的调试实录
5.1 求解失败和矩阵装配错误
我调试过程中遇到最多的报错是维度和边界不匹配。scipy的milp函数要求约束矩阵的维度、lb/ub数组的长度、bounds的个数必须严格一致，任何一处少了一个变量，都会抛出ValueError。排查方法很简单：在求解前打印NVAR、矩阵行数、lb长度这三个数值，核对是否满足“行数等于lb长度”这个基本关系。
更难定位的是逻辑错误而不是语法错误。比如SOC递推约束中，如果t>0时忘记把上一时段的e_bat[t-1]加入矩阵，约束就变成每一时段只与初始值相关，储能能量会在一天结束后出现明显漂移。这类错误不会导致求解失败，但结果会非常离谱，需要靠绘图才能发现。
5.2 结果“看起来不对”的几个典型症状
症状一：电池和抽蓄同时大功率充放电。这几乎可以肯定是互斥约束写错了，比如互斥方向弄反，或者用了近似松弛却没加功率合限。症状二：SOC长期贴着上限或下限。这说明惩罚权重或电价结构有问题——如果目标函数里缺少对“储能被用完”的约束，优化器会在最后几个时段把储能全部放空，以便少购电，单日调度勉强可以，但滚动调度次日就无电可用。我之前在单日测试里没加末态约束，结果最后两小时抽蓄强制满功率发电，电量全部清零。加一条E(T) ≥ E_init的约束，这个问题立刻解决。
症状三：弃风弃光量为全零，但购电成本异常高。这通常是λ取值过大，优化器宁可高价购电也不弃新能源，逻辑上没毛病，但违背了经济调度的常识。λ需要标定，我建议和最低电价值对比：λ略高于低谷电价，系统才会优先消纳新能源。
5.3 参数敏感性：效率和容量的影响比想象中大
最后说一个经验：在做参数分析时，储能效率对调度结果的影响远大于容量。我把电池效率从0.95降到0.85，整个系统的购电成本立刻上升十几个百分点，因为每一次储能循环都在“亏电”。而把电池容量翻倍，对成本的影响反而不明显，因为成本瓶颈在系统高峰时段的功率支撑能力，而不是总能量。
这引出一个更深的体会：调度优化的重点不只是“让模型跑通”，而是理解参数变化时系统行为的迁移规律。比如废弃矿井抽蓄的综合效率如果实测下来只有65%而不是模型的76%，它的调度优先级会大幅下降，很多原本分配给它的能量搬移任务会被电池接管。建模时留出参数接口，运行后做敏感性分析，这个习惯能帮你避开很多“模型做完才发现参数离谱”的坑。
这套代码我目前还在持续迭代，最近在尝试把预测模块从统计分布换成长短期记忆网络，并引入滚动时域控制框架，让每一天的调度决策都能根据最新的预测信息修正。从个人体验来看，先把线性规划和混合整数规划版本的调度闭环跑通，把约束条件、数据流和参数标定这些基础打牢，再逐步加复杂度，是最高效的路径。希望这篇分享能帮你少走几步弯路。
