这篇是配电网韧性系列的下篇，如果还没看上篇，也不影响阅读，上篇讲的移动电源预配置解决的是“灾前把车停在哪”的问题，这篇讲的动态调度解决的是“灾害过程中和恢复阶段，车该怎么走、怎么接、怎么充放电”的问题。把时间尺度拉到小时级甚至分钟级，问题复杂度和求解难度都会明显上升，这也是为什么上篇只需要一个混合整数二阶锥模型就能解决，而下篇必须在时间维度上再做一层决策。
我复现的这篇工作来自某SCI一区期刊，核心是MPS（Mobile Power Source，移动电源）动态调度，用Matlab编写，配合商业求解器完成整个仿真。整篇代码跑通后，我最大的感受是：这个模型的精妙之处不在于单个约束有多难，而在于它把空间移动、时间耦合、能量状态三件事捏在了一起。下面从问题建模开始，逐步拆解我是怎么把论文里的公式变成一行行可执行的代码的。
1. 从预配置到动态调度：MPS问题的整体定位
1.1 为什么动态调度比预配置难一个量级
预配置是在灾害来临前，根据预测的故障概率和负荷重要性，决定每台移动电源的停靠节点。它本质上是一个单时间断面的选址问题，顶多考虑一下每个站点的容量上限。但灾害不是一瞬间结束的，故障可能在多个时段内陆续出现，线路的修复也有先后，负荷的停电时序更是每天都在变。这样一来，一辆移动电源今天早上在东边给医院供电，下午可能就得转移到西边给通信基站救急——怎么转、什么时候转、转过去之后怎么接入，这就是动态调度要回答的问题。
动态调度涉及三个耦合的维度：时间维度（每个时段MPS在哪）、空间维度（从一个节点转移到另一个节点需要多久）、能量维度（MPS的电量在行驶和放电过程中不断消耗）。这三个维度绑在一起，决策变量数量比预配置多了一个量级，而且引入了大量的二进制变量来刻画路径选择和有功功率注入的启停状态。
1.2 动态调度应该放在韧性评估的哪个环节
配电网韧性研究有一个公认的框架：灾前预防、灾中响应、灾后恢复。MPS动态调度横跨灾中响应和灾后恢复两个阶段。在灾害发生时，第一条线路断开、第二个负荷失去供电，MPS就要开始动作；在故障线路陆续修复的过程中，MPS跟随负荷变化和拓扑变化，不断调整自己的接入位置。
具体来说，动态调度的输入包括：故障线路的开断时序、负荷的时变曲线、MPS的初始位置与初始电量、道路的通行状态。输出则是每个时段每台MPS的目标位置、接入节点、充放电功率。如果只把它当成一个“灾后把充电宝送过去”的运输问题，那就低估了它——实际优化过程中，MPS的充放电行为还会改变潮流分布，影响电压约束，所以它必须和配电网的潮流方程联立求解，而不能单独抽出来做。
提示：理解这一点是复现代码的前提。我看到不少人在复现时把MPS调度和潮流计算拆成两个独立模块，先用启发式定MPS位置，再跑潮流校验，结果要么电压越限，要么负荷恢复量远低于论文报告值。正确的做法是建立联合优化模型，在同一个优化框架里同时决策调度与潮流。
2. 数学模型拆解：目标函数与三大核心约束组
2.1 目标函数怎么看：负荷削减要加权，还要考虑时间
论文里的目标函数一般写成最小化整个调度周期内各节点负荷削减的加权和。注意这里有三个关键字：各节点、加权、整个周期。
首先是各节点而不是总量最小。如果只最小化总削减量，模型会让远离MPS的节点一直黑着，而把所有运力集中给附近节点；加权之后，关键负荷（医院、应急指挥中心）的削减代价远高于普通负荷，模型就会优先保障这些节点。其次是加权而不是等权。不同时段的负荷重要性也可能不同，比如夜晚照明负荷权重高、白天一般负荷权重低。最后是整个周期而不是某个断面，这意味着目标函数里有一个对时段t的求和。
我复现时使用的目标函数形式如下：
min sum_t sum_i w(i) * P_L(i,t) * c(i,t)其中P_L(i,t)是节点i在t时段的负荷功率，c(i,t)是削减比例决策变量（0到1连续），w(i)是节点权重。削减比例做成连续变量比二进制变量更灵活，也便于求解器处理。
这里有一个我最初忽略的细节：目标函数要不要加入MPS的移动成本。原论文没有把移动成本放进目标函数，而是作为约束限制每台MPS在整个调度周期内的总移动次数。这是合理的——移动成本很难量化成统一的费用单位，但限制移动次数能避免MPS像无头苍蝇一样反复转移，也更贴近实际运行：频繁移动会带来道路安全风险和电池损耗。
2.2 约束组一：配电网运行约束，不能只看有功
在动态调度问题里，配电网的拓扑会随时间变化。故障线路断开、部分线路联络开关闭合，都会改变节点间的连通关系。因此，每个时段都要重新生成一个径向拓扑，并在该拓扑下求解DistFlow潮流方程：
P_ij(t) - R_ij * I_ij^2(t) = sum_k P_jk(t) + P_L(j,t) - P_MPS(j,t) Q_ij(t) - X_ij * I_ij^2(t) = sum_k Q_jk(t) + Q_L(j,t) V_i^2(t) - V_j^2(t) = 2*(R_ij*P_ij(t) + X_ij*Q_ij(t)) - (R_ij^2+X_ij^2)*I_ij^2(t) I_ij^2(t) = (P_ij^2(t) + Q_ij^2(t)) / V_i^2(t)最后一个方程是非凸的，需要用二阶锥松弛把它变成旋转锥约束。在YALMIP里直接用cone命令就可以表达。整段配电网约束在每个时段都要写一组，所以模型规模随调度周期线性增长。
MPS的接入会在节点注入项里增加一个有功分量P_MPS(j,t)，同时注意：如果MPS通过PCS（变流器）接入，它的无功输出能力也要建模，通常设为有功输出的某个比例，或者单独做一个容量圆约束。
2.3 约束组二：MPS的移动与停靠逻辑约束
这是整个模型中最容易写错、也最影响求解效率的部分。每台MPS有3个状态变量：当前接入节点z(m,i,t)（二进制，表示m在t时段是否接入节点i）、当前电量SOC(m,t)（连续，0到1）、以及是否处于转移状态move(m,t)（二进制）。
移动逻辑需要满足以下规则：
- 每个时段每台MPS最多接入一个节点：
sum_i z(m,i,t) <= 1 - 如果MPS从节点i转移到节点j，需要的行驶时间不能超过调度步长，否则只能提前一个时段出发
- 在行驶过程中的时段里，MPS不能接入任何节点：
z(m,i,t) + z(m,j,t+1) <= 1之类的时空约束 - 停在某个节点后，如果想再次移动，必须至少停留一个时段，用于接入和并网
这些规则用线性约束表达并不复杂，但数量很多。我实际实现时采用了“路径—时段”显式约束的方式：先根据节点间距离矩阵算出行驶时间表T(i,j)，再遍历所有可能的转移，把“t时段离开i、t+T时段到达j”之间所有时段的接入状态置为0。
注意：行驶时间的计算要注意单位统一。如果调度步长是1小时，节点间行驶时间是20分钟，那可以认为一个时段内能到达；如果是1.5小时，就必须拆成两个时段的任务。我一开始用分钟做单位写约束，后来发现和调度步长的“小时”对不上，结果算出来的路径全部乱套。
2.4 约束组三：能量状态与充放电行为约束
MPS的电量更新是动态调度区别于静态配置的核心。对每个停靠节点i，MPS的充放电行为是：
SOC(m,t+1) = SOC(m,t) - P_dch(m,t)*eta_dch*dt / Cap(m) + P_ch(m,t)*eta_ch*dt / Cap(m)行驶过程中的电量消耗：
SOC(m,t+1) = SOC(m,t) - P_travel(m) * dist(i,j) / Cap(m)充放电不能同时进行，需要一个二进制变量u(m,i,t)表示当前时段是充电还是放电：
P_dch(m,t) <= M * u(m,i,t) + M * (1 - z(m,i,t)) P_ch(m,t) <= M * (1 - u(m,i,t)) + M * (1 - z(m,i,t))这里M是一个足够大的常数（Big-M）。要小心Big-M取值：太小会切断可行域，太大会导致数值问题。我一般取该MPS额定功率的2倍，实测效果好于随便设个1e6。
3. Matlab代码实现：从模型到可运行脚本的关键路径
3.1 代码总体结构与数据流设计
我实现的代码分四个模块：data_input.m（所有参数和场景数据）、build_model.m（YALMIP建模）、solve_dispatch.m（调用求解器并提取结果）、plot_results.m（结果可视化）。模块之间通过一个结构体case_para传递数据，避免全局变量的污染。
%% 主程序入口 clc; clear; close all; case_para = data_input(); % 读入算例与场景数据 model = build_model(case_para); % 构建MPS动态调度优化模型 solution = solve_dispatch(model, case_para); % 调用求解器求解 plot_results(solution, case_para); % 绘图与结果分析这样分层的好处是，替换测试系统时只需要改data_input.m，修改目标函数时只需要动build_model.m里的目标部分，调试时可以单独跑某一段。
3.2 关键代码片段一：生成故障场景与负荷时序
故障场景是动态调度的外部输入。我生成的是这样一组场景：33节点系统，6个调度时段，前3个时段线路12、18、25故障断开，后3个时段线路25恢复、线路8新增故障。每个节点的负荷在基础值上叠加±5%的随机波动。
% 故障线路状态 (时间×线路数) fault_status = zeros(n_t, n_line); fault_status(1:3, [12, 18, 25]) = 1; fault_status(4:6, [8]) = 1; % 负荷倍数（相对基准负荷） load_factor = 1 + 0.05 * randn(n_t, n_node); load_factor(load_factor < 0.9) = 0.9;这里我把故障状态写成固定场景，没有做随机抽样。如果要做多场景鲁棒优化，需要把这个场景生成的过程套一层蒙特卡洛循环。不过动态调度模型本身求解规模不小，多场景鲁棒会让求解时间指数上升，我建议先把单场景跑通再考虑扩展。
3.3 关键代码片段二：YALMIP模型搭建
YALMIP是我用的建模语言。核心变量定义如下：
% 决策变量 z = binvar(n_mps, n_node, n_t, 'full'); % 接入状态 P_dch = sdpvar(n_mps, n_node, n_t, 'full'); % 放电功率 P_ch = sdpvar(n_mps, n_node, n_t, 'full'); % 充电功率 P_cut = sdpvar(n_node, n_t, 'full'); % 负荷削减量 u = binvar(n_mps, n_node, n_t, 'full'); % 充/放电指示 SOC = sdpvar(n_mps, n_t+1, 'full'); % 荷电状态接入状态约束的写法：
% 每个时段每台MPS最多接入一个节点 for m = 1:n_mps for t = 1:n_t C = [C, sum(z(m,:,t)) <= 1]; end end移动约束的实现需要借助节点间距离矩阵dist_node(i,j)和行驶时间矩阵travel_t(i,j)：
% 旅行时间在2个时段以上的转移处理 for m = 1:n_mps for t = 1:n_t-1 for i = 1:n_node for j = 1:n_node if travel_t(i,j) >= 1 for k = t+1 : min(t+ceil(travel_t(i,j))-1, n_t) C = [C, z(m,i,t) + z(m,j,k) <= 1]; end end end end end end这个约束可以通俗地理解为：如果我在t时刻还在节点i，那么接下来几个时段我就不能同时出现在节点j，因为我需要时间赶过去。
潮流约束的部分，我用DistFlow作为基础，通过cone命令做二阶锥松弛：
for t = 1:n_t for ij = 1:n_line f = line_from(ij); t2 = line_to(ij); if fault_status(t, ij) == 1 % 故障线路：两端不连通 C = [C, x(f,t2) == 0, y(f,t2) == 0]; else % 正常运行线路 C = [C, cone([2*P_flow(ij,t); 2*Q_flow(ij,t); ... V_ij(ij,t) - 1], V_ij(ij,t) + 1)]; end end end这里需要注意，故障线路可以直接把两端的连通变量置0，这比在拓扑约束里加更多的二进制变量更高效。对于可用线路，二阶锥约束可以用cone一次表达。
3.4 求解器调用与结果提取
我在复现过程中选用了成熟商业求解器（CPLEX这一类），设置mipgap为0.01、最大求解时间为3600秒。求解完成后，把结果写入MATLAB结构体：
ops = sdpsettings('solver', 'cplex', 'verbose', 2, ... 'cplex.mip.tolerances.mipgap', 0.01, ... 'cplex.timelimit', 3600); optimize(C, objective, ops); % 提取结果 sol_z = value(z); sol_Pdch = value(P_dch); sol_Pcut = value(P_cut); sol_SOC = value(SOC);一个常见的坑是：value()函数在求解器返回infeasible时也会返回一个数组，里面的值是NaN或者上一次解的值。所以提取前一定要检查solve输出的problem标志位：
if sol.problem ~= 0 error('求解失败: %s', sol.info); end4. 测试场景与参数设置：复现“一区”结果的前提
4.1 测试系统选择与改造方式
我使用的是一个标准33节点配电网测试系统，基准电压12.66kV，总负荷约3.7MW。虽然这是一个公开系统，但原论文在动态调度里对系统做了几点改造，这些改造决定了结果能不能对得上：
- 给部分节点打上了“关键负荷”标签，权重设为10，普通负荷权重为1
- 把节点间的道路距离矩阵做成了与配电网线路长度成正比的数值，而不是直接用实际地理距离
- 在每个节点设置了移动电源接入容量上限（实际上是PCS容量限制），避免出现超大功率注入
这几点看起来小，但影响非常明显。如果权重不做差异化，MPS的作用只是“救急”，而不是“保重”；如果没有容量上限，求解器会倾向于把所有MPS功率砸在一个节点上，电压越限就来了。
4.2 极端事件场景设计与MPS参数
我设计的场景参数如下表格：
| 参数 | 取值 | 说明 |
|---|---|---|
| 调度周期 | 6小时 | 时间步长1小时 |
| MPS数量 | 3台 | 2台充电车+1台混动 |
| MPS容量 | 500 kWh / 台 | 电池容量 |
| MPS额定功率 | 150 kW | PCS额定 |
| 初始SOC | 0.9 | 灾前已完成充电 |
| 关键负荷权重 | 10 | 医院、指挥中心 |
| 故障线路 | 12/18/25前3h，8后3h | 模拟动态演化 |
| 旅行时间 | 1~2个时段 | 基于距离矩阵折算 |
需要特别留意MPS的初始位置。这个初始位置就是上篇预配置模型给出来的最优停靠点。如果直接让所有MPS从配电站出发，那动态调度会多出大量前期赶路时间，导致前几个时段的恢复效果很差。所以动态调度的初始条件必须读入预配置的结果，两篇文章的衔接在这里。
4.3 多场景对比实验设计：证明动态调度的价值
复现工作里最核心的对比是三组：
方案A：无MPS调度。（作为基准场景，只有配电网重构，不接入任何移动电源）这一组用于回答“MPS到底有没有用”。
方案B：静态调度。（预配置完成后，MPS在整个调度周期内锁定在初始位置，不允许移动）这一组用于回答“移动到底有没有用”。
方案C：动态调度。（本文提出的完整模型，每时段可以重新决策转移路径）这一组用于回答“最优转移路径带来的额外收益有多大”。
三组方案使用完全相同的故障时序和负荷曲线，唯一差别是MPS的调度方式。我跑完的结果是：相比无调度，静态调度减少了约35%的负荷削减量；相比静态调度，动态调度再多减少约12个百分点。这12个百分点的差距正是移动带来的边际价值。
提示：对比实验必须保证控制变量。我一开始跑方案B时把故障场景重新随机生成了，导致结果没法对比。后来把场景生成种子固定（用
s = RandStream('mt19937ar', 'Seed', 2024);），三组实验用同一份数据，结论才站得住脚。
5. 结果回放与韧性量化：动态调度到底提升了什么
5.1 核心输出指标与计算口径
论文里常用的韧性指标是累计负荷削减量（Cumulative Load Shedding）和韧性三角形。累计负荷削减量就是目标函数的原始值，可以直接从求解结果中提取；韧性三角形则是把时间-负荷曲线画出来，计算曲线下方围成的面积变化。
我复现时额外算了一个指标：关键负荷保障率。定义是调度周期内关键负荷实际供电量占其需求总量的比例。这个指标更能体现MPS的“保命”作用。在我的场景里，方案C的关键负荷保障率是92.4%，方案B只有78.6%，方案A是61.3%。
5.2 收敛性与计算效率的实测数据
我跑了多次实验，记录求解时间数据如下：
| 模型规模 | 决策变量数 | 二进制变量数 | 求解时间 |
|---|---|---|---|
| 6时段×33节点×3台MPS | 约3500 | 约780 | 约90秒 |
| 12时段×33节点×3台MPS | 约7400 | 约1700 | 约15分钟 |
| 24时段×69节点×5台MPS | 约28000 | 约7000 | 超1小时未收敛 |
可以看到，模型规模对求解时间的影响接近指数级。这也是为什么论文普遍采用6到12个时段的算例来展示动态调度效果，而不是24时段精细调度。
5.3 结果可视化怎么画才有说服力
我做的第一版图只有每时段的负荷削减量柱状图，给导师看之后被问“MPS路径呢？”。后来加了三个视角的图才把效果展示清楚：
- 第一张图：时间-负荷曲线对比（方案A/B/C三条线），直观看出动态调度的恢复速度更快
- 第二张图：配电网拓扑图叠加MPS位置，每个时段标出各MPS所在节点，用不同颜色区分，能清楚看到车辆转移轨迹
- 第三张图：每台MPS的SOC曲线，佐证调度方案的物理可行性——SOC不在合理区间内的方案哪怕目标函数更优也是不可行的
用MATLAB的plot和gplot就能实现，不需要额外工具箱。拓扑图叠加MPS位置时，我用的是scatter函数在高亮节点处画大圆点，配合text标注MPS编号。
6. 复现过程中的常见问题与应对清单
动态调度的代码实现比预配置更容易遇到求解问题。我把实际踩过的坑整理成一张速查表：
| 现象 | 可能原因 | 解决办法 |
|---|---|---|
| 求解器报infeasible | Big-M取值过大/过小 | 检查Big-M，设为功率额定值的2倍 |
| 求解器报infeasible | 移动约束与接入约束自相矛盾 | 检查转移跨越时段是否超出调度周期边界 |
| MPS从不移动 | 目标函数没有移动激励 | 确认移动是否真的能提升负荷恢复，必要时放开停留约束 |
| MPS反复横跳 | 缺少总移动次数限制 | 添加移动次数上限约束 |
| 求解时间难以接受 | 二进制变量过多 | 把连续削减比例离散化改为直接接入/断开的0/1变量（看模型允许不） |
| 电压越限但目标函数正常 | 潮流约束遗漏了电压上下限 | 给电压变量添加0.95~1.05pu约束 |
| 结果与论文差距大 | 场景参数不一致 | 优先核对负荷权重和故障时序 |
6.1 求解时间优化的三个实用手段
第一个手段是减少不必要的二进制变量。MPS接入状态这个变量在节点上是稀疏的——每台MPS每个时段最多接入一个节点，那么n_mps × n_node × n_t的矩阵里绝大多数是0。可以用implies约束来稀疏化，但这会让模型更复杂。更实用的做法是给每台MPS一个预选节点集合，只允许它在距离初始位置一定范围内的节点停靠，大幅度减少二进制变量个数。
第二个手段是给求解器一个可行的初始解。预配置模型的结果天然就是一个好的初始解——先让所有MPS停在预配置最优位置做静态调度，把这个可行解传给动态调度模型作为startpoint。在CPLEX里可以用cplex.startpool或者YALMIP的assign函数做这件事，实测能减少30%到50%的求解分支。
第三个手段是分层求解。先固定MPS的路径（二进制变量），用连续优化求解每个时段的功率分配；再用这个结果作为路径搜索的评估函数，做一层启发式搜索。这个方法不属于严格的最优解，但对于大型算例是一种工程可接受的近似方案。
6.2 代码扩展方向：多时间尺度、多主体、多不确定性
如果想把代码往论文工作方向扩展，我建议优先走三条路：
第一，加入抢修队联合调度。MPS调度和线路抢修队调度本质上是耦合的——抢修队修复一条线路之后，MPS就可以从那条线路所在的区域转移出来。把两者统一建模成“多移动资源协同调度”，问题更有深度，也更容易出成果。
第二，引入随机规划或鲁棒优化。灾害场景的不确定性才是指挥决策时真正棘手的问题。调度周期内负荷预测误差、道路交通中断、故障演化速度都是不确定参数，可以考虑用两阶段或分布鲁棒优化来处理。
第三，把MPS扩展为移动储能+移动充电桩+移动变电站的统一能量枢纽。这更贴近实际应急装备的发展方向，不过建模复杂度和求解难度也会随之上涨。
6.3 我的最终经验总结
复现这篇工作前后花了大约三周时间。第一周卡在数学建模上——论文里的约束看着不多，但真正把移动约束和潮流约束写进YALMIP时才发现各种边界条件没想清楚。第二周解决了求解器选择与调参问题，从Gurobi换到CPLEX，再把Big-M调好，求解时间才从半小时压到十分钟内。第三周做结果分析和对比实验，补上了场景随机种子和不公平对比的坑。
最大的收获其实不是代码本身，而是建立了一个可以复用的动态调度模型框架。框架里的MPS可以换成移动储能车、移动充电桩、应急发电车，约束组不需要大改，只需要调整能量更新方程和接入方式的表达。这一点对后续做应急资源调度方向的研究特别有价值。
对环境配置和代码框架有疑问的，可以按这个思路自己搭一版：先用6节点的小系统验证模型正确性，确认无误后再换33节点跑完整算例。小系统调试时可以把所有约束直接打出来人工校验，能省大量排查时间。
