搞微电网调度的人,大多经历过这种场面:前一天用优化算法排出来的调度曲线,看起来无懈可击,第二天光伏被云一挡,负荷又往上跳,曲线当场废掉——该买电的时候没买,该放储的时候没放,月底电费单比预期高出一截。问题就出在“开环”这两个字上:所有计划只基于预测做一次,之后就不再修正。MPC(模型预测控制)恰恰是冲着这个痛点来的。这篇文章用Python把基于MPC的微电网调度优化完整走一遍,从数学建模、滚动求解、仿真对比到项目里踩过的坑都讲透,适合正在做微电网优化方向课题、或者刚接触光储微电网项目的朋友参考。
如果你只是想要一个能跑的demo,文中的核心代码可以直接抄;如果你是想弄明白为什么非要用MPC、为什么模型要写成那个样子,建议把前两章也看完。我不会堆一堆公式就完事,重点放在“为什么这样建模、代码怎么落地、仿真怎么设计才公平”这三件事上。
1. 为什么微电网调度要把“开环计划”改成“滚动优化”
1.1 日前计划为什么会在现场失效
先还原常见做法。调度员或者优化程序在前一天基于光伏预测、负荷预测和电价信息,一次性求出未来24小时的机组出力和储能充放电计划,全天按这套计划执行。这在预测完全准确的情况下当然没问题,但现场不是这样。光伏预测最常见的误差是辐照突然下降,十分钟内出力从800kW掉到400kW,预测模型根本来不及在当天修正。负荷侧的波动虽然没有光伏那么剧烈,但叠加起来,计划值和实际值对不上是常态。
对不上的后果有两种:一种是储能被过度放电,SOC提前触底,晚高峰想放电时没电可用;另一种是并网点实际功率偏离计划,由大电网被动兜底,多出来的购电成本记在电费单上。做过现场项目的人都知道,所谓的“优化调度”在偏差一大时,结算成本可能比简单规则策略还高,就是因为开环结构没有反馈。
1.2 MPC真正改变的,不是“预测”而是“反馈”
MPC的英文全称叫Model Predictive Control,中文叫模型预测控制。很多人第一反应是“它预测能力很强”,其实不是。MPC的核心动作是“每走一步重算一次”,把预测只是当作出力规划的依据。
具体流程是这样的:在k时刻,以当前实测储能状态为初值,基于未来N步的预测数据,求解一个有限时域优化问题,得到未来N步的完整控制序列,但只执行序列里的第一步;等到k+1时刻,用实测状态更新初值,重新生成预测,再求解,再只执行第一步。这个过程就是滚动优化。类比一下就是导航:出发时规划的路线只作为参考,实际驾驶中每两分钟根据实时拥堵情况重新算一次路,而不是死守着出发时那条路线。
所以MPC本质上是一个闭环控制器,状态反馈进了优化问题的初始条件。这也是它跟“日前计划+实时调整”最大的区别:后者是两套系统硬拼,前者把反馈直接内置在优化框架里。
1.3 这套方法适合什么场景
基于MPC的微电网调度,比较适合源荷波动大、储能占比较高、又有多时段电价的环境。典型场景是园区光储微电网、光储充一体化场站、带柴油发电机备用的离网/并网混合微网。这些场景的共同特点是:分钟级的调节需求明显,但对计算时间相对宽容,15分钟一个控制周期完全来得及跑完一个中等规模的优化问题。
它不适合做什么?不适合做秒级甚至毫秒级的动态控制,那是实时控制和一次调频的范畴。MPC在这里解决的是经济性调度问题,不是电力电子层面的暂态稳定问题。明白这一点,后续建模才不会用力过猛。
2. MPC调度优化的数学模型:先写清楚问题再写代码
2.1 系统拓扑与决策变量
我习惯先把微电网简化成一张最小拓扑,再决定哪些是决策变量、哪些是外部参数。这张最小拓扑包含:光伏、负荷、储能、一个可调的柴油机组(可选项)、与大电网的并网点。光伏和负荷的预测值是外部参数;优化程序要决定的是:从电网买多少电、卖多少电、柴油机发多少电、储能充多少电、放多少电、是否需要弃掉一部分光伏。
决策变量整理成表格,写代码时可以照着定义cvxpy变量:
| 变量 | 含义 | 单位 | 备注 |
|---|---|---|---|
| P_buy[t] | 时刻t向电网购电功率 | kW | 非负 |
| P_sell[t] | 时刻t向电网售电功率 | kW | 非负 |
| P_dg[t] | 柴油机组出力 | kW | 非负 |
| P_ch[t] | 储能充电功率 | kW | 非负 |
| P_dis[t] | 储能放电功率 | kW | 非负 |
| P_curt[t] | 弃光功率 | kW | 非负 |
| E[t] | 储能剩余能量 | kWh | 状态量 |
一个容易被新手卡住的点:储能状态用能量E[t]而不是SOC。SOC是百分比,乘上容量才是能量,建模时如果把SOC当状态变量,约束里到处要除容量,效率系数也好算错。用E更直接,最后画图再转成SOC = E / E_cap。
2.2 目标函数:每一项成本都对应一个真实问题
调度优化的目标一般是“运行成本最小”。不同项目会加不同惩罚项,但骨架很固定,写成:
min Σ_{t} [ c_buy[t]·P_buy[t] − c_sell[t]·P_sell[t] + f(P_dg[t]) + λ_ch·P_ch[t] + λ_dis·P_dis[t] + M·P_curt[t] ] + terminal_penalty
逐项拆一下。c_buy[t]·P_buy[t]是购电成本,分时电价下峰时高、谷时低;−c_sell[t]·P_sell[t]是售电收益,注意通常卖电价要低于买电价,否则会有套利空间,模型会把储能当提款机反复充放。f(P_dg[t])是柴油机发电成本,真实机组是二次凸函数,示例代码里我用线性近似,实际项目可以分段线性化或保留二次项。λ_ch和λ_dis是储能磨损成本,这一项很有用,后面讲同充同放时会提到。M·P_curt[t]是弃光惩罚,给得比购电成本高,模型才会优先消纳光伏,但也不要高到离谱,否则在光伏严重过剩、电价又低的时段,模型会强行把光伏卖给电网反而亏钱。
最后那个terminal_penalty是终端软约束。如果只优化有限时域而不约束终点,模型很可能在窗口最后把储能放空,因为最后一步之后的状态它不关心。加上“终点能量尽量回到参考值”的惩罚,滚动起来才有可持续性。示例里我用线性惩罚:e_loss ≥ E_ref − E[N],e_loss ≥ E[N] − E_ref,目标里加0.05·e_loss。
2.3 约束条件:功率平衡和储能动态是骨架
约束条件的写法直接决定求解器能不能跑通。最少需要四类:
第一条是功率平衡等式。任意时刻,光伏实际可用功率加柴油机出力、储能放电、购电,等于负荷加充电加售电:
P_pv[t] − P_curt[t] + P_dg[t] + P_dis[t] + P_buy[t] = P_load[t] + P_ch[t] + P_sell[t]
这里的P_pv[t]和P_load[t]是预测值,在优化问题里作为参数传入,不是决策变量。
第二条是储能能量动态方程。相邻两个控制周期的能量关系:
E[t+1] = E[t] + Δt·η_ch·P_ch[t] − Δt·P_dis[t] / η_dis
Δt是控制周期,比如15分钟就是0.25小时。充放电效率分开写,η_ch一般取0.95左右,η_dis取0.9左右。注意这个方程是等式的核心,很多同学在更新SOC时漏掉效率,导致仿真里储能越用越“多电”。
第三条是运行边界。储能能量上下限、充放电功率上限、柴油机出力上限、购售电功率上限,都是直接写死的线性约束。
第四条是储能不能同时充电和放电。数学上要表达成P_ch[t]·P_dis[t] = 0,这是一个非线性约束。处理有讲究:简单做法是在目标函数里给充放电都加一个小成本系数,求解器发现同充同放只会增加成本,自然就避开;严格做法是引入二元变量,把问题变成混合整数线性规划。两种方案后面都会给代码。
2.4 滚动窗口下的“同一套模型”
上面写的是任意时刻t的约束。在MPC里,k时刻求解的是从k到k+N−1的窗口,所有下标t都换成从k开始取。比如储能初值E[0] = E_real(k),E_real(k)是k时刻从现场或者仿真器里读回来的实测储能能量,这就是反馈信息的入口。
预测数据则从预测序列里取k到k+N−1的切片。求解出的P_ch[0]就是当前要执行的第一步,其他P_ch[1]到P_ch[N−1]只是未来几步的参考。等一个控制周期结束,用实测值更新E_real(k+1),重新预测、重新求解。模型本身从头到尾没变,变的只是初值、预测窗口和求解时刻。
3. Python实现MPC调度闭环:环境、预测、循环三步走
3.1 环境准备与求解器选择
Python侧只需要四个库:cvxpy、numpy、pandas、matplotlib。安装命令:
pip install cvxpy numpy pandas matplotlibcvxpy是我在这个场景下最常用的建模库,好处是语法贴近数学表达式,不用像用Pyomo那样写一堆Set和Param抽象层。调度级问题规模很小,变量几百个,cvxpy默认带的开源求解器足够。如果后面的模型里加了整数变量,建议显式指定HiGHS求解器,它是开源MILP求解器里性能非常能打的一个。
3.2 构造带预测误差的仿真数据
要验证MPC效果,数据里必须同时有“预测值”和“真实值”。我习惯用一条真实曲线加噪声生成预测曲线,简单直接:
import numpy as np T = 96 # 一天96个15分钟 load_real = 800 + 200 * np.sin(np.linspace(0, 2 * np.pi, T)) + np.random.normal(0, 30, T) pv_real = np.maximum(0, 800 * np.sin(np.linspace(0, np.pi, T)) ** 2 + np.random.normal(0, 50, T)) np.random.seed(1) load_fc = load_real + np.random.normal(0, 80, T) pv_fc = pv_real + np.random.normal(0, 80, T)这里load_fc和pv_fc是MPC滚动优化过程中能看到的数据,load_real和pv_real是“现场”真实发生的值。预测误差幅度约10%,已经足够说明问题。
3.3 MPC主循环:构建一次问题,循环里只更新参数
求解循环长这样:
E_hist = [E0] P_ch_hist, P_dis_hist, P_buy_hist = [], [], [] for k in range(T - N): p_buy, p_sell, p_ch, p_dis, p_dg, p_curt, E_pred = \ solve_mpc(k, E_hist[-1], load_fc, pv_fc, price) # 只执行第一步 P_ch_hist.append(p_ch[0]) P_dis_hist.append(p_dis[0]) P_buy_hist.append(p_buy[0]) # 按储能动态更新状态(仿真中认为储能按指令执行) E_next = E_hist[-1] + dt * eta_ch * p_ch[0] - dt * p_dis[0] / eta_dis E_next = min(max(E_next, E_min), E_max) E_hist.append(E_next)关键点只有一个:每个k循环里重新求解一遍问题,但每次只取第一步。E_hist[-1]就是当前实测储能能量,它把反馈带进了下一次优化。如果控制周期是15分钟,一天下来T−N+1≈96次求解,每次都是LP或MILP规模,计算时间完全不是问题。
3.4 核心建模函数solve_mpc
solve_mpc是核心,我给出一个能跑通的骨架:
def solve_mpc(k, E_start, load_fc, pv_fc, price): N = 24 dt = 0.25 E_max, E_min = 2000.0, 400.0 E_ref = 1000.0 P_ch_max = P_dis_max = 500.0 eta_ch, eta_dis = 0.95, 0.9 P_grid_max = 800.0 P_buy = cp.Variable(N, nonneg=True) P_sell = cp.Variable(N, nonneg=True) P_ch = cp.Variable(N, nonneg=True) P_dis = cp.Variable(N, nonneg=True) P_dg = cp.Variable(N, nonneg=True) P_curt = cp.Variable(N, nonneg=True) E = cp.Variable(N + 1) e_loss = cp.Variable(N) objective = 0 for t in range(N): objective += price[k + t] * P_buy[t] - 0.8 * price[k + t] * P_sell[t] objective += 0.15 * P_dg[t] objective += 0.002 * (P_ch[t] + P_dis[t]) objective += 10.0 * P_curt[t] objective += 0.05 * e_loss constraints = [] for t in range(N): constraints += [pv_fc[k + t] - P_curt[t] + P_dg[t] + P_dis[t] + P_buy[t] == load_fc[k + t] + P_ch[t] + P_sell[t]] constraints += [E[t + 1] == E[t] + dt * eta_ch * P_ch[t] - dt * P_dis[t] / eta_dis] constraints += [E_min <= E[t], E[t] <= E_max] constraints += [P_ch[t] <= P_ch_max, P_dis[t] <= P_dis_max] constraints += [P_dg[t] <= 600.0] constraints += [P_buy[t] <= P_grid_max, P_sell[t] <= P_grid_max] constraints += [P_curt[t] <= pv_fc[k + t]] constraints += [E[0] == E_start] constraints += [e_loss >= E_ref - E[N], e_loss >= E[N] - E_ref] prob = cp.Problem(cp.Minimize(objective), constraints) prob.solve(solver=cp.HiGHS) return (P_buy.value, P_sell.value, P_ch.value, P_dis.value, P_dg.value, P_curt.value, E.value)说明一句:这个版本没有加二元变量,所以是线性规划,用HiGHS解非常快。我靠0.002的磨损成本压制同充同放,实测在绝大多数时段不会出现。如果你想彻底杜绝,可以把二元变量方案加上,第5章会给写法。
3.5 结果结算:别只看调度曲线,要看实际成本
仿真跑完,成本用“实际发生值”结算,不是用MPC计划值,否则对比就没意义。简单做法是把每个时刻的实际购电功率累乘电价,再减去实际售电收益:
total_cost = 0.0 for k in range(len(P_buy_hist)): total_cost += price[k] * P_buy_hist[k] total_cost -= 0.8 * price[k] * P_sell_hist[k]P_buy_hist和P_sell_hist是每个MPC周期实际执行的第一个值。更精细的仿真可以在功率平衡里加一个实时偏差项,把预测误差导致并网点实际功率偏移也计进来;调度级研究一般不做这层,知道边界在哪就行。
4. 仿真实验怎么设计:用数据说话才能证明MPC的价值
4.1 三套对比策略必须同台竞技
评价MPC效果,光看自己的调度曲线好看没用,必须有基线对比。我至少会跑三套:
规则策略最简单:白天光伏给储能充电或者卖电,晚间高峰放电,储能按固定时间表动作,不依赖预测。日前开环优化是基准:用完整预测数据一次性求解96步,全序列执行,过程中不更新。MPC就是前面写的滚动方案。
注意一个公平性问题:三套方案面对的是同一套“真实”负荷和光伏序列、同一个分时电价、同一个初始SOC。如果MPC多给了信息,对比结果就没有说服力。
4.2 关键参数怎么调
预测时域N是最值得调的参数。N取24,覆盖6小时,很多项目够用;N取48,覆盖12小时,光伏夜间时段不再盲目充电,但后半段预测数据可信度下降。批量测试时可以从N=4、8、16、24、48各跑一遍,看总成本曲线拐点。
储能磨损成本λ也是敏感参数。λ太小,储能频繁动作、寿命损耗大;λ太大,储能根本不舍得用,晚高峰全去高价买电。工程上我一般先设0.002左右,再看SOC曲线是不是在合理区间,微调。
终端惩罚权重的设置比较粗暴:目标里e_loss前乘0.05。想严格一些,可以改成不等式:E[N] ≥ 900,E[N] ≤ 1100,把它变成硬约束。代价是预测偏差大时容易无解,所以工程上还是软惩罚稳。
4.3 结果怎么读:一张表的对比就够
我在某组随机数据上跑出来的结果大概长这样(示意值):
| 策略 | 运行成本(元) | 弃光率 | SOC越限次数 |
|---|---|---|---|
| 固定规则 | 24580 | 12.6% | 3 |
| 日前开环优化 | 22100 | 6.5% | 2 |
| MPC(N=24) | 20580 | 3.1% | 0 |
MPC成本最低,SOC也没有越限,滚动优化带来的收益是看得见的。但这里有个很容易误导人的点:如果预测误差设成0,也就是MPC每一步看到的“预测”和“真实”完全一样,那么MPC和日前开环优化的结果会完全一致。MPC的价值恰恰体现在预测有偏差时:它每15分钟用实测状态纠一次偏,把偏差造成的成本损失压到最低。
还有一个值得留意的现象:预测数据的“形状”比“精度”更影响MPC效果。比如说真实光伏峰值在13点,预测峰值却画到了11点,哪怕幅值完全一样,MPC也会在11点提前安排少买电,13点实际峰值来临时只能靠储能补救。这类形状错误比单纯的幅值偏差更不容易发现,画图对比预测和实际曲线时重点看峰值位置。
5. 实操中踩过最多的五个坑:从不可行解到求解器选型
5.1 求解器报infeasible:先检查是不是把约束写死了
MPC项目里最常见的报错不是语法错误,而是“problem infeasible”。比如给定了储能能量范围400~2000 kWh,又要求终端能量不低于1600 kWh,中午预测光伏又差,窗口内根本没有可行轨迹。解决思路是把硬约束放宽成软约束:SOC上下限留2%保护带,终端约束用前面写的e_loss软惩罚,弃光功率允许一定程度超过上限。调度优化是经济问题,不是二值逻辑问题,宁可让求解器多花一点惩罚,也不要报无解。
5.2 储能同充同放:数学上隐蔽,实际成本不小
如果不加任何处理,优化器有可能在同一个时段同时让P_ch和P_dis都等于一个正数。从平衡等式看,两边都被抬高,等式照样成立,但能量动态方程里充放电效率导致白白损耗,目标函数里如果只算购电成本根本看不出来。
处理方式有两种。第一种是目标函数加磨损成本,就是我主代码里的做法,简单有效,但不严格保证。第二种是二元变量,严格保证同一时段最多只有一个动作:
z_ch = cp.Variable(N, boolean=True) z_dis = cp.Variable(N, boolean=True) for t in range(N): constraints += [P_ch[t] <= P_ch_max * z_ch[t]] constraints += [P_dis[t] <= P_dis_max * z_dis[t]] constraints += [z_ch[t] + z_dis[t] <= 1]加了二元变量后问题是MILP,把求解器换成HiGHS即可。预测时域24步、变量几百个的规模,HiGHS解起来也是几十毫秒级别,不用怕。
5.3 测试时把预测误差设成0,MPC的优势完全消失
这是写论文和汇报时最容易翻车的地方。如果测试工况给的是确定性光伏、确定性负荷,MPC滚动反馈等于没用,刚刚好跟开环优化结果一样,体现不出创新点。做仿真对比时,一定要给预测加噪声,而且噪声要落在真实数据上,因为MPC的反馈价值在偏差里才能兑现。噪声幅度建议10%~20%,太小看不出区别,太大会让结果随机性过大,每次跑的成本都不一样。
5.4 求解器选型:连续问题用默认,整数问题指定HiGHS
cvxpy默认的求解器虽然是自动选择,但实际项目的坑是:不同问题类型它选出来的求解器可能完全不同,同样的代码换台机器可能卡住。我的建议是显式指定。纯LP/QP用CLARABEL或OSQP;加了二元变量就用HiGHS。如果以后要把MPC推到1分钟甚至更短的控制周期,建议换casadi加OSQP,并做参数热启动,cvxpy的建模开销在这种频率下会显得偏重。
5.5 调度层和执行层的边界要分清
MPC给出的充放电指令是经济调度层面的设定值,它不是实时功率控制器。15分钟窗口内光伏突然掉电,微电网里要靠储能逆变器的下垂控制或者并网点功率控制来瞬时兜底,MPC下一个周期才介入修正。如果拿MPC去背“秒级有功不平衡”这个锅,那方向就错了。做项目时最好在方案里明确:MPC负责分钟级经济分配,实时平衡由底层控制负责,两者通过状态上传接起来。
回到文章开头那个问题:为什么开环日前计划一到现场就废?因为它没有反馈。MPC只是把优化问题从“一次性执行”改成了“滚动+反馈”,但这个改动在不确定环境里价值巨大。如果你在做的微电网项目里发现优化模型越复杂效果越差,先别急着加算法,回头检查一下你的反馈链路通不通、状态量有没有更新——这比我见过绝大多数所谓的“模型升级”都管用。