企业官网建设流程全解析

先把一个最容易被问住的点说在前面：这篇所谓的“SCI一区复现”，真正的难点并不是把论文里的公式敲进Matlab，而是怎么把“配电网韧性提升”这个大词，落到“应急移动电源（MPS）到底该在哪个时刻、哪个节点、以多少出力接入”这样一个可计算的问题上。上篇讲的是预配置——也就是灾前把这些移动电源放在哪里蹲点；这篇讲的是动态调度——也就是灾难发生之后，怎么指挥这些MPS在各节点之间转移、接入、出力和退出。后者才真正让人头疼，因为这里不仅有空间维度，还有时间维度，变量之间互相牵制，优化模型一不小心就无解或者解出来一堆违反物理直觉的“神仙调度”。

如果你正在复现这类工作，或者准备用Matlab搭一个类似的两阶段韧性优化框架，那这一篇应该能帮你少走不少弯路。下面我把MPS动态调度的问题定义、数学建模、求解器选型、Matlab代码骨架，以及我在实操中踩过的坑完整过一遍。这些经验都是我拿着实际算例跑出来的，不是课本上的标准答案。

1. 问题定位：MPS动态调度到底在解决什么

1.1 配电网韧性提升的语境

先说“韧性”（ resilience ）。电力系统里常说的“可靠性”，指的系统在常态扰动下维持供电的能力，关注的是年平均可用度；而“韧性”面向的则是小概率、高影响事件，比如台风、洪涝、极端冰雪天气。这类事件会在短时间内造成多线路同时断开，常规的N-1校验根本不覆盖这种场景。此时配电网很可能出现大面积失电，而抢修队伍不可能瞬时恢复所有线路，于是就有了一个时间窗口，在这个窗口内，必须靠现有资源尽量保住重要负荷。

移动应急电源MPS就是在这种背景下进入模型的。它本质上就是一个可以跟着卡车跑的柴油发电机或储能单元，有容量上限，有额定出力上限，从节点A转移到节点B需要时间，且转移期间无法供电。这些特性决定了它的调度不是简单的“把电源放哪”，而是一个带时序的路径—出力联合优化问题。

我复现时选的是IEEE 33节点配电系统，把其中几条线路设为故障，主网失电，形成一个时间段内的孤岛或局部失电区。MPS的任务就是在这个失电窗口内，通过移动和接入，动态改变局部供电拓扑，把重要负荷的断电时间压到最低。

1.2 预配置与动态调度的衔接逻辑

标题带有“(下)”，说明这篇讨论的是两阶段决策的后半段。上篇的预配置解决的是“灾前”：在灾害预测信息的基础上，决定每台MPS的初始停放节点和预备容量。预配置的结果不是MPS最终位置，而是一个起点；动态调度从这个起点出发，拿到灾害实际造成的故障场景，再决定每台MPS在每一时段去哪里、发多少电。

两阶段之间不是割裂的，而是递进关系。预配置模型的目标函数里通常已经包含了“预期恢复效果”的影子，比如最小化期望失负荷量；动态调度则是把这个预期落在具体场景里。所以复现时不能只看下篇的公式，要把上篇给的初始停放位置当成固定输入写进动态调度模型里。

如果你没有上篇的代码，也可以自己随便给一个合理的初始位置，比如根据故障预测概率把MPS停放负荷中心附近。这样动态调度模型照样能跑通，只不过曲线会差一些。

1.3 为什么用Matlab做这件事

我知道很多电力方向的优化项目现在都流行用Python+Gurobi，但Matlab在配电网这块依然有不可替代的地方。一是配电网的潮流计算、拓扑分析有很多现成工具包，数据处理随手就能画图；二是Yalmip这个建模层让复杂约束的表达非常直观，尤其适合二值变量多、时序关系复杂的模型，代码写出来基本就是数学公式的翻译。

另外，评审或课程组里总有那么一批人的主干代码是Matlab，采用Matlab也方便后续交接和扩展。如果你手上正好有Matlab R2021b以上的版本，再装上Yalmip和一款MILP求解器（Gurobi或CPLEX），就完全够用。

2. 核心模型拆解：动态调度问题的数学表达

2.1 目标函数怎么定才合理

动态调度的目标函数并不只有一种写法，我在复现过程中先后试过三种，结果差异很大。最常见的写法是最小化总失负荷量（ENS）或者说最大化累计供电量：

