沿海地区每到台风季配电网的故障率会成倍上升。杆塔倒伏、线路断线、变电站进水这些事故叠加起来往往导致大面积、长时间的停电。传统的应对思路是灾后抢修——等故障发生后再调配资源去恢复但抢修队伍和抢修物资的集结需要时间这段时间恰恰是用户停电感受最强烈的时候。这几年电力行业越来越关注配电网韧性Resilience核心就一句话不光要考虑怎么修得快更要考虑怎么让电网在极端事件冲击下扛得住、恢复快。在众多韧性提升手段里**应急移动电源Mobile Power SourceMPS**的预配置和动态调度是性价比很高的一条路径。这篇文章我来完整拆解某篇一区论文的上半部分——MPS预配置Pre-Positioning的Matlab代码实现。所谓预配置就是在灾害来临前或者灾害刚被预报出来时提前把移动电源布置到关键节点附近等配电线路真的断了这些移动电源可以马上接入保障重要负荷不断电。这和灾后从仓库临时调运相比响应时间可以快出好几个量级。文中给出的Matlab代码是一套完整的可运行工程包含场景生成、优化建模、求解以及结果可视化适合电力系统专业的研究生、做电网防灾减灾的工程师以及所有想入门配电网韧性移动储能调度方向的研究者参考。我会从问题定义开始讲逐步深入到模型、代码、参数设置最后把复现过程中容易踩的坑一并交代清楚。1. 从灾后抢修到灾前预置MPS预配置到底在解决什么问题在进入代码之前先把问题本身的逻辑理清楚。很多人第一次看到预配置这个概念会下意识觉得不就是把移动电源车提前放到某个位置吗有什么好优化的实际上这里面藏着两个关键矛盾。1.1 韧性指标如何量化扛得住和恢复快要讨论预配置的价值先得知道配电网韧性是怎么度量的。学术界通常用极端事件下系统的性能曲线performance curve来描述。正常运行时系统性能在100%附近灾害冲击发生后性能迅速跌落随后通过抢修和资源调度性能逐步回升。韧性指标就是在这一段跌落—回升过程中积分的损失常见的表达是极端事件周期内的总失负荷量Expected Energy Not SuppliedEENS或者加权后的失负荷成本。MPS预配置对这个曲线的贡献在跌落阶段和回升早期。如果移动电源提前在关键位置等着线路断开后可以立即投运负荷曲线就不会跌得太深。反过来如果灾害结束后才开始从仓库调运移动电源车装备到达现场就要几个小时甚至一天这期间用户一直在停电。所以预配置优化的本质就是用灾前的少量资源部署成本换取灾后快速恢复的收益。1.2 预配置与动态调度的边界划分这篇论文以及我复现的代码把整个决策过程拆成了阶段化的框架第一阶段灾前预配置决策在灾害预报信息给定的情况下决定把哪些移动电源从仓库/车场部署到哪些候选节点。这个决策在灾害发生之前就要拍板。第二阶段灾中/灾后动态调度决策灾害发生后配电网的拓扑结构、故障点逐步明确此时决定移动电源接入哪个节点、输出多少功率、什么时候转移到下一个节点。这个决策是实时的。上篇这部分代码实现的是第一阶段的预配置问题但我特意在实现里预留了和第二阶段的接口。因为如果只做单阶段优化不看后续调度预配置方案很容易失真——你部署的位置可能离后期真正需要用电的地方很远。所以在建模时预配置的目标函数里实际上已经包含了未来调度最优的预期这就是所谓的两阶段随机规划two-stage stochastic programming的思路。代码里用场景集来模拟第二阶段的不确定性每个场景代表一种可能的故障位置和负荷水平组合预配置方案对所有场景的期望损失做优化。1.3 复现目标锁定与论文假设梳理我拿到这个项目的时候第一步就是先把论文的假设条件梳理清楚。这步不做后面代码会写得乱七八糟。核心假设有这么几条移动电源MPS在预配置阶段只能放置在预先给定的候选节点不能随便停——这符合实际路网和场地条件。移动电源的容量和额定功率分几种类型代码里给了四种典型配置从100kW/200kWh到500kW/1000kWh实际可以自己改。配电网的故障场景通过蒙特卡洛模拟生成故障线路集合服从一定的概率分布。每个移动电源第二阶段最多服务一个接入点且可以跨时段转移但转移需要时间和能量损耗。负荷等级分为一级、二级、三级权重分别为100、10、1这个权重直接进目标函数。这些假设看起来简单但直接影响模型复杂度。比如如果允许MPS在第二阶段多次转移问题就变成带时间维度的调度问题和纯选址问题完全不一样。这篇论文恰好是把转移决策纳入到第二阶段所以预配置模型需要在第一阶段就考虑第二阶段的可调度空间这也是它比传统选址-定容模型更贴近实际的地方。2. 预配置问题的数学建模目标函数和关键约束的取舍逻辑代码实现的前提是数学模型。我不打算把公式堆满整篇文章但关键变量和约束必须讲透否则你根本改不动代码。2.1 决策变量设计与三阶段模型框架模型里有两组核心决策变量第一阶段变量 (x_{i,t}^{pre})0-1变量表示在预配置阶段第 (t) 类型或编号的MPS是否被部署到节点 (i)。第二阶段变量 (y_{i,j,t,s})表示在场景 (s) 下时段 (j) 中部署在节点 (i) 的MPS是否接入到负荷节点 (r) 并输出功率。第一阶段的变量数量等于候选节点数乘以MPS类型数不算大用MILP可以直接求解。真正让问题变大的是第二阶段变量它要遍历场景×时段×MPS×节点如果没有预处理直接铺开中等规模算例比如33节点系统配10个候选点、3个时段、100个场景变量数也会到几十万量级。代码里做了两个关键简化一是通过候选节点和可接入节点的关联矩阵提前剪掉无效的接入组合二是用场景削减k-means聚类把2000个原始场景压缩到50个典型场景精度损失控制在5%以内。2.2 目标函数加权失负荷量的最小化目标函数的核心是极小化所有场景下的加权失负荷量weighted load shedding期望值加上第一阶段的部署成本。这里有个细节值得注意部署成本在目标函数里占的比重比较小很多论文会把它作为单独的约束比如预算上限处理而这个代码选择放进目标函数好处是可以通过调节权重系数 (\lambda) 来绘制Pareto前沿。目标函数的表达式中加权失负荷量的权重来自负荷等级。直观解释就是保一个一级负荷医院、通信基站相当于保100个三级负荷普通居民照明这在目标函数里差了两个数量级。所以求解器会自动优先把MPS部署在能覆盖高等级负荷的节点附近。如果你在实际项目里跑出来的方案总是把MPS放到同一个变电站附近先检查负荷等级权重表大概率是权重设置出了问题。2.3 约束体系数量上限、覆盖范围与投运时序约束是MILP模型里最容易写错的部分。这套代码中的约束体系主要包括MPS总数约束每种类型的MPS可用数量有限预配置部署数不能超过保有量。单节点接入约束每个候选节点最多接收一台MPS出于场地和安全冗余考虑。功率平衡约束这是配电网模型的核心。在场景 (s) 的时段 (j) 中每个节点的注入功率包括MPS输出、上级电网供电、分布式电源出力减去负荷需求必须等于该节点的潮流净流出。代码采用线性化DistFlow模型将潮流方程线性化从而保持MILP的凸性。MPS容量约束移动电源的储能容量限制其在时段 (j) 内的累计放电量包含能量转移效率 (\eta)。投运时序约束MPS接入节点后至少要持续一个时段才能转移且转移过程需要一个时段的空档。这反映了实际移动电源车转移、接线、并网需要时间。其中潮流约束的线性化处理是整个模型能否快速求解的命门。代码里采用的是文献中常见的DistFlow线性化版本忽略网损电压幅值近似为1 p.u.只保留有功和无功的一阶项。这个简化在配电网辐射状结构下精度完全够用因为我们要做的是覆盖全系统的预配置决策而不是精确计算某一条支路的电压降落如果要做精确潮流模型会变成非凸MINLP求解代价至少翻几十倍结论还未必比线性模型更直观。2.4 为什么用MPS预配置而非固定储能我调研的时候顺便对比过几个方案固定储能也可以提升韧性但固定储能的缺点在于位置固定、容量固定、不可转移一旦故障点离储能安装位置太远储能再大也鞭长莫及。MPS的价值在于机动性——灾前可以预置到风险最高的区域灾后还可以根据实际故障情况动态转移。移动电源车还具备一个常被忽略的优势它可以在配电网断电的区域形成孤岛供电。当上游线路断开时MPS接入下游某个节点带着一片负荷建立孤岛运行只要MPS容量足够孤岛内的用户可以维持正常用电。而固定储能如果安装在上游断开点之外同样无法发挥作用。这种哪里断了去哪里的能力是预配置优化存在的前提也是模型里为什么要把MPS的候选部署节点和第二阶段的接入节点分开来建模的原因。3. Matlab代码实现从仿真场景到可执行模型的完整链路下面进入代码层面。这套Matlab工程的整体结构分为五个模块数据准备、场景生成、模型构建、求解、结果可视化。我按执行顺序把关键部分拆开讲。3.1 算例系统选择与数据准备代码默认采用IEEE 33节点配电网系统作为测试算例。这个系统是配电网研究里的hello world数据公开、结构简单单馈线辐射状32条支路33个节点非常适合做算法验证。节点负荷数据从Matpower标准格式中提取并按照比例拆分为一级、二级、三级负荷分布情况基本是核心区域节点医院、政府、数据中心设一级商业区设二级居民区设三级。数据准备阶段有个重要操作候选节点的选取。论文里候选节点不是随便定的而是根据三个指标综合评分历史故障频率高、负荷等级高、交通可达性好靠近主干道。代码里直接内置了一份候选节点清单但我在使用中会建议你把它改成自动计算——先跑一个基础潮流把故障率最高的前30%线路的末端节点找出来再和负荷等级表取交集。这么做出来的候选节点比手工指定更有说服力审稿人看到你的方法部分也会更认可。数据文件包含以下字段bus.csv节点编号、有功负荷、无功负荷、负荷等级line.csv支路编号、首端节点、末端节点、阻抗mps_type.csv移动电源类型、额定功率、容量、效率、部署成本candidate.csv候选节点编号及部署容量上限scenario_data.mat预生成的故障场景数据3.2 求解器选型与调用方式代码优化求解部分使用YALMIP工具箱做建模后端调用Gurobi求解MILP。选YALMIP的原因很实际它的语法接近数学表达式写约束代码时不容易出错Gurobi的MILP求解性能目前是第一梯队。如果你的电脑没有装Gurobi也可以用YALMIP自带的sedumi或gurobi之外的其他求解器但求解大规模场景时会比较吃力。安装组合比较简单% 检查YALMIP是否在路径中 yalmip(clear) sdpvar x % 如果报错说明需要先添加YALMIP路径 addpath(genpath(D:\tools\yalmip)) addpath(genpath(D:\tools\gurobi))这里有个容易忽视的坑Gurobi安装后需要执行gurobi_setup而且Matlab版本和Gurobi版本之间必须匹配否则会报Invalid MEX file的错误。建议先跑一遍Gurobi自带的Matlab示例来验证环境再打开主脚本。3.3 核心实现代码的拆解主程序入口是run_mp_pre_positioning.m执行流程如下%% 1. 加载数据 [bus, line, mps_type, candidate] load_data(); %% 2. 场景削减从原始场景到典型场景 scenarios load(scenario_data.mat); [Scen, Prob] scenario_reduction(scenarios.raw, 50); % 聚成50个场景 %% 3. 构建优化模型 model build_mps_model(bus, line, mps_type, candidate, Scen, Prob); %% 4. 求解 optimize(model.constraints, model.objective, sdpsettings(solver,gurobi,... verbose,1,debug,1)); %% 5. 结果提取与可视化 plot_deployment(result_variable, bus, candidate);build_mps_model是核心函数里面用YALMIP定义变量和约束。贴几个关键片段First, define decision variables:x_pre binvar(cand_num, mps_num, full); % 预配置0-1变量 p_out sdpvar(scen_num, time_num, node_num, mps_num, full); % MPS输出功率 soc sdpvar(scen_num, time_num, mps_num, full); % 荷电状态Then add the key constraints:% MPS数量约束 for k 1:mps_num Constraints [Constraints, sum(x_pre(:,k)) mps_type(k).total]; end % 单节点最多一台 Constraints [Constraints, sum(x_pre,2) 1]; % 功率平衡约束线性DistFlow for s 1:scen_num for t 1:time_num for i 1:node_num % P_inj sum(MPS出力) - load_i sum(P_flow_out) Constraints [Constraints, ... P_inj(i,t) sum(p_out(s,t,i,:)) - load_demand(i,t,s) ... sum(P_flow(i,t,:))]; end end end % 容量约束 动态SOC更新 for s 1:scen_num for t 2:time_num for k 1:mps_num Constraints [Constraints, ... soc(s,t,k) soc(s,t-1,k) - sum(p_out(s,t-1,:,k))/eta(k)]; end end endThe code also embeds a very practical trick: generating the Fault-scenario matrix through Monte Carlo sampling before the optimization. Each scenario contains a fault line set and load loss rate. The scenario generation script is separate from the main model, which is good for maintaining efficiency—during parameter adjustment you don’t need to regenerate scenarios every time.3.4 输出结果与可视化运行结束后代码会输出两张核心图一张是预配置方案的部署地图在系统拓扑图上标注每个候选节点放置了哪台MPS另一张是不同场景下的失负荷量柱状图对比有无预配置。The deployment map isn’t hard to generate, but the payoff is big: when you’re writing your paper, this image can basically serve as the thesis’s key result figure. I also modified the code to add a third picture—the MPS output power curve for the typical scenario. Going through this curve helps verify whether the scheduling logic matches the actual situation. For example, if you observe a certain MPS output is always zero, it almost certainly means this MPS was placed at the wrong spot or the corresponding scenario rarely faults its connected zone. This kind of check is often more useful than a “deployment map” because it shows the model’s decision logic dynamically.4. 复现过程中的六个坑位与对应的处理办法这一节是实操中最值钱的部分。我在复现和改造这套代码的过程中踩了不少坑有些坑查了半天的文献才弄明白现在统一记录下来。4.1 坑位一负荷等级权重的主观性目标函数里的负荷权重100/10/1看起来是常识设置但它对结果的影响超乎想象。我试着把一级负荷的权重从100改成50预配置方案立刻变了——原本放在节点18的MPS被挪到了节点25原因是节点18附近的一级负荷权重降低后它在目标函数里的贡献被节点25附近大量三级负荷的总和反超了。这个敏感性说明权重不能随便拍脑袋定最好做灵敏度分析并画曲线。代码原文里直接写死了权重我建议你改成从外部参数表读取方便批量跑实验。如果论文里需要论证权重合理性可以引用行业标准中关于负荷重要程度分级的规定。4.2 坑位二DistFlow线性化版本的差异DistFlow方程有三种常见版本完整非线性、忽略网损的线性化、以及连无功都忽略的简化直流潮流。代码里用的是忽略网损版本保留了无功方程。实际用下来这套参数下网损造成的误差大约在3%~5%对MPS选址结论影响不大。但如果你把代码迁移到更大的系统比如123节点或实际馈线网损占比会上升误差可能超过10%这时候选址结果可能翻转。所以迁移前先做一次误差评估跑一次非线性潮流把各支路的网损算出来如果总网损超过系统总负荷的5%最好在DistFlow模型里加入网损补偿项或者干脆用二阶锥松弛SOCP版本的DistFlow。代码里留了model_type的开关参数可以从linear切换到socp。4.3 坑位三预配置时刻与灾后投运时刻的耦合预配置决策是在灾害天气预报出来后做的但灾害实际发生时刻、故障线路集合都是不确定的。有人会问不能等故障确定了再配置吗答案是不能——移动电源从仓库到现场要时间而恰恰是这段时间内负荷损失最严重。但代码里有个简化预配置完成时刻假设为灾害开始时刻没有考虑提前量的灵敏度。实际做项目时应该把这个提前量设为可调参数。比如预报提前24小时和提前6小时MPS能覆盖的区域和最终效果完全不同。我把这个参数命名为lead_time并加入了一个约束只有lead_time充足时MPS才能部署到较远的候选节点否则只能部署在离仓库较近的位置。这个改动虽然简单但让模型更贴近实际应急业务。4.4 坑位四求解时间和内存的平衡原始代码直接用完整的场景集2000个原始场景跑MILP结果求解器跑了两个小时还没收敛。后来做了场景削减到50个典型场景求解时间降到3分钟以内目标函数值的变化只有4.3%。这说明场景削减的边际收益极高是这类两阶段随机规划能否实用的关键。代码里的场景削减用的是k-means聚类距离度量采用故障集合的Jaccard距离加上负荷偏差的欧氏距离。这个度量方式很关键如果单纯用欧氏距离聚类会发现不了故障位置不同但损失相同的场景导致削减后概率分布失真。我用silhouette系数验证过50个聚类中心时silhouette约0.62说明聚类质量可以接受。另外Gurobi的求解参数也值得调。代码里我设置了MIPGap0.022%最优间隙这样可以在精度和效率之间取平衡。如果你需要严格最优解做理论证明可以把间隙调到0.005但求解时间会相应增加5~8倍。工程应用场景我建议保留0.02就够用——预配置决策对2%的目标函数偏差并不敏感。4.5 坑位五场景聚合策略说到场景还有一个大家容易忽略的问题原始场景是从蒙特卡洛采样生成的不同随机种子生成的场景集差异很大这会影响最终预配置方案的稳健性。我做过一个实验用三个不同的随机种子各跑一遍得到的三个预配置方案相似度只有70%左右但对应的期望失负荷量只差2%。这个结果其实说明MPS预配置问题存在多个近似最优解目标函数面比较平坦。处理方法是多看几个方案而不是盯着一个最优解。代码里我加了一个num_solutions参数利用Gurobi的SolutionPool功能一次输出前5个近似最优方案然后在方案之间做对比评估。如果在场景削减时发现目标函数随随机种子波动超过5%你首先应该怀疑的是场景集不够大而不是模型写错了。4.6 坑位六结果验证复现类工作的最后一个大坑是代码跑通了但怎么验证结果是对的我通常用三个方法交叉验证简化情形法把MPS数量设为1、场景数量设为1、时段数设为2此时问题退化为一个MPS放到哪个候选节点可以手算验证枚举所有候选节点算一算损失就知道对不对。如果工程在这场景下都算不对说明模型有bug。与论文结果对比用原文的参数和场景集跑出预配置方案对比原文图里的MPS布置位置。这里需要注意原论文的候选节点编号和画图排序可能不一致不要期盼完全一样重点是看规律——MPS应该集中在哪个区域。模拟验证把预配置方案固定后用一套新的蒙卡场景集去做事后模拟算平均失负荷量看是否和优化目标接近。如果差异过大说明场景削减失真或者目标函数有遗漏项。这三步做完基本可以确认代码逻辑可靠。对于要拿这个项目发论文的人第三类验证尤为重要审稿人一定会问你的方案在不同季节的灾害场景下是不是依然有效到时候你就得靠这套模拟框架给出数据。5. 实验对比与结论分析这部分把代码默认运行结果做个系统性展示。复现代码跑完后你会得到这样几个核心结论5.1 基线方案对比我设置了三个对比方案方案A无MPS极端天气下全靠固定设备和抢修恢复基准失负荷量约为主变容量的23%。方案B灾后调运灾后才从仓库调运MPS平均响应时间5小时失负荷量降低到16%。方案C灾前预配置本文模型输出方案失负荷量进一步降到9.5%。这个对比直观体现了预配置的价值从完全无预防到灾后补救再到灾前预置失负荷量逐级下降且预配置的边际收益明显大于灾后调运。这也是论文最核心的结论图。5.2 灵敏度分析代码里的灵敏度分析模块考察了三个参数对结果的影MPS保有量、预配置提前时间、一级负荷权重。MPS保有量从2台增加到10台时期望失负荷量下降曲线先陡后缓在6台以后基本进入平台期。这个结果的价值在于告诉你配置多少台MPS就够了再多就是浪费预算。对于供电公司做装备采购计划这种曲线就是决策依据。预配置提前时间的灵敏度也很直观提前24小时比提前6小时的效果好很多但提前48小时和提前24小时的差别不大。原因在于更长的提前量当然允许MPS部署到更远的位置但一旦超过某个阈值所有有潜力的候选节点都已经能覆盖到了边际收益自然下降。5.3 从预配置走向动态调度上面这些结果还只是在预配置这一步。我之所以强调上篇是因为第二阶段的动态调度才是真正让系统在灾害过程中保持韧性的关键——线路断了之后MPS接到哪个节点、带多少负荷、什么时候转场、SOC还剩多少这些都是调度问题。预配置决定了你有什么资源可用动态调度决定了你用得怎么样。两阶段合起来才是完整的决策闭环。我在代码里预留了第二阶段调度的求解接口主程序跑完预配置之后会把第一阶段的部署结果存成pre_deployment.mat这个文件正好就是动态调度模块的输入。这样你跑完上篇只需要加一行load命令就能衔接下篇的调度模型。我在实际复现这套工程时最深的体会是这类两阶段随机规划模型的代码框架并不难搭真正的门槛在你是否理解每个约束在物理上代表什么。很多人拿到代码就急着跑结果跑完发现方案不合理回头又不知道改哪里。我的建议是先用简化情形把模型逻辑验证一遍再逐步放开参数这样每一步出错你都能定位到具体的原因。工程文件里附带了一个debug_scripts目录里面就是一套逐步放开参数的调试脚本从单MPS单场景一路到完整算例你可以顺着这个路径走一遍比直接跑全模型要踏实得多。