1. 从“预配置”到“动态调度”为什么韧性提升必须走完这最后一步前一阵子把一篇配电网韧性提升方向的SCI一区论文完整复现了一遍上篇我们已经聊过应急移动电源Mobile Power SourceMPS的预配置也就是说在灾害来临前基于气象预报和历史数据把有限的移动电源车提前部署到最可能受灾、负荷又重要的节点附近。当时我留了一个尾巴预配置只是“开局布阵”灾害真正发生之后电网的结构、故障位置、负荷损失情况都在不断变化如果移动电源只会停在预配置的位置不动那韧性提升的效果至少要打对折。这篇就接着把这个下半场讲透也就是MPS的动态调度并且在Matlab里给出可跑的代码实现思路。先说清楚什么是“动态调度”。灾害发生之后配电网的故障信息是逐步明确的不是一次性全部暴露。比如台风过境可能先是一回馈线跳闸抢修过程中又出现第二处故障原计划接入MPS的节点可能本身已经被水淹了或者负荷转移后其他节点的缺电更严重。这时候就需要一个决策机制在每一个“信息更新时刻”重新回答三个问题我手头这几台移动电源车现在应该接到哪个节点充了多少电、还剩多少能量接下来往哪个方向移动最划算这套机制就是MPS动态调度。它是典型的“信息不完全、滚动决策”问题跟预配置那种“基于场景集合的一次性优化”有本质区别。我之前做的时候最大的感受是预配置更像下棋开局摆子动态调度则是棋到中盘、每步都要重新算。两者配合起来才是完整的韧性提升调度框架。这篇文章适合正在做配电网韧性、极端灾害下供电恢复、移动储能调度这类方向的硕士博士和工程师阅读。如果你手里恰好有Matlab和Gurobi或者CPLEX这类求解器基本上可以照着我给出的模型和代码结构把动态调度模块搭起来再跟你自己的预配置模型对接。文章中会包含数学模型拆解、滚动时域框架设计、关键代码实现、以及我复现过程中踩过的坑尤其是那些论文里不会写的细节。2. 动态调度的数学建模目标函数和约束条件怎么定2.1 动态调度的目标函数先说建模目标。预配置阶段的目标往往是在灾害发生前的多个预测场景下最大化全时段的关键负荷恢复期望值动态调度阶段不一样它面对的是一个正在演化的确定性或部分可观测状态目标函数通常写成“从当前时刻到调度结束累计恢复的失电负荷最大化”或者“累计切负荷量最小化”。我复现的这篇论文采用的是后者也就是最小化加权失负荷量。写成数学表达式大致是这样的[ \min \sum_{t \in \mathcal{T}} \sum_{i \in \mathcal{N}{load}} w_i \cdot P{i}^{shed}(t) \cdot \Delta t ]其中 (w_i) 是负荷节点 (i) 的权重医院、应急指挥中心这类关键负荷权重高普通商业和居民负荷权重低(P_{i}^{shed}(t)) 是时段 (t) 在节点 (i) 的切负荷量(\Delta t) 是每个时段的时长。目标函数很简单但动态调度真正的难度在约束条件。2.2 约束条件里的“硬骨头”MPS动态调度的约束可以从四个维度来理解每一个都对应一个工程现实。第一个是MPS自身的能量和功率约束。每台移动电源车有额定容量比如500kWh最大放电功率比如150kW。它在某个节点接入时输出功率不能超过最大放电功率同时在整个调度周期内累计放出的能量不能超过初始电量加如果能充电的话补充电量。这个约束决定了MPS不能“永动机式”地到处救火跑地方越多、出力越大能量下降越快后面时段可能就无电可放了。第二个是MPS的空间移动约束。这是动态调度区别于预配置的核心所在。MPS在当前时段接入节点 (n_1)如果下一时段要接到节点 (n_2)必须要“花时间移动”。移动时间取决于两个节点之间的地理距离和移动速度。在离散时段模型里这个约束通常写成“若MPS在时段 (t) 和 (t1) 接入不同节点则它们之间的距离必须小于等于移动速度乘以时段时长”。如果配电网的网络拓扑出现断线还要考虑哪些节点之间通过路网仍然可达。现实中移动电源车走的是公路不一定沿着电网走所以这里需要引入路网距离矩阵我早期图省事直接用电气距离结果计算出来的调度路径在现实中根本不可行。第三个是配电网络的运行约束。MPS接入到某个节点后它本质上是一个分布式电源向该节点及其下游或相邻节点供电。因此需要满足潮流方程。对于辐射状配电网常用的是DistFlow线性化模型[ P_{ij}(t) \sum_{k: (j,k) \in \mathcal{E}} P_{jk}(t) P_{j}^{load}(t) - P_{j}^{MPS}(t) - P_{j}^{DG}(t) ][ Q_{ij}(t) \sum_{k: (j,k) \in \mathcal{E}} Q_{jk}(t) Q_{j}^{load}(t) - Q_{j}^{MPS}(t) ]这套方程通过线性化把支路潮流简化成节点功率平衡加线路容量约束避开了二次项使得整个模型可以交给混合整数线性规划求解器处理。此外还要考虑节点电压幅值约束、线路电流容量约束。在动态调度里网络拓扑还可能因为故障和抢修而改变所以拓扑参数要做成时变参数输入而不是静态定值。第四个是决策变量的时序一致性约束。比如MPS不能凭空消失再凭空出现也就是每个时段每台MPS必须而且只能在一个节点两台MPS不能同时接同一个节点如果节点容量不够。这些在模型里都是线性约束写起来不难但容易漏。我把这几类约束用一张表整理一下方便对照约束类别代表约束物理意义易错点能量/功率约束累计放电量 ≤ 初始电量输出功率 ≤ 最大放电功率电池不是无限能量漏掉时段累计约束导致MPS反复“满血复活”移动约束(\sum_{n} x_{m,n,t}1)跨节点移动需满足时间窗移动需要时间不能瞬移忽略路网距离直接用电气距离网络运行约束DistFlow功率平衡、电压上下限、支路容量接入后不会过载、压降不越限线性化时忽略二次项误差累积时序一致性(x_{m,n,t}) 为0-1变量每台每个时段位置唯一决策逻辑自洽缺少状态转移约束位置跳变无轨迹2.3 动态调度与预配置在模型层面的根本差异预配置和动态调度虽然目标一致但模型结构有两处根本不同。第一不确定性处理方式不同。预配置通常使用场景法设定多个灾害场景每个场景对应一组故障和负荷值对场景求期望或最坏情况。动态调度用的是滚动时域加实时信息更新每一轮优化只针对当前已知信息和未来短时预测不追求一次性覆盖全场景。这样做的原因是灾害演化过程很难精确预测过早绑定某个场景反而会让决策失真。第二决策变量类型不同。预配置的决策变量主要是“哪台MPS停到哪个节点”的离散变量一旦确定基本不再改变动态调度则要在整个灾后恢复期内反复调整位置和出力每一个时段都有对应的位置变量和功率变量所以决策变量维度至少比预配置多一个时间维度。这也直接导致动态调度的求解难度显著上升。3. 滚动时域调度框架为什么不能“一次算到底”拿到完整的数学模型后我一开始的想法非常简单粗暴直接对整个灾后恢复期比如24小时按1小时间隔分成24个时段建立一个大混合整数线性规划模型一次性求出所有时段MPS的位置和出力。当时我还以为这样更“全局最优”。但跑下来结果很不理想不是求解器算不动就是给出的方案在实际场景里根本没法执行。原因在于信息是逐步到达的。论文里的动态调度本质上是一个闭环反馈决策流程每个决策周期开始时收集当前配电网的实时状态包括哪些线路仍然故障、哪些负荷恢复、MPS当前电量然后在这个信息集上求解从当前时段到预测时域的优化问题但只执行第一个时段的决策。等到下一个周期信息更新再重新求解如此反复。这就是模型预测控制Model Predictive ControlMPC的基本思想也是滚动时域调度框架。用伪代码把这个流程表达出来是这样的初始化读取预配置阶段的MPS初始位置和电量 for k 1:K 更新当前网络拓扑、故障信息、负荷曲线、MPS状态 求解MPC优化模型决策时域为[k, kH] 执行第k时段的调度指令 更新MPS位置和电量状态 end这里的 (H) 是预测时域长度取太小会让决策只顾眼前取太大会引入过多不确定信息导致模型保守或求解困难。我实测下来当时间粒度为1小时、预测时域取4到6小时、MPS数量不超过5台、节点规模在100个左右时Gurobi求解一轮大概需要几十秒到几分钟基本能满足在线调度的需求。如果节点规模超过200耗时就会明显上升此时需要考虑有效不等式或者分解算法。滚动时域的另一个好处是天然支持“预配置—动态调度”两阶段衔接。预配置给出的MPS初始位置只要在动态调度第一轮作为初始状态传入即可两者不需要强行耦合在同一套模型里求解这大大降低了建模和求解的复杂度。4. Matlab代码实现从模型到可运行的程序4.1 数据结构设计讲完框架进入代码环节。我基于YALMIP加Gurobi的搭配在Matlab里实现了完整的滚动时域动态调度。代码结构上我建议按三层来组织数据层、模型层、调度层。数据层负责将配电网参数和灾害场景封装成统一的struct。我用的核心数据结构如下%% 配电网基础数据 net.N 33; % 节点数 net.E 32; % 支路数 net.from [1 2 3 ...]; % 各支路首端节点 net.to [2 3 4 ...]; % 各支路末端节点 net.R [0.0922 0.493 ...]; % 支路电阻(p.u.) net.X [0.0470 0.251 ...]; % 支路电抗(p.u.) net.rater ones(32,1)*100; % 支路容量 %% 负荷数据 loads.active [100 90 120 ...]; % 各节点有功负荷(kW) loads.reactive [60 40 80 ...]; % 各节点无功负荷(kvar) loads.weight ones(33,1); % 节点负荷权重关键负荷置为10或100 %% MPS参数 mps.num 3; % MPS台数 mps.capacity 500 * ones(3,1); % 每台容量(kWh) mps.pmax 150 * ones(3,1); % 最大放电功率(kW) mps.init_pos [12 25 7]; % 初始接入节点来自预配置阶段 mps.init_soc [0.8 0.8 0.8]; % 初始SOC %% 路网距离矩阵 mps.road_dist roadDistMatrix; % 任意两节点间公路距离(km)从地图数据生成 mps.speed 30; % 移动速度(km/h)这里特别提醒一下roadDistMatrix一定要基于实际路网生成最简单的办法是用经纬度坐标通过地图API算路网距离或者用各节点之间的直线距离乘以绕路系数。论文里一般不会讲这个细节但实际复现的时候这一笔数据直接决定调度路径是否可执行。4.2 核心决策变量定义模型层是代码的核心。按动态调度的要求决策变量需要分三类定义%% 决策变量定义 T 6; % 预测时域长度当前MILP的时段数 N net.N; M mps.num; % 1. 位置变量MPS在时段t是否接入节点n x binvar(T, M, N, full); % x(t,m,n) % 2. MPS出力变量 p_mps sdpvar(T, M, N, full); % p_mps(t,m,n) 每台MPS在各时段的输出功率 % 3. 切负荷变量 p_shed sdpvar(T, N, full); % p_shed(t,n) 各节点各时段的切负荷量 % 4. 线路潮流变量 p_ij sdpvar(T, E, full); % 支路有功潮流 q_ij sdpvar(T, E, full); % 支路无功潮流位置变量 (x(t,m,n)) 是三维0-1变量这是模型里规模最大的部分。如果预测时域取6、MPS有3台、系统有33个节点位置变量就有594个0-1变量加上连续变量和约束这个MILP规模对求解器来说不算轻松但对Gurobi完全可解。如果节点数到100以上建议把MPS可接入节点限制在候选节点集合内而不是所有节点能大幅减少变量数量。4.3 约束条件的代码表达约束条件是最容易出现bug的部分我把动态调度特有的几类约束代码贴出来%% 约束每台MPS每个时段必须且只能在一个位置 Constraints []; for t 1:T for m 1:M Constraints [Constraints, sum(x(t,m,:)) 1]; end end %% 约束MPS累计放电量不能超过初始电量 % soc(t,m)表示第m台MPS在时段t结束时的剩余电量比例 soc sdpvar(T, M, full); Constraints [Constraints, soc(1,m) mps.init_soc(m) - sum(p_mps(1,m,:))*dt / mps.capacity(m)]; for t 2:T for m 1:M Constraints [Constraints, soc(t,m) soc(t-1,m) - sum(p_mps(t,m,:))*dt / mps.capacity(m)]; Constraints [Constraints, soc(t,m) 0.1]; % SOC下限保护 end end %% 约束MPS移动时间窗核心动态调度约束 % 从节点n1到节点n2移动时间 road_dist(n1,n2)/speed % 如果MPS在t时段在n1t1时段在n2且n1~n2则需要满足移动时间 dt for t 1:T-1 for m 1:M for n1 1:N for n2 1:N if n1 ~ n2 travel_time mps.road_dist(n1,n2) / mps.speed; Constraints [Constraints, x(t,m,n1)x(t1,m,n2)-1 travel_time/dt]; end end end end end这段移动时间约束写得朴素但有效。逻辑上是这样的如果 (x(t,m,n1)1) 且 (x(t1,m,n2)1)那么等式左边等于1约束强制 (travel_time/dt \geq 1)也就是说MPS在两时段之间必须拥有至少一个时段的移动时间。如果 (travel_time/dt) 大于1意味着从一个节点到另一个节点需要超过一个时段那模型会在下一时段继续约束。这个写法属于最直接的Big-M思想但胜在直观、不容易出错。效率上如果候选节点数量大这个双重循环会产生大量约束建议先预计算所有节点对的旅行时间并只保留 (travel_time dt) 的组合来生成约束。网络运行约束的代码相对常规%% 约束DistFlow潮流线性化 for t 1:T for j 2:N % 功率平衡 parent_idx find(net.to j); % 这里需要根据你的数据结构匹配 in_flow 0; for k 1:length(parent_idx) in_flow in_flow p_ij(t, parent_idx(k)); end Constraints [Constraints, in_flow - sum(p_mps(t,:,j)) loads.active(j)*load_factor(t,j) - p_shed(t,j)]; end end注意负荷曲线是随时间变化的load_factor(t,j)表示节点 (j) 在时段 (t) 的负荷系数这一步容易漏掉。如果直接把静态最大负荷代入MPS出力计划会偏保守在低负荷时段浪费宝贵的移动电源能量。4.4 求解器调用与结果提取模型组装完之后调用求解器就比较机械了%% 求解 ops sdpsettings(solver, gurobi, verbose, 0, showprogress, 0); ops.gurobi.MIPGap 0.01; % 允许1%次优解大幅提升求解速度 ops.gurobi.TimeLimit 300; % 单轮求解时限5分钟 optimize(Constraints, Objective, ops); %% 结果提取 x_val value(x); p_mps_val value(p_mps); p_shed_val value(p_shed); soc_val value(soc);这里有个经验MIPGap设为0.01跟设为0相比求解速度往往能快几倍甚至一个量级而工程决策上1%的次优性完全可以接受。尤其是滚动时域调度本身就在不断更新信息没必要在单轮上追求绝对最优否则下一轮信息一更新上一轮算的“最优”就失效了。这一点是好多初学者没想通的地方。4.5 滚动时域主循环把上面的模型封装成函数solveMPSDynamicScheduling(mpc_data, mps_state, current_time)然后在主循环里反复调用%% 滚动时域主循环 horizon 24; % 总共24个时段即24小时 pred_len 6; % 滚动时域长度 dt 1; % 时段长度(h) % 当前MPS状态位置、SOC来自预配置结果或上一轮执行结果 mps_state.position mps.init_pos; mps_state.soc mps.init_soc; for t 1:horizon-pred_len1 % 更新当前故障信息哪些线路断开、哪些节点失电 current_status updateNetworkStatus(t); % 更新负荷预测曲线 load_forecast updateLoadForecast(t, pred_len); % 求解这一轮的动态调度优化问题 decision solveMPSDynamicScheduling(net, loads, mps, ... mps_state, current_status, load_forecast, dt); % 执行当前时段的决策 executeDispatching(decision(:, 1)); % 只看第一个时段 mps_state updateMPSState(mps_state, decision(:, 1)); end这里面的工程细节在于预测时域内的故障状态和负荷曲线是“预测”的而当前时段的故障状态是“实测”的。求解时当前时段用实测值、未来时段用预测值是一种比较务实的处理方式。我第一次做的时候把未来时段的故障当作全未知处理结果调度方案非常短视MPS总是原地不动后来改成已知预测值性能才明显提升。5. 复现过程中踩过的坑常见的求解失败与逻辑陷阱5.1 求解器报“infeasible”不可行这是我遇到最多的问题。明明MPS数量和容量都足够但模型说无解。排查后发现大部分情况是约束之间互相矛盾。典型矛盾场景一MPS初始位置与移动时间窗约束冲突。比如MPS初始在节点12但第一时段要求它必须出现在节点25而这两个节点之间的路网距离需要2小时、时段时长只有1小时模型当然无解。解决方法是把初始位置作为变量传入约束而不是固定值或者放宽第一时段的移动约束允许MPS从初始位置出发后在到达目标节点前处于“移动中”状态。典型矛盾场景二SOC约束与负荷需求的冲突。当MPS的剩余电量不足以支撑某个节点的关键负荷时模型会选择切负荷而不是无解但如果我把SOC下限设得过高比如0.3MPS可用的能量太少导致所有可行方案都满足不了负荷需求。这时候需要把SOC下限调低到0.1或者允许MPS在指定节点充电。提示遇到infeasible时不要一上来就改约束先把YALMIP的check结果打印出来看看哪条约束违反最严重。配合assign一个已知可行解用check逐条检查能快速定位问题。5.2 数值病态Big-M取值不当移动时间约束其实是一种Big-M约束还有其他类似 (x \cdot P) 的约束也用到Big-M。Gurobi对Big-M数值特别敏感。取值过大会导致求解精度下降出现“决策变量明明为1结果却违反约束”的诡异现象取值过小则可能会剪掉可行解。我的处理经验是Big-M尽量取紧。比如MPS最大放电功率已知为150kW那约束里用到Big-M的位置就取2001.3倍而不是图省事取10000。位置变量的Big-M取候选节点数量和时段数的上界再加余量。这个看起来是小细节但对求解稳定性影响非常明显。5.3 拓扑与供电逻辑陷阱还有一个特别隐蔽的逻辑陷阱MPS接入某个节点后它到底能给哪些负荷供电这取决于配电网当前的拓扑结构。如果MPS接入的节点与某个关键负荷节点之间的线路处于断开状态MPS出力再大负荷也恢复不了。而论文中的模型通常通过二进制变量描述“某节点是否恢复”内部隐含了电源与负荷必须处于同一连通区域这一条件。我在复现的时候曾经因为简化这个问题把MPS出力直接加到节点功率平衡里完全不考虑连通行导致结果中MPS在孤岛里“发电”负荷侧却毫无改善目标函数值虚高。后来老老实实把连通性约束加进去计算量变大了但结果才真正对得上论文的图和表。5.4 求解时间爆炸如果MPS数量多、候选节点全放开、预测时域长MILP规模会迅速爆炸。我的实测数据可以给个参考场景参数0-1变量数Gurobi求解时间MIPGap0.0133节点3台MPS预测时域6h59415~40秒33节点5台MPS预测时域8h13202~5分钟118节点5台MPS预测时域6h354010分钟以上需加有效约束118节点5台MPS候选节点限制为30个9001~3分钟我建议当候选节点数量明显大于MPS数量时必须引入候选节点选择逻辑只把预配置阶段选出的高潜力节点和当前失电节点附近的关键恢复节点作为候选其余节点不许接入。这样的话0-1变量数量会急剧下降求解速度提升非常显著而调度效果几乎不受影响。6. 可视化与结果分析怎么判断调度效果好坏算完调度方案不能只看一个目标函数数值还得从结果里看出门道。我会把MPS的位置轨迹和负荷恢复曲线画出来检查几点MPS是否出现了频繁无意义的移动比如来回折返、关键负荷是否在灾害初期就得到优先恢复、MPS电量是否均匀消耗、切负荷总量是否随着时间单调下降。%% 绘制MPS位置轨迹 figure; plot_traj squeeze(x_val(:,1,:)); hold on; for m 1:M pos_seq zeros(T,1); for t 1:T [~, pos_seq(t)] find(squeeze(x_val(t,m,:))); % 提取每个时段位置 end plot(1:T, pos_seq, -o, LineWidth, 1.5); end xlabel(时段); ylabel(接入节点编号); legend(MPS1,MPS2,MPS3); grid on;还有一个非常实用的分析角度对比有动态调度和无动态调度即MPS固定不动的负荷恢复曲线。用上一轮预配置的结果做初始值再固定位置不移动重新求解一遍得到的就是“静态方案”的恢复效果。把两条曲线画在同一张图上动态调度的收益一目了然。我在复现过程中会把这种对比图放进论文中作为佐证审稿人非常吃这一套因为直接量化了动态调度相对静态调度的增益。7. 最后一个实操建议参数敏感性分析怎么做模型本身跑通之后如果你想用它支撑一篇论文或者一个实际工程方案千万不要停在“能出结果”这一步。建议做几组敏感性分析一是MPS数量对恢复效果的影响通常曲线会先快速上升后饱和帮你确定配置多少台MPS最经济二是预测时域长度对调度效果的影响这个参数跟MPS移动速度强相关移动速度越快、时域可以越短三是SOC下限对方案可行性的影响调高SOC下限意味着给系统留更多的备用但会牺牲当前时段的供电能力两者之间有个权衡。我在这几组敏感性分析里发现一个有意思的规律MPS动态调度的收益在灾害初期最大。原因很简单故障刚发生时故障信息不完整预配置的位置很可能不再合适这时候动态调整的边际收益极高随着抢修推进、网络逐步恢复正常MPS的位置调整空间变小动态调度的收益也自然下降。这个结论反过来指导实际应用移动电源车的通信和远程调度系统一定要在灾害发生后的第一个小时内保持高可用别等到灾后三小时再搭调度链路那时候移动电源车的价值已经衰减一半了。复现这套论文代码的过程最值钱的不是跑通的那一刻而是中间反复理解模型假设、调试约束、跟论文图表核对结果的过程。每一步都在加深对配电网韧性调度的理解尤其是“信息—决策—执行”这条闭环逻辑放到真实的应急抢修场景里完全成立。如果你也在复现这篇或者类似方向的论文建议先把自己手里的Matlab环境梳理干净再把预配置模块和本文的动态调度模块按接口对接剩下的就是耐心调参了。