企业官网建设流程全解析

做一个配电网韧性的SCI复现项目，最大的感受就是：论文里一句“采用两阶段方法”带过的内容，落到代码里是两种完全不同的问题。上篇我们拆过灾前预配置——在灾害来临前，把有限的应急移动电源（MPS）放到哪些关键节点待命。这篇继续往下走，解决灾后“怎么动”的问题：故障已经发生，负荷一片片地失电，手头这几台移动电源车，每个时段应该开到哪个节点、给谁供电、充多少放多少，才能让整个配电网的恢复效果最好。

这篇复现笔记聚焦的是MPS动态调度这部分，对应Matlab代码的完整实现逻辑。内容适合电力系统方向准备复现论文的研究生，也适合刚接触配电网韧性、移动储能调度、MILP建模的工程师。我会把建模思路、代码框架、求解器配置，以及我在实际复现时踩过的高频坑都整理出来。代码不是单纯地把公式翻译成脚本，关键是理解为什么这么建模，约束项之间是怎么互相牵制的。

1. 先把问题定义清楚：动态调度到底在解什么

1.1 为什么“预配置”和“动态调度”必须分两个阶段做

配电网韧性问题里，应急移动电源的部署天然就是一个两阶段决策过程。原因很简单：灾害发生前，你只知道台风要来了，但具体哪条馈线会断、哪些杆塔会倒，是不确定的。这时候做的决策是“把MPS放在哪儿”——这是典型的先验决策，必须在不确定信息尚未揭晓时落地，对应的就是预配置阶段。

一旦灾害过去，故障位置、故障范围、失电负荷基本确定了，这时候MPS可以从预置点出发，根据实际受灾情况动态调整自己的位置和出力。这就是动态调度阶段，属于后验决策。两个阶段的目标虽然一致——都是最小化失电损失、提升恢复水平——但信息结构和约束差别很大，放在一个模型里求解会非常复杂，解空间爆炸，所以文献里几乎都采用两阶段或序贯框架。

理解这个拆分逻辑很重要。我复现的时候一开始贪心，想用一个统一的模型把“布点+调运”全做了，结果变量量级直接失控，求解器跑到天亮都出不来可行解。后来老老实实按论文思路拆成两段：上篇解决“从哪里出发”，下篇解决“往哪里去、供多少电”，问题的颗粒度一下子清晰了。

1.2 动态调度阶段的输入输出到底有哪些

在做代码之前，必须先把数据流理清。动态调度模块的输入不是随便拍的，它来自两路：一路是预配置阶段输出的MPS初始驻点，另一路是灾后的故障场景和负荷信息。

具体到我的Matlab实现里，输入包括：

• 配电网拓扑和线路参数（节点数、支路阻抗、线路容量）；
• 失电区域的节点编号和负荷曲线（有功/无功需求，分时段）；
• 故障支路集合（哪些线路断开，系统被划分为多少个孤岛）；
• MPS参数：数量M、单台容量、最大输出功率、移动速度、初始预置节点；
• 调度周期T和时段步长，比如把24小时按1小时分成24个时段。

输出则是决策变量的最优取值：每台MPS在每个时段移动到了哪个节点，在哪个节点以多大功率放电，以及整个T个时段内系统的负荷恢复情况。

这里要特别提醒一个前置检查：预配置阶段的输出文件里，每台MPS的初始位置必须是有效的节点编号，并且在动态调度阶段要保持相同的数据格式。我遇到过两次，上篇代码里节点编号是从0开始存成数组索引，下篇读取时却按1开始的Matlab习惯去索引，结果全是维度报错。建议在主程序里加一段断言，检查初始驻点编号是否在[1, N]范围内，并且和拓扑矩阵的列数能对上。

2. MPS动态调度建模：目标函数与约束体系

2.1 时间、空间、状态三要素怎么落到变量上

动态调度模型的核心是刻画MPS在时间和空间上的运动行为。MPS是一辆移动的车，它在某个空间节点待着，并且可以移动到别的节点。建模的第一步就是把“位置”状态表达出来。

