我最初接触基于模型预测控制(MPC)的微电网调度优化是帮朋友复现一篇论文。项目标题里写着“Python代码实现”我当时以为主要工作是调库真正动手才发现MPC本身并不复杂难的是把微电网里的设备模型、约束条件、滚动时域逻辑串成一条能闭环跑通的链路。这篇文章就围绕这条链路展开把我从建模到Python实现再到调参排错的经验全部摊开讲。想复现课题的学生或者在做微电网能源管理系统EMS的工程师都可以直接拿这套思路当脚手架。MPC在微电网调度里的定位其实很直白每个调度时刻读取当前实测状态拿到未来一段时间的预测序列在线求解一个有限时域优化问题然后把序列里的第一步动作下发执行等下一个调度时刻再重复“预测—优化—执行—滚动”。这个闭环正是它比传统规则调度和一次性开环优化更稳的根本原因。下面我把建模细节、代码实现和踩坑记录按顺序展开希望能帮大家把这条链路一次走通。1. 为什么是MPC微电网调度绕不开的三个坑1.1 传统调度方法在微电网面前的三道坎先说传统调度为什么会在微电网这里吃瘪。最直观的要数“一刀切”的开环优化按预测的光伏和负荷曲线一次性求解未来一整天的调度计划并全部执行。如果预测偏差很小这套方案看起来非常完美但光伏出力受云层影响几分钟之内就可能掉一半等运行到那一刻计划里的功率平衡就兜不住了储能SOC可能越界柴油机也可能被迫爬坡。问题不一定出在优化算法本身而是开环结构没有“回头路”后面的误差会一步步累积。第二种常见做法是规则调度比如储能SOC低于30%就充电、高于80%就放电柴油机作为备用电源在负荷高时启动。这类规则实现简单、调试直观很多现场系统还在这么用但规则写得简单就难以同时兼顾经济性优化和多种安全约束。光伏大发时要不要充电、电价高时要不要让储能放电参与峰谷套利这些逻辑一条条堆进规则表最终会变成一个庞大的if-else王国。微电网设备一多维护成本会高到让人想重写。再往上一步有人会自然想到PID负反馈。但微电网本质上是带约束的多输入多输出系统PID设计时依赖固定工作点无法天然处理储能SOC上下限、柴油机爬坡限制这类不等式约束强行串级整定参数会非常痛苦。所以当问题同时具备随机扰动大、约束条件多、需要经济优化这三个特征时MPC几乎是最顺手的工具。1.2 滚动优化把“一次性计划”改成“步步修正”你可以把传统开环调度想象成闭着眼睛开车出发前按地图规划好全程路线中间无论路况怎么变你都照着原计划转向。MPC则像正常开车时不断看路每走一小段就重新评估前方路况只把最近一段的路线定下来等开到下一个路口再重新规划。这个比喻里预测时域就是你“看多远”的距离控制时域就是你“先走多少”这一段。MPC每个调度周期做的事情可以拆成三步读取当前实测状态比如储能当前SOC拿到未来N步的光伏、负荷预测序列然后在线求解一个有限时域最优化问题把未来每步的设备出力算出来最后只执行当前这一步动作。等到下一个15分钟再重新读取状态、更新预测、重新求解。因为每次决策都用最新状态做了修正预测误差带来的连锁反应会被一次一次打断。在微电网调度里MPC的核心优势有两点。第一目标函数可以直接写成经济运行成本各类运行约束也能原样写进优化模型不需要像传统控制那样绕来绕去第二滚动时域带来天然的反馈校正能力光伏预测错了没关系下一轮重新算一遍就行。这也是为什么这个研究方向最常见的实验设计就是拿MPC和开环优化、规则调度各跑一遍仿真对比总成本和约束违规次数。从工程角度看MPC的计算量完全够实时一个96时段的优化问题在普通电脑上毫秒级就能解完。2. 微电网建模把蓄电池和柴油机写成MPC能懂的数学表达式2.1 先搭一个典型微电网微电网不是把一堆设备简单接在一起建模需要明确每类设备的出力范围、状态变化和相互之间的功率平衡。本文用一个最经典的光储柴并网型微电网作例子光伏作为随机电源储能负责平移电量柴油机作为可控备电配网联络线允许与外部电网双向买卖电。这个结构覆盖了大多数研究场景后续扩展风电、电动汽车充电桩也比较方便。系统内部的功率关系可以用一条等式概括任意时刻光伏出力、柴油机出力、储能放电功率、电网交互功率之和等于当前负荷功率。如果允许甩负荷就在等式右侧留一个松弛变量。这条等式是微电网调度模型的骨架所有设备模型和约束都挂在它上面。2.2 设备模型与关键变量每个设备都需要用变量和约束描述清楚。下面是我常用的最小模型集合设备变量模型与约束光伏p_pv0 ≤ p_pv ≤ p_pv_pred可弃光储能p_ess、soc-P_ess_max ≤ p_ess ≤ P_ess_maxSOC_MIN ≤ soc ≤ SOC_MAXsoc[k1] soc[k] - p_ess[k] * Δt / E_ess柴油机p_dg0 ≤ p_dg ≤ P_dg_max爬坡约束 |p_dg[k1] - p_dg[k]| ≤ R_dg联络线p_gridp_grid 0 表示购电p_grid 0 表示售电负荷p_load预测序列p_unmet 表示失负荷量这里最需要注意的是储能模型的符号约定我习惯把 p_ess 定义为“放电为正”所以SOC差分方程里是减号。如果符号定义反了整个闭环跑起来SOC会朝反方向飞代码越调越乱。很多初学者第一版代码跑不通最后查出来都是这个小地方。光伏模型用的是“预测序列作为上限”的思路意思是不要求光伏必须全额消纳优化器可以主动弃光。这在光伏渗透率较高的微电网里很有用避免为了消纳光伏而让储能过充或让柴油机停机带来更大的代价。负荷侧我保留了失负荷变量正常情况下通过高额惩罚把它压在零附近只在极端场景下才能放开这比硬性要求功率平衡更容易保证模型有解。2.3 目标函数把“省钱”变成数学表达式确定了设备模型和平衡约束接下来就是让优化器知道“什么方案更好”。目标函数就是经济成本最小化我通常写成minimize Σ (C_dg * p_dg C_buy * max(p_grid, 0) - C_sell * max(-p_grid, 0)) 惩罚项其中 C_dg 是柴油机度电成本C_buy 是购电价C_sell 是售电价。既然允许跟电网交互卖电收入就作为负成本放进目标购电价通常高于售电价所以优化器不会傻到同时买和卖。惩罚项一般有两类。一类是失负荷惩罚系数要设得比任何正常调度成本高比如50元/kWh让优化器只有在极极端情况下才会切负荷另一类是储能磨损惩罚我常用 p_ess 的平方项比如 C_ess * Σ p_ess²。平方项的作用是抑制储能功率频繁反转和大幅波动防止MPC为了几毛钱价差让电池剧烈充放。不加这个项目标优化会表现得过于“贪”调度曲线锯齿状跳变实际系统根本不敢这么用。2.4 一套可直接用的仿真参数给出一套我常用的参数新手可以直接照抄先跑通再根据场景调整参数数值含义Δt0.25 h调度步长15分钟N24预测时域共6小时光伏最大出力80 kW预测曲线峰值负荷均值60 kW负荷曲线基值储能容量 E_ess200 kWh可用容量储能功率 P_ess_max±80 kW充放电最大功率SOC范围0.2 ~ 0.9上下限柴油机最大出力120 kW单机容量柴油机爬坡限制30 kW/15min每步爬坡能力选15分钟步长是一个比较常规的折衷时间粒度太小计算量和数据量都会成倍增加粒度太大又难以反映光伏的分钟级波动。另外这组参数在物理上是自洽的储能容量200kWh功率限值80kW意味着一个调度步长最多能充放20kWh从SOC下限到上限需要多个步长不会出现一步充饱这种不合理现象。开始建模前先手动检查这类数值一致性能省掉后面排查不可行解的大量时间。3. Python实战把MPC闭环完整跑起来3.1 工具链选型cvxpy就够了Python环境建议用3.10以上的干净虚拟环境避免系统Python里各种包的版本互相打架。依赖其实很少主要就是numpy、cvxpy和matplotlibpip install numpy cvxpy matplotlibcvxpy是声明式优化建模库写出来的代码和数学表达式几乎一一对应不需要手写KKT条件也不需要自己实现梯度下降。相比之下scipy.optimize.minimize也能做小规模问题但处理成千上万个变量和约束时的代码量和调试成本会高很多。求解器方面cvxpy默认带OSQP、ECOS等我们这里的问题本质是带线性约束的二次规划直接用OSQP就足够。如果以后加入柴油机启停这类整数变量需要额外安装HiGHS但第一版先别急着上整数连续MPC跑通后再扩展。3.2 参数定义与预测序列生成先定义系统参数并模拟一组光伏和负荷预测数据。这里的预测序列就当作外部预测模块的输出实际工程中它可以来自数值天气预报加历史负荷模型。import numpy as np import cvxpy as cp import matplotlib.pyplot as plt DT 0.25 # 调度步长 15min单位小时 N 24 # 预测时域 TOTAL 96 # 预测数据总长度 SIM_STEPS 72 # 实际滚动调度的步数 PV_MAX 80.0 # 光伏最大出力 kW ESS_E 200.0 # 储能容量 kWh ESS_P_MAX 80.0 # 储能最大功率 kW SOC_MIN, SOC_MAX 0.2, 0.9 DIESEL_MAX 120.0 # 柴油机最大出力 kW RAMP_LIMIT 30.0 # 柴油机每15min爬坡限制 kW C_DG 0.8 # 柴油机成本 元/kWh C_BUY 1.2 # 购电价 元/kWh C_SELL 0.5 # 售电价 元/kWh C_LOSS 50.0 # 失负荷惩罚 元/kWh C_ESS 0.02 # 储能磨损惩罚系数 np.random.seed(42) t_hour np.arange(TOTAL) * DT # 简单模拟白天6点到18点有光照夜间为0 pv_base PV_MAX * np.clip(np.sin(np.pi * (t_hour - 6) / 12), 0, None) pv_forecast pv_base * (0.85 0.15 * np.random.rand(TOTAL)) # 负荷预测白天偏高叠加小幅随机波动 load_forecast 60 20 * np.sin(np.pi * (t_hour - 8) / 12) 5 * np.random.randn(TOTAL)PV预测和负荷预测我都加了随机噪声用来模拟真实预测不可能完全正确的情况。没有这个噪声闭环MPC和开环优化会得到几乎一样的结果体现不出滚动优化的优势。3.3 核心MPC求解函数下面这步是整个项目的核心把第2章的数学表达式翻译成cvxpy代码。这个函数返回未来N步的最优功率序列但我们最终只会取第一步动作。def solve_mpc(soc0, pv_pred, load_pred): # 决策变量 p_pv cp.Variable(N) # 光伏实际出力 p_ess cp.Variable(N) # 储能功率正为放电 p_dg cp.Variable(N) # 柴油机功率 p_grid cp.Variable(N) # 电网交互正为购电 p_loss cp.Variable(N) # 失负荷量 soc cp.Variable(N 1) # SOC轨迹 constraints [soc[0] soc0] for k in range(N): constraints [ 0 p_pv[k] pv_pred[k], -ESS_P_MAX p_ess[k] ESS_P_MAX, 0 p_dg[k] DIESEL_MAX, SOC_MIN soc[k 1] SOC_MAX, p_loss[k] 0, # 功率平衡光伏柴油储能电网 负荷 - 失负荷 p_pv[k] p_dg[k] p_ess[k] p_grid[k] load_pred[k] - p_loss[k], # 储能SOC差分方程 soc[k 1] soc[k] - p_ess[k] * DT / ESS_E, ] # 柴油机爬坡约束 for k in range(N - 1): constraints [cp.abs(p_dg[k 1] - p_dg[k]) RAMP_LIMIT] objective cp.sum(C_DG * p_dg) \ cp.sum(C_BUY * cp.pos(p_grid)) \ - cp.sum(C_SELL * cp.neg(p_grid)) \ cp.sum(C_LOSS * p_loss) \ cp.sum(C_ESS * cp.square(p_ess)) prob cp.Problem(cp.Minimize(objective), constraints) try: prob.solve(solvercp.OSQP, verboseFalse, eps_abs1e-5, eps_rel1e-5, max_iter5000) except Exception: return None if prob.status ! optimal: return None return p_pv.value, p_ess.value, p_dg.value, p_grid.value, p_loss.value, soc.value有几个细节值得展开说一下。一是cp.pos和cp.neg它们分别提取正部和负部天然适合处理买卖电价格不对等的情况自己用np.maximum去截断在numpy里可以但在cvxpy的建模阶段会破坏DCP规则。二是储能磨损惩罚里的cp.square这个平方项让目标函数变成二次规划仍然可以交给OSQP高效求解。三是我对决策变量直接用了0 p_pv[k] ...这种区间写法cvxpy会自动扩展成两条不等式约束代码可读性高很多。函数返回前检查了prob.status不是optimal就返回None。这个判断非常重要因为在滚动时域里某个时刻的预测数据特别离谱时硬约束的模型确实可能无解调用方有了None这个信号才能外面做容错处理不至于让整个仿真崩溃。3.4 滚动主循环只执行第一步有了单步求解函数主循环的逻辑就简单了每一步用当前真实SOC作为初始状态拿从当前时刻开始的未来N步预测调用solve_mpc得到最优调度序列然后只取序列里的第一个值执行再按真实模型更新SOC。soc_now 0.5 soc_history [soc_now] pv_log, ess_log, dg_log, grid_log, loss_log [], [], [], [], [] for k in range(SIM_STEPS): pv_pred pv_forecast[k:k N] load_pred load_forecast[k:k N] res solve_mpc(soc_now, pv_pred, load_pred) if res is None: print(fStep {k}: MPC infeasible, break.) break p_pv, p_ess, p_dg, p_grid, p_loss, _ res # 只执行第一步 pv_log.append(p_pv[0]) ess_log.append(p_ess[0]) dg_log.append(p_dg[0]) grid_log.append(p_grid[0]) loss_log.append(p_loss[0]) # 状态更新用真实SOC方程递推 soc_now soc_now - p_ess[0] * DT / ESS_E soc_history.append(soc_now)这里我先把模型自身当成了“真实系统”假设预测误差为零所以状态更新跟优化器里的SOC方程完全一致。实际工程中真实系统会受到各种扰动这时状态更新要改用现场采集的SOC测量值这也正是滚动时域的意义所在每一步都拿实测SOC重新优化而不会傻傻沿用开环计划里的SOC轨迹。循环结束后可以打印累计运行成本total_cost (C_DG * np.sum(dg_log) C_BUY * np.sum(np.clip(grid_log, 0, None)) - C_SELL * np.sum(np.clip(-grid_log, 0, None))) * DT print(fTotal operation cost: {total_cost:.2f} CNY)注意乘上DT因为所有功率变量的单位是kW一个步长是0.25小时功率乘时间才是电量再乘价格才是钱。这一步单位换算是新手最容易漏的地方漏掉后成本计算会大四倍。3.5 结果可视化与解读画图建议直接看四类信息各电源出力曲线、储能SOC轨迹、负荷平衡情况和电网交互功率。一个简单版本是这样fig, axes plt.subplots(3, 1, figsize(11, 8), sharexTrue) k np.arange(len(pv_log)) axes[0].plot(k, pv_log, labelPV) axes[0].plot(k, dg_log, labelDiesel) axes[0].plot(k, np.maximum(ess_log, 0), labelESS discharge) axes[0].plot(k, np.maximum(-ess_log, 0), labelESS charge) axes[0].plot(k, grid_log, labelGrid) axes[0].set_ylabel(Power (kW)) axes[0].legend(locbest) axes[1].plot(k, soc_history[:-1], labelSOC) axes[1].axhline(SOC_MIN, linestyle--, colorgray) axes[1].axhline(SOC_MAX, linestyle--, colorgray) axes[1].set_ylabel(SOC) axes[2].plot(k, np.array(pv_log) np.array(dg_log) np.array(ess_log) np.array(grid_log) np.array(loss_log), labelGenerated - Load) axes[2].set_ylabel(Balance (kW)) axes[2].set_xlabel(Time step (15min)) plt.tight_layout() plt.show()跑出来的结果会出现几个典型现象。白天光伏出力上来时MPC会让储能充电把多余电量存起来傍晚光伏归零后储能开始放电避峰柴油机只在负荷高峰和SOC较低的时候启动如果当天购电价高于柴油成本优化器会优先用柴油机而不是电网买电。SOC曲线会一直在0.2~0.9之间游走不会越界这就是约束起作用的表现。这些行为不需要很复杂的分析一眼就能看出MPC是在“提前充、高峰放”而不是被动响应。4. 调参避坑求解失败、SOC越界和预测误差4.1 求解器报“infeasible”从哪里查滚动仿真跑到一半突然break大概率就是某个时刻的优化问题无解。infeasible的原因通常就三类单位不一致、约束互相矛盾、初始值不在可行域内。单位不一致是最隐蔽的。比如储能容量用了kWh但功率用了kWΔt却忘了乘以小时数那么SOC差分方程和功率平衡约束就会差一个数量级模型在数学上根本不可能满足。我建议在建模前把所有单位统一写成kW、kWh、h三个基本单位的组合然后手动检查一个步长内储能能量变化是否合理。约束互相矛盾也很好理解。比如SOC被限制在0.2~0.9但初始SOC只有0.2而柴油机爬坡又慢光伏夜间为零负荷又高优化器可能找不到一组设备出力让功率平衡。遇到这种情况先把个别硬约束换成软约束比如给SOC边界加上松弛变量或者允许失负荷变量非零且付出惩罚至少让每个时刻都有可行解。之后再根据求解结果逐步收紧。还有一个小技巧调试时把prob.status打印出来再用类似constraints[5].violation()的手段看哪些约束在求解结束时仍然不满足。cvxpy在求解器返回non-optimal时也提供了残差信息这个比人肉盯代码快得多。4.2 cvxpy报“DCP”错误怎么破cvxpy会强制校验问题是否满足DCP规则目的是保证问题可解且解是全局最优。报这个错误时第一反应是检查目标函数和约束里是否出现了“变量乘变量”“max(a,b)里放变量”这类写法。前面我特意用cp.pos和cp.neg处理买卖电而不是写np.maximum(p_grid, 0)就是为了避免在建模阶段混入numpy的不透明操作。如果只是简单的分段成本cp.maximum本身也可以用它是凸函数符合DCP规则但如果目标和约束里混了cp.maximum(变量, 0)后再取负号整个表达式可能就不再凸。遇到这类报错最稳妥的办法是把目标拆成几个子表达式一个个加进去看到底是哪一项触发报错。另一个容易遇到的坑是版本兼容性。numpy版本过新和cvxpy版本不匹配有时候导入cvxpy就直接报错或者求解时出现诡异的C回溯。我踩过几次后习惯在虚拟环境里固定安装numpy 1.26.x和cvxpy 1.4.x工程稳定性优先于版本追新。4.3 预测误差不是MPC的万能药MPC能抗预测误差不代表它不怕误差。我做过一组简单实验把pv_forecast叠加不同方差的高斯噪声分别运行闭环MPC和一次性开环优化。结果很典型噪声小时两者成本接近噪声大时开环优化开始出现SOC越界甚至失负荷而MPC因为有滚动修正成本上升幅度明显更小违规次数也更少。这说明滚动时域面对不确定性确实有优势。但噪声大到一定程度MPC也会无能为力因为预测窗口内整段数据都严重偏了再修正也只是“信息滞后”。这时候就需要更高级的处理了比如第5章会提到的鲁棒MPC和随机MPC。如果你做的课题是仿真对比用不同噪声水平做敏感性分析是非常加分的内容也能帮你真正理解MPC的边界在哪。4.4 实时性优化与热启动连续MPC在普通PC上一般毫秒级就能解完但如果把预测时域拉长到96步实时求解次数成倍增加调度周期又短还是会紧张。能做的优化有几个一是减小预测时域比如从24步减到12步代价是对长期约束的预判变弱二是给求解器设置更宽松的精度OSQP的eps_abs调到1e-4就够工程使用了三是复用上一时刻的解作为热启动。cvxpy的热启动写法很简单在solve里传入warm_startTrue并把上一轮得到的变量value赋值给新的变量。这一招在连续滚动仿真里尤其有用因为相邻两个调度周期的预测数据只滑动了一格最优解变化不大热启动能显著减少迭代次数。4.5 常见问题速查表现象可能原因处理方案求解器返回infeasible单位不一致、约束过紧、初值不在边界内统一单位降低ΔtSOC边界加软约束DCP规则报错目标里有变量乘变量或numpy函数改用cp.pos、cp.neg、cp.squareSOC曲线反向飞p_ess正负号定义反了检查储能差分方程放电为正时SOC递推用减号调度曲线锯齿严重缺少储能磨损惩罚目标里加平方项惩罚系数从小到大试仿真结果MPC与开环几乎一样预测无误差或扰动太小给预测序列叠加噪声再对比成本与违规次数求解太慢N太大或有整数变量减小N连续松弛设置求解器eps和max_iter5. 从论文到工程MPC调度还能往哪走5.1 鲁棒MPC把预测误差放进优化模型如果预测误差不可忽略并且你希望方案在“最坏情况”下也不越界可以走鲁棒MPC路线。核心思路是把光伏和负荷预测建模成不确定集比如“预测值上下浮动15%”然后在优化约束中让所有可能的不确定场景都满足安全边界。代价是模型规模增大计算量上升但换来的是更强的保守性和安全性。5.2 随机MPC用场景树描述不确定性随机MPC是另一种处理不确定性的方式。它用一批采样场景代表光伏和负荷可能的未来轨迹目标函数改为期望成本。相比鲁棒MPC随机MPC不会过于保守代价是要生成足够多的场景问题规模成倍膨胀。在微电网调度研究里这种方法和蒙特卡洛模拟经常结合在一起非常适合做对比实验。5.3 分布式MPC多微电网协同当研究对象从单微电网变成微电网群集中式优化在通信和计算上的压力都会变大。分布式MPC按区域把问题拆开每个子微电网自己跑一个MPC再通过协调层交换边界功率信息迭代若干轮后得到全局近似最优解。常见迭代手段是ADMM工程实现上要处理通信延迟和异步问题但可扩展性比集中式强很多。5.4 与强化学习结合MPC和强化学习的结合是最近较火的方向思路也很自然用强化学习学习复杂的不确定性预测或策略再让MPC充当“安全护栏”保证输出不会违反调度约束。对这种混合方案Python生态很友好cvxpy可以当作安全优化层RL框架负责策略学习两者通过接口对接。如果你已经有了一套能跑的MPC代码往RL方向扩展其实比从零开始要省力很多。我个人在实际操作中的体会是MPC这个方向真正难的不是最优控制理论而是把状态更新、符号约定、单位换算这些“细枝末节”一次弄对。建议先跑通一个最简单的小模型确认SOC轨迹和调度行为在常识上是合理的再逐步加复杂度。如果发现SOC在跳变或求解频繁失败不要急着怀疑求解器先回去检查p_ess的正负号和Δt的单位。把这些细节理顺之后剩下的就是不断调权重系数看成本曲线和约束违规情况之间的此消彼长。这篇代码和调试经验就是顺着这条思路沉淀下来的希望对正在复现这个课题的人有所帮助。