简介本资源面向电力系统韧性研究与配电网优化调度方向的研究生、科研人员及工程技术人员聚焦IEEE Trans. Smart Grid 2019年文献中MPS动态调度阶段的完整复现。针对灾后应急移动电源MPS派遣与配电网运行在时间尺度上的耦合、道路网络与电力网络的双重耦合难题资源给出了基于YALMIP建模的混合整数规划求解方案并在IEEE 33节点与123节点测试系统上验证了路由与调度策略对系统韧性的提升效果。压缩包共5个文件约49KB包含3个m脚本文件、1个p文件与1个vsdx拓扑图分别承担主调度程序、结果可视化、参数数据与节点接线示意等功能结构紧凑、便于直接运行与二次开发。目前已有322人学习下载适合希望快速掌握MPS动态调度建模思路、复现顶刊算例并迁移至自身研究场景的读者参考。1. 从预配置到动态调度MPS 在韧性配电网里到底扮演什么角色台风季一到沿海某 10kV 馈线因为主干线路倒杆下游三个台区同时失电。固定式储能还没充满柴油发电车从城区调过来要 40 分钟。这 40 分钟里如果有一台移动电源Mobile Power SourceMPS提前停在关键节点附近失电范围能压掉一大半。问题在于MPS 数量有限、容量有限、路网通行时间不确定停在哪、什么时候动、动几台全是耦合决策。上篇讲的是灾前预配置本质是一个「选址定容」的静态优化这一篇落到灾中的动态调度。MPS 动态调度要解决的是在故障集合已知或滚动更新的前提下如何分配每台移动电源的接入节点、接入时段和出力曲线使加权失负荷量最小。它和经典机组组合的区别在于——决策变量里多了「空间位置」这一维时间尺度也更短通常按 15 分钟一个时段滚动。适合谁看做过配电网优化、用过 YALMIP 建模、想把这套「预配置调度」两阶段框架跑通的人。下面所有模型都在 MATLAB YALMIP 里实现求解器用 Gurobi 或 CPLEX 都行代码结构不依赖具体求解器。2. MPS 动态调度的数学模型与 YALMIP 建模2.1 决策变量与目标函数怎么定动态调度的核心变量分三类。第一类是空间决策$x_{k,i,t}$ 表示第 $k$ 台 MPS 在时段 $t$ 是否接入节点 $i$第二类是时间耦合变量$y_{k,i,t}$ 表示第 $k$ 台 MPS 是否在时段 $t$ 从节点 $i$ 移动到节点 $j$这里通常用 $y_{k,i,j,t}$ 表达第三类是功率变量$p_{k,t}$ 为第 $k$ 台 MPS 在时段 $t$ 的放电功率$s_{k,t}$ 为其剩余电量。目标函数是加权失负荷最小加上 MPS 移动成本$$\min \sum_{t}\sum_{i} w_i \cdot LS_{i,t} \sum_{k}\sum_{t} c^{move} \cdot z_{k,t}$$其中 $w_i$ 是节点 $i$ 的负荷权重医院、通信基站权重高$LS_{i,t}$ 是节点 $i$ 在时段 $t$ 的失负荷量$z_{k,t}$ 是第 $k$ 台 MPS 在时段 $t$ 是否发生移动的 0-1 变量。移动成本这一项不能省——不加的话求解器会让 MPS 每个时段都换位置物理上不可行。2.2 用 YALMIP 把约束写成可求解形式约束分四组功率平衡、MPS 电量递推、移动逻辑、接入唯一性。下面这段是可直接跑的骨架代码。% MPS 动态调度主模型骨架 % K: MPS 数量, N: 节点数, T: 时段数 % w: 节点权重, Pload: 各节点各时段负荷, dt: 时段长度(h) % Tmove: 节点间移动时间矩阵, Emax: MPS 容量, Pmax: MPS 最大出力 x binvar(K, N, T, full); % 接入决策 y binvar(K, N, N, T, full); % 移动决策 z binvar(K, T, full); % 是否移动 p sdpvar(K, T, full); % 放电功率 s sdpvar(K, T1, full); % 电量状态 LS sdpvar(N, T, full); % 失负荷量 Constraints []; % 1) 每个时段每台 MPS 只能接入一个节点 for t 1:T Constraints [Constraints, sum(x(:,:,t), 2) 1]; end % 2) 电量递推s(k,t1) s(k,t) - p(k,t)*dt for k 1:K Constraints [Constraints, s(k,1) Emax]; % 初始满电 for t 1:T Constraints [Constraints, ... s(k,t1) s(k,t) - p(k,t)*dt, ... 0 p(k,t) Pmax, ... 0 s(k,t1) Emax]; end end % 3) 移动逻辑从 i 到 j 需要 Tmove(i,j) 个时段 for k 1:K for t 1:T for i 1:N for j 1:N if i ~ j t Tmove(i,j) T Constraints [Constraints, ... y(k,i,j,t) x(k,i,t), ... y(k,i,j,t) x(k,j,tTmove(i,j))]; end end end end end % 4) 功率平衡MPS 出力 剩余供电 负荷 for t 1:T for i 1:N Constraints [Constraints, ... sum(p(:,t) .* x(:,i,t)) LS(i,t) Pload(i,t), ... LS(i,t) 0]; end end Objective sum(sum(w .* LS)) 0.5 * sum(sum(z)); ops sdpsettings(solver,gurobi,verbose,1); diagnostics optimize(Constraints, Objective, ops);逻辑说明第 1 组约束保证「一台 MPS 同一时刻只在一个节点」这是接入唯一性第 2 组是电量递推注意 $s(k,1)Emax$ 表示灾前已充满如果预配置阶段有部分充电这里改成对应初值第 3 组是移动逻辑用 $y(k,i,j,t) \le x(k,i,t)$ 和 $y(k,i,j,t) \le x(k,j,tT_{move})$ 把「移动」和「两端接入」绑死避免求解器凭空瞬移第 4 组是功率平衡$LS$ 是松弛变量代表没被满足的负荷。参数说明Tmove矩阵由路网距离除以平均车速得到单位是时段数向上取整w建议按负荷等级设一级负荷取 10二级取 3三级取 1dt取 0.25 对应 15 分钟。求解器选 Gurobi 时MIPGap 默认 1e-4节点数超过 30 建议放宽到 1e-2否则求解时间会爆炸。2.3 线性化处理把乘积项拆开上面代码里p(:,t) .* x(:,i,t)是连续变量乘 0-1 变量属于双线性项YALMIP 不会自动线性化。常见做法是引入辅助变量 $q_{k,i,t} p_{k,t} \cdot x_{k,i,t}$用 Big-M 法约束% 线性化 p*x 乘积项 q sdpvar(K, N, T, full); M Pmax; % Big-M 取 MPS 最大出力 for k 1:K for i 1:N for t 1:T Constraints [Constraints, ... q(k,i,t) M * x(k,i,t), ... q(k,i,t) p(k,t), ... q(k,i,t) p(k,t) - M * (1 - x(k,i,t)), ... q(k,i,t) 0]; end end end这样功率平衡约束改成sum(q(:,i,t)) LS(i,t) Pload(i,t)。Big-M 取 Pmax 是最紧的界取大了会拖慢分支定界。如果 MPS 数量超过 5 台、节点超过 20 个建议改用 Gurobi 的 indicator constraintsYALMIP 里用implies函数写求解效率更高。3. 滚动调度与故障场景更新的实现细节3.1 滚动时域怎么切灾中调度不是一次算完而是每 15 分钟滚动一次。每次滚动只优化未来 4 个时段1 小时但只执行第一个时段的决策下一轮再重新优化。这样做的原因是故障集合会变——可能又有新线路跳闸也可能某条线路抢修恢复。实现上把上面的模型包成一个函数solve_mps_dispatch(fault_set, t_start, horizon)每轮传入当前故障集合和起始时段。关键点是电量状态要跨轮传递上一轮结束时的 $s(k,T1)$ 作为下一轮的 $s(k,1)$。% 滚动调度主循环 horizon 4; % 每次优化 4 个时段 total_T 24; % 全天 24 个时段 s_current Emax * ones(K,1); % 初始电量 for t0 1:horizon:total_T t_end min(t0 horizon - 1, total_T); [x_opt, p_opt, s_next] solve_mps_dispatch(... fault_set, t0, t_end, s_current); % 只执行第一个时段 apply_dispatch(x_opt(:,:,1), p_opt(:,1)); s_current s_next; % 电量传递到下一轮 end逻辑说明horizon取 4 是经验值太短会频繁切换 MPS 位置太长则对故障变化响应迟钝。s_current跨轮传递是滚动调度的核心漏掉这一步会导致电量凭空重置。apply_dispatch是实际下发指令的接口仿真里可以只记录不执行。3.2 故障场景生成与权重设置故障集合通常来自蒙特卡洛抽样或历史故障库。每个场景是一组「线路编号故障时段」的列表。权重 $w_i$ 的设置直接影响调度结果下面这张表是常见负荷等级的参考值。负荷类型权重 $w_i$典型节点说明一级10医院、数据中心失电后果严重二级3商业区、水厂可短时中断三级1居民区可较长时间中断特殊20通信基站影响面广权重不是越大越好。如果一级负荷权重设成 100求解器会优先保它但可能导致 MPS 长距离移动反而增加总失负荷。我一般先用 10:3:1 跑一版看结果再调。3.3 求解失败的排查顺序YALMIP 报 infeasible 时按这个顺序查第一看diagnostics.info如果是 1 说明模型不可行用diagnostics里的infeasible字段定位冲突约束第二检查Tmove矩阵如果某两个节点间移动时间超过 horizon移动约束会直接冲突第三检查电量递推如果Pmax * dt * horizon EmaxMPS 会在 horizon 内放空后续时段无电可放功率平衡就崩了。常见做法是加一个「松弛变量」到功率平衡里允许少量失负荷先让模型可行再逐步收紧。YALMIP 里用optimize的usex0选项给初值也能帮求解器更快找到可行解。4. 从调度结果反推预配置两阶段闭环与验证4.1 调度结果怎么反馈给预配置上篇的预配置模型假设 MPS 停在某节点就不动了但动态调度会移动它。如果预配置阶段没考虑移动可能出现「预配置选了 A 点但调度发现 B 点更优MPS 全跑 B 点去了」的情况。闭环做法是用动态调度的平均移动距离和平均接入时长修正预配置模型里的「服务半径」参数迭代 2-3 轮。具体操作跑完 24 小时滚动调度后统计每台 MPS 的实际移动次数和总移动距离如果平均移动距离超过预配置假设的 30%说明预配置的选址偏了把预配置模型里的覆盖半径缩小重新求解。4.2 验证调度方案是否合理的三个指标第一个指标是加权失负荷率即 $\sum w_i LS_{i,t} / \sum w_i P_{load,i,t}$低于 5% 算合格。第二个指标是 MPS 利用率即实际放电量除以容量低于 40% 说明 MPS 配多了或位置不对。第三个指标是移动次数如果每台 MPS 一天移动超过 6 次说明调度太频繁实际执行会有困难。% 计算三个验证指标 WLS sum(sum(w .* value(LS))); TotalLoad sum(sum(w .* Pload)); loss_rate WLS / TotalLoad; utilization sum(sum(value(p))) / (K * Emax); move_count sum(sum(value(z))); fprintf(加权失负荷率: %.2f%%\n, loss_rate*100); fprintf(MPS 利用率: %.2f%%\n, utilization*100); fprintf(总移动次数: %d\n, move_count);逻辑说明value()是 YALMIP 提取求解结果值的函数LS、p、z都是优化后的变量。loss_rate越低越好但低于 2% 可能是权重设置不合理导致求解器过度优化。utilization在 40%-70% 之间比较合理太低浪费太高说明 MPS 不够用。4.3 一个容易忽略的坑时段粒度与移动时间如果时段粒度是 15 分钟但两个节点间移动需要 20 分钟Tmove向上取整为 2 个时段30 分钟这会导致 MPS 实际到达时间比需要的时间晚 10 分钟。常见做法是把时段粒度设成 5 分钟或者把移动时间矩阵按实际路网算得更细。我一般用 5 分钟粒度跑调度输出结果再聚合到 15 分钟上报这样精度和计算量都能接受。另一个坑是初始电量。如果预配置阶段 MPS 只充到 80%动态调度里s(k,1)要改成0.8*Emax否则结果会偏乐观。这个初值最好从预配置阶段的充电计划里读不要硬编码。本文还有配套的精品资源点击获取