用二进制变量x(m, n, t)表示第m台MPS在时段t是否位于节点n，取1表示在，取0表示不在。因为一台MPS在一个时段内只能在一个地方，所以每个时段对每台MPS要加约束：所有节点上的x(m, n, t)之和等于1。

同时还需要一个连续变量E(m, t)表示第m台MPS在时段t结束时剩余的电量。MPS的移动会消耗能量，对外放电也会消耗能量，这些都需要在能量平衡约束里反映出来。如果你还考虑MPS在某个节点可以给负荷供电，还需要P_out(m, n, t)这样的放电功率变量。

这里我习惯把所有变量按维度组织好再开始写约束，不然代码容易乱。比如用三维变量存放位置状态时，第一个维度是MPS编号、第二个是节点编号、第三个是时段编号。在YALMIP里直接用binvar(M, N, T, 'full')定义即可，但要注意YALMIP默认维度顺序是全索引的，写入约束循环时别搞错顺序。

2.2 目标函数怎么选：加权失电量、恢复收益还是综合成本

动态调度阶段的目标函数，文献里常见的有三类：最小化系统加权失电量、最大化恢复的负荷量、最小化MPS调度成本加上失电惩罚的综合目标。虽然在数学上等价的，但实现难度和解释性不一样。

我复现时采用的是“加权失电量最小化”目标。具体写法：把每个节点在时段t的基础有功负荷记为P_load(n, t)，实际恢复供电量记为P_res(n, t)，两者之差就是失电量。根据节点的重要程度不同，给每个节点一个权重w(n)，目标就是使所有时段的失电量加权总和最小。

用YALMIP表达目标函数，核心一句话：

objective = sum(sum(sum(w .* (P_load_3d - P_res_3d))));

其中w是N维权重向量，需要先扩展成N x T的矩阵再参与三维乘法。很多新手在这里栽跟头：w是节点权重，不是时段权重，直接写成w' * ones(1, T)再逐元素乘，就能避开维度不匹配的警告。

为什么不直接用“恢复负荷最大”这个目标？因为两者实质等价，但加权失电量在表达式上更直观，而且后续如果要加惩罚成本（比如MPS行驶费用），加权失电量目标更容易和惩罚项线性叠加，不会出现两项量纲不一致的问题。

2.3 关键约束：位置唯一性、能量迁移、功率容量、移动耗时

动态调度建模最核心的是约束体系。任何一个约束出错，结果要么是求解器报错，要么是输出一个看似合理但物理上根本不可行的方案。

位置唯一性约束：每个时段每台MPS只能在一个节点，这约束很简单但必不可少。一旦缺失，求解器会给同一台MPS分配多个节点，模型输出会完全失真。

MPS供电能力约束：MPS在节点n的放电功率不能超过其最大输出功率，同时放电功率必须小于该节点负荷需求，因为MPS只能恢复负荷，不能反向倒送功率制造新的安全隐患。

能量平衡约束：这是移动储能类问题的灵魂。MPS在t时段结束时的剩余能量，等于上一时段剩余能量，减去本时段输出电能，再减去移动消耗的能量。写成计算公式：

E(m, t) = E(m, t-1) - P_out(m, t) * delta_t - P_move(m, t) * delta_move

但这里有个坑：P_move只在移动发生时起作用，而移动事件本身需要用x(m, n, t)的转移来刻画。如果直接写成总移动耗能常量，就会导致MPS在“原地待命”时也平白扣掉移动能量，结果严重偏离实际。

解决的办法是为每台MPS引入一个移动状态标志变量move_flag(m, t)，它不等于1意味着该时段有位置变化，而这个变化是和x变量的转移关系耦合的。通常的做法是加一条约束：如果x(m, n, t)和x(m, n, t-1)在同一个节点n上都是1，则MPS没有移动；否则视为移动。在混合整数线性规划（MILP）里，这种“如果-那么”逻辑需要用大M法转成线性不等式约束，不能直接写if语句。

移动时间约束：移动不是瞬时的。MPS从节点i移动到节点j需要travel_time(i, j)个时段。这个约束写起来很考验建模功力。简化的做法是假设一个时段内MPS最多移动一次，且不能在同一时段内移动并给负荷供电——移动中和供电是互斥的。

