★ FEATURED ARTICLE 上篇把MPS预配置模型聊完之后很多做配电网韧性方向的朋友都来问同一个问题预配置只是决定了电源车的“初始车位”灾害过程中这辆移动电源车到底该怎么跑什么时候去哪个节点、接入之后供多久、何时撤出来奔向下一个关键负荷这些就是今天要聊的MPS动态调度。把这篇论文的算例复现完我最大的感触是动态调度真正的难点不在“优化”本身而在“时空耦合”四个字。电网状态随时间变、交通路网约束电源车行动、电源车容量又限制供电时长三者拧在一起之后模型规模一下子就上来了。这篇文章我把整个复现过程做一次完整梳理先从数学建模讲清楚每个约束的物理含义再给出MatlabYALMIPGurobi的落地写法最后把我调试时踩过的几个坑和排查思路一并分享出来。如果你是准备复现SCI一区电力调度论文、或者在做配电网韧性相关研究这份笔记应该能帮你少走不少弯路。1. 问题背景与建模动机1.1 为什么韧性课题里一定要有“动态调度”先把这个概念捋清楚。配电网韧性Resilience关注的是系统在极端灾害下的表现灾害来了电网扛不扛得住、受损之后能不能快速恢复关键供电。传统可靠性Reliability讨论的是常态运行下的供电连续两者不在一个维度上。打个比方可靠性是“平时身体好”韧性是“大病一场后能不能挺过来、恢复快不快”。极端灾害下配电网可能出现多条线路同时故障上级变电站也可能失压局部区域变成一个个“电力孤岛”。这时候想让关键负荷恢复供电一个快速有效的办法就是投入移动电源。但是移动电源车数量有限、容量有限、行动还受路网约束所以必须靠优化来决定怎么用。MPS动态调度解决的正是灾中阶段的问题已知每台移动电源车的初始位置预配置阶段的输出、路网通行时间、电网故障场景和负荷预测曲线优化给出电源车在整个灾害时段内的移动路径、接入节点、各时段出力计划目标是让加权恢复供电量最大。这里值得强调一下“动态”二字。灾害过程不是静止的负荷曲线在变、故障线路可能修复、电源车容量在消耗。如果只是灾前定一个接入节点就把电源车固定在那里本质上就是静态调度根本没有发挥移动电源“移动”的价值。动态调度的意义在于让有限的电源车跟着关键负荷和系统状态走一个节点救完之后再去救下一个。1.2 与预配置阶段的分工与衔接整个项目分上下两篇上篇做预配置下篇做动态调度两者不是孤立关系。预配置阶段的输出是每台MPS的初始驻留节点。这个节点通常选取在那些“灾后最可能成为重要供电枢纽”的位置同时要兼顾移动电源车从驻地出发的可达性。到了动态调度阶段这些初始位置就作为 t0 时刻的边界条件输入模型。换句话说预配置决定了故事的起点动态调度决定了后续整个旅程。我复现时在数据准备里做了一个衔接处理把预配置算出的最优车位直接写成一个向量 init_pos作为动态调度模型的初始状态。这个衔接看起来简单实际很关键。很多初学者在复现时喜欢把两个阶段完全分开结果预配置给出的初始位置在动态阶段根本不通比如到关键负荷的通行时间过长最后算出来的“最优解”完全不符合实际。所以复现时一定要把两个模块放在同一个数据环境下跑网格拓扑、路网矩阵、故障场景都要对齐。2. 数学模型拆解从物理过程到约束方程2.1 目标函数恢复多少才算“最优”动态调度的目标函数在论文里通常写成最小化未恢复的加权负荷或者等价地最大化加权恢复负荷。两种写法在求解上没有本质区别我用的是最大化恢复加权负荷max ∑_{t∈T} ∑_{i∈N} α(i,t) · w_i · D(i,t)其中 α(i,t) 是节点 i 在 t 时段的恢复比例取0到1之间的连续值w_i 是节点负荷权重医院、通信、供水这些关键负荷权重高D(i,t) 是该时段节点预测负荷。为什么要用“加权”而不是单纯加总因为极端灾害下电源容量远小于需求必须保重点。医院断电和普通居民楼断电的损失差一个数量级如果目标函数不给权重差异优化结果很可能把电源车派到一个负荷大但不重要的节点这在应急场景里是不能接受的。实际复现时权重不是随便拍的可以参考论文里对负荷等级的划分一级负荷权重10、二级负荷5、三级负荷1也可以根据实际场景调整。这里有一个容易忽略的细节α(i,t) 如果是连续变量模型允许节点部分恢复比如恢复了50%的负荷这在低压配网中其实不太符合实际。很多节点的负荷是一个整体要么恢复要么不恢复。因此我在复现时对关键负荷节点采用了二进制恢复变量对非关键负荷用连续恢复比例两者混用。这种处理既能保持模型可解性又更贴近应急场景下的实际调度逻辑。2.2 时空耦合约束MPS的移动轨迹怎么建模这是整个模型中最核心、也最容易写错的部分。动态调度需要一个“时空网络”的概念把电网节点沿时间轴复制T层形成一张“时间-空间”扩展图。MPS在这张图上的路径就是它的完整调度轨迹。如果 t 时刻位于节点 it1 时刻要么仍在 i继续发电要么已经在前往其他节点的路上要么已经到达新节点开始供电。我采用的变量定义如下z(i,t,k)0-1变量表示第 k 台MPS在 t 时段是否位于节点 ix(i,j,t,k)0-1变量表示第 k 台MPS在 t 时段开始从节点 i 前往节点 jTT(i,j)路网最短通行时间从电网节点 i 对应位置到节点 j 对应位置单位取调度时段整数倍第一组约束是位置唯一性与移动唯一性∑_i z(i,t,k) ≤ 1每个时段MPS至多在一个节点允许在路上 ∑_{i,j} x(i,j,t,k) ≤ 1每个时段至多发起一次移动第二组是流守恒约束也就是MPS位置随时间的演变关系z(i,t,k) - z(i,t-1,k) -∑_j x(i,j,t,k) ∑_j x(j,i,t-TT(j,i),k)这个式子的物理含义是如果 t-1 时刻在第 i 节点要么 t 时刻还在这里等式右边两个求和都为零要么本时段从 i 出发离开出现 −1同时如果之前某个时段从 j 出发来 i、旅途时间刚好等于 TT(j,i)那么 t 时刻应该到达 i出现 1。第三组是移动前提约束x(i,j,t,k) ≤ z(i,t-1,k)这条保证MPS不可能从一个它根本不在的节点出发。这套约束写起来直觉上很顺但实际编程时有几个边界坑t-TT(j,i) 可能小于1要对索引做保护另外 z 变量在“途中”时段哪都不在∑z ≤ 1 而不是 1我见过不少人写成 1然后把“在路上”的时段强行分配一个位置结果移动时间约束直接崩溃。这个细节在后面代码部分会再强调。2.3 电网运行约束孤岛场景下的潮流与二阶锥MPS接入配电网节点后会和该节点附近未故障的区域形成一个或多个孤岛。孤岛内的潮流需要满足物理约束这里用的是DistFlow支路潮流模型。对每条支路 (i,j)记 P(i,j,t)、Q(i,j,t) 为有功、无功潮流U(i,t) 为节点电压幅值的平方I(i,j,t) 为支路电流幅值的平方。DistFlow的约束包括节点功率平衡注入功率 上游流入 负荷 下游流出 电压降落方程U(i,t)-U(j,t) 2(rPxQ) - (r²x²)I 二阶锥条件‖2P; 2Q; U_i-zI‖₂ ≤ U_izIZ 是线路阻抗模的平方相关量这条二阶锥约束是模型能否高效求解的关键。DistFlow原始方程是非凸的直接求解很费劲通过引入 UV²、I|I|² 的变量替换并把等号松弛为锥条件问题就变成了一个可高效求解的混合整数二阶锥规划。更重要的是在辐射状配电网中这个松弛在很多情况下是精确的不会因为松弛而丢失可行解的物理意义。这个性质是论文选SOCP松弛而不是直接上交流潮流迭代的本质原因。复现时我遇到一个和拓扑相关的问题故障线路是时变的。灾害初期损坏的线路可能在中后期修复所以每条线路的投运状态需要做成时变参数。具体做法是给每条线路单独设一个二进制参数 line_status(l,t)故障时段为0修复时段为1。在潮流约束里支路潮流和电压降落约束都要乘以这个状态参数故障线路相当于从网络中摘除。这个处理听起来简单但做不好会出很多莫名其妙的不可行解尤其是支路功率没被强制清零时潮流会“穿”过故障线路。2.4 MPS运行约束容量、出力与能量连续性MPS本质上是带容量限制的移动电源运行约束主要有三块。第一块是接入出力上限。第 k 台MPS在节点 i 的注入功率 P_inj(i,t,k)不能超过其额定出力 P_max(k)且必须满足“只有接入才能出力”的耦合关系0 ≤ P_inj(i,t,k) ≤ z(i,t,k) · P_max(k)这条约束在代码里用大M或直接乘积实现YALMIP里可以直接写 P_inj z * P_max。第二块是能量连续性。MPS在 t 时段的剩余能量 E(t) 等于上一时段剩余能量减去本时段放电量E(t,k) E(t-1,k) - ∑_i P_inj(i,t,k) · Δt这里 Δt 是单个调度时段的时长单位要和P_inj对应起来。比如P_inj单位是kW、Δt是小时那么能量单位就是kWh。如果不考虑充电假设电源车在灾害期间没有充电条件约束就是能量单向递减。第三块是能量下限。为了避免把电池完全放空通常要留一定裕度E(t,k) ≥ SOC_min · E_cap(k)这个下限的设置直接影响调度结果下限设高了MPS不敢长时间供电会频繁转移设低了虽然单节点恢复时间长但后期可能无电可用。论文里那个30%的SOC下限是有讲究的复现时不要随意调整成一个很低的数否则最优解会变得非常激进现场调度时根本执行不了。3. 求解方案与工具选型3.1 为什么选MISOCP而不是启发式算法整个优化模型包含二进制变量位置、移动、恢复状态和二阶锥潮流约束数学上属于混合整数二阶锥规划也就是MISOCP。我见过不少人在复现这类论文时第一反应是写遗传算法或粒子群。不是说启发式不能解但在论文复现这个场景下用精确求解器有三个明显优势第一能拿到带最优性间隙的全局解和论文结果对比时有底气第二MISOCP模型本身是凸松弛之后的结构求解器可以利用分支定界高效搜索中小规模算例的计算时间完全可以接受第三换场景换参数时不需要重新调试算法超参数省去大量玄学调参。以IEEE 33节点系统、3台MPS、24个调度时段、5条故障线路为例我实测下来YALMIPGurobi求解时间大约在1到5分钟量级这个速度对离线灾中调度决策来说完全够用。如果算例规模再大一倍可以设置MIPGap提前终止获取近似最优解依然比大多数启发式算法可控。3.2 求解器配置与参数选择我用的求解链路是Matlab YALMIP Gurobi。YALMIP负责建模把优化问题写成约束和目标函数Gurobi负责底层求解。Gurobi里我推荐设置三组参数ops sdpsettings(solver, gurobi, ... verbose, 2, ... gurobi.MIPGap, 0.01, ... gurobi.TimeLimit, 600, ... gurobi.NumericFocus, 1);MIPGap设成0.01意味着允许1%的次优性求解速度会明显快于默认的1e-4。对动态调度这种场景1%的误差完全不影响决策质量。TimeLimit设上界防止极端算例把时间耗完。NumericFocus设成1缓解由于Big-M取值偏大导致的数值问题。这里提醒一下YALMIP版本和Gurobi版本有兼容性问题。我一开始用较新版本的Gurobi搭配老版本YALMIP结果YALMIP识别不了求解器折腾了很久。建议直接用官方文档里标注兼容的版本组合不要两个都是“最新版”最新版之间往往接口对不上。4. 关键代码实现与复现笔记4.1 数据准备从电网拓扑到路网通行时间矩阵动态调度需要两类基础数据电网拓扑数据和路网移动时间数据。电网拓扑我用的是IEEE 33节点测试系统节点数不多但足以验证模型正确性。基础数据包括节点坐标、支路电阻电抗、额定负荷、关键负荷权重。在Matlab里我习惯把支路数据整理成两个数组branch_from 和 branch_to分别记录每条支路的首端和末端节点编号。路网通行时间矩阵TT是动态调度的灵魂。如果研究场景里没有现成的路网数据可以按“欧氏距离/移动速度”估算再加一个系数模拟道路曲折。比如节点i和节点j的坐标距离是2公里电源车平均行驶速度40公里/小时曲折系数1.5那么通行时间就是 2×1.5/40×604.5分钟取整到小时时段后就是1个时段。% 路网时间矩阵生成简化为基于坐标的估算 dist_mat pdist2(node_coord, node_coord); TT ceil(dist_mat / (speed_kmh / 60) * detour_factor); % 单位分钟 TT max(1, round(TT / time_step_min)); % 对齐到调度时段 % 场景参数 T 24; % 调度时段数 nmps 3; % MPS数量 init_pos [6, 14, 25]; % 预配置阶段给出的初始位置 E_cap [500, 500, 800]; % kWh容量 P_max [100, 100, 150]; % kW额定出力注意TT矩阵的对角线是0自己到自己通行时间为0不可达的节点对可以设置一个很大的数但在建模时建议直接把这些“路径”从移动变量里排除而不是给一个大成本值否则求解器会在无意义的道路上浪费分支。4.2 变量定义与核心约束的YALMIP写法变量定义是整个建模过程里最需要条理的一步。我把变量分成三组位置与移动变量、潮流变量、MPS出力变量。% 位置与移动变量 z binvar(nb, T, nmps, full); % 位置指示 x binvar(nb, nb, T, nmps, full); % 移动决策出发 % 电网运行变量 P_line sdpvar(nl, T, full); % 支路有功 Q_line sdpvar(nl, T, full); % 支路无功 U sdpvar(nb, T, full); % 节点电压平方 I_line sdpvar(nl, T, full); % 支路电流平方 % MPS注入功率 P_inj sdpvar(nb, T, nmps, full); E sdpvar(nmps, T, full); % 剩余能量约束拼装的核心是前面第2节讲的那几组。这里展示两个关键代码片段。第一段是位置流守恒约束。注意我用了一个辅助函数来对到达时刻做边界保护Constraints []; for k 1:nmps % 初始位置 Constraints [Constraints, z(init_pos(k), 1, k) 1]; % 位置唯一性 for t 1:T Constraints [Constraints, sum(z(:, t, k)) 1]; end % 流守恒 for t 2:T for i 1:nb depart sum(x(i, :, t, k)); arrive 0; for j 1:nb tt TT(j, i); if (t - tt) 1 arrive arrive x(j, i, t - tt, k); end end Constraints [Constraints, ... z(i, t, k) - z(i, t - 1, k) -depart arrive]; end end % 移动只能从当前所在节点出发 for t 2:T for i 1:nb for j 1:nb if i ~ j Constraints [Constraints, x(i, j, t, k) z(i, t - 1, k)]; end end end end end这段代码我实际跑的时候犯过一个隐蔽的错把 x(i,j,t) 的出发时刻和 z 的时段对齐搞错了。一开始我把 x 定义成“t 时段结束时出发”结果流守恒里 depart 和 arrive 全部错位轨迹乱成一团MPS甚至能同时出现在两个位置。后来用一个极简的三节点算例逐步打印每个时段的 z 值才把这个语义对齐的问题揪出来。建议所有复现动态调度的朋友先拿一个 3 节点、5 时段的小算例验证轨迹方程确认符合直觉之后再上完整算例这一步能省下大量debug时间。第二段是DistFlow和二阶锥约束。在线路故障时段直接强制支路功率和电流为零for t 1:T for l 1:nl i branch_from(l); j branch_to(l); % DistFlow 电压降落 Constraints [Constraints, ... U(i, t) - U(j, t) ... 2 * (r(l) * P_line(l, t) x(l) * Q_line(l, t)) - ... (r(l)^2 x(l)^2) * I_line(l, t)]; % 二阶锥松弛 Constraints [Constraints, ... norm([2 * P_line(l, t); 2 * Q_line(l, t); ... U(i, t) - U(j, t)]) U(i, t) U(j, t)]; % 故障线路强制开断 if line_status(l, t) 0 Constraints [Constraints, ... P_line(l, t) 0, Q_line(l, t) 0, I_line(l, t) 0]; end end end这里有个性能优化点当 line_status 全为0时支路潮流变量仍然存在只是被约束为0。这会产生大量冗余的等式约束拖慢求解。更高效的做法是直接收缩变量只建立“有效线路有效时段”的潮流变量不过那会让索引管理复杂很多。我实际测下来33节点规模下冗余变量影响不大不需要做这个优化但如果扩展到上百节点量级建议做收缩处理。4.3 目标函数与求解调用目标函数我按最小化未恢复加权负荷来写load_shed sdpvar(nb, T, full); % 各节点各时段未恢复比例 alpha 1 - load_shed; % 恢复比例 Objective sum(sum(W .* (load_shed .* D_pred))); % 节点功率平衡电源出力 上级注入 负荷恢复量 for t 1:T for i 1:nb Constraints [Constraints, ... sum(P_line(branch_to i, t)) - sum(P_line(branch_from i, t)) ... sum(P_inj(i, t, :)) ... alpha(i, t) * D_pred(i, t)]; end end % MPS容量约束 for k 1:nmps for t 2:T Constraints [Constraints, ... E(k, t) E(k, t - 1) - sum(P_inj(:, t, k)) * dt]; Constraints [Constraints, E(k, t) 0.3 * E_cap(k)]; end Constraints [Constraints, E(k, 1) E_cap(k)]; end % 求解 ops sdpsettings(solver, gurobi, verbose, 1); sol optimize(Constraints, Objective, ops);节点功率平衡这一条是最容易和DistFlow重复的。如果前面已经写了完整的支路潮流约束这里的节点平衡本质上由DistFlow的第一个方程隐含了。但我在建模时习惯显式写出来宁可冗余不要遗漏。对YALMIP来说冗余等式通常不会带来严重性能问题。5. 调参与避坑实录5.1 Big-M取值解决了模型一大半的数值问题模型里有两个地方隐式用了大M思想一个是 P_inj ≤ z * P_max一个是节点注入功率与恢复状态之间的耦合。YALMIP会自动引入Big-M但Big-M的具体值取决于变量的上下界而不是用户显式设置。电压幅值平方 U 的处理要特别小心。如果系统基准电压是12.66kV有名值下的 U 可以高达1.6e8和1这种小量级的功率值放在同一个约束里数值尺度严重失衡很容易把求解器的数值精度拖垮。解法是全部采用标幺值电压 U 在0.9到1.1之间功率在几十到几百之间二阶梯锥的三个分量数量级接近求解器处理起来非常顺。我一开始没有用标幺值结果Gurobi报Numerical Error模型偶发不可行。把电压、电流、功率全部转成标幺值之后同样的约束逻辑一次就跑通了。5.2 初始位置不合法最常见的不可行来源动态调度的初始位置来自预配置阶段但预配置通常不保证这个位置在所有场景下都“走得通”。假设预配置给出的初始位置在第6节点但灾害发生后从第6节点到任何关键负荷节点的通行时间都超过调度时段总数的一半那这台MPS的前半程基本废掉了。排查这种问题时我建议先把移动约束和电网约束分开验证把电网潮流约束全部注释掉只保留MPS的移动约束和目标函数比如最大化总移动距离看有没有可行解。如果移动约束单独跑都可能不可行那基本可以断定是初始位置或路网时间矩阵的问题。如果移动约束可行但加上潮流约束后不可行那问题出在电网侧。5.3 同质MPS对称性爆炸如果两台MPS的容量和出力上限完全一样它们在模型里就是“同质”的。同质变量会让分支定界算法在对称解空间里反复打转白白浪费时间。3台同质MPS的对称性会让分支树膨胀好几倍。破解方法很简单给MPS编号加对称性破缺约束。比如强制第一台MPS的初始位置编号小于第二台第二台小于第三台或者按“同一时段位于较小节点编号的MPS优先”来写一个偏序约束。加了这条约束之后求解时间经常能缩短30%以上。5.4 可移植结构用场景数据表驱动模型动态调度模型涉及的场景很多不同故障线路组合、不同负荷水平、不同MPS数量。我在刚开始复现时把场景写死在代码里每次换场景都要改好几处数据非常容易出错。后来痛定思痛把所有场景参数集中到一个结构体里用数据表驱动整个模型scenario.nb 33; scenario.nl 32; scenario.branch_from [...]; scenario.branch_to [...]; scenario.r [...]; scenario.x [...]; scenario.T 24; scenario.TT TT; scenario.init_pos [6, 14, 25]; scenario.line_status line_status; % nb x T 或者 nl x T scenario.D_pred D_pred; % nb x T scenario.W W; % nb 权重向量这样在批量跑多个故障场景时只需要循环更换 scenario 的部分字段模型代码完全不用动。这是整个复现工程里性价比最高的一次重构。6. 算例效果与结果分析6.1 算例设置我用IEEE 33节点系统做了验证。设置5条线路同时在初期故障故障修复时间分别为8、12、16、20、24时段。3台MPS的初始位置分别在第6、14、25节点容量为500kWh/500kWh/800kWh额定出力100kW/100kW/150kW。调度周期24小时时间步长1小时。关键负荷集中在医院节点8、通信基站节点12、水厂节点20、社区应急中心节点29。6.2 动态调度 vs 静态接入对比为了验证动态调度的价值我设置了两个对照组。对照组A是无移动电源接入对照组B是静态接入也就是电源车在预配置位置原地发电、不移动只通过配电网络向周边供电。场景加权负荷恢复量关键负荷恢复率平均供电时间无MPS000静态接入2140 kWh52%5.3 h动态调度2875 kWh78%8.1 h动态调度相比静态接入加权恢复量提高了约34%。这个增益主要来自两方面第一个增益来源是MPS的逐点救援。静态接入的MPS只能给初始位置周边的负荷供电距离远的关键负荷即使权重再高也鞭长莫及。动态调度允许MPS在完成一个节点的主要供电任务后转移到下一个高权重节点把每度电都用在最需要的地方。第二个增益来源估计很多人会忽略动态调度带来的负荷恢复“时间前移”。由于MPS可以主动向着关键负荷移动接入节点的选择不再受限于初始位置关键负荷往往能提前几个时段恢复供电。在应急场景里早一小时恢复医院供电价值远不止多出来的那几度电。6.3 MPS轨迹的两点观察观察MPS的时序轨迹有几个有趣的规律。第一台MPS初始位置6在前8个时段围绕节点6、8区域供电服务医院负荷第9时段开始向节点12转移因为此时通信负荷权重在灾后上升。第二台MPS初始位置14前半程在节点20水厂附近恢复供电后半程转移至节点29应急中心。第三台大容量MPS初始位置25几乎没有移动因为它所在片区本身有高权重负荷且容量大原地供电比移动更划算。这个结果说明一个问题动态调度不意味着每台MPS都要频繁移动移动是有代价的移动期间不发电只有移动带来的负荷恢复增益大于“途中停电”损失才值得动。7. 从复现到工程落地的几点心得这篇代码跑通之后我尝试往更实际的方向扩展发现几个值得深入的点。第一是把不确定性纳入调度。灾害过程中线路故障不是完全已知的实际调度往往采用“滚动优化”的思想每2到4小时根据最新灾情信息重新求解一次动态调度更新MPS的后续轨迹。我在复现时尝试过在Matlab里套一层循环每个滚动周期内调用一次优化模型效果比一次性的全局调度更贴近实际不过计算开销会成倍增加。第二是多类型移动电源的混合调度。论文标题里说的是“应急移动电源”实际项目里往往同时存在移动储能车和移动柴油发电车。两者的特性差异很大储能车响应快但容量有限柴油车能持续供电但燃料供应受路网影响。把两者放在同一个优化框架里约束会多一组燃料补给的逻辑问题规模又上一个台阶但这才是工程实际。第三是MPS的就地充电问题。灾害中后期如果配电网部分区域恢复正常MPS可以接入电网补充电量形成“电网给MPS充电、MPS再机动转移”的模式。这个逻辑在模型里就是给能量连续性约束增加一个充电项求解框架不需要大改但调度结果会发生质变——MPS不再是单次放电的消耗品而变成可持续参与恢复的机动力量。从个人实际复现经验看这套代码最大的价值不在于精确复现某篇论文的某个算例而在于把“配电网韧性提升”这个宏大概念落成了一个可以反复实验、换场景、加约束的量化工具。故障场景一换、负荷权重一调、MPS参数一改马上就能看到调度策略和恢复指标的联动变化。这种把抽象概念变成可计算模型的能力才是复现一篇论文最值得带走的东西。
阅读完成 · 觉得有帮助？
