简介面向电力系统韧性研究者和配电网调度算法开发人员的 MATLAB 实现包复现顶刊文献中应急移动电源MPS事件后动态调度策略。该策略基于混合整数规划建立优化模型解决道路网络与电力网络耦合、MPS 派遣与配电网运行多时间尺度协调等难点并在 IEEE 33/123 节点系统上验证了有效性。压缩包共5个文件含3个 .m 脚本、1个版权说明文件和1个 vsdx 拓扑示意图整体仅49KB轻量便于运行。目前已有322人学习下载。可直接运行 main_mps_Scheduling.m 查看 MPS 动态调度结果配合 show_result.m 绘制恢复效果图通过 IEEE33.vsdx 可快速理解测试系统拓扑及接入位置适合正在研究韧性配电网、应急电源调度或 YALMIP 建模的读者作为复现与扩展参考。1. 配电网韧性场景里MPS 动态调度为什么比“预先放好”更棘手台风刚过境配电网里几条馈线同时断开抢修队还没进场最有效的临时供电手段往往是应急移动电源。MPS 动态调度要回答的是“车从哪个待命点出发、什么时间发车、走哪条可达道路、接入哪个节点、给哪些负荷充多久”这一整组时空耦合问题。相比预配置阶段的站点选址动态调度在灾情逐步揭晓的过程中做多阶段决策既要考虑道路是否中断又要盯住剩余电量这个状态量的消耗是配电网韧性计算里公认最难复现的部分。这篇笔记适合正在复现两阶段配电网韧性论文、或设计应急电源调度方案的人按模型、求解、滚动时域、踩坑四条线把它讲透。2. 把动态调度写成优化问题MPS 状态机、目标函数与 MISOCP 选型先别急着写代码。动态调度和预配置最大的差别在于时间维度被真正拉出来了。预配置只回答“灾害前把 MPS 放在哪”是一组空间决策动态调度要把“发车时刻、到达时刻、接入时段、转场路径”全放进模型里稍不留神就会得到一个变量成倍膨胀、求解器直接卡死的混合整数非凸问题。这个章节先把模型骨架立住。2.1 预配置环节传下来什么动态调度的输入数据清单复现这类两阶段框架时常见做法是让预配置模型输出四类东西动态调度模型直接继承。第一类是 MPS 初始待命站点集合也就是每台移动电源在 t0 时刻所在的节点第二类是每台 MPS 的额定功率和初始能量它决定了后续供电能力上限第三类是候选接入点集合通常是灾后仍具备接入条件的配电变压器节点第四类是道路可达性矩阵这个矩阵在动态调度阶段会反复更新。预配置环节输出动态调度环节的用途初始待命站点决定 t0 时刻 MPS 的初始位置状态额定功率 / 初始 SOC作为接入后的放电功率上限和总能量边界候选接入节点集合削减二进制变量的枚举范围避免求解器白算道路可达性矩阵发车路径选择直接读取不可达路段从候选中删除我在复现时经常看到有人把预配置和动态调度耦合成一个大模型一次性求解。那在数学上没错但计算规模很容易失控而且失去了“动态”的本意。真正的动态调度应该是把预配置结果当作固定参数读入只对灾后恢复过程做滚动优化这也是标题里“下”这个系列的合理切分方式。2.2 MPS 只有三种状态状态机决定决策变量怎么设MPS 在任意时刻只有三种状态待命、在途、接入供电。待命是指停在某个安全节点等指令在途是指已经发车、正在沿道路移动接入供电是指已经到达目标节点并以额定功率向负荷放电。动态调度模型的所有决策变量本质上是围绕这三种状态的切换来设计的。待命状态不需要消耗能量但这台 MPS 也暂时不能为负荷出力在途状态要占用道路资源同时有行程时间时间窗口内到不了就接不上接入供电状态产生放电功率但它把 MPS 钉死在了某个节点上短期内不能再移动。因此建模时最常见的一套变量是二进制变量 x[m][k][t] 表示 MPS m 在时段 t 是否接入节点 k另一个二进制变量 y[m][k][t] 表示 MPS m 在时段 t 是否正在从 k 节点发车转移连续变量 p[m][t] 表示第 m 台 MPS 在时段 t 的放电功率连续变量 soc[m][t] 表示第 m 台 MPS 在时段 t 结束时的剩余能量。这套变量组合能充分表达“一台 MPS 上午在 A 节点供电、中午转场、下午给 B 节点供电”的动作序列。状态机对约束设计的直接影响是一台 MPS 同一时刻只能处于一种状态。用约束表达就是 x 和 y 不能同时为 1一台 MPS 同一时段最多接入一个节点。很多复现方案翻车不是输在目标函数上而是忘了给状态机补互斥约束。2.3 目标函数和约束怎么选为什么用 DistFlow 加二阶锥松弛目标函数优先选择加权恢复负荷的最大化形式是“每个时段被恢复负荷的状态指示变量乘以权重系数后求和”再乘上时段长度折算成电量。权重系数怎么定医院、养老院、通信基站这类关键负荷权重高普通居民负荷权重低具体数值可以根据灾后重要性排序来标定。也可以退一步用“失负荷量最小化”两者在数学上等价但加权恢复形式更容易解释韧性提升的效果。等式约束里最核心的是潮流方程。配电网跟输电网不一样电阻电抗比大直流潮流那种忽略无功和电压的近似误差太大不能直接用。业界和顶刊里最常用的是 DistFlow 模型它把支路有功、无功、电流平方和节点电压平方写成递推关系。直接求解 DistFlow 是非凸问题为了能用商用求解器保证收敛通常对支路电流项做二阶锥松弛也就是把原来等号连接的电流平方项放成不等式让模型变成混合整数二阶锥规划简写就是 MISOCP。选 MISOCP 而不是其它形式理由很实际第一Gurobi、SCIP 这类求解器对二阶锥约束支持成熟不需要自己写外逼近算法第二在辐射状配电网的标准算例上二阶锥松弛的误差通常可以控制在可接受范围内第三二进制变量和锥约束可以放在同一个模型框架里处理不需要做两层迭代。此外还要补充节点功率平衡、MPS 接入容量、能量守恒、移动行程时间约束这些约束会在下一章直接落到代码里。3. 复现最小可算模型Pyomo 建模代码与求解器参数设置模型结构定下来之后下一步就是把它翻译成求解器能吃的代码。我不建议一上来就堆完整算例那会让排查问题变得很痛苦。正确做法是先搭一个三到五个节点的最小网络跑通一个 MPS 的接入和放电再逐步放大。这个章节给出一套可以直接改的 Pyomo 骨架求解器用 Gurobi 或 SCIP 都行。3.1 先按 p.u. 整理数据再写核心约束代码复现这类模型最容易忽视的就是单位。我一般先把所有功率、电压、能量统一到标幺值 p.u. 体系基准功率取系统总容量负荷功率和 MPS 额定功率都除以基准值。这一步做完数值尺度问题就会少一大半。下面是最小算例的核心建模代码import pyomo.environ as pyo # 最小链式馈线: 节点0是上级电源, 节点1/2带负荷 NODES [0, 1, 2] BRANCHES [(0, 1), (1, 2)] # (from, to) 支路列表 LOAD {1: 0.6, 2: 0.4} # 有功负荷, 单位 p.u. TIME [0, 1, 2, 3] # 4 个时段, 每时段 15 分钟 DT 0.25 # 时段长度 0.25 h # MPS 参数: 两台移动电源, 额定放电功率 0.3 p.u., 初始能量 4.0 p.u.*h MPS_LIST [MPS_1, MPS_2] MPS_RATED_P 0.3 MPS_INIT_ENERGY 4.0 model pyo.ConcreteModel() model.T pyo.Set(initializeTIME) model.M pyo.Set(initializeMPS_LIST) # 决策变量: 接入指示、负荷恢复指示、放电功率、剩余能量 model.x pyo.Var(model.M, NODES, model.T, withinpyo.Binary) model.s pyo.Var(NODES, model.T, withinpyo.Binary) model.p pyo.Var(model.M, model.T, bounds(0, MPS_RATED_P)) # 约束1: 负荷恢复必须由至少一台 MPS 接入该节点支撑 def load_rule(mdl, k, t): return mdl.s[k, t] sum(mdl.x[m, k, t] for m in MPS_LIST) model.load_cons pyo.Constraint(NODES, model.T, ruleload_rule) # 约束2: 每台 MPS 同一时段只能接一个节点 def one_node_rule(mdl, m, t): return sum(mdl.x[m, k, t] for k in NODES) 1 model.one_node_cons pyo.Constraint(model.M, model.T, ruleone_node_rule) # 约束3: 放电功率只有在接入时才能为正 def discharge_rule(mdl, m, t): return mdl.p[m, t] MPS_RATED_P * sum(mdl.x[m, k, t] for k in NODES) model.discharge_cons pyo.Constraint(model.M, model.T, ruledischarge_rule)这段代码里有个关键逻辑负荷恢复约束用的是“小于等于”而不是“等于”目的是允许优化器选择不恢复某个非关键负荷把有限的 MPS 能量留给权重更高的负荷。放电功率约束则把 p[m][t] 和接入指示变量 x 绑定在一起MPS 没接入就不可能出力。如果去掉这条约束求解器会在 x 全为 0 的时候还给负荷送电得到一个物理上不可能的解。3.2 补全潮流方程与能量约束让模型真正可解上面的代码只是 MPS 侧约束还差两块潮流方程和能量守恒。潮流方程在 DistFlow 框架下需要定义支路功率变量 P[j][t]、支路电流平方变量 I[j][t] 和节点电压平方变量 V[i][t]然后对每条支路写电压递推关系再对非根节点写功率平衡方程。二阶锥松弛落在电流平方约束上形式是“支路电流平方大于等于支路有功平方加无功平方除以首端电压平方”。# 支路功率与电压变量 model.P pyo.Var(BRANCHES, model.T) model.Q pyo.Var(BRANCHES, model.T) model.I pyo.Var(BRANCHES, model.T, withinpyo.NonNegativeReals) model.V pyo.Var(NODES, model.T, bounds(0.8**2, 1.2**2)) # 二阶锥松弛: I (P^2 Q^2) / V def soc_rule(mdl, b, t): i, j b return mdl.I[b, t] * mdl.V[i, t] mdl.P[b, t]**2 mdl.Q[b, t]**2 model.soc_cons pyo.Constraint(BRANCHES, model.T, rulesoc_rule) # 能量守恒: 累计放电量不能超过初始能量, 每时段末更新剩余能量 def energy_rule(mdl, m): return sum(mdl.p[m, t] * DT for t in TIME) MPS_INIT_ENERGY model.energy_cons pyo.Constraint(model.M, ruleenergy_rule)能量约束用“累计放电量小于等于初始能量”是一种保守写法它假设所有放电都发生在调度时域内。如果 MPS 在灾前已经消耗过一部分电量就把 MPS_INIT_ENERGY 换成预配置阶段输出的剩余电量值。还有一个细节distflow 的电压递推关系我没在这里展开因为完整实现还要定义线路电阻电抗数据并逐支路写两行递推约束但核心的非线性难点已经在 soc_cons 里体现出来了。把这个约束交给求解器后模型就从原来的非凸问题变成了 MISOCP这是能稳定求解的关键一步。3.3 Gurobi 参数、MIPGap 与求解超时的判断标准模型建好之后求解器参数设置决定了你是十分钟拿到结果还是跑一夜不出解。我复现时常用的求解配置是TimeLimit 设置为 300 秒MIPGap 设置为 0.01Threads 设置为 4MIPFocus 设置为 1前两个参数是为了在可接受时间内拿到次优解后两个参数是为了让求解器在分支定界时更积极。opt pyo.SolverFactory(gurobi) opt.options[TimeLimit] 300 opt.options[MIPGap] 0.01 opt.options[Threads] 4 opt.options[MIPFocus] 1 res opt.solve(model, teeTrue)求解结束之后绝对不能只看求解器状态就收工。Gurobi 返回 Optimal 不代表解可靠要同时检查两个指标最终 MIPGap 是否真的小于 0.01以及目标函数值和最优界 ObjBound 之间的相对偏差。如果 TimeLimit 到了还没收敛求解器会返回 feasible 解经常被误当成最优解用这个问题我会在避坑章节专门展开。4. 动态体现在哪里滚动时域调度主循环与状态传递到这里模型已经能求解“给定完整灾情信息下MPS 应该怎么动”。但真实场景里灾情信息不是一次性到位的而是随着巡检、用户报修和故障定位逐步确认。所以复现顶刊里的“动态调度”核心不是把单次优化模型写好而是要把时间轴滚动起来。4.1 为什么不用“全知求解”而用滚动时域如果假设灾情信息全局已知把所有时段的所有线路故障一次性告诉求解器得到的自然是全局最优方案。这在算例演示里很好看但工程上不成立。真实场景是台风过境后 2 小时调度员只掌握部分馈线跳闸信息再过 4 小时又有新的线路被确认断线。如果一开始就按“未来不会新增故障”来调度前面把 MPS 派到了某个节点后边真正重要的负荷反而接不上。滚动时域调度的核心思想是每个时段只做未来一个有限窗口的优化执行窗口内的第一步决策等新的灾情信息到达后再重新优化。窗口长度可以取 4 到 8 个时段每个时段 15 分钟也就是覆盖未来 1 到 2 小时的调度计划。这样做的代价是放弃了全局最优换取了决策对最新信息的响应能力。4.2 滚动主循环代码新故障进入、固定已执行决策、只下发第一步滚动主循环代码的逻辑比模型代码更值得仔细看。每一步接收新的故障状态、更新候选接入节点、固定上一轮已经发车的移动决策、求解、只取第一步执行。H 4 # 滚动窗口长度 horizon list(range(H)) # MPS 运行状态字典, 在循环外部维护 mps_state { pos: {MPS_1: 0, MPS_2: 0}, # 当前所在节点 soc: {MPS_1: 4.0, MPS_2: 4.0}, # 剩余能量 tasks: [] # 已下发且正在执行的任务 } for t in range(TOTAL_STEPS): # 1) 读取本时段新确认的断线信息和新增故障节点 new_faults read_fault_status(t) blocked_nodes update_blocked_nodes(new_faults) # 2) 构建当前窗口的模型, 并把已执行任务固定到模型里 mdl build_rolling_model(mps_state, blocked_nodes, horizon) # 3) 求解当前窗口, 只取第一步动作 solve_and_log(mdl, labelfstep_{t}) action extract_first_step(mdl) # 4) 执行第一步, 推进 MPS 状态: 移动的车按路程前进, 接入的车放电 mps_state advance_state(mps_state, action, dtDT) # 5) 窗口向前滑动 horizon [h 1 for h in horizon]构建滚动窗口模型时有一个强制操作不能省对 mps_state[tasks] 里已经下发但还没完成的任务要把对应的接入或发车二进制变量固定为 1然后传给求解器。这样本轮优化器不能反悔上一轮的决定。如果不加这个固定操作优化器会经常“把已经开到半路的 MPS 又调回原节点”表面上看目标函数更好实际方案根本不可执行。按照我复现这类算法的经验这一步是动态调度和普通多时段静态调度的本质区别。4.3 状态推进里最容易被忽略的细节剩余能量与在途位置状态推进函数 advance_state 看起来很机械但有两个细节特别容易翻车。第一移动中的 MPS 也要消耗时间不能假设它瞬间到达。如果两台节点之间的行车时间是 30 分钟而每个调度时段是 15 分钟那这辆 MPS 应该占用两个时段处于“在途”状态第二个时段结束时才能更新到目标节点。我一般会给每台在途 MPS 维护一个剩余行程时间字段每滚动一个时段就减 1减到 0 才更新位置。第二剩余能量的更新不能只算放电还要考虑移动途中是否耗电。绝大多数论文模型假设 MPS 移动作业不消耗车载能量但实际移动电源车在长途转场时会用自身电池驱动辅助系统这部分耗电虽然小却会让 soc 字段出现微小偏差。复现时可以在状态推进函数里加一个可配置的运输能耗系数默认设 0需要精细化时再打开。这种细节不影响主结果但会让你的复现代码在面对审稿人或者部门评审时更有说服力。def advance_state(mps_state, action, dt): for m in mps_state[pos].keys(): # 更新剩余行程时间: 只对在途车辆生效 if m in mps_state[traveling]: mps_state[traveling][m] - dt if mps_state[traveling][m] 0: mps_state[pos][m] mps_state[dest][m] del mps_state[traveling][m] # 接入供电时扣除放电能量, 按 15 分钟时段累计 if action.get(dispatch, {}).get(m): mps_state[soc][m] - action[dispatch][m] * dt return mps_state这里 action 是从求解结果里解析出来的“第一步动作”包括哪些车发车、哪些车在哪个节点放电。注意 soc 的扣除必须只针对当前时段不能把整个窗口的放电量一次性扣掉否则下一轮滚动窗口会用错误的电量状态重新规划误差越滚越大。5. 复现 MPS 动态调度绕不开的 5 个坑这个部分的每一条都是我实打实踩过的有的坑甚至让我在无解的模型里耗了两天。按“现象 → 原因 → 解决”的顺序写希望你能直接绕过。5.1 模型报 infeasible先查数值尺度而不是约束逻辑现象模型写好后一求解Gurobi 直接报 Model is infeasible而且 Irreducible Inconsistent Subsystem 指向的约束看起来完全合理。原因大部分情况不是约束逻辑错了而是数值尺度跨了太多量级。真实电量单位里负荷功率是几百千瓦线路电阻是零点几欧能量却是几千千瓦时这些数值混在一起后求解器内部的容差处理会出问题明明可行也判定为不可行。解决把所有功率、电压、能量全部转换到标幺值体系基准功率选系统容量基准电压选额定电压。转换完之后负荷、MPS 功率、线路参数都集中在 0.001 到 10 之间数值尺度问题基本消失。这也是我为什么在第 3 章要求先做数据预处理再写模型。5.2 SOCR 锥约束不紧解出来电压虚高现象模型能求解目标函数也很漂亮但把解回带到原始 DistFlow 方程里一验算支路电流平方和功率不满足原等式约束电压曲线异常平滑。原因二阶锥松弛把等号变成了不等式如果目标函数对这个锥约束没有“向边界拉”的作用松弛后的锥体内部点会被选中导致结果偏离真实物理可行域。这在某些目标函数和参数组合下特别明显。解决检查每个支路每个时段的锥约束松弛间隙也就是 I 和 P、Q、V 之间的差值。如果间隙超过 1%就在目标函数里加一个极小的惩罚项让求解器倾向于把锥约束压紧。惩罚系数从 1e-6 试起每次放大十倍观察松弛间隙变化直到间隙回到可接受范围。5.3 滚动时域里的“后悔药”已发车的移动被回退现象滚动调度结果里某台 MPS 第一轮从节点 A 发往节点 B第二轮又重新从节点 A 发往节点 C仿佛它从来没出发过。原因构建新一轮滚动窗口时没有把上一轮已经确定的发车变量和接入变量固定下来。求解器看到当前窗口的故障信息更好就把未来时段的决策改掉了这在数学上是合法的在物理上毫无意义。解决在 build_rolling_model 函数里强制加入固定变量逻辑对所有已经下发且未完成的任务把对应的 x 或 y 变量的上下界同时设为 1。另外建议在日志里记录每轮下发任务的时间戳和变量名排查问题时能直接对账。5.4 道路通行时间固定导致方案不可执行现象调度方案显示两台 MPS 都规划了完美的路径但现场反馈实际到达时间比计划晚了一个半小时窗口内的供电计划全部错位。原因模型里把任意两个节点间的行车时间设成了固定值没有考虑灾后道路拥堵、绕行、桥梁封闭造成的时间增加。这种情况在预设算例里不明显一旦接入真实路网数据误差会直接影响 MPS 到达时序。解决给行程时间加上场景系数。常用做法是维护一个道路可用性矩阵把断交道路对应的时间设成无穷大并在模型中删除该路径同时对剩余可用道路统一乘一个 1.2 到 1.5 的保守系数还可以把行程时间设成离散化的小时数而不是连续值这样更贴合实际调度习惯。5.5 求解时间耗尽把 suboptimal 误当最优解现象TimeLimit 设了 300 秒求解器在 300 秒时返回了一个可行解状态显示 feasible 或 suboptimal但日志里没细看直接拿这个解去画恢复率曲线结果下一轮实验换了参数后目标值突然提升了 8%。原因时间耗尽后返回的解只是当前最好可行解离真实最优还有一段距离。滚动时域每轮都产生一个新模型后一轮可能碰巧找到一个更好的解导致结果前后不一致复现结论站不住脚。解决每个滚动窗口求解结束后强制打印并记录 MIPGap 和目标函数值。MIPGap 超过 5% 的窗口需要延长 TimeLimit或者把窗口长度缩短。更重要的一点状态推进时用的 soc 和位置必须以“已经找到的可行解”为准不能假设后续窗口会补出一个更优解来修正之前的动作。6. 用 IEEE 33 节点算例验证复现结果三个进阶技巧模型和算法都跑通之后还有一个躲不掉的步骤怎么让结果可信。这个章节不讲新模型只讲三个我每次复现都会用的验证手段做完这三件事结果才敢拿出去见人。6.1 先搭“不移动”的对照组再用 MPS 移动性算增益动态调度的价值必须通过对照实验体现。我一般先在 IEEE 33 节点系统上跑两套方案第一套把 MPS 当成固定式储能落地后整个调度时域内不许移动第二套启用完整的移动和转场能力。两套方案的输入数据完全一致只差“能否移动”一个变量。对比两张恢复率曲线就能量化 MPS 移动性到底带来了多少韧性增益。方案恢复电量 (p.u.*h)MIPGap (%)单轮求解时长 (s)固定式储能由实验记录由求解器输出由求解器输出MPS 动态调度由实验记录由求解器输出由求解器输出这样填出来的表格才是你自己的复现结果而不是从哪篇论文里抄来的数字。如果算下来的增益低于 10%先别怀疑算法回去检查行程时间是不是设得过大、负荷权重是不是拉不开差距这两项直接决定了移动性的价值空间。6.2 用恢复率曲线下面积量化韧性提升配电网韧性最直观的指标是恢复率曲线横轴是时间纵轴是“当前时刻已恢复的关键负荷 / 总关键负荷”。动态调度响应越快、恢复率爬升越陡曲线下面积就越大。计算时把恢复率曲线按时间积分再除以总调度时长得到一个 0 到 1 之间的归一化韧性指标方便和不同方案直接比较。我在复现时会额外把曲线覆盖在同一个坐标系里和固定式储能方案做对比同时把每个时段的 MPS 位置标注在节点上。这样不仅能看到恢复率数值差异还能解释差异来源——是 MPS 转场更早还是能量分配更优。单纯看一个目标函数数字很难定位问题但曲线一画漏洞马上暴露。6.3 从确定性动态调度走向随机场景的改造方向如果目标是继续往上做学术延伸确定性滚动时域的下一步就是随机规划。常见改造是在每个滚动窗口内部挂多个故障场景把目标函数改成所有场景下恢复电量的期望值。这个改造的工作量不大但模型规模会乘以场景数需要把 TimeLimit 放宽松或者对场景做缩减。最后说一个我的习惯每次跑实验都会把每轮滚动窗口的 MIPGap、求解时长、下发任务明细存成日志文件。这个习惯帮我多次避开了“结果前后矛盾”的尴尬局面尤其是在调参数卡住的时候翻日志比对就能定位是哪一轮决策出了问题。动态调度复现拼的不只是模型功底更是对求解器行为细节的把控。希望这篇笔记能帮你在复现 MPS 动态调度这条路上少走几步弯路。本文还有配套的精品资源点击获取
