主动配电网故障恢复：网络重构与孤岛划分的统一MISOCP模型

主动配电网故障恢复这个方向做配电自动化的人应该都不陌生。尤其是分布式电源DG大范围接入之后配电网从原来“变电站单向送电到负荷”的无源网络慢慢变成了一个多电源、可重构、局部还能独立运行的有源网络。故障恢复里的两个核心手段——网络重构和孤岛划分以前大多数是分开做的先重构再划岛或者先划岛再重构。这个“升级版本”做的最重要的一件事就是不再把这两个问题拆开而是统一放进一个混合整数二阶锥规划MISOCP模型里让开关状态和孤岛范围同时被优化然后用Matlab加YALMIP加求解器直接求解。我完整复现了这套代码在IEEE 33节点系统上做了多组故障场景验证也踩了不少坑。这篇博客把模型设计思路、Matlab实现细节、算例结果和调试经验一次性讲清楚适合正在做配电网故障恢复、主动配电网运行优化相关课题的研究生以及想把这套方法落到实际馈线自动化里的工程师参考。1. 项目概述与核心需求解析1.1 为什么要把重构和孤岛划分放进同一个模型先说说传统做法的问题。早期没有DG或者DG渗透率很低的时候故障恢复基本就一件事通过改变联络开关和分段开关的状态把失电区域转移到相邻的健康馈线上这叫网络重构。做法通常是先定位并隔离故障然后从剩余网络里找一个满足辐射状运行要求的开关组合让失电负荷尽可能多地恢复。这个过程本质上是个大规模组合优化问题过去大家常用遗传算法、粒子群算法或者专家规则去搜。DG大面积接入后情况变了。故障点上下游被隔离开之后如果某个区域里挂着DG和负荷在满足孤岛运行条件的情况下这个区域可以不依赖主网自己形成一个微型供电系统继续运行这就是主动孤岛。主动孤岛能显著提升恢复率所以后来很多人单独研究“孤岛划分”做法一般是在重构方案确定之后再去判断哪些DG能带多大的负荷、孤岛边界画在哪里。听起来没问题但实际一算就会发现两个问题是强耦合的。重构改变了网络的拓扑结构拓扑一变孤岛能够覆盖的负荷范围就变了反过来如果先孤岛划分再去重构重构的可行域又被孤岛边界死死限制住。两阶段法在这种强耦合场景下经常得到次优解甚至在某些故障场景下直接求不出可行解。这个升级版本的出发点很直接把开关状态和孤岛决策放在同一个优化模型里由求解器同时决定哪些开关断开、哪些闭合、哪些DG进入孤岛运行、孤岛涵盖哪些负荷。模型虽然比两阶段法复杂但恢复率能做到全局最优而不是“先重构再划岛”的局部最优。而且所有约束放在一起天然不会自相矛盾后续扩展电压约束、DG出力约束、不确定性建模都方便很多。1.2 这个版本解决的实际问题这个“升级版本”不是简单把两个模型拼在一起而是对常见的统一建模方案做了几处关键改进。第一拓扑约束的改进。很多早期统一模型直接用环路枚举或者深度优先搜索来保证辐射状算例一大就失效。这个版本采用基于虚拟流的生成树约束把“网络必须连通且无环”这个图论要求转化成一组线性不等式求解器可以直接处理123节点系统也能跑。第二潮流约束的改进。配电网潮流方程是非凸的直接塞进混合整数模型解不了。升级版本采用DistFlow二阶锥松弛把非凸潮流方程松弛成二阶锥约束配合支路开关状态变量形成一个MISOCP模型。这种模型好处是有严格的数学最优性保证不像启发式算法那样只能给出“还不错”的解。第三孤岛平衡约束的改进。常见模型里孤岛运行时的频率和电压支撑往往被忽略或者简单假设只要DG容量大于负荷就行。这个版本显式考虑了孤岛内部的功率平衡、电压上下限约束并且在目标函数里加入可中断负荷的松弛项避免因为约束太强导致模型无解。1.3 适用人群与使用场景这套代码和思路覆盖的场景很典型配电网发生永久性故障之后调度系统需要在秒级到分钟级内给出一个恢复方案既要保证负荷恢复率又要保证操作安全还要避免过多操作开关。适合的研究方向包括主动配电网故障恢复、配电网自愈控制、DG接入后的运行优化、以及源网荷储协调下的配电网韧性提升。对初学者来说这个项目也是一个很好的混合整数二阶锥规划入门案例能直观看到图论约束、潮流约束、离散变量是怎么在一个模型里协同工作的。如果是现场工程师可能关心的不是模型到底多精细而是能不能快算、结果可不可执行。我的建议是33节点甚至135节点系统用这个模型直接求解没问题但如果是数百节点真实馈线需要考虑模型降维或分解算法后面第5部分会展开聊求解时间的问题。2. 统一模型的核心设计与建模思路2.1 决策变量和网络拓扑表达方式配电网重构的经典数学表达是支路开关状态变量。在这里用z_ij表示支路是否闭合z_ij 1表示闭合z_ij 0表示断开。配电网正常运行时要求辐射状也就是不能有环路。注意配电网的“辐射状”和一般网络图里的“树”不完全一样它允许有一个主电源根节点所有负荷节点必须能从根节点通过闭合支路到达但整体不能形成闭合环。这个“无环且连通”的约束很难用一条简单式子表达常见做法有三种生成树约束、单商品流约束、以及基于可达性的约束。这个版本用的是生成树约束思路非常巧妙给每个非根节点定义一个虚拟需求假设根节点向整张网络注入虚拟流量每条闭合支路上虚拟流量有一个下限根节点总注入等于所有非根节点虚拟需求之和。这样如果网络中存在环虚拟流量就没有办法在环上同时满足有向下限和守恒关系于是环被自动禁止。决策变量方面除了支路开关z_ij还有节点电压平方变量u_i、支路电流平方变量l_ij、支路有功和无功功率P_ij、Q_ij以及DG的孤岛运行状态和出力变量。把这些变量放一起模型规模看起来很大但实际上YALMIP处理这种规模的符号建模完全没问题。2.2 孤岛划分约束如何与重构约束统一这是全篇最容易绕晕的地方。孤岛划分本质上是要把节点集合分成若干子集每个子集至少有一个DG并且子集内部能够实现功率平衡。在统一模型里我没有给每个节点显式打一个“属于几号孤岛”的标签那样会引入太多离散变量。更聪明的做法是通过功率平衡方程和连通性约束隐式实现。具体思路是这样主网恢复区域和孤岛区域之间的断点其实就是一些断开的支路一旦某些支路断开原本从主网流入的功率路径就断了。如果这个断开的区域恰好有一个DG在供电那么该区域内部的节点电压就不再由主网平衡节点决定而是由这个DG决定。模型里把这个逻辑表达成系统总功率平衡方程中主网注入功率只供给与根节点连通的区域孤岛内部由DG出力与负荷匹配。为了让孤岛内部能运行需要确保每个孤岛都有一个“松弛节点”。传统配电网潮流只有一个平衡节点也就是变电站母线出现孤岛后每个孤岛内至少要有一个DG承担电压支撑角色。升级版本里引入了一个虚拟母线把孤岛内所有在线DG的公共连接点连接起来并对同一孤岛内的节点施加相角一致性约束。从数学上看这些约束都是线性的或者二阶锥可表达的不会破坏整个模型的MISOCP结构。这里还要提醒一点孤岛划分的可行性非常依赖DG的容量和调节能力。如果孤岛里只有一个DG那它的有功和无功出力范围必须覆盖孤岛内负荷的变化范围如果有多个DG可以通过下垂控制或协调控制共同支撑。如果模型中不考虑这些算出来的孤岛方案在实际中很可能无法稳定运行。2.3 目标函数的设计与权重处理目标函数是模型“指挥棒”直接决定了恢复方案的倾向。这个版本的目标函数由三部分组成失电损失最小、开关动作次数最少、网损尽量小。但三者不是简单相加必须处理量纲和优先级。第一部分是恢复供电的负荷量最大。实际工程里负荷不是平等的医院、通信基站、重要工业负荷的优先级远高于普通居民负荷。所以这一项要写成加权失电损失即[ \min \sum_{i \in N} w_i \cdot P_i^{\text{load}} \cdot (1 - u_i) ]其中u_i表示节点i是否恢复供电w_i是负荷等级权重。这样求解器会自动优先恢复重要负荷。第二部分是开关动作次数最少。每次都用一个绝对值项|z_ij - z_ij0|其中z_ij0是故障前的开关状态。绝对值在数学上不利于求解但可以通过引入辅助变量进行线性化不影响整体模型。第三部分是网损。注意故障恢复的核心矛盾是尽可能多送电而不是让网损最小。如果把网损权重设得太大求解器可能为了降低网损而牺牲恢复率。我实测下来的经验是网损项权重取失电惩罚权重的1%到5%比较合适纯粹作为一个“打破对称性”的辅助目标存在。升级版本里把目标函数最终写成了“加权失电损失最小 开关动作次数惩罚 小权重网损”这样既保证主目标明确又不会让结果出现一堆无意义的开关抖动方案。3. Matlab代码实现的关键模块拆解3.1 数据组织与算例初始化Matlab实现里最影响效率的往往不是求解器而是数据组织方式。我建议所有算例数据用结构体统一封装不要散落成一堆全局变量。对IEEE 33节点系统核心数据包括支路参数、负荷数据、DG参数、联络开关位置。推荐用边列表格式存储因为后续生成约束和遍历都很方便。%% 基础数据定义 % 支路数据编号, 起点, 终点, 电阻(Ω), 电抗(Ω), 容量(kVA) branch [ 1 1 2 0.0922 0.0470 300 2 2 3 0.4930 0.2511 300 % ... 33节点系统完整支路数据见文献 ]; % 负荷数据节点, 有功(kW), 无功(kVar), 负荷等级(1重要/2一般/3可中断) load_data [ 2 100 60 2 % ... ]; % DG数据节点, 有功上限(kW), 无功上限(kVar), 是否可孤岛运行 dg_data [ 18 400 200 1 % ... ]; % 联络开关和分段开关的初始状态 s0 ones(size(branch,1),1); % 1表示闭合0表示断开这里有一个非常容易踩的坑IEEE 33节点原始数据里电阻电抗的单位是欧姆负荷单位是kVA和kV但有的算例数据会用标幺值。建议统一转换成有名值并由变量名明确单位。我在调试时因为单位混用导致结果完全不对的情况出现过不止一次。3.2 生成树约束的YALMIP实现生成树约束是这个模型的核心。使用YALMIP定义决策变量和约束的代码比较简洁。以虚拟流方式构造生成树约束%% 决策变量 z binvar(nbranch, 1); % 支路状态1闭合 v sdpvar(nnode, 1); % 节点电压幅值平方 l sdpvar(nbranch, 1); % 支路电流幅值平方 P sdpvar(nbranch, 1); % 支路有功 Q sdpvar(nbranch, 1); % 支路无功 w sdpvar(nnode, 1); % 虚拟流变量 %% 生成树约束 Constraints []; M 1000; for k 1:nbranch i branch(k,2); j branch(k,3); % 支路闭合时虚拟流可以传输支路断开时虚拟流为0 Constraints [Constraints, -M*z(k) w(i) - w(j) M*z(k)]; Constraints [Constraints, -M*z(k) w(j) - w(i) M*z(k)]; end % 根节点的虚拟注入为简单起见设根节点虚拟需求为节点数-1 Constraints [Constraints, sum(w(2:end)) nnode - 1]; Constraints [Constraints, w(1) nnode - 1]; % 每个非根节点必须有“源” for i 2:nnode Constraints [Constraints, w(i) 1 - M*(1 - sum(z(find(branch(:,2)i | branch(:,3)i))))]; end上面这段只是示意完整实现时还要考虑支路方向。重点是大M不能随便取。M太大会导致求解器数值病态M太小又可能把可行解切掉。对33节点系统M取1000没有问题但在更大系统上要根据节点数量自适应调整或者干脆用指标约束的方式避免大M。升级版本里我改成了带implies的逻辑约束写法求解稳定性会更好但约束数量略增。3.3 潮流约束与二阶锥松弛DistFlow方程是配电网潮流计算中最常用的支路潮流模型。对于辐射状网络中的一条支路(i,j)其精确形式为[ P_{ij} - \sum_{k:(j,k)} P_{jk} - r_{ij} l_{ij} p_j ][ Q_{ij} - \sum_{k:(j,k)} Q_{jk} - x_{ij} l_{ij} q_j ][ v_i - v_j 2(r_{ij}P_{ij} x_{ij}Q_{ij}) - (r_{ij}^2 x_{ij}^2) l_{ij} ]再加上二阶锥约束[ (2P_{ij})^2 (2Q_{ij})^2 (l_{ij} - v_i)^2 \le (l_{ij} v_i)^2 ]这个二阶锥约束就是整个模型能被称为MISOCP的原因。YALMIP里可以直接用cone函数写也可以用( ... ) norm([...])的形式。对于断开状态的支路需要把该项松弛也就是当z0时支路功率和电流都强制为0但电压关系不再约束。这部分处理不好会让求解器报“约束不一致”排查起来非常让人头疼。for k 1:nbranch i branch(k,2); j branch(k,3); r branch(k,4); x branch(k,5); % 支路闭合时DistFlow方程 Constraints [Constraints, P(k) r*l(k) - ... sum(P(branch(:,2)j)) -load_data(j,2) dg_output(j)]; % 类似写出Q方程和电压方程 % 二阶锥约束 Constraints [Constraints, cone([2*P(k); 2*Q(k); l(k)-v(i)], l(k)v(i))]; % 支路断开时功率为0 Constraints [Constraints, -M*(1-z(k)) P(k) M*(1-z(k))]; Constraints [Constraints, -M*(1-z(k)) Q(k) M*(1-z(k))]; Constraints [Constraints, l(k) M*(1-z(k))]; end注意这里为了表达简单做了简化实际模型里负荷节点的注入和DG出力要分别定义且要乘上恢复状态变量u_i。如果负荷没有恢复即使支路闭合该节点也不应该有功率注入这个耦合关系特别容易漏。3.4 求解器配置与结果可视化Matlab环境里我建议用YALMIP作为建模层求解器选Gurobi或者Mosek。YALMIP安装非常简单下载后添加路径即可。Gurobi需要申请学术许可证但性能比开源求解器好太多。如果暂时没有Gurobi用CBC或SCIP也能跑33节点小算例但求解时间会明显增加而且MISOCP支持不稳定。求解命令如下%% 求解配置 ops sdpsettings(solver,gurobi,verbose,2, ... gurobi.TimeLimit,3600, ... gurobi.MIPGap,0.001); optimize(Constraints, Objective, ops); %% 输出结果 z_opt round(value(z)); P_opt value(P); Q_opt value(Q); v_opt value(v);结果可视化用Matlab的graph和highlight比较直观。把恢复后的闭合支路画出来故障支路用红色标记孤岛区域用不同底色填充。平时做汇报或者写论文的时候一张清晰的网络拓扑图比几百个字都有说服力。%% 网络拓扑可视化 G graph(s_edge(:,1), s_edge(:,2)); p plot(G, Layout, force); highlight(p, fault_branch, EdgeColor, r, LineWidth, 2); highlight(p, open_branch, EdgeColor, k, LineStyle, --);4. 典型故障场景与算例验证4.1 单点故障与多点故障场景构造故障场景设计是算例验证里最容易水掉但其实最重要的一环。最基本的场景是单条支路发生永久性三相短路故障。以IEEE 33节点系统为例假设支路5-6发生故障变电站侧的断路器跳闸故障点下游节点6到18、以及部分联络线路附近节点失去主网供电。此时如果节点8附近有DG统一模型会自动判断这个DG能不能单独带起一部分负荷如果不能全部带起是切掉一部分负荷还是通过闭合联络开关转移到其他馈线具体场景我设置了三种场景1网络中无DG只能靠纯重构恢复作为基准场景。场景2在一个失电区域内配置单个DGDG容量略大于该区域内的重要负荷。场景3配置多个DG并设置DG出力存在10%到20%的波动范围测试模型在不确定性边界下的鲁棒性。多点故障场景则是在两条馈线上同时设置故障比如支路5-6和支路9-10同时断开。这种情况下失电区域更大、可选转供路径更少对模型的全局寻优能力要求更高。两阶段法在这种场景里经常出现“为了恢复一个区域而放弃另一个区域”的问题而统一模型能够同时权衡两个失电区域的恢复优先级。4.2 算例结果怎么看每次求解完成后我建议固定记录五类指标恢复率、重要负荷恢复率、开关动作次数、目标函数值、求解时间。恢复率是所有人第一眼看的指标但重要负荷恢复率更能反映方案质量。有时候整体恢复率差不多但重要负荷恢复率差很多方案的实际价值就差很多。以下是一组基于IEEE 33节点算例的典型结果。注意不同版本代码和不同参数下数字会变但趋势比较稳定场景整体恢复率(%)重要负荷恢复率(%)开关动作次数求解时间(s)是否有孤岛场景1无DG78.480.031.2否场景2单DG91.295.022.8是场景3多DG波动93.597.045.6是可以看到场景2比场景1恢复率高出12.8个百分点这主要就是靠孤岛把原本无法通过重构恢复的远端负荷带了起来。场景3的恢复率更高但求解时间明显增加主要是因为鲁棒约束让二阶锥可行域变小求解器需要更多节点来证明最优性。4.3 与两阶段法的对比为了说明统一模型的优势我特意复现了一个两阶段法作为对比基准第一阶段只做网络重构用遗传算法搜索转供路径第二阶段在剩余失电区域里做孤岛划分判断哪些节点可以由DG带起来。两个阶段各自求解最后合成一个恢复方案。对比结果很能说明问题场景2中两阶段法的恢复率约为86.4%统一模型是91.2%场景3中两阶段法恢复率约为86.8%统一模型是93.5%。差距在6到7个百分点左右。这个差距主要来自两阶段法的“信息损失”第一阶段重构时不知道后面孤岛能带多少负荷可能会把本该留给孤岛的有利路径提前用掉第二阶段划岛时又受限于第一阶段的开关状态无法调整。统一模型没有这个问题代价是求解时间多一倍但对故障恢复这种离线优化场景来说这个代价完全值得。5. 调试经验与常见问题排查5.1 找不到可行解的根源多半是约束打架我调试这个模型时遇到最多的报错是“infeasible problem”也就是求解器说这个模型无解。一开始我还觉得是潮流约束写错了后来逐步排查发现真正的原因经常出在目标函数和约束的边界条件上。举个典型例子如果某个孤岛内的总负荷大于DG容量模型就不可能有解。但实际运行里孤岛可以通过切除部分非重要负荷来维持平衡所以不应该直接把负荷全部作为必须恢复的刚性约束。升级版本的处理办法是给每个负荷节点增加一个恢复状态变量u_i当负荷不恢复时该节点不消耗功率。这样“孤岛内总负荷大于DG容量”就不会导致无解而是让求解器自己决定切掉哪些负荷。这个修改非常关键算是判断一个模型“能不能用于实际”的重要分水岭。另一个导致无解的原因是节点电压约束设置太紧。如果要求所有节点电压都在0.95到1.05标幺值之间而孤岛内DG到末端负荷的电气距离太长末端电压很可能越限。这种情况下要么提高DG无功出力要么切除部分远端负荷。模型里需要把DG功率因数约束放开允许无功输出在一定范围内调节。5.2 求解时间爆炸怎么办33节点系统用Gurobi求解通常几秒就结束但换成135节点或者更多节点时求解时间可能从几秒膨胀到几小时。我踩过的坑主要有三个。第一个坑是大M参数选得太大。所有含大M的约束都会让混合整数规划的线性松弛非常松分支定界树变得很大。解决方法是尽量用指标约束替代大M约束或者在求解前用预求解信息把M收紧。第二个坑是生成树约束写得太冗余。我一开始加了一大堆“任意两个节点之间路径长度”的约束想让模型更精确结果约束数量爆炸预求解就花了半天。删掉冗余约束后求解时间降了一半以上。第三个坑是目标函数里开关动作惩罚的权重太小导致求解器在“多动作一次但网损好一点点”的方案之间反复试探MIPGap很难收敛。把权重上调后求解时间显著减少。5.3 开关动作次数与恢复率如何平衡开关动作次数是现场非常关心的指标。配电开关的机械寿命有限而且每一次操作都有失败风险所以调度员不喜欢那种为了多恢复几十千瓦负荷而连环动作好几个开关的方案。这里有个调参的实用技巧先把恢复负荷的目标函数权重设为1开关动作惩罚设为0求解一次得到最大恢复率下的最小开关动作数k_min。然后以k_min为基准把开关惩罚权重设置为每动作一个开关相当于损失0.1个标幺负荷再重新求解。这样得到的解既不会为了省一次操作而牺牲大量负荷也不会为了多恢复一点点负荷而乱动开关。这个思路比凭经验调权重靠谱得多。5.4 常见问题速查表问题现象可能原因处理办法求解器直接报不可行孤岛内负荷大于DG总容量且负荷被当成刚性约束给每个负荷增加恢复状态变量允许切负荷求解时间几小时不收敛大M参数太大松弛过松收紧M或改用逻辑约束给求解器设置MIPGap结果出现环网生成树约束只加了一半检查所有非根节点是否都有虚拟流需求根节点虚拟注入是否等于总需求孤岛内节点电压越限DG无功调节能力没放开修改DG无功上下限增加电压约束松弛结果开关动作太多开关惩罚权重太小采用“先求最小动作数、再定权重”的方法调整联络开关闭合后线路过载支路容量约束缺失或容量值过大逐一检查场景下的潮流结果确保所有闭合支路容量约束生效6. 个人使用体会与后续扩展方向6.1 做这个项目最大的体会如果只让我说一点体会那就是统一模型不是简单把两个问题塞进一个文件而是要从数学结构上理解它们为什么耦合。网络重构影响的是拓扑连通性孤岛划分影响的是局部功率平衡和电压支撑能力。这两者通过“哪些支路断开”这个共同决策发生联系。建模的时候抓住这根主线后面的约束都是围绕它展开的。另一个体会是Matlab加YALMIP这种组合非常适合快速验证理论模型。YALMIP让符号化建模变得非常顺手不用手搓约束矩阵修改约束时逻辑特别清楚。但代价是模型规模大了之后YALMIP的符号开销会拖慢预求解。如果以后要做数百节点级别的算法研究建议转到Python里的Pyomo或JuMP这类更高效的建模工具模型思路可以完全平移。6.2 后续可以怎么扩展这个模型框架的扩展价值很大。第一个方向是多时段协调。故障恢复不是一瞬间完成的从故障隔离到开关动作再到DG出力调整每段时间尺度都不同。把单时段模型扩展成多时段MISOCP可以合理安排各阶段的操作顺序。第二个方向是考虑DG出力的不确定性采用鲁棒优化或分布鲁棒优化让恢复方案在DG出力波动时依然可行。第三个方向是把模型从“故障恢复”迁移到“正常运行方式的网络重构优化”只需要把故障隔离约束删掉目标函数改成网损最小模型框架可以直接复用。另外如果需要对实际馈线系统建模建议把节点和支路数据自动化导入而不是手动维护。我在做某地区实际馈线算例时看得到最花时间的不是求解而是把GIS导出的CIM模型转成Matlab能用的边列表和节点表。这一步做扎实了整个模型的应用价值能提升一大截。最后分享一个小技巧调试混合整数模型时一定要给决策变量命名加前缀比如z_switch、u_load、dg_on并且把所有负荷的单位统一成kW。不要问我为什么强调这个——我因为单位混用花了整整一个晚上排查一个根本不存在的潮流约束错误。数据和代码的整洁比求解器性能更重要。