做配电网韧性提升方向的复现,最常被问到的就是:应急移动电源(MPS)到底该怎么动态调度?这个问题看起来简单,真上手做Matlab实现时,你会发现它同时牵扯潮流、时序、移动路径和抢修进度四条线,哪一条没处理好,结果都会被带偏。
先说个实际场景。台风过境后的凌晨,调度大屏上还有七八条线路处于故障状态,库房里的应急移动电源车早就到位了,但停在哪、什么时候出发、先去哪个负荷点,都不是一眼能看出来的。如果只做"预配置",那就相当于把移动电源当成固定储能,放在一个猜好的位置等着——这在灾后修复时间不确定的现实里几乎不可用。这篇博文专门聊动态调度部分,也就是MPS在灾中如何根据实时故障状态、抢修进度和自身电量,一步步滚动作出"下一站去哪、出多少功率"的决策。
文章里所有代码框架都以Matlab实现为主线,参考了某篇SCI一区论文中关于两阶段韧性能源调度的思路,但实现细节是我自己重构过的,算例用的是公开的33节点配电网测试系统。读完之后,你应该能直接照着搭出一套可运行的MPS动态调度仿真。
1. 从"放在哪"到"走去哪":动态调度到底在解什么问题
1.1 极端事件下的三个时间窗口
配电网韧性研究和传统可靠性研究最大的区别,在于它面对的是小概率、高影响的极端事件。台风、冰灾、地震这类事件会造成多条线路同时段故障,而且故障是陆续被发现、陆续被修复的。整个过程大致分成三段:灾前预配置阶段、灾中持续响应阶段、灾后逐步恢复阶段。
灾前阶段能做的事情有限,主要是根据气象预报和灾害模型,把移动电源、抢修队伍、应急物资预先布置到风险最高的区域。这也就是"上篇"讲的预配置问题:从哪里调移动电源、放到哪个变电站、放几台。到了灾中,情况开始变得不可控,实际受损位置和预测往往对不上,某条线路可能提前修复,也可能比预计多停一天。这时候,如果移动电源还傻乎乎地待在预配置的位置,很可能错过真正需要供电的负荷点。
所以动态调度的价值,是在灾害展开之后,根据"已经确认的故障信息 + 短期修复预测 + 自身剩余电量",不断重新规划移动电源接下来的行动轨迹。它把"静态放置"变成了"动态路由",这也是配电网韧性和传统配电网调度之间最本质的差异。
1.2 移动电源与固定储能的核心区别:时空强耦合
固定储能接入电网后,位置不变,优化的只是充放电时序。移动电源多了一个自由度——位置决策,但它不是白来的:
- 移动需要时间,在途中的时段无法提供供电;
- 转移需要消耗自身电量,赶到新地点之后SOC可能已经下降;
- 移动到某个节点之后,该节点能否形成有效供电孤岛,取决于周边配电网的拓扑连通性;
- MPS被一个节点"锁住"的时间越长,它支援其他节点的机会成本就越高。
我复现过程中最直观的感受是:移动电源调度本质上是把"空间路由问题"和"能量管理问题"叠在一起,再加上配电网潮流的约束。如果只用经验规则,比如"哪里负荷重就往哪走",大概率会在某个时段把MPS卡在一个既没负荷、又离故障点很远的中间位置。
1.3 为什么要动态调度,而不是一次性算出全局方案
理论上,如果灾后故障集合、修复时间全部已知,完全可以用一个全时段的大模型一次性算出MPS全生命周期的最优路径。但现实是,故障信息是逐步确认的,修复时间也经常变化。全时段模型一旦信息更新,就必须重新求解;而窗口一大,变量数和整数变量规模会急剧膨胀,常规商业求解器在合理时间内根本解不动。
滚动时域(rolling horizon)是对这个问题最自然的妥协:每一步只求解未来若干个时段的小规模子问题,执行当前时段的决策,然后滚动推进。论文里面它经常被包装成MPC(模型预测控制)或者滚动优化,底层逻辑都一样。我的实现也采用这个框架,代码结构非常清晰:更新状态、构建窗口模型、求解、执行、再更新。
2. 动态调度的数学内核:时序耦合与移动约束建模
2.1 目标函数:加权失负荷量是最直观的韧性度量
我在这套实现里采用的目标函数是"最小化加权失负荷量"。写成文字化的公式就是:
min ∑(t=1..T) ∑(i=1..N) ω_i · P_i^L · (1 - u_i,t) · Δt其中:ω_i是节点i的负荷权重,P_i^L是该节点的有功负荷需求,u_i,t表示该节点在t时段是否恢复供电,Δt是时段长度(我取1小时)。
权重怎么设非常关键。一级负荷比如医院、通信基站、供水系统的权重可以给到50以上,普通居民负荷给5,一般工商业给1。权重设平了,MPS会变得"就近原则",哪里近就先去哪,失去了负荷重要性维度;权重差距拉太大,模型又会过度偏向少数几个节点,导致部分区域完全不管。
算例中我会固定三组对比场景:无MPS、静态预置、动态调度。目标函数值越小,说明加权失负荷量越小,韧性表现越好。
2.2 关键约束组:一个都不能少
光有目标函数还不够,动态调度的核心难度在约束。我复现时整理的约束组大致分五类:
- 潮流约束:采用DistFlow线性化。节点功率平衡、支路功率、电压递推关系都用线性等式约束描述,电压幅值限制在0.95到1.05 pu之间。
- MPS位置唯一性约束:每一个时段,MPS只能待在一个节点,所以∑ x_i,t = 1,x_i,t是二进制变量。
- 移动时间约束:MPS从节点i转移到节点j需要一定的旅行时间 τ_ij,这段旅行时间内它既不能在i出力,也不能在j出力。实际建模时,我会预先算好所有节点对之间的旅行时间,然后用一个"禁止快速跨越"的时序约束表达。
- SOC耦合约束:SOC(t+1) = SOC(t) − (p_t^dis / η_dis − p_t^ch · η_ch) · Δt,外加SOC上下限(一般0.1~0.95),充放电功率不能超过额定值,且同一时段不能同时充放电。
- 供电路径与恢复状态约束:某节点在t时段若由主网供电,前提是上游配电路径已修复;若由孤岛供电,前提是MPS当前接入该孤岛内的一个节点,且孤岛内部连通。
2.3 一个最小示例看懂约束耦合
拿一个三节点系统举例:节点1连变电站,节点2和节点3各带一个医院负荷。故障发生后,线路1-2断开,节点2、3全部失电。MPS在节点3,容量能支撑两个负荷半小时,但转移去节点2需要1小时。
这时候模型会怎么决策?如果只看当前时段,MPS留在节点3给医院3供电是"最优"的,因为目标函数立刻下降。但如果再往未来一个时段看,抢修队预计会在2小时之后修复线路1-2,那么节点3很快可以恢复主网供电,而节点2还需要再等更久。滚动窗口看到这一点后,就会倾向于让MPS尽早转移到节点2,哪怕当前时段节点3短暂失电也在所不惜。
这就是动态调度"牺牲局部、换取全局"的典型场景。如果只用单时段贪婪策略,这种跨时段的决策根本做不出来。
3. Matlab代码实现:滚动时域求解框架的完整链路
3.1 数据准备:33节点测试系统与故障场景生成
我用的是经典33节点配电网测试系统的公开数据,支路阻抗、节点负荷都是现成的。导入时为了避免脚本到处是魔法数,我封装了一个函数:
% load_case33.m function CASE = load_case33() CASE.N = 33; % 节点数 CASE.NB = 32; % 支路数 CASE.refbus = 1; % 变电站根节点 BASE.MVA = 10; BASE.KV = 12.66; BASE.ZBASE = (BASE.KV^2) / BASE.MVA; % 阻抗标幺化基准 % 支路数据:首端、末端、电阻(ohm)、电抗(ohm) br = [... 1 2 0.0922 0.0470; 2 3 0.4930 0.2511; 3 4 0.3660 0.1864; ... 31 32 0.3720 0.5740; 32 33 0.2090 0.3280]; CASE.br = br; CASE.r_pu = br(:,3) / BASE.ZBASE; CASE.x_pu = br(:,4) / BASE.ZBASE; % 节点有功、无功负荷(kW/kVar,转标幺) Pload = [... 0.1; 0.09; ... ]; % 实际按测试系统数据填 CASE.Pload = Pload / BASE.MVA; CASE.Qload = Qload / BASE.MVA; CASE.w = ones(CASE.N, 1); % 权重后续按负荷等级覆盖 end故障场景生成我建议做成独立的结构体,方便做多场景对比。包括:故障线路编号和发生时段、每条故障线路的修复时间分布、负荷波动系数。我的做法是写一个随机生成器,种子固定,保证每次复现结果一致:
% generate_scenario.m function SC = generate_scenario(rng_seed) rng(rng_seed); SC.fault_branch = [4; 9; 18]; % 三条线路同时故障 SC.fault_time = [1; 1; 1]; % 从第1时段开始 SC.repair_time = [5; 8; 11]; % 分别在第5、8、11时段修复 SC.load_factor = 0.85 + 0.3 * rand(CASE.N, 1); % 节点负荷波动 end3.2 YALMIP模型构建:变量定义与约束组装要点
Matlab端我采用的是 YALMIP 建模,求解器用商业MILP求解器。变量定义如下:
% build_rolling_model.m function [Constraints, Objective, vars] = build_rolling_model(CASE, MPS, SC, t_cur, H) N = CASE.N; T = H; u = binvar(N, T, 'full'); % 节点供电恢复状态 x_mps = binvar(N, T, 'full'); % MPS所在节点位置 p_dch = sdpvar(1, T, 'full'); % 放电功率 p_ch = sdpvar(1, T, 'full'); % 充电功率 soc = sdpvar(1, T+1, 'full'); % SOC轨迹 v = sdpvar(1, N, T); % 节点电压(可选向量化写法) Constraints = []; % 初始SOC由外部传入 Constraints = [Constraints, soc(1) == MPS.soc_init]; for t = 1:T Constraints = [Constraints, 0 <= p_dch(t) <= MPS.Pmax]; Constraints = [Constraints, 0 <= p_ch(t) <= MPS.Pmax]; Constraints = [Constraints, soc_min <= soc(t) <= soc_max]; end % SOC递推 for t = 1:T Constraints = [Constraints, ... soc(t+1) == soc(t) - (p_dch(t)/MPS.eta_dis - p_ch(t)*MPS.eta_ch) * CASE.dt]; end % 位置唯一性 Constraints = [Constraints, sum(x_mps, 1) == 1]; % 移动时间约束(示意版) % travel(i,j) 是MPS从i到j需要的时段数 for t = 1:T for i = 1:N for j = 1:N if travel(i,j) > 0 && t > travel(i,j) Constraints = [Constraints, ... x_mps(j,t) <= 1 - x_mps(i,t-travel(i,j))]; end end end end % 潮流约束、电压约束、恢复状态约束此处省略, % 完整代码里用 DistFlow 线性化逐时段的等式约束。 ... end这里有一个非常关键的点:移动时间约束的时序对齐。很多第一次复现的人会直接写成"当前位置变量在连续两个时段不能跨越",但忽略了一个事实——如果MPS从节点5移动到节点9需要2个时段,那么它在中途的第1个时段实际上出现在"系统的任何节点之外",可如果强制∑ x = 1,这个时段就不可能合法。我简化处理的办法是:引入一个虚拟节点表示"在途",MPS在途中时段就落在虚拟节点上,这样位置唯一性约束始终满足,同时虚拟节点不允许出力。
3.3 主循环求解流程
滚动时域的主循环是所有流程的中枢。我把它拆得很清楚:
% run_mps_dynamic_dispatch.m CASE = load_case33(); SC = generate_scenario(20240831); MPS = init_mps(); % capacity=1.0MWh, Pmax=0.5MW, soc0=1.0 H = 4; % 滚动预测窗长度 T = 24; % 仿真总时段 result.model_traj = struct(); result.soc_hist = zeros(1, T+1); result.served_hist = zeros(1, T); for t = 1:T % 1. 根据当前实际发生的事件更新场景信息 SC_now = update_scenario(SC, t); % 2. 生成未来H个时段的预测场景(可以带误差) SC_pred = forecast_scenario(SC_now, t, H); % 3. 构建并求解当前窗口的MILP [Constraints, Objective, vars] = build_rolling_model(... CASE, MPS, SC_pred, t, H); ops = sdpsettings('solver', 'gurobi', 'verbose', 0, 'mipgap', 1e-3); diagnostics = optimize(Constraints, Objective, ops); if diagnostics.problem ~= 0 warning('时段%d求解失败: %s', t, diagnostics.info); break; end % 4. 只执行当前时段决策 x_now = value(vars.x_mps(:,1)); p_now = value(vars.p_dch(1) - vars.p_ch(1)); % 5. 更新MPS真实状态 MPS = execute_mps_decision(MPS, x_now, p_now, CASE.dt); result.soc_hist(t+1) = MPS.soc; % 6. 记录结果 result.model_traj(t).x = x_now; result.model_traj(t).u = value(vars.u(:,1)); result.model_traj(t).p = p_now; end这套流程你跑通之后会发现:真正让代码变复杂的是第1步和第2步的状态更新,而不是优化模型本身。因为每个时段都要确认"哪些线路在这个时段突然恢复了、哪些负荷点需要重新从MPS转向主网、MPS当前到底在哪个节点",这些状态如果不同步,滚动优化就变成了"睁着眼睛跳崖"。
4. 动态调度与抢修进度的联动逻辑
4.1 抢修进度如何改变MPS下一步去向
灾后恢复过程中,MPS的动态调度和抢修队进度是强耦合的。线路修复后,原本靠MPS支撑的孤岛负荷会重新切回主网供电路径,这时候MPS就多了一个"可选行动":要么原地降低出力待命,要么立刻转移到另一个仍然失电的区域。
我在代码里用一个接口函数解决这个联动:
% update_scenario.m function SC_now = update_scenario(SC, t) SC_now = SC; % 根据修复时间表,更新当前已恢复线路集合 SC_now.repaired_lines = SC.repair_branch(SC.repair_time <= t); % 重新计算当前失电节点集合:基于已修复线路做一次尾部可达判断 SC_now.outage_nodes = compute_outage_nodes(SC_now); end这里最值得注意的细节是:节点是否恢复供电,不要只在模型里用二进制变量"硬推"。配电网是辐射状结构,某个节点能否恢复,取决于它到变电站的路径上所有故障线路是否都已经修复。我最初直接在YALMIP里用大量连通性变量建模,模型规模爆炸,求解器跑得非常吃力。后来改成在Matlab侧用可达集预计算,把每个时段的"可恢复节点集合"算出来,再以参数形式传入优化模型,速度和稳定性都好了很多。
4.2 线路修复与MPS转移时机:谁该等谁
抢修进度和MPS转移这两个动作,在时间上经常存在一个"岔路口":
- 如果MPS先转移到某孤岛,但故障线路迟迟修不好,MPS电量耗尽后就只能陪跑;
- 如果MPS在一个即将修复的区域死守,等线路一修好它就立刻想去下一站,转移时间又可能来不及赶上当晚的负荷高峰。
这种问题在纯优化模型里靠目标函数自动权衡,但工程上我会加一个额外约束:MPS在某个节点接入后,至少停留两个时段才能撤离。这个"最小停留时长"约束虽然会牺牲一点点目标函数值,但能避免滚动优化因为预测误差导致MPS频繁移动,实际运行中完全不可执行。
实测下来,加了这个约束之后,SOC曲线更平稳,MPS路径也更像一个"有经验的人"会做的选择:去一个地方,稳定供一段电,然后才转移。如果完全没有这个约束,MPS在极端情况下甚至会出现"每过一个时段就换一个节点"的抖动,目标函数数字倒是好看,工程上一看就不现实。
5. 结果可视化与韧性指标定量评价
5.1 三条曲线:负荷恢复率、MPS位置、SOC
复现论文代码时,可视化不是最后画个图那么简单,它其实是你校验模型行为是否合理的最直接手段。我每次跑完都强制自己看三张子图:
% plot_dispatch_result.m figure('Color','w','Position',[100 100 800 700]); subplot(3,1,1); stairs(0:T, [1, served_percent], 'LineWidth', 1.5); xlabel('时段'); ylabel('负荷恢复比例'); title('系统加权负荷恢复率时序'); ylim([0 1]); grid on; subplot(3,1,2); stairs(0:T, [MPS.soc0, result.soc_hist(2:end)], 'LineWidth', 1.5); xlabel('时段'); ylabel('SOC'); title('MPS 荷电状态轨迹'); ylim([0 1]); grid on; subplot(3,1,3); stem(0:T, [MPS.node0, node_hist], 'LineWidth', 1.2); xlabel('时段'); ylabel('节点编号'); title('MPS 位置轨迹'); grid on;负荷恢复率曲线看的是宏观韧性表现,SOC曲线看的是能量管理是否合理,MPS位置轨迹看的是空间决策是不是有解释性。我见过不少复现版本结果数字很漂亮,但画出来位置轨迹乱跳,说明约束漏了或者滚动窗口设置有误。可视化在复现过程中的作用,是作为"合理性检验器",而不是论文配图工具。
5.2 韧性提升怎么量化:曲线下面积对比
SCI一区论文里常见的韧性指标之一是系统恢复过程中的负荷时间曲线下面积,反映整个事件过程的总供电水平。我的实现中归一化处理方式是:
% compute_resilience_index.m function RI = compute_resilience_index(served_percent) RI = sum(served_percent) / length(served_percent); end也就是把曲线下面积除以总时段数,得到一个0到1之间的韧性指数。对比时固定同一故障场景跑三组实验:
| 方案 | 加权失负荷量 | 韧性指数 | 说明 |
|---|---|---|---|
| 无MPS | 基准值 | 较低 | 纯粹靠抢修恢复 |
| 静态预置MPS | 略有改善 | 中等 | MPS全程固定在预配置节点 |
| 动态调度MPS | 明显降低 | 最高 | MPS跟随抢修进度不断转移 |
从我的算例看,动态调度相对静态预置的加权失负荷量降幅一般在10%到20%之间,具体取决于故障点的空间分布和修复时间差。如果故障点非常集中,两者差距会被压缩;如果故障点分散且修复时间错开,动态调度的优势会非常显著。
6. 复现过程中踩过的四个坑与处理经验
6.1 大M参数过大导致的决策错位
很多经典教材喜欢用"大M法"把逻辑约束转成线性约束,但大M取了1e6这种量级,求解器数值上会非常难受。我遇到过一次很典型的怪象:模型求解结果里,MPS明明在某个时段报告"已接入节点3",但节点3的功率平衡等式却完全不包含MPS注入功率。排查到最后发现就是大M取值过大,导致二进制变量和连续功率之间的逻辑约束数值病态,求解器在默认数值容差下"认为"约束已经满足了。
我的处理原则:大M值永远按实际物理上界取,且不超过该上界的1.2倍。比如放电功率上界0.5MW,大M就取0.6;电压偏差上界0.05pu,大M就取0.06。这样既保留了逻辑约束的效力,又不给求解器埋数值炸弹。
6.2 初始SOC和预配置阶段脱节
这是做"上篇+下篇"两阶段复现时最容易犯的错误。如果上篇预配置阶段算完,MPS到达某个节点时已经用掉了30%的电量,下篇动态调度开始时就必须把它接过来,而不是默认满电。我写代码时专门给init_mps()加了一个soc_init字段,由预配置阶段的终态SOC传入。
如果忽略这一步,动态调度会把本不该存在的电量白白释放,第一时段的"最优决策"完全失真。这个坑的隐蔽性在于:单跑动态调度看不出问题,只有把两阶段串起来做整体仿真时,结果才对不上。
6.3 滚动窗口太短导致MPS反复横跳
我最早把滚动窗口H设成1,纯粹是"在线决策",结果MPS路径像喝醉了一样在两个失电区域之间来回跑。原因很直观:窗口太短时,模型看不到未来区域2的修复时间更长,只觉得"现在区域1有负荷、去区域1有收益",等下个时段又觉得"区域2也有负荷、返回去也有收益"。
把H扩到4个小时后,路径才稳定下来。这个经验也直接呼应了4.2节说到的"最小停留时间约束"——两件事配合起来,动态调度结果才具备工程可行性。
6.4 对比实验的求解精度不一致
做有MPS、无MPS、静态预置三种方案的对比时,如果每个场景都用默认的求解器收敛设置,某些场景会停在不同的MIP gap上,结果差异可能来自求解精度而不是调度策略本身。我把所有场景统一设置为mipgap=1e-3,同时加一个'maxtime', 300的兜底限制,任何一个场景超时都会在日志里明确标注,而不是默默输出一个次优解。
这样处理之后,至少保证对比结果在数值上是一致的。学术复现尤其要注意这个细节,否则审稿人或者同行拿你的代码重跑,很可能得到完全不同的趋势。
6.5 补充提醒:孤岛电压校验别只信线性化结果
最后再补一个不算坑但容易忽略的提醒。DistFlow线性化在小幅电压偏差下比较可靠,但MPS在某个孤岛内满载支撑时,节点电压可能逼近甚至低于0.95pu,这时候线性化结果会偏乐观。我建议在动态调度仿真结束后,对每个时段的断面对MPS所在孤岛做一次前推回代精确校验,看是否出现低电压。如果出现,就把该时段的MPS出力上限稍微调低,或者把电压下界在优化模型里收紧一点,例如设成0.96pu,留出线性化误差的余量。
整套代码跑下来,我对滚动时域调度感触最深的一句话是:优化模型给的是决策逻辑的上限,而真正让系统落地的,是状态更新、时序对齐和数值处理这些"看不见的细节"。你现在拿到的这个骨架,换掉数据加载函数、改一下MPS参数,就可以去接自己的算例;如果再往上加入抢修队路径协同优化,这套框架也完全撑得住。之后如果有时间,我打算把预配置阶段和抢修队模型的实现也像这样整理一遍,把两阶段韧性调度完整串起来。