这里把时段按小时离散，T是调度周期（比如24个时段），ω_n表示节点n的重要程度权重。这个目标简单直接，求解器也友好，因为它是求和，不会有太奇怪的数值行为。但是纯用这个目标有个问题：它不区分“早恢复”和“晚恢复”，只要总供电量一样，早恢复晚恢复在目标值上是完全等价的。

所以更贴合韧性语义的做法是引入恢复性指标。很多论文用韧性三角形（Robustness, Resourcefulness, Recovery）下的时变韧性函数来做，最常见的是某个时刻的供电恢复率与理想恢复曲线围成的面积。实际建模时，可以通过在目标函数中给每个时段的失负荷量乘一个随时间递增的惩罚系数来实现，越晚恢复损失越大，这样求解器才会倾向于把MPS调度到“先恢复重要负荷”的路径上。

我最终用的目标是“最小化加权失负荷量之和+时间惩罚项”。具体地，把目标函数写成：

min Σ_t Σ_n ω_n · α_t · (P_n^load - P_n^supply(t)) · Δt

这里的α_t是随t递增的权重因子，比如α_t = (t/T)^0.5，让恢复速度对目标的影响更贴合韧性语义。当然这只是个人处理的技巧，不同论文对韧性指标的定义不同，但你只要有办法把“韧性”翻译成可计算的目标，整个模型就成立了。

2.2 三类核心约束

动态调度模型的难点主要在约束，而不是目标函数。我把约束归纳成三类，对应三种物理世界的限制。

第一类是MPS自身的时序物理约束。每台MPS同一时刻只能接入一个节点，这是一个二值变量x(k,t,n)的唯一性约束：对任意k、t，Σ_n x(k,t,n) ≤ 1。MPS的输出出力只有在接入状态为1时才允许非零，否则出力必须为0，这需要Big-M约束：0 ≤ P(k,t) ≤ M · Σ_n x(k,t,n)，这个M取MPS的额定容量。同时还有电量状态约束，比如S(k,t) = S(k,t-1) - P(k,t)·Δt，并且S有上下限，对应MPS的储能容量和最低荷电比例。