也正因为这个约束，模型中的时间粒度不能太大。比如把一天切成24个1小时时段，MPS跨节点移动要2小时，那从n1到n2的过程就要跨越两个时段，必须用辅助变量把“途中状态”表达出来，否则求解器会直接让MPS“瞬移”，结果就非常离谱。

我复现的这篇论文里，移动时间是通过预先计算好的travel_time矩阵给出的。所以在代码里我直接用了一个二维参数矩阵T_travel(i, j)表示从节点i到节点j所需的时段数，并用线性不等式将x(m, j, t)与x(m, i, t-T_travel(i,j))关联起来，保证了路径时序上的可行性。

拓扑与潮流约束：配电网韧性研究中，动态调度通常要和“恢复孤岛”的拓扑约束联合建模。简化做法是用DistFlow线性化表达节点功率平衡，同时考虑故障线路的隔离。复现时，如果只关注MPS本身，可以把潮流简化为有功功率平衡加线路容量不等式；如果想要更贴近原文，就要保留含无功和电压幅值的线性化潮流约束。

这里我提一个建议：第一版复现先跑通“位置调度+有功功率平衡+容量约束”的简化版，确定代码能出结果、能复现韧性曲线的走势，再逐步把潮流、无功、电压约束加回去。一上来就上完整约束，很容易因为模型太紧导致求解器找不到可行解，你根本分不清是模型问题还是代码bug。

3. Matlab代码实现：从数据准备到求解器的完整框架

3.1 输入数据准备与参数初始化

代码实现的第一步是数据准备。我通常用一个结构体Config把整篇的参数集中管理，比如：

Config.N = 33; % 节点数 Config.T = 24; % 调度时段数 Config.M = 3; % MPS数量 Config.dt = 1; % 每个时段长度（小时） Config.MPS_capacity = 500; % 单台MPS容量（kWh） Config.MPS_Pmax = 100; % 单台最大放电功率（kW） Config.P_base = 100; % 基准功率（kVA）

节点负荷数据通常是从IEEE标准算例改出来的。我用的案例是简化后的33节点系统，在每个节点乘以一个随时间变化的负荷系数，模拟一天内的负荷波动曲线。故障场景用一个数组记录，例如faulted_lines = [5, 18, 25]，表示这几条线路断开。

预配置结果读入也很关键。我写了一个load_pre_position.m脚本，把上篇输出的MPS初始节点读成向量init_pos，长度为M。在动态调度主程序里，这个向量用来给位置变量x在t=1时段赋初值——注意这里只赋初值，不锁死，否则动态调度就成了纯预配置的搬运工，没有意义。

3.2 用YALMIP构建MILP模型的核心骨架

我复现用的是YALMIP+外部求解器Gurobi的组合。YALMIP的好处是建模语言接近数学表达式，调试起来可视性强；求解器负责真正的求解计算。核心骨架如下：

% 定义变量 x = binvar(Config.M, Config.N, Config.T, 'full'); % 位置变量 E = sdpvar(Config.M, Config.T, 'full'); % 剩余能量 P_out = sdpvar(Config.M, Config.N, Config.T, 'full'); % 各节点放电功率 P_res = sdpvar(Config.N, Config.T, 'full'); % 各节点恢复负荷 % 目标 objective = sum(sum(weight .* (P_load - P_res))); % 约束集 Constraints = []; for t = 1:Config.T for m = 1:Config.M % 位置唯一性 Constraints = [Constraints, sum(x(m, :, t)) == 1]; % 放电功率不超过最大功率 Constraints = [Constraints, sum(P_out(m, :, t)) <= Config.MPS_Pmax]; end for n = 1:Config.N % 恢复功率不超过该节点负荷 Constraints = [Constraints, P_res(n, t) <= P_load(n, t)]; % 恢复功率由各MPS出力之和提供 Constraints = [Constraints, P_res(n, t) == sum(P_out(:, n, t))]; end end % 能量更新约束 for t = 2:Config.T for m = 1:Config.M Constraints = [Constraints, E(m, t) == E(m, t-1) ... - sum(P_out(m, :, t)) * Config.dt - move_energy(m, t)]; end end % 求解 optimize(Constraints, objective, sdpsettings('solver', 'gurobi', 'verbose', 2));

