配电网韧性研究这几年最热门的方向之一就是应急移动电源MPS在灾害前后如何发挥最大价值。而MPS预配置作为整个“预配置动态调度”体系的第一步恰恰是最容易被初学者轻视、却最影响复现结果的一环。我翻了不少SCI一区论文也自己动手用Matlab复现过类似的模型发现很多人卡住的地方并不是动态调度而是第一步预配置就没弄明白MPS应该提前停在哪、为什么停在那、这个决策和灾后的调度是什么关系。这篇先把(上)写透——只讲MPS预配置的建模思路和Matlab实现动态调度细节留到(下)再展开这样两篇合起来正好构成一个完整的两阶段随机规划复现闭环。先说这系列文章的适用对象正在复现配电网韧性相关论文的研究生、想用Matlab做灾害场景模拟的工程师、以及需要给移动储能电站做前期规划方案的同行。这篇文章不会只丢一个“论文模型”给你而是会拆清楚模型每一层约束背后的物理含义再给出可以落地运行的Matlab代码框架。看完之后你至少能回答三个问题预配置为什么是两阶段问题里的第一阶段预配置目标函数到底在优化什么怎么用Matlab把预配置模型搭出来并求出合理的方案。1. 研究问题拆解配电网韧性提升与MPS预配置到底解决什么问题1.1 配电网韧性从“停电”到“快速止血”传统配电网可靠性研究关注的是一般性故障某条线路跳闸、某台配变过载这类故障持续时间短、影响范围可控常规保护与调度就能应付。但配电网韧性Resilience面对的完全是另一类问题——极端灾害场景下多条馈线可能同时跳闸变电站出线甚至母线都会失压整个配电网可能陷入大面积停电。台风的狂风暴雨、冬季冰灾的线路覆冰、地震造成的杆塔倒伏都属于典型的“高影响-低概率”事件不可能用平时那套可靠性指标去评估。学界描述配电网韧性时最常用的工具是韧性曲线横轴是时间纵轴是系统剩余供电能力或恢复负荷比例。一场极端事件发生后曲线会先快速跌落然后随着故障隔离、应急供电、抢修恢复逐步回升。韧性研究的核心目标就是要把这条曲线的“跌得慢一些、谷底浅一些、回升早一些”给做出来。MPS预配置就是典型的事前措施在灾害发生前的预测预警期把有限的移动电源布置到合理节点为灾后快速建立应急供电走廊做好物理准备。这里有必要定义一下MPS。在不少论文里MPS写作Mobile Power Source也可以理解为移动储能车MES、移动柴油发电车或者广义的“可移动能量源”。它可以自己带电量也可以从变电站或健康馈线取电后转运到失电区域。MPS的数量、容量、移动速度都是有限资源所以预配置不可能覆盖所有负荷只能针对重点区域做布局。1.2 MPS预配置与动态调度两阶段问题的时间分工从时间维度看一场典型极端天气过程可以粗略分成三个阶段预测预警期灾害发生前几小时到几天、灾害发展期几小时、恢复期几小时到几天。MPS预配置发生在预测预警期动态调度发生在灾害发生与恢复期。预配置解决的是“车提前停在哪”的问题。动态调度解决的是“灾后车怎么走、接在哪、功率怎么出”的问题。这两个决策在数学上正好构成了两阶段随机规划的骨架预配置位置是第一阶段决策变量也叫here-and-now变量必须在灾害场景真实发生之前确定移动路径、接入位置、充放电功率、负荷恢复量都是第二阶段决策变量也叫wait-and-see变量需要等具体的故障场景信息暴露后再决定。这里有一个初学者极容易踩的坑既然灾后还要动态调度那预配置是不是随便选个点就行答案是完全可以反过来理解——预配置位置决定了动态调度的初始终端集合也就决定了灾后调度的可行域上限。如果预配置点选在一段极易被台风刮断的馈线末端灾后这条路一断车就被困在孤岛里连出车的机会都没有。如果预配置点太分散灾后移动距离拉长车辆在路上白白消耗恢复窗口。所以预配置是在“拓扑位置可达性”和“负荷需求空间分布”之间做风险对冲跟“就近原则”完全是两回事。1.3 为什么“预配置”值得单独开一篇从代码复现的角度讲很多论文把预配置和动态调度一起塞进一个两阶段模型拿到手的第一步就是变量爆炸——根本分不清哪个变量是一阶段的哪个变量是每个场景一套的二阶段变量。我最初复现时也花了很久才理顺。这篇先把预配置单独拆出来目标很明确把阶段一彻底吃透。目标函数长什么样、约束怎么建、为什么这么建、Matlab里对应代码如何组织都在这一篇里交代清楚。把预配置这个模块整明白再看动态调度那一篇就是水到渠成的事。2. 数学建模精读预配置阶段的目标函数与约束体系2.1 目标函数在不确定性场景下最小化加权负荷损失预配置优化的标准数学框架是两阶段随机规划。记第一阶段预配置决策为向量x第二阶段在场景s下的调度决策为y_s目标函数可以写成min C_place(x) Σ_s ω_s · Q(x, s)其中C_place(x)是预配置本身的动作代价比如车辆从驻地出发到达预配置点的移动成本如果为便于简化也可以直接设为常数或0ω_s是场景s的发生概率Q(x, s)是在预配置方案x已经确定、场景s的故障线路已知的前提下灾后动态调度所能实现的最小加权失负荷量。Q(x, s)的内部形式是Q(x, s) min Σ_t Σ_n w_n · L_curt(n, t, s)这里的w_n是节点负荷重要性权重L_curt是未恢复被削减的负荷量。权重设置非常关键医院、通信基站、消防设施这类敏感负荷权重可以设为10甚至更高普通居民负荷权重就是1。韧性提升和传统经济调度的目标差异就在这里——它追求的不是电价最低或者网损最小而是关键负荷的连续供电时间最长、累计失电量最小。理解这个目标函数需要抓住一点预配置不是针对某一个确定的灾难场景做优化而是对所有可能场景的期望表现负责。这个“期望值”恰恰是两阶段随机规划区别于确定性选址的核心。比如某节点在单一场景下可能是最优位置但它在其他几个高概率场景里完全不可达那么它的期望得分就很低大概率不会被选中。2.2 第一阶段约束预配置的唯一性与离散性预配置阶段的约束看起来很简单只有三类每个MPS只能预配置在一个节点Σ_n x_v,n 1对每辆车v成立。每个节点最多停靠一辆MPSΣ_v x_v,n ≤ 1对每个节点n成立。MPS可用总数限制Σ_v Σ_n x_v,n N_mps。这三条约束保证了预配置方案是一个“一对一”的离散匹配关系本质上是组合优化里的匹配问题。第二类约束每个节点最多一辆很容易被初学者漏掉。漏掉之后模型可能会把两辆车都放在同一个变电站节点上虽然目标函数值可能好看但物理上完全不合理——同一个节点多放一辆车并不会带来更多供电能力只会浪费资源。为了在空间上形成覆盖这个约束必须加上。这里我提一个经验性的补充解释为什么预配置阶段不做移动时间约束因为预配置发生在灾害之前车辆从驻地到预配置点有充足的时间窗口通常以小时计完全可以在目标函数里用固定成本处理或者干脆令其为0只保留位置决策。把移动时间约束放在第二阶段才是符合物理规律的。2.3 场景内约束移动、容量、潮流三类缺一不可预配置方案的质量必须通过第二阶段的场景内约束来检验。因为历史原因配电网韧性论文里动态调度的约束体系相对统一主要分成三类。第一类是移动约束。车辆v在时段t的接入位置Z_v,t必须满足空间可达性约束。设预配置节点为n0那么灾后初始时刻t1车辆必须出现在n0。如果后续要转移到另一个节点n1必须满足从n0到n1的移动时间要小于等于t-1个时段的累计时长。实际操作中一般会提前计算一个可达性矩阵reach_time(i,j)表示从节点i到节点j的最短移动时间然后写成整数约束。第二类是容量与功率约束。每辆MPS的电量SOC_v,t要在上下限之间充放电功率不能超过额定功率且充电与放电不能同时进行。这个“不能同时”在数学上是非凸约束通常用大M法线性化p_ch ≤ M·γ、p_dch ≤ M·(1-γ)其中γ是0-1辅助变量。还需要注意充放电效率实际SOC变化是SOC_v,t1 SOC_v,t η_ch·p_ch - p_dch / η_dch效率系数如果取0.95左右仿真结果才能贴近真实储能设备的特性。第三类是潮流约束。配电网在恢复过程中通常要把网络盘活所以用的是支路潮流模型。为了可计算性绝大多数论文都用线性化的DistFlow方程V_j(t) V_i(t) - 2·(r_ij·P_ij(t) x_ij·Q_ij(t))P_ij(t) Σ_k P_jk(t) P_j_load(t) - P_j_MPS(t) - P_j_DG(t)其中P_j_MPS(t)是MPS注入该节点的功率P_j_DG是分布式电源出力。这里需要特别提醒如果只是做选址预配置可以先用直流潮流或忽略电压约束的版本做一个快速判断但如果要把结果写进论文必须回到带电压约束的DistFlow模型否则优化方案可能在弱网场景下过于乐观——因为实际配电网节点电压偏差是会限制MPS实际可输出功率的。2.4 线性化处理与模型简化能省则省两阶段随机规划天然是混合整数线性规划MILP求解规模本身就很大。建模时要注意把一切可以线性化的环节都做掉充放电互斥条件用大M法。目标函数里min-max结构通过把min放到外层、内层Q(x,s)作为各场景最小失负荷量不需要额外线性化。潮流方程直接在DistFlow线性化版本下构建不做非线性交流潮流。如果有联络线开关变量也需要加0-1整数变量但这类变量数量过多会严重影响求解速度。实操经验是能用连续变量表达的绝不用整数变量能用线性约束表达的绝不用非线性约束。预配置本身是整数变量避不开但第二阶段的移动路径选择、MPS接入位置也是整数变量这一块如果能压缩就需要尽量压缩。比如移动路径中关于“车辆不得越过故障线路”的空间拓扑约束可以通过预计算可达集合并用0-1系数矩阵参与约束避免引入路径级联变量。3. Matlab代码实现核心模块搭建与求解策略3.1 整体代码架构数据、场景、模型、结果四层分离复现这类模型我最推荐的方式是把代码拆成几个独立脚本不要全都堆在一个main文件里。原因很简单两阶段随机规划的调试周期长中间任何一层出错都会导致结果莫名其妙模块化之后才能逐层定位问题。一个典型的目录结构是main.m算例配置和运行骨架。case_ieee33.m拓扑、线路参数、节点负荷、MPS参数。scenario_generation.m蒙特卡洛生成故障场景并做场景削减。mps_pre_config_model.mYALMIP建模、求解、返回预配置位置和期望指标。plot_resilience.m输出韧性曲线和对比图。main.m只需做三件事加载数据、生成场景、调用模型脚本。具体业务逻辑全部放到子函数里。3.2 场景生成模块蒙特卡洛抽样与场景削减配电网韧性研究里灾害场景的生成是决定结果可信度的底座。最常用的方法是按线路故障概率做蒙特卡洛抽样。比如给每条线路设定一个故障概率p_break然后对每条线路独立采样一个0-1变量1代表故障。这样做S次就得到S个场景。场景数据结构建议用struct数组每个元素存故障线路索引、故障时段范围、以及发生概率。% 场景生成按独立故障概率抽样 S 200; % 初始场景数 K 12; % 削减后的代表场景数 prob_line 0.06; % 每条线路故障概率示例值 scenario_list struct(broken, cell(1,S), prob, 0); for s 1:S is_break rand(num_line, 1) prob_line; scenario_list(s).broken find(is_break 1); scenario_list(s).prob 1 / S; end % 用k-means做场景削减减少代表性场景数量 [scenario_red, ~] reduce_scenario(scenario_list, K);我习惯先抽200个场景然后削减到10~15个代表场景再用削减后的场景做两阶段优化。削减算法不复杂把每个场景的故障线路集合编码成一个one-hot向量做k-means聚类取每个聚类中心最接近的真实场景作为代表该代表的概率等于簇内所有场景概率之和。这步做完求解规模能降一个数量级而且对预配置结论影响很小。3.3 YALMIP建模变量定义、约束组装、目标函数Matlab里做这类优化我基本离不开YALMIP工具箱。用YALMIP有两个好处一是变量定义直观二是不用关心底层求解器的接口细节。核心建模逻辑如下% 载入数据 data load_case_ieee33(); N data.num_node; M data.num_line; V 3; % 3辆MPS T 24; % 24时段 % ---------- 第一阶段变量 ---------- x_place binvar(V, N); % 每辆车是否停在节点n % ---------- 第二阶段变量每个代表场景一套 ---------- for k 1:K z{k} binvar(V, N, T); % 时段t车辆接入节点 p_dch{k} sdpvar(V, T); % 放电功率 p_ch{k} sdpvar(V, T); % 充电功率 soc{k} sdpvar(V, T); % 荷电状态 L_rec{k} sdpvar(N, T); % 恢复的负荷量 end Constraints []; % 每辆车只能选一个预配置节点 Constraints [Constraints, sum(x_place, 2) ones(V,1)]; % 每个节点最多停一辆 Constraints [Constraints, sum(x_place, 1) ones(1,N)];第二阶段约束的关键是“预配置位置与灾后初始接入位置”的耦合。用x_place直接约束z的第1个时段% 场景k中车辆在t1必须在预配置节点上 for k 1:K Constraints [Constraints, sum(z{k}(:, :, 1), 3) x_place]; end这一段是整篇代码最容易出错的地方。必须注意维度匹配z是一个(V,N,T)的三维数组sum(z(:,:,1),3)会得到(V,N)矩阵让它等于x_place就保证了“车停在哪就必须从哪出发”。如果漏掉这条约束预配置变量和场景内调度变量就完全脱节了算出来的结果毫无意义。目标函数要按场景概率加权Objective 0; for k 1:K L_curt_k data.L_demand - L_rec{k}; % 失负荷矩阵 (N,T) loss_k sum(sum(data.w .* L_curt_k)); % 加权失负荷 Objective Objective scenario_red(k).prob * loss_k; end ops sdpsettings(solver, cplex, verbose, 2, miprelgap, 0.01); result optimize(Constraints, Objective, ops); x_best value(x_place);如果Cplex或Gurobi没装YALMIP自带的内置求解器也能跑小规模问题但要求解MILP还是建议装一个商用求解器学术许可证通常是免费申请的。我实测过同一模型在Cplex和Gurobi下求解速度差别明显预配置问题在IEEE 33节点系统上Gurobi往往比Cplex快20%到50%。3.4 求解策略当大规模MILP跑不动时的替代方案两阶段随机规划全变量展开后变量数量会非常大。如果不想一上来就碰大规模优化有一个非常实用的替代方案枚举预配置组合 固定第一阶段后的子问题求解。当节点数较少、MPS数量也不多时比如IEEE 33节点系统、3辆MPS预配置的组合数量是C(33,3)5456。这个数字虽然看起来不小但每组合只需要解一个第二阶段的确定性等价问题完全可以用parfor并行跑完。这样得到的预配置方案是全局最优的因为枚举覆盖了所有可行组合。代码框架大致是combos nchoosek(1:N, V); % 5456行每行一个三节点组合 n_combo size(combos, 1); obj_each zeros(n_combo, 1); parfor c 1:n_combo x_fixed zeros(V, N); for v 1:V x_fixed(v, combos(c, v)) 1; end % 用写的子问题函数求期望失负荷 [loss, ~] eval_subproblem(x_fixed, scenario_red, data, params); obj_each(c) loss; end [min_loss, idx_best] min(obj_each); best_combo combos(idx_best, :);这个“枚举子问题验证”的思路本质上就是把两阶段问题的第一阶段摘出来做粗暴搜索第二阶段用线性规划少数整数变量快速求出期望最优值。它最大的好处是不需要处理复杂的分解算法也不用担心大规模MILP求解器跑一天还没收敛代码调试极其友好。当节点规模升级到100以上、MPS数量到5以上时再用Benders分解或L-shaped方法也不迟。4. 实验设计与结果可视化4.1 算例配置以IEEE 33节点系统为默认测试床做配电网韧性仿真的默认算例几乎绕不开IEEE 33节点系统。它是8节点系统的经典扩展版12.66kV33个节点32条支路5个联络开关总负荷约3.7MW。在做应急移动电源研究时我习惯在原始拓扑基础上做如下配置节点负荷采用一天24小时的归一化曲线峰值负荷按各节点原始额定值乘以1.2倍放大给“削峰填谷”留空间。3辆MPS每辆额定功率100kW容量300kWhSOC运行范围0.1~0.9初始SOC设为0.5。车辆平均移动速度取20km/h节点间移动时间按地理距离折算向上取整到1小时。场景的故障线路集从10条易损线路中抽取每条线路故障概率设为0.1模拟台风路径覆盖区域。联络开关在故障发生后允许闭合以形成额外转供路径。这些参数是我复现同类论文时整理的一组“接近实际、又方便计算”的值。如果你只是为了对论文结果先按原论文参数跑如果是为了验证自己代码逻辑用上面这组默认值就足够。4.2 对比方案设计无MPS、随机预配置、优化预配置、事后最优有了预配置模型之后如何证明它的价值光看优化结果不够必须做对比实验。我常用的对比方案有四种无MPS完全不做应急移动电源参与反映系统基础韧性水平。随机预配置在33个节点里随机选3个位置放MPS重复做20次取平均值代表“没有规划、凭经验放”的结果。优化预配置用两阶段随机规划求出的方案。事后最优假设故障场景已知再回代求解最优预配置位置这是理论上的“上帝视角”上界只用来衡量预配置优化的极限增益。评价指标用两个就够了累计加权失负荷量相当于失电能量损失单位kWh和恢复比例曲线下面积恢复比例对时间积分越接近1说明韧性越好。4.3 韧性曲线绘制与解读绘制韧性曲线时我习惯将纵轴定义为“系统恢复供电比例”即当前时段恢复负荷除以全部负荷。代码很简单figure; plot(t_hour, ratio_opt, b-, LineWidth, 1.6); hold on; plot(t_hour, ratio_random, g--, LineWidth, 1.4); plot(t_hour, ratio_no_mps, r:, LineWidth, 1.6); xlabel(时间/h); ylabel(恢复供电比例); legend(优化预配置, 随机预配置, 无MPS, Location, southeast); grid on;真正有价值的是读曲线时的判断逻辑。通常优化预配置方案会在灾后早期第2到第4小时窗口就把恢复比例拉起来而随机预配置可能会在前期“干瞪眼”一两小时因为车辆要么离故障点太远要么所在区域本身就失电。这种差异反映的不只是调度能力更是“预配置位置”决定的时间窗口优势。所以预配置优化的核心收益往往不体现在最终稳态稳态时大家都恢复而是体现在恢复过程的“响应速度”上。5. 常见问题与踩坑清单5.1 建模阶段的高频错误两阶段变量层次混淆是最常见的问题。新手容易把预配置变量写成每个场景一套导致第一阶段变量变成了场景相关模型结果不可解释。记住所有与场景下标k相关的变量都是第二阶段变量第一阶段变量必须不带k。目标函数漏乘场景概率也很常见。如果只对每个场景求和而不加权高概率场景和低概率场景的影响就一样了预配置方案会严重偏向那些概率极低但破坏极大的场景得到的方案“看起来很保守实际很浪费”。负荷权重w设置不统一是复现论文对不上号的一大原因。很多论文对权重默认全部取1但实际应用里一定是有差别的。我在复现时建议在说明文档里明确写出哪些节点是医院/应急负荷权重多少哪些节点是普通负荷权重多少。否则你跑出的韧性与原论文有差异先不怀疑代码先检查权重表。5.2 求解性能优化与数值问题如果初始场景数取300个又做了k-means削减到15个整体规模合理Cplex求解通常几十秒到几分钟不等。但如果场景削减没做好或者MPS数量和时段数偏大MILP求解时间就可能呈指数增长。遇到这种情况先检查是不是整数变量太多了。大M系数也容易踩坑。取值过大比如10^8会造成数值病态求解器可能给出满足容差但物理意义错误的解取值过小比如只是功率上限的2倍又可能约束不够松。我习惯把M设为“相关变量理论最大值的10倍”比如放电功率上限100kWM取1000就够再往上没必要。还有一个容易被忽视的点是求解器的整数容忍度设置。Cplex默认的整型变量容忍度有时候太紧导致求解很慢。我一般会把miprelgap设为0.01把IntFeasTol放到1e-5级别求解速度能提升不少而结果差距通常小于千分之一。做科研复现完全可以接受。5.3 复现时容易忽略的细节移动时间约束的参数设定直接影响动态调度的可行性。我最初复现时忽略了这个约束直接把车辆从一个节点瞬间转移到另一个节点结果预配置方案特别激进——车辆可以任意飞到任何地方导致模型给出“某车初始配置在A节点灾后瞬间出现在B节点”这种物理不合理的调度。后来加上reach_time矩阵约束后预配置方案立刻就变得更保守、更贴近实际。联络开关的状态设置也是一个隐藏变量。IEEE 33系统的5个联络开关默认是断开的但如果允许故障后闭合联络开关相当于给恢复过程增加了很多转供路径系统的韧性基线会大幅抬高。这也意味着MPS预配置的“增量价值”会变小因为网架自己就能依靠联络线恢复不少负荷。论文里两种假设都存在你复现时必须看清原论文是“固定开环”还是“允许动态重构”否则结果差异会非常大。场景故障线路的设定方式也需要看仔细。有的论文假设所有线路都有相等故障概率有的论文只对特定区域线路设定故障概率。我在复现时更倾向于用“故障场景由线路故障概率抽样生成”的方式因为它可以灵活模拟台风范围也让预配置模型真正面对不确定性。如果想要复现一篇具体论文建议先核对它用的场景生成方式再决定自己怎么抽。我自己的体会是MPS预配置这个阶段看似只是“放几个点”但它几乎决定了后面所有内容的上限。预配置放得准灾后调度哪怕用最朴素的最短路优先策略韧性指标也不会差预配置放得偏后面调度模型再漂亮也只能在给定的初始点上小修小补。所以复现这类论文我强烈建议你在预配置阶段多花时间先吃透两阶段结构再理解约束背后的物理动机最后才动手写Matlab代码。顺序反了就会在一个个数值怪癖里反复横跳。接下来准备写(下)的时候我会把动态调度部分遇到的两个有意思的问题也整理出来算是给这个系列的收尾。