第二类是转移行为约束。MPS在时段t接入节点n，在时段t+1接入节点n'，如果n和n'不一样，就必须保证“转移时间≤单位时段时长”。这是动态调度区别于静态选址的关键。在实际建模中，需要把“节点发生切换”这件事用二值变量显式表达出来，然后给一个带Big-M的约束：travel_time(n, n') ≤ Δt + M(1 - move_flag)。这里move_flag就是表征“这台MPS在t到t+1之间是否移动”的变量，而travel_time矩阵可以事先通过Floyd最短路算法算出来。

第三类是配电网潮流约束。这是最费精力的部分。要描述MPS接入节点后对整个配网电压和功率分布的影响，严格做法是用DistFlow模型，这会引入非线性的P²+Q²项，常见处理是用二阶锥松弛（SOCP）把它变成可解的凸问题。但如果你只是做24个时段的MILP，不想引入过多非线性，还有一个更工程化的替代方案：做一次线性化DistFlow，假设电压在1.0 p.u.附近波动，忽略网损，用线性方程描述支路潮流和节点功率平衡。在故障范围有限、MPS容量不大的场景下，这种线性化误差对调度决策的影响通常可以接受。

节点功率平衡可以用一个线性等式来描述：从上游注入的功率+本地MPS出力+其他分布式电源出力=该节点负荷+流向子节点的功率。每个节点的电压幅值偏移量要限制在±5%以内，线路容量限制用不等式约束。

2.3 为什么要用Big-M法

我在建模时大量使用Big-M技巧，这里多说两句。很多新手在写“出力只有接入时才允许非零”时，会试图直接用一个if判断或者逻辑表达式建模，这在Yalmip里是跑不通的，必须把它拆成线性不等式。用Big-M的本质是“用一个足够大的常数把约束‘松弛’掉”：当接入状态为0时，右边变成M·0=0，约束退化为P≤0，配合P≥0就把出力强制为0；当接入状态为1时，约束变成P≤M，而M只要大于容量上限，约束实际上不起作用。

M的选择要谨慎。太大会导致数值病态，求解器报数值警告甚至错误；太小又会把可行域砍掉。经验做法是根据该MPS的最大出力加上一个小裕度来取，比如额定容量1.2MW，就取M=1.5，而不要为了“绝对保险”取到1000。

3. Matlab代码实现细节：从框架到每个矩阵

3.1 数据组织和预处理

拿到代码的第一步不是写约束，而是把数据组织好。我习惯用一个结构体或MAT文件把所有场景参数集中管理。

以IEEE 33节点为例，节点基础数据包括：节点编号、有功负荷（变成时段负荷矩阵后是T×N）、负荷权重（一级负荷权重高，二级次之，三级最低）。线路数据包括：起始节点、终止节点、电阻（标幺值或欧姆值）、电抗、容量上限。故障信息是一个N×N的线路断开矩阵，把故障线路直接置0，代表断开。

MPS数据包括：每台MPS的编号、额定容量、初始位置（来自预配置结果）、当前电量、移动速度（km/h）、节点间最短路径时间矩阵travel_time。travel_time的计算不要自己手写，直接用图论最短路函数，把配电网的拓扑当成无向加权图处理，权重取线路长度/移动速度。

代码大概长这样：

% 读取基础数据 mpc = loadcase('case33.m'); % 用Matpower格式读取IEEE 33节点系统 [N, branch] = size(mpc.branch); % 定义MPS参数 K = 3; % MPS台数 MPS_capacity = [0.5, 1.0, 0.8]; % 单位MW MPS_energy = [2, 4, 3]; % 单位MWh MPS_pos0 = [12, 18, 22]; % 初始位置，来自预配置结果 % 用最短路算法计算转移时间矩阵 travel_time = zeros(N, N); for i = 1:N for j = 1:N if i ~= j travel_time(i,j) = shortest_path_time(i,j, mpc, speed); end end end

这一步看似基础，但其实我在这里浪费过不少时间。因为IEEE 33节点系统的节点编号和线路起始节点不一定严格对应0/1索引，必须统一坐标，不然后面构造Yalmip变量的维度会错位，求解器报“Inconsistent dimensions”能让人排查到崩溃。

3.2 变量定义与目标函数实现

Yalmip的关键优势在于，你可以直接用sdpvar和binvar定义变量，然后把约束写成一个cell数组。动态调度的核心变量有三个。

% 定义决策变量 x = binvar(K, T, N, 'full'); % x(k,t,n)=1表示第k台MPS在t时段接入节点n P = sdpvar(K, T, 'full'); % MPS输出的有功功率 S = sdpvar(K, T+1, 'full'); % 荷电状态/剩余电量 % 定义目标函数 Objective = 0; alpha = sqrt((1:T) / T); % 时间惩罚因子 for t = 1:T for n = 1:N load_supply = sum(P(:,t) .* squeeze(x(:,t,n))); % 节点n上接入的MPS出力之和 Objective = Objective + weight(n) * alpha(t) * (load_profile(t,n) - load_supply); end end

这里有个非常容易踩的坑：如果直接像上面这样把P和x相乘，就会引入双线性项，模型变成非凸MINLP，Gurobi直接拒绝求解。正确做法是引入每台MPS在每个节点上的辅助变量P_n(k,t,n)，让它满足P_n ≤ M·x(k,t,n)和P_n ≤ P(k,t)，然后节点功率平衡用的是P_n而不是P乘x。

这也是我在复现时踩过最大的坑之一，后面在问题章节会详细说。

修正后，目标函数应这样写：

% 引入辅助变量: P_node(k,t,n) P_node = sdpvar(K, T, N, 'full'); Constraints = []; % 添加约束: 每台MPS同一时刻最多接入一个节点 for k = 1:K for t = 1:T Constraints = [Constraints, sum(x(k,t,:)) <= 1]; end end % 添加约束: P_node与x、P的关系 for k = 1:K for t = 1:T for n = 1:N Constraints = [Constraints, 0 <= P_node(k,t,n) <= M * x(k,t,n)]; Constraints = [Constraints, P_node(k,t,n) <= P(k,t)]; Constraints = [Constraints, P(k,t) <= P_node(k,t,n) + M * (1 - x(k,t,n))]; end end end

这样做的逻辑是：如果x=0，那么P_node被压到0；如果x=1，P_node可以等于P(k,t)，从而在节点功率平衡中使用P_node。这是一个纯线性表达。

3.3 转移约束和电量约束怎么落地

转移约束是MPS动态调度里最微妙的部分。核心逻辑是：如果MPS在t时段接入节点n，在t+1时段接入节点n'，且n≠n'，那么从n到n'的最短转移时间必须不超过一个时段长度Δt。

在Yalmip里，这个约束不能直接用if写，要用“切换指示变量”来显式定义。代码如下：

% 定义切换变量: move(k,t)=1表示第k台MPS在t到t+1时段发生移动 move = binvar(K, T-1, 'full'); for k = 1:K for t = 1:T-1 % 移动指示变量与接入状态的关系 for n = 1:N for np = 1:N if n ~= np Constraints = [Constraints, x(k,t,n) + x(k,t+1,np) - 1 <= move(k,t)]; end end end % 如果移动，则必须满足时间约束 Constraints = [Constraints, MPS_time_penalty * move(k,t) <= 1]; % 上面这行是示意，真正的约束是用travel_time矩阵和Δt比较： % travel_time(n,np) <= dt + M*(1 - move(k,t)) end end

真正的移动时间约束需要结合接入状态与travel_time矩阵写，建议单独封装成一个函数build_move_constraints(x, travel_time, dt, M)，返回约束集。这里的大M不能取太大，否则转移时间约束形同虚设；我调试时取M为最大travel_time的三倍就足够。

电量约束也需要小心处理。S(k,t)表示第k台MPS在时段t结束时的剩余电量，初值来自预配置阶段的S(k,0)。电量的动态过程是：

for k = 1:K for t = 1:T Constraints = [Constraints, S(k,t+1) == S(k,t) - P(k,t) * dt]; Constraints = [Constraints, S_min(k) <= S(k,t+1) <= S_max(k)]; end end

需要注意，MPS在移动过程中自身也会消耗少量能量，但很多论文会忽略这点。如果你想把模型做细，可以在移动时电量减少一个固定值energy_move · move(k,t)。我一开始没加这一项，导致调度结果在最后一刻出现“满电移动”的不合理情况，加上移动耗电项之后才正常。

3.4 潮流约束的简化处理与实现

如果不引入完整潮流，用线性化DistFlow的话，需要定义两个主要变量：每个节点t时段的电压偏移量delta_V(n,t)和每条支路的有功潮流P_branch(l,t)。

节点功率平衡约束对每个节点、每个时段写：

Σ_{l∈in(n)} P_branch(l,t) - Σ_{l∈out(n)} P_branch(l,t) + Σ_k P_node(k,t,n) = load_profile(t,n)

这个约束保证了每个节点注入和流出的功率与负荷平衡。电压偏移约束则采用如下形式：

delta_V(child,t) = delta_V(parent,t) - (r_l * P_branch(l,t)) / V_base

这里V_base取1.0 p.u.，r_l是支路电阻。之所以叫线性化，是因为把DistFlow中的非线性项(P²+Q²)省略掉，只用实功率近似。在故障恢复场景下，MPS容量通常比较小（几百kW到几MW），这种近似不会导致决策上的重大偏差。

3.5 求解器配置与计算规模控制

把所有约束组装好之后，用一行代码交给求解器：

ops = sdpsettings('solver', 'gurobi', 'verbose', 2, 'mipgap', 0.01); result = optimize(Constraints, Objective, ops);

我建议不要一上来就求解24时段、3台MPS的完整模型。第一次先跑4时段、1台MPS的模型，把逻辑验证通，再逐步扩大。24时段×33节点×3台MPS的完整MILP，变量规模大概是K×T×N = 2376个二值变量，加连续变量后大概5000个变量，Gurobi默认参数下通常几分钟到十几分钟能收敛到1%的MIPGap。如果超过30分钟没动静，多半是约束写错了，而不是算力不够。

4. 复现过程中的坑与经验清单

4.1 求解器无解时的排查思路

这是我最想分享的实战经验。当你写了完整的模型、约束也完全照着论文来，结果优化返回“Infeasible problem”，第一反应不应该是怀疑求解器，而是回头检查约束之间是否互相冲突。最常见的无解原因是“转移时间约束和唯一性约束冲突”：MPS在t时段接入节点n，但t+1时段必须接入节点n'，然而n到n'的转移时间大于1个时段，于是约束同时要求“必须移动”和“移动时间不足”，直接无解。

解决办法有两个：一是把时间步长放大（比如从0.5小时改成1小时）；二是允许MPS在无法按时到达下一节点时保持原地待命，也就是加一个“允许不移动”的逻辑，把移动判定做成≤而不是==。如果论文里是硬性约束，那你至少要保证初始位置的选择本身满足“到任何目标节点的转移时间都能在单时段内完成”，这又会反向约束预配置的结果。

另一个无解常客是“负荷权重矩阵里有零值”。如果某些节点负荷为0但root节点权重也是0，那么优化器为了加快求解可能直接砍掉这些节点，结果某个强加的等式约束让某个节点必须满足功率平衡，而它又没有负荷又没有电源，自然无解。干脆把所有负荷为0的节点权重设成一个小正数，比如1e-3，能有效避免这类问题。

4.2 求解结果不合理：MPS到处乱跳

我跑出的第一版结果，MPS在t=1时在节点12，t=2时在节点25，t=3时又回到节点12，来回转移，看起来非常“智能”，但实际上完全违反物理。原因出在我遗漏了“接入节点切换后，转移期间无法供电”的约束。

如果节点本身没有供电能力，MPS从一个节点移到另一个节点，至少需要跨越转移时间，而转移时间内不能接入任何节点。如果模型里没表达这个“转移空窗”，优化器就会让MPS瞬移。严格做法是：对于每台MPS，如果发生了移动，则需要在转移时段内强制所有接入状态为0。这等价于在t时段结束后进入转移，在t+1时段必须还在路上，直到t+2时段才能接入新节点。

我最后采用的简化处理是设定一个“切换冷却时段c”：如果MPS在t时段执行了移动，则t+1时段必须保持所有x都等于0，t+2才能重新接入。这个冷却时段直接等于ceil(最大转移时间/Δt)。加了这个约束之后，调度轨迹变得稳妥很多。

4.3 求解速度优化：Big-M别乱设

M≤1.5倍额定容量，转移时间约束的M≤最大travel_time的3倍。按这个原则设，数值稳定性好很多。25倍以上就等着看Gurobi输出一堆numerical trouble吧。

另外，变量对称性也是拖慢求解的一大杀手。3台MPS如果容量都相同，模型会认为把它们交换一下还是同一个解，分支定界里就会出现大量对称分支。消除对称性很简单：给MPS加编号约束，比如规定MPS1的初始电量≥MPS2的初始电量≥MPS3的初始电量，或者强制MPS1的接入位置字典序小于MPS2的接入位置。这个技巧能把求解时间压缩到原来的三分之一到五分之一。

4.4 算例对比与结果展示的技巧

复现论文时，一定要把结果的可视化做到位，这是展示工作时最有说服力的部分。我通常会输出三类图：第一类是各时段的负荷恢复率曲线，横轴是时段，纵轴是恢复比例，把“无MPS”“静态MPS”“动态调度MPS”三条线放一起对比，动态调度的优势一目了然。第二类是每台MPS的时序位置热力图：横轴是时段，纵轴是节点，某台MPS在哪个时段在哪个节点，用颜色块标出来。第三类是每台MPS的出力和电量曲线，看它是不是合理利用了存储容量。

在表格里，我习惯把关键指标整理成：累计失电量（kWh）、重要负荷失电时间（h）、电压最低点（p.u.）、MPS总移动次数。动态调度相对静态方案通常能降低15%到30%的失电量，具体数值取决于故障场景的严重程度和MPS的数量。

5. 一些补充分享

做这类复现工作时，我逐渐意识到一个很现实的问题：论文里写得无比严密的模型，落实到代码里总会有一堆“论文没写但你绕不开”的细节。比如预配置阶段的初始电量如何进入动态调度模型的初始状态，比如MPS移动时的自身能耗怎么处理，比如负荷曲线的时变特性该怎么从T时段的时间序列映射到每个节点的负荷剖面。这些细节往往才是复现质量的真正分水岭。

我的建议是，每拿到一篇论文，先不要急着抄公式，而是把它的因果链路捋一遍：预配置阶段输出的到底是什么变量，动态调度阶段输入的初始条件又有哪些，两阶段之间通过哪些参数耦合。把这个链条理顺了，Matlab代码不过是体力活。

最后再分享一个小技巧：如果你在调试时总也找不到模型无解的原因，可以先固定一部分决策变量，比如固定所有MPS的接入位置为预配置阶段的初始位置，只优化出力，看子问题是否可行。如果子问题也不可行，就说明问题不在移动约束上，而是在潮流或电量上。这种“二分式排查”在复杂优化模型里非常实用。