这只是一个框架性代码，真正的论文模型里还需要加入移动状态约束和移动时间约束。YALMIP的优势在于，即使你加了几十条约束，写起来依然接近自然语言，出bug后能比较快地定位是哪条约束导致不可行。

3.3 与预配置阶段结果衔接：如何读入初始驻点

很多复现者会忽略“接口”这个环节，造成两段代码各自为政。我在工程上做了这么一件事：预配置阶段结束后，不仅保存结果变量，还保存了一份result_summary.mat文件，里面包含MPS数量、每台初始驻点、优化目标值、求解耗时。

动态调度主程序在启动时执行：

load('result_summary.mat', 'MPS_init_pos', 'Config'); assert(length(MPS_init_pos) == Config.M, '预配置MPS数量不一致');

然后把MPS_init_pos转成第一个时段的位置约束：

for m = 1:Config.M Constraints = [Constraints, x(m, MPS_init_pos(m), 1) == 1]; end

这一步相当于告诉求解器：第m台车在t=1必须出现在预配置指定的节点。至于之后怎么动，完全交给动态调度优化器自己去权衡。如果你在预配置阶段同时得到了每个时段的预调度计划，也没关系，动态调度时可以只利用初始位置，丢弃后面时段的固定计划——因为动态调度本来就是要根据实时故障信息重新优化。

3.4 求解器选择、求解流程和常见设置

MATLAB自带的intlinprog也能解MILP，但实际算下来，遇到节点数几十个、MPS台数10台以内、时段24个这样的规模，intlinprog勉强能跑；一旦MPS数量增加到十几台、时段细化到48个，求解时间会指数级上升，性能和稳定性都跟不上一线商用求解器。

我用的是Gurobi，通过YALMIP的sdpsettings('solver', 'gurobi')接入。关键参数设置我一般这样调：

options = sdpsettings('solver', 'gurobi', 'verbose', 2); options.gurobi.MIPGap = 0.05; % 5%的gap就停止，速度快很多 options.gurobi.TimeLimit = 600; % 单次求解时间上限（秒） options.gurobi.NumThreads = 8; % 多线程

MIPGap设为0.05意思是：求解器找到的最优解和上界之间的相对差距在5%以内就停止。对韧性研究这种工程问题，5%的误差完全可以接受，但求解时间可能能省下一大半。如果你追求论文里的精确最优值，可以把gap调到0.01，但要做好等很久的心理准备。

求解流程上，我强烈建议按“简化模型跑通—完整模型验证—参数敏感性分析”三步走。第一次先不加入移动能量消耗，把MPS当作“瞬移”的电源来解，快速得到一版韧性曲线。确认曲线形态正确、恢复趋势符合常识后，再逐步加上移动耗时、能量消耗、潮流约束。我就是这么干的，第一次加移动约束后结果全乱，但因为有简化版的对照，很快发现是能量更新约束里把移动耗能重复计算了。

4. 复现过程中的高频问题和排查实录

4.1 维度对不上：节点编号、时间索引和矩阵形状的错误

这个问题出现的频率最高。Matlab的数组索引从1开始，而有些论文的节点编号从0开始。如果你从论文公式直接翻译代码，几乎100%会遇到“Index exceeds array bounds”或“Subscript indices must either be real positive integers or logicals”的报错。

我吃过一次大亏：把节点编号当成数组索引直接用，结果节点0对应Matlab中实际不存在的下标。后来我写了一个映射工具函数：

function idx = node2idx(node_num) idx = node_num + 1; end

所有拓扑矩阵、负荷矩阵统一用idx访问，公式推导时保持论文的0编号习惯，代码内部全部转成1编号。这样既能追溯论文公式，又不会让Matlab崩溃。

4.2 移动约束里的“瞬移”坑：辅助变量与逻辑约束

动态调度最容易出现的“伪合理”结果就是：MPS在t时段还在节点a，t+1时段突然出现在距离老远的节点b，但模型里并没有约束去限制这种瞬移。

