“预配置是解决‘电站在哪里等着’动态调度解决‘电站怎么跑过去’。真正让一区论文落地后者才是硬骨头。”上篇我们把MPS应急移动电源的预配置模型讲完了在台风灾害来临前基于场景生成把有限数量的移动电源提前布置到配电网的关键节点本质上是一个两阶段随机规划的第一阶段决策。但如果你只做预配置文章最多发到二区。一区论文真正的戏肉在于第二阶段——灾害发生后MPS怎么根据故障演化的实际情况动态调整位置和出力。这篇我们专门拆MPS动态调度也就是论文(下)对应Matlab代码实现的核心部分。1. MPS动态调度要解决的核心问题先说清楚为什么“动态”这两个字值钱。配电网在极端灾害下的故障不是一次性全暴露的而是随时间逐步演化的。台风过境时某条馈线可能先跳闸半小时后另一条线路才因为树障短路再过一小时某个开关站才出现设备损毁。如果我们在预配置阶段把所有MPS一次性锁定在初始位置不动就会出现一个很尴尬的局面某个节点明明开始的故障很严重但修复到一半已经不再需要MPS支援而另一个节点故障才开始爆发却分不到资源。MPSMobile Power Source移动电源本身是两个特征的复合体它既是分布式电源有发电能力又是移动资源有运输属性。这两个特征叠加起来数学建模的复杂度不是加法是乘法。你既要刻画它接入电网后的潮流约束还要刻画它在路网中的移动路径和时序约束这两套系统还在时间维度上耦合——移动期间MPS无法发电接入之后才能出力而出力大小又和当前配电网的拓扑、负荷水平、故障位置强相关。动态调度的本质是一个多时段耦合的混合整数规划问题而且是在不确定性环境中做的决策。论文里通常的处理思路有两种一种是滚动时域控制RHC在每个决策时刻基于当前观测到的故障状态求解未来一段时间的最优调度另一种是全局优化假设整个故障演化过程是已知轨迹一次性求一个开环最优解。前者更接近实际工程后者数学上更干净、代码实现更简单也是大多数复现代码采用的方式。这篇帖子里我们先讲清全局优化版本的实现再讨论怎么改造成滚动时域。与预配置模型相比动态调度的变量维度和约束复杂度都上了一个台阶。预配置阶段每台MPS只有一个位置变量动态调度阶段每台MPS在每个时段的每个可选节点都可能有接入/断开动作同时还要处理移动路径上“从哪到哪、花了多久”的时空一致性问题。代码里最容易写崩的就是这里人们经常在写目标函数时只考虑运行成本结果解出一个MPS在几个孤岛间反复横跳的荒谬调度方案或者在写移动约束时忘了给移动过程加上时间门槛导致模型认为MPS可以瞬间位移。这些都是我们在复现过程中第一次跑代码时真的会遇到的坑。2. 动态调度的数学模型长什么样2.1 目标函数到底在优化什么SCI一区论文里动态调度的目标函数通常不是简单的成本最小化而是“韧性导向”的多目标。所谓韧性Resilience在电力系统语境里一般用“灾害期间系统损失的负荷电量”来度量同样一场台风损失的电量越少系统韧性越高。因此动态调度目标函数的第一项是最大化或最小化其负值网络总供电恢复量也可以写成每个时段每个节点负荷削减量的加权求和。第二项是MPS的移动成本。注意这里的移动成本不只是油费还包括一套时间代价系数。为什么要用代价系数而不是直接把移动时间写进约束因为从数学上说如果你只约束移动时间、不在目标里给移动惩罚求解器会倾向于把移动步骤安排到负荷损失已经发生之后的位置——尤其是当恢复收益在时间上分布不均时模型可能为了一个远期的小收益放弃当前时段的快速恢复这不符合调度直觉。典型的目标函数可以写成min Σ_t Σ_b [ λ_b · L_b^cut(t) · Δt ] Σ_m Σ_t [ c_m · y_m(t) π_m · travel_time(m,t) ]其中λ_b是节点b的单位失负荷价值系数L_b^cut是t时段节点b的切负荷量y_m(t)是MPS m在t时刻的移动状态指示变量π_m是移动时间惩罚系数。这个目标函数里隐含了“削峰填谷”的调度逻辑把失负荷当成一种“负收益”MPS的调度本质上是拿有限的移动资源去填失负荷的谷底。谁的λ大、谁当前失负荷最严重谁就有优先获得MPS支援的资格。2.2 约束条件里的时间-空间耦合动态调度的约束远比预配置复杂核心是三个第一配电网潮流约束和辐射状运行约束。这部分和静态优化一致采用DistFlow线性化模型P_j P_i - r_ij·I_ij^2线性化后I_ij^2用常数替代或分段线性加上每个节点的有功/无功平衡方程。需要注意的是故障线路在动态调度时段内必须从拓扑中摘除而且摘除后不能恢复——除非论文考虑了修复过程这时需要用额外的二进制变量表示“该线路是否已修复”。强烈建议第一版代码先忽略修复过程把故障线路全部作为“永久断线”处理跑通后再逐步加复杂度。第二MPS的时空一致性约束。这是动态调度模型最核心、也最容易写错的部分。每一台MPS在任意时刻只处于两种状态之一接入某节点并网运行或在路网中移动。可以用一个二维二进制变量 u_{m,i,t} ∈ {0,1} 表示“MPS m在t时段接入节点i”。当 u_{m,i,t}1 时MPS出力必须满足其功率容量上限当所有 u_{m,i,t}0 且MPS未接入任何节点则意味着它在移动。移动有速度限制因此从节点i到节点j的行程时间必须满足捷运距离除以通行速度再除以步长Δt向上取整。这个约束在YALMIP里可以用一个“移动完成时间标记”变量来实现否则会出错。第三MPS出力与网络状态的耦合。同一时段内接入节点的MPS出力可以注入本节点但若本节点处于孤岛状态与主网失去联络MPS只能带岛内负荷潮流约束会被“孤岛识别”逻辑限制。更严谨的做法是在模型中显式引入“孤岛索引”变量把MPS接入节点和岛内负荷之和建立联系。很多复现代码在这里做了简化通过设置一个很大的故障线路变量摘除后的节点分块利用潮流方程自动限制孤岛内的发电与负荷平衡。2.3 场景削减与不确定性处理如果你是严格复现一区论文大概率还会遇到一个模块场景削减。动态调度的输入是预配置阶段生成的多个典型故障场景而在动态调度阶段论文往往采用“时间序列上的故障演化场景”而不是静态的单时刻场景。这就意味着每个场景本身是一个二维数组——故障线路ID × 故障时刻。代码里常用同步回代消除法Backward Reduction对原始场景树做削减比如从500个原始场景削减到20个等概率场景。手写这个算法大概需要80行左右Matlab代码核心是不断计算场景两两之间的Kantorovich距离矩阵然后把与其它场景距离之和最小的那个场景合并到最近邻场景合并时把概率加起来。这一步如果求快可以直接调用MATLAB的probability distribution object或者自己写两层for循环实测100×100的场景矩阵根本不需要向量化优化直接硬算。3. Matlab代码实现从数学到可运行3.1 代码整体框架设计拿到一区论文的复现任务建议代码不要一上来就写主脚本。先把整体框架按“数据层—模型层—求解层—分析层”拆开。数据层负责生成IEEE 33节点系统参数、负荷曲线、故障场景模型层用YALMIP定义决策变量和目标、约束求解层调用Cplex或Gurobi求解MILP分析层把最优解落回13维数组里画出曲线。主脚本的组织逻辑大概是初始化参数节点数、时段数、MPS数量、移动速度、最大行程时间生成基础配电网拓扑从Matpower格式导入或直接在代码里定义支路表生成故障场景矩阵每种场景包含故障线路编号-断线时刻-持续时间用同步回代消除法削减场景定义YALMIP优化变量二进制变量为binvar连续变量为sdpvar循环添加约束潮流约束、MPS移动约束、出力约束设置Cplex/Gurobi求解参数MIP GAP设为0.5%~1%左右求解并把结果存入结构体绘图与分析建议从“全场景同时优化”起步。所谓全场景同时优化Sample Average Approximation, SAA指的是每个场景都对应一组动态调度决策目标函数取所有场景的期望值。这样做的好处是变量多但求解结构规整Gurobi处理这种带块结构的MILP效率很高坏处是求解时间可能长达几十分钟。如果你想做滚动时域那就要把模型参数化成一个带调度时刻输入的函数然后在主循环里反复调用——这一步我们后面单独讲。3.2 关键参数的选取与计算过程跑动态调度之前有几个参数必须自己先算清楚否则结果必定离谱。时间步长Δt。论文里一般取30分钟或1小时。如果你取1小时那么一台移动速度为40 km/h的MPS在配电网尺度节点间路网距离一般取1~5 km内几乎可以每时段迁移一个节点调度自由度极高解出来的方案会偏向“哪里故障去哪里”如果取15分钟MPS移动被显著限制你会看到调度方案更保守倾向于提前顶在关键节点不动。这是一个很微妙的敏感性参数建议多跑几组对比。MPS的容量配置。预配置阶段已经定了MPS的数量和类型动态调度阶段需要把它们当作参数传入。这些参数包括每台MPS的额定功率常见取100 kW~1 MW区间段、最大电量/油量对应持续时间、移动速度、接入切负荷响应时间。论文里经常会做“同构MPS”假设也就是认为所有MPS参数一致但真实现场大概率是异构的。异构化之后代码量会膨胀不少而且需要额外增加一组类型索引变量。负荷时序曲线。动态调度必须绑定时变负荷曲线。IEEE 33节点系统每节点都有基准负荷建议用日负荷率线性缩放。这里有个经验如果负荷曲线是24点调度周期拉长到48或72小时那么MPS电量约束会显著起作用——不能只在单一时刻考虑有功出力上限还要考虑累计放出电量。这个约束很多复现代码会漏漏掉之后MPS会变成“永动机”在仿真过程中反复满功率出力和真实物理不符。3.3 YALMIP核心建模代码详解先把最容易劝退的YALMIP代码部分拆解一遍。第一步定义变量u binvar(nMPS, nBus, NT, full); % MPS接入状态 P_mps sdpvar(nMPS, nBus, NT, full); % MPS有功注入 gamma binvar(nMPS, NT, full); % MPS是否处于移动状态 theta sdpvar(nMPS, NT, full); % 累计移动距离用于行程约束注意二维/三维变量的索引顺序。我建议统一用(mps, bus, time)的顺序避免最后reshape时把行、列搞混。YALMIP对nMPS×nBus×NT的三维binvar支持没问题但Cplex的输入接口在个别版本中需要把三维变量展开成向量这时用reshape(binvar(...), [], 1)转换。如果先定义二维binvar再循环赋值约束拼接会慢到让你怀疑人生。接着是核心约束MPS同一时刻只能接入一个节点或处于移动状态Constraints []; for t 1:NT Constraints [Constraints, sum(u(:,:,t), 2) gamma(:, t) 1]; end然后是移动时间一致性约束。假设MPS从节点i到节点j需要的最少时间步数为D_ij由距离矩阵除以速度并向上取整得到那么下面这个约束限制了“在t时段接入节点j的MPS必然在t-D_ij之前离开了某个节点i”for t 1:NT for m 1:nMPS Constraints [Constraints, sum(u(m,:,t), 2) ... 1 - gamma(m, max(1, t - D_min(m)) : t) ... ]; end end这个写法需要小心边界t-D_ij小于1时要从1开始。但更精细的写法需要引入“到达时刻”变量比如定义arrive(m,j,t)表示MPS m在t时刻到达节点j然后用u的差分约束来表示接入和到达的关系。潮流约束用DistFlow逐时段添加for t 1:NT for i 1:nBus % 从父节点注入的有功 - 负荷有功 MPS有功 流入子节点的功率和 Constraints [Constraints, ... P_inj(i,t) - P_load(i,t) sum(P_mps(:,i,t), 1) ... sum(P_line(:, i, t), 2)]; end end这个写法是示意实际实现中需要先构建支路-节点关联矩阵。求解设置方面Cplex参数一般这样设置options sdpsettings(solver, cplex, ... cplex.mip.tolerances.mipgap, 0.005, ... cplex.mip.tolerances.integrality, 1e-5, ... verbose, 2, ... savesolveroutput, 1, ... savesolveroutput, 1); optimize(Constraints, Objective, options);Gurobi对应参数名是gurobi.MIPGap。MIP Gap设0.5%是个不错的平衡点低于0.1%会把求解时间拉长数倍而高于2%可能让结果出现明显次优。复现论文时建议先保持默认Gap跑一次再收紧到0.5%验证结果是否一致确认模型没有隐藏的数值问题。4. 动态调度案例设计与结果分析4.1 IEEE 33节点测试系统搭建测试系统采用标准IEEE 33节点辐射状配电网基准电压12.66 kV总有功负荷约3.72 MW无功负荷约2.3 Mvar。支路参数直接用经典数据表填进去就行网上很好找。建议在代码里用两个矩阵存储bus(:, [编号, 有功负荷, 无功负荷]) 和 branch(:, [首端, 末端, 电阻, 电抗])。故障场景怎么设我这边用的是“台风灾害导致3条线路断线断线时刻依次是t4、t6、t8单位时段”。比如线路16、21、25断线会造成若干个负荷节点脱离主干网形成孤岛。MPS数量设为2台额定功率各300 kW移动速度取20 km/h考虑台风天路面通行能力下降调度时域取24个时段每时段1小时。预配置阶段的结果可能是在节点25和节点8预置了2台MPS。动态调度的任务就是在故障演化过程中决定这2台MPS是留在原地支撑还是移动到其他孤岛节点。留原地的好处是当前已有接入关系无需额外移动时间坏处是如果当前节点的负荷恢复优先级不高MPS可能处于“大马拉小车”的状态。动态调度就是在权衡这种“恢复收益”和“移动时间成本”之间的替代关系。4.2 最优调度结果解读跑通模型后输出目标函数值、各时段各节点的切负荷量、MPS接入位置和移动轨迹。你大概率会看到这样一个规律故障初期t1~4MPS按预配置位置原地待命或就近接入高优先级负荷节点中部时段t5~12当某些孤岛节点因故障失电且负荷量大、单位失负荷价值高时MPS会移动过去通常移动1~2个节点距离后期t13~24故障逐步排除、主网恢复供电MPS回到初始充电位置或继续补充供电。从曲线图上看有MPS动态调度与无MPS的“系统恢复面积”差异非常明显无MPS情况下孤岛负荷在故障期间基本全损失预配置固定接入情况下只有预配置位置的负荷被恢复动态调度情况下恢复的负荷量显著提高还能看到MPS移动过程中短暂出力下降的“缺口”——这是移动时间造成的必然代价。韧性指标建议计算两个负荷损失期望值EPNS单位MW·h和系统韧性三角面积。前者直接反映调度效果后者用来做论文对比图非常合适。所谓韧性三角横轴是时间纵轴是系统可用负荷比例灾害发生时曲线下降恢复过程中曲线回升下降-低点-回升形成的三角形面积越小系统韧性越好。4.3 场景削减对调度结果的影响同步回代消除法在代码里的实现值得单独提醒削减后的场景如果数量太少比如只剩3~5个动态调度结果会被少数极端场景主导解出的MPS调度策略可能过于冒险——集中资源去赌某个小概率严重故障场景数量太多比如50个以上求解时间指数级增长但目标函数值不再明显改善。我在复现时用20个削减场景Gurobi求解时间大约7~12分钟结果和500个原始场景的期望值差距在2.3%以内这个误差水平对工程应用完全可以接受。你可以跑一个“场景数量-目标函数值”的收敛曲线作为论文里的敏感性分析图期刊审稿人很吃这一套。5. 复现过程中的常见问题与排查技巧5.1 模型求解时间爆炸怎么办动态调度的MILP模型变量数很容易突破六位数。如果解到两三个小时还出不来先检查下面几项第一是否有多余的大M约束。YALMIP建模时如果用了M法处理逻辑约束M的值不要取“超大数”比如1e9取1e3或1e2就够过大的M会让Cplex/Gurobi的松弛解空间出现严重的数值病态分支定界效率暴跌。第二是否可以把一些连续变量固定上下界。MPS有功出力P_mps虽然是非负变量但它的上界受额定功率限制主动设置ub后能让求解器做更好的边界传播。第三对称性问题。多台MPS如果参数完全一致模型天然存在置换对称性求解器会重复探索等价分支。解决办法是加对称破缺约束比如按编号规定“MPS 1接入节点的编号不小于MPS 2的接入节点编号”前提是节点按顺序合理编号。这一步往往能节省30%以上的求解时间。5.2 解出来出现MPS“瞬移”或“穿墙”现象这个问题最典型几乎每个人第一次跑动态调度都会遇到。现象是明明MPS从t2在节点12到了t3出现在节点25中间路程需要4小时模型却认为1小时就到了。原因不外乎三个忘记加移动状态变量gamma只定义接入变量u。模型可以在任意两个节点间“跳变”因为它只限制了每个MPS同时最多接入一个节点没有限制“从接入到再次接入必须经过移动时间”。距离矩阵D算错了。D(i,j)存的是节点编号差而不是路网距离比如节点12到节点25的曼哈顿距离可能只有13但路网长度可能是8 km除以速度后是8个时段——你存的却是13个时段导致移动过于宽松。动态约束的边界条件没处理好。如果在t1之前没有给MPS赋“初始接入状态”模型可能让MPS从“虚空”中凭空出现在任意节点。正确做法是加入初始接入约束sum(u(m,:,1)) 1且初始位置必须等于预配置方案给出的节点。排查方法是把最优解的u和gamma矩阵打印出来逐时段看你很快能定位到问题出现在时段衔接还是距离矩阵。5.3 场景削减算法里的概率归一化问题用同步回代消除法削减场景之后如果忘了把剩余场景概率重新归一化到1目标函数值会出现系统性偏移——偏大或偏小取决于削减过程中丢掉的场景总概率偏大还是偏小。设原始场景概率为p_s削减后保留集合S每个保留场景的新概率为p_s p_s / Σ_{k∈S} p_k。不加归一化这行代码后面的动态调度目标函数期望值可能误差5%甚至更高。6. 从复现到改进还能怎么扩展如果你已经完整跑通上篇的预配置和下篇的动态调度恭喜你已经拿到了一区论文复现的整套流程。但让我说句实在话复现只是第一步。真正想用这个模型发自己的文章至少还有三个值得扩展的方向。一是把全局优化改成滚动时域控制。全局优化的假设是“我们知道整个故障演化轨迹”现实中这不成立。滚动时域控制把动态调度拆成多个决策窗口每个窗口只对未来4~6个时段做优化执行一个时段后滚动更新故障信息。这个改动的工程意义非常大代码量增加也不多——只需要把模型封装成函数在循环里反复调用实测从全局优化改到RHC大约是半天工作量但审稿人眼里这就是“从理论到应用”的一步。二是加入多类型移动资源协同。MPS之外很多论文还会考虑维修队Crew的调度。维修队修复故障线路MPS承担临时供电两者之间存在时序配合——维修队还没到现场MPS先顶上维修完成MPS撤出前往下一处。多资源联合调度能让MPS的利用率提升一个台阶代价是模型里需要额外引入“资源到达时间”耦合约束MILP的块结构会被破坏求解难度陡增。三是计及交通路网不确定性。MPS移动速度不是常数台风天气下路面积水、树木倒伏会影响路网通行能力。更精细的模型可以把交通路网作为独立的图模型MPS的移动时间和路况强耦合——这个方向更像运筹学里的“灾难物流”问题跨学科属性更强发文章的高地也更宽阔。我自己在实际复现中的体会是动态调度模型的风险不在数学推导而在编码习惯。YALMIP的语法确实友好但三维变量的索引、约束循环的边界、求解器参数的单位Gurobi的MIPGap要小数、Cplex的mipgap也要小数但有的版本要整数百分比这些细节稍不注意就会浪费一整天。如果让我给一条最实用的建议先把目标函数里MPS的相关项全部注释掉只求解“纯切负荷最小化”问题验证潮流约束写对了再放开MPS变量验证移动约束写对了最后加上惩罚系数调参调权重。分步调试虽蠢但稳能省下至少一个通宵的排查时间。完整Matlab代码实现框架、参数文件与测试算例数据已按模块拆好可以直接对照本文流程逐段跑通。