原因在于，移动时间矩阵虽然存在，但位置变量x的相邻时段约束没有真正关联。如果你只写了位置唯一性约束，没有写“移动必须耗时”的逻辑，求解器当然会走捷径，让MPS直接跳过去。

要堵住这个bug，我用的是“路径中间状态”思路。具体做法：引入中间变量transit(m, i, j, t)表示第m台MPS在时段t从节点i向节点j移动；通过一系列线性不等式将transit与x关联，并确保移动期间MPS不能同时对外供电。

我建议你在建模阶段先在白纸上画一个时空坐标系：横轴是时段，纵轴是节点。把MPS的轨迹画成一条折线。每一条水平线段表示待命/供电，每一条斜线段表示移动。这样画完之后，每条约束对应哪条线段、哪个交点，一目了然。我最初建模就是靠这种“时空轨迹图”自查，才发现自己漏了移动期间的供电互斥约束。

4.3 预配置结果接口的变量类型与数值精度问题

YALMIP的binvar默认是0/1的二进制变量，但优化完成后，result中拿到的值有可能是double类型的0或1，也可能因为求解器数值精度问题出现1e-6这种微小的非整数值。如果你把这些值直接传递到下一个阶段，或者用来做节点索引，就会报错或让结果出现莫名其妙的偏移。

我实际采用的处理是，在接口处统一做一个舍入：

init_pos = round(value(MPS_init_pos_bin)); init_pos(init_pos < 1) = 1;

这样既保证了节点编号落在有效范围内，又保护了后续索引的安全。还想多说一句：在预配置和动态调度之间，最好把两个阶段的Config都对一遍，MPS台数、总容量、节点数任何一处不一致，下面所有结果都会是空中楼阁。

4.4 求解时间爆炸：big-M、变量初始化与MIPGap调参

变量规模一大，MILP求解时间就很容易失控。我在复现时碰到过一次，把故障场景从单一故障扩展到多个故障后，求解时间从十几秒直接涨到几十分钟，gap还降不下去。

排查下来主要有三个原因。

第一，big-M值设得太大。大M约束是线性化的常用手段，但M太大会让松弛后的线性规划可行域过于宽松，分支定界效率极低。建议M尽量取和问题量级匹配的值。比如供电功率的M，取“所有节点最大负荷之和+一点裕量”而不是随便给个10^6。

第二，没有给求解器提供好的初始解。YALMIP的assign和initial_guess功能可以传一个初始可行解。我的办法是先用一个贪心算法生成MPS最初的分配方案，把它作为MILP的start vector。Gurobi会从这个可行解出发做分支，速度快得很明显。

第三，MIPGap和TimeLimit设置不合理。我之前追求gap=0.01，结果一跑就是半小时。后来调整成0.05，再加上600秒的时限，大多数场景都能在两三分钟内收敛到一个工程可用的解，曲线形态和gap=0.01的结果差别很小。做学术复现不是做竞赛调优，性价比很重要。

4.5 不可行问题：模型约束过紧时怎么定位

另一个常见问题是求解器直接提示模型不可行（infeasible）。这种情况通常是约束之间互相矛盾，但代码里又没有直接提示是哪条约束导致的。

我推荐一个方法：用YALMIP的“约束松弛法”逐个解放约束。比如先注释掉供电容量约束，跑一遍看能否可行；再加回来，注释掉移动时间约束，跑一遍。靠二分法定位到具体是哪类约束导致矛盾。这个过程快的时候几分钟就能定位。

我印象很深的一次：加移动时间约束后一直不可行，排查发现是因为移动时间矩阵中，某些节点对所需的移动时段数大于总调度时段数。MPS从该节点出发，根本不可能在24时段内到达目标节点，模型自然不可行。修正方法是限制MPS初始驻点必须落在合理范围内，或者在调度周期外增加“前置移动时段”，让MPS在灾前就已经开始朝目标方向移动。

5. 仿真结果长什么样：从曲线看图说话

5.1 韧性曲线与恢复进程的判读

跑完代码，第一步不是看数据表格，而是画韧性曲线。典型画法是横轴是时间，纵轴是系统失电负荷总量或累计失电量。曲线会呈现一个“断崖—低谷—爬坡—恢复”的形态：故障发生瞬间失电量飙升，随后MPS逐步调度到位，失电量逐渐回落直到负荷大部分恢复。

我建议同时画三张图：单次场景下的失电曲线、多个故障场景的失电曲线族、以及有无MPS调度的对比曲线。第三条最关键——如果没有MPS或者MPS不动，恢复曲线基本上是一条平线，失电量不会下降；有了动态调度，曲线才有明显的“回勾”。这个视觉差异就是整个模型价值的直接体现。

我在复现时，还特别喜欢把每台MPS的时空轨迹单独画出来。横轴时段、纵轴节点，MPS的轨迹就显示为一条跨节点移动的折线。看到轨迹合理，比如MPS优先向负荷权重高的大用户节点移动，说明模型行为符合预期，比单看目标数值可靠得多。

5.2 MPS数量、容量、速度对恢复效果的敏感性

参数敏感性分析是复现论文精度的关键步骤，也是最容易扩成自己成果的地方。我做了三组分析：

• MPS数量的影响：把M从2依次增加到6，看累计失电量随数量的下降曲线。逻辑上，M越多恢复越快，但边际效果递减。如果你发现MPS增加到某个数量后目标基本不变，说明电网的瓶颈已经不是移动电源不足，而是配电网拓扑限制或线路容量限制了功率输送。
• 单台容量的影响：容量从200 kWh增加到800 kWh，主要影响的是持续供电时长。容量虽然大，但如果MPS不移动，它只能守护一个节点，整体恢复效果有限。这组实验能直观证明“移动”本身的价值，也是论文里比较出彩的一张图。
• 移动速度的影响：把travel_time矩阵整体缩小一半，观察韧性曲线变化。速度提升带来的收益，本质上体现的是MPS的网络响应速度。如果你发现这个指标很敏感，说明论文场景里“抢时间”是关键矛盾。

5.3 从复现到自己的课题：扩展方向

论文复现不是终点，而是起点。我把这套代码跑通后，自己做了两个实际扩展：

第一个是把确定性动态调度改成两阶段鲁棒或分布鲁棒优化——把故障位置、故障持续时间作为不确定参数集合，MPS调度计划在不确定集下依然保证一定的恢复水平。这需要在动态调度模型外面再套一层“最大最小”结构，用列与约束生成算法迭代求解。

第二个是考虑MPS的充电行为。现实中的移动储能车在放电之后需要回到充电桩补能，不能无限循环。我在模型里增加了一个充电状态变量，让MPS在低电量时自动决策“去充电节点补能”还是“继续支撑负荷”。虽然模型规模变大，MIPGap难压下去，但结果更贴近工程实际。

如果你准备往韧性方向发论文，我建议把重点放在不确定性建模和MPS与其他柔性资源（比如配电网重构、需求响应、分布式电源）的协同调度上。这套复现代码里，位置变量、能量变量、功率变量本身就是一个很好地承载协同机制的框架，改起来并不难。

复现心得与最后一点建议

把这篇一区论文的下半部分复现完，我最大的体会是：动态调度模型的难点不在“优化”本身，而在“时空关联约束”能不能建得干净准确。位置唯一、移动耗时、能量变化这三类约束环环相扣，任何一个细节错位，结果都会差出十万八千里。论文里一句“MPS在各节点间的移动时间由道路距离决定”，落到代码里就是一个travel_time矩阵和一连串带大M的不等式。

如果你正准备复现类似论文，我给你几个实操建议：第一，先在一张纸上画时空轨迹图，把MPS运行逻辑彻底想清楚再动手写代码；第二，第一版一定要跑简化模型，快速看曲线形态，确认方向对了再往上堆复杂度；第三，Gurobi配YALMIP是当前性价比最高的组合，别在自带intlinprog上死磕大模型；第四，每一版代码都用系统性方法去排查不可行问题，不要靠运气试。

最后分享一个小技巧：每次跑完优化，把目标值、求解时间、Gap、每台MPS的轨迹长度都打印到日志文件里。参数多轮调整后，回看这些日志能一下看出哪个版本的行为发生了突变。这比调参时凭感觉试要靠谱太多，也省下了大量重复实验的时间。