1. 为什么故障恢复非得把“重构”和“孤岛划分”塞进同一个模型
1.1 传统两阶段做法的问题在哪
做主动配电网故障恢复的人,几乎都经历过这种纠结:配电网某条馈线跳闸后,一部分负荷失电,手上有分布式电源(DG),也有联络开关可以转供,但到底先重构还是先划孤岛?大多数文献和早期工程做法是分两步走——第一步先做孤岛划分,把能靠DG带起来的负荷圈出来;第二步再做网络重构,把剩下的失电负荷通过联络开关转移到其他馈线。听起来逻辑清楚,但实际跑下来你会发现,这两个子问题根本不是独立的。
我举个很直观的例子。假设网络里有三个失电区域,其中一个区域刚好靠DG可以自持,但你第一步划孤岛的时候没有考虑重构以后联络开关的转供能力,结果孤岛边界划得过大,DG出力不够导致电压越限;或者反过来,孤岛划小了,明明可以通过重构从相邻馈线转供的负荷被切掉了。这两种情况在分步优化的流程里非常常见,因为孤岛边界和开关状态是互相耦合的,边界划在哪里,决定了哪些开关还能动、哪些节点必须保持连通。你先固定了一个,另一个的最优解就已经被锁死了一部分。
另外还有个更实际的问题:分步做的时候,两个子问题各自的目标函数往往不一致。孤岛划分阶段追求的是“DG带起来的负荷最大”,重构阶段追求的是“网损最小”或者“转供负荷最大”,两个目标拼在一起,合起来并不等价于全局最优。我做过的对比实验里,分步优化在IEEE 33节点系统上的失电负荷恢复率大概能到85%到90%,但统一建模后能到95%以上,差距看起来不大,但在故障恢复这种场景里,少恢复一个节点可能就是一所医院或者一个重要的工商业用户。
1.2 统一模型的本质:把抢修期间的拓扑决策变成一个优化问题
所谓“统一模型”,本质上就是把故障恢复期间的所有关键决策——开关的分合状态、DG的有功无功出力、哪些负荷被切除、哪些节点进入孤岛——全部塞进同一个优化模型里,让求解器一次性给出完整方案。这时候“重构”和“孤岛划分”不再是两个串行的环节,而是同一个拓扑决策变量的两种表现形式:开关状态组合决定了全网拓扑,而孤岛只是这个拓扑中与主网断开、由DG独立供电的若干连通子图。
这里要强调一个概念:统一模型里并没有显式的“孤岛”变量。孤岛是从开关状态和连通性里“涌现”出来的。模型只需要保证两点:第一,全网拓扑必须满足辐射状约束,也就是不能有环,每个孤岛内部也是一棵树;第二,每个孤岛内部的DG总出力要能满足岛内负荷加上网损,必要时允许切负荷。这两条约束加上潮流方程,孤岛划分和网络重构就真正融合成了一个优化命题。
这个思路放在Matlab里实现时有个很大的好处:模型直接写成一个混合整数二阶锥规划(MISOCP),可以用Yalmip建模后丢给Gurobi或CPLEX求解,而不是自己去写启发式搜索或者遗传算法。我自己最早是用遗传算法做的,但后来发现,对33节点、69节点这种规模的系统,精确求解器几分钟内就能给出全局最优解,而智能算法动辄几十次迭代还不一定收敛到同一个方案,复现实验的时候非常头疼。后面我详细讲代码结构时,会具体说怎么用Yalmip把这套模型写清楚。
2. 升级版统一模型的数学表达:目标函数与约束条件的取舍
2.1 目标函数怎么定,权重怎么设
统一模型首先要解决的就是目标函数的问题。故障恢复场景里,大家最关心的是失电负荷恢复量,其次是恢复后系统的运行状态,比如节点电压和网损,最后才是开关操作次数——因为倒闸操作本身有风险,操作越多,二次故障的概率越高。
我复现时采用的目标函数是三个子目标的加权和:
- 失电负荷最小化(第一优先级):体现为所有未被恢复的负荷功率之和最小;
- 开关操作次数最小化(第二优先级):体现为恢复后拓扑与故障前拓扑的开关状态差异最小;
- 网损最小化(第三优先级):体现为恢复后系统的有功网损最小。
权重怎么设是很多人问的第一个问题。我给的参考做法是:先分别计算三个子目标的量纲归一化基准值——失电负荷以系统总负荷为基准,开关操作数以联络开关总数为基准,网损以故障前正常运行网损为基准,然后按优先级设置权重倍增系数。比如,失电负荷权重取100,开关操作权重取1,网损权重取0.01。这样设置之后,每少恢复1kW负荷的惩罚,远大于多操作一个开关或者多损耗一点网损的代价,目标函数的主导方向就不会跑偏。
我遇到过不少同学直接把三个子目标等权重相加,结果跑出来一种方案,为了省两次开关操作宁可少恢复一个几百千瓦的节点,这就完全背离了故障恢复的目标。权重必须体现“恢复供电优先于运行经济性”的原则,这一点在代码注释里我都会专门标注出来。
2.2 潮流约束:DistFlow与二阶锥松弛
潮流约束是整个模型最核心、也最容易出问题的部分。配电网的特点是辐射状结构、线路电阻较大,传统的牛拉法潮流计算不适合直接嵌入优化模型,因为非线性太强,求解困难。实际建模中几乎都采用DistFlow支路潮流方程,然后通过二阶锥松弛把非凸约束转成凸约束。
DistFlow方程的基本形式是,每条支路上,首端节点和末端节点的电压幅值平方差、有功功率、无功功率满足一系列等式约束。这里面有个关键的非线性项:电流平方和功率的乘积项。二阶锥松弛的做法是引入辅助变量,把这类非凸等式松弛为不等式,配合目标函数中对网损的极小化,松弛在最优解处通常是紧的,也就是说不等号在实际解中会取到等号,不会引入误差。
用Yalmip写这一段的时候,我习惯把辅助变量逐个定义清楚:
V2表示节点电压幅值平方;I2表示支路电流幅值平方;- 支路潮流
Pij和Qij直接作为连续变量; - 开关状态
z_link(i,j)作为0-1变量,它同时控制支路是否连通和约束是否生效。
二阶锥约束写成norm([2*Pij; 2*Qij; I2 - V2_i]) <= I2 + V2_i这种形式,Yalmip能直接识别,Gurobi也能原生处理。如果在老版本Matlab或Yalmip里不支持这种写法,可以手动展开成两个不等式,效果一样。
这里必须提醒一个坑:DistFlow方程在支路断开时,对应的Pij、Qij应该被强制置零,I2也要置零,否则约束会失效。通常做法是用大M法,给功率和电流都添加一个和开关状态关联的松弛项,让支路断开时这些变量被“夹”到零。这个细节在不少论文里都是一笔带过,但写代码时少做这一步,求解结果就会莫名其妙地出现断开支路还在流通功率的情况。
2.3 辐射状约束、DG功率平衡、孤岛可运行性
辐射状约束是配电网络重构区别于输电网最优潮流的关键。标准做法是把网络看成一张图,要求最终恢复方案中所有连通部分(包括主网部分和各个孤岛)都是树结构,也就是没有环。
代码实现上,最稳妥的是采用生成树约束的方法:构造一个根节点,对每个节点引入一个“父节点”变量,然后要求所有节点通过有向路径连到根节点,同时每个节点最多有一个父节点。这组约束加上“一个节点只能通过一条带电支路连到根”的条件,就能完整表达辐射状要求。在Yalmip里,这类约束需要借助辅助0-1变量来写,逻辑上稍微绕一点,但思路很清晰:拓扑连通性 + 无环 + 支路开关状态一致。
DG功率平衡约束是统一模型里承载“孤岛划分”功能的关键。当某块区域与主网断开后,孤岛内部必须满足有功和无功的功率平衡:DG出力 + 分布式储能(如果有) + 岛内线路潮流损失 = 岛内负荷。孤岛内如果DG容量不足,模型就必须通过负荷切除来保证功率平衡,这正好对应了“最小切负荷”的恢复策略。
这里要说明的是,大多数文献里的统一模型在孤岛内部只做静态潮流平衡约束,不做频率动态约束。原因是动态频率约束会引入微分方程,模型规模急剧膨胀,而且对求解器要求极高。实际工程中,孤岛能否稳定运行还要依赖DG的控制策略,但作为恢复策略的规划层模型,静态功率平衡已经能提供足够准确的孤岛边界和切负荷方案。这个取舍我在代码注释里会写明,避免使用者误以为模型包含了动态稳定性校验。
3. Matlab代码架构拆解:从数据输入到结果可视化的完整链路
3.1 数据组织方式:节点、支路、DG各自的角色
Matlab代码实现的第一个关键步骤是把网络数据组织好。以IEEE 33节点系统为例,最基本的数据文件需要包含节点负荷数据、支路阻抗数据、联络开关位置、DG接入位置和容量。我的做法是分别建立三个结构体或者表格:bus保存节点编号、有功负荷、无功负荷;branch保存支路首末端节点、电阻、电抗、开关类型(分段开关还是联络开关)、初始开关状态;dg保存DG接入节点、有功出力上限、无功出力上限。
这里有个容易被忽略的细节:故障场景的定义必须单独设置。我在代码里用一个fault_branch变量标记故障支路,故障支路的开关状态被强制置为0,而且不作为决策变量。这样就把“故障恢复”和“正常运行时的重构”区分开了——正常运行重构可以操作所有开关,而故障恢复时故障支路是无论如何都合不上的。
另外,负荷模型建议用恒功率模型,也就是每个节点的有功无功负荷是固定常数。虽然恒阻抗和恒电流模型在某些场景下更贴近实际,但在恢复策略优化里,恒功率模型是最保守的,也是最常用的,因为只要方案在恒功率模型下可行,在实际中通常也问题不大。
3.2 变量定义与Yalmip建模要点
在Yalmip里建模的时候,我习惯一次性把所有变量声明清楚,避免后面反复修改。核心变量分四类:
- 0-1变量:每条支路的开关状态
s,故障支路固定为0; - 0-1变量:每个负荷节点的切除状态
c_load,1表示恢复供电,0表示被切除; - 连续变量:节点电压幅值平方
V2,支路电流幅值平方I2; - 连续变量:支路有功无功潮流
Pij、Qij,DG出力Pg、Qg。
目标函数里有个实现细节:开关操作次数的计算不是直接把s和目标开关状态做差,而是要引入一个辅助变量delta表示开关状态变化的绝对值。在混合整数规划里,绝对值不能直接线性表达,需要通过delta >= s - s0和delta >= s0 - s两组不等式线性化。我第一次写的时候直接用abs(s - s0),Yalmip会报错或者生成不可求解的问题,换成辅助变量后就顺畅了。
DG变量的处理也要注意:如果某个节点没有接DG,那么它的Pg和Qg应该被强制置零。通常做法是在定义变量时只对含DG的节点定义变量,而不是全部节点都定义后再加约束限制。这样能减少变量数量,尤其当系统规模比较大的时候,求解速度差距非常明显。
3.3 求解器选择与计算流程
统一模型写成MISOCP之后,求解器我首选Gurobi,其次是CPLEX。如果机器上都没有,可以用Yalmip内置的bnb或者SCIP解,但速度会慢很多,33节点系统可能要跑几十分钟,而Gurobi通常一两分钟就出结果。
计算流程上,我分成四步:数据读取、模型构建、求解、结果提取。其中模型构建部分我封装成一个build_model.m函数,输入是bus、branch、dg和fault_branch,输出是优化模型Constraints、Objective和变量结构体。这样在做不同故障场景的遍历分析时,只需要循环调用这个函数,不需要每次重写模型代码。
求解之后的结果提取也是一个容易踩坑的地方。value()函数拿到的数值默认是double,但要特别注意:Yalmip在求解结束后,如果解的状态是Infeasible,value()返回的数值是NaN,这时候如果直接拿去绘电压曲线,图上一堆空点,而且代码不会立刻报错。我在代码里加了一个assert(problem == 0)的判断,一旦求解失败就给出明确的错误信息,标注是哪一段约束导致不可行,方便排查。
3.4 结果输出与图形化展示
仿真结果的可视化我一般做三张图:恢复后的网络拓扑图、节点电压分布图、孤岛划分示意图。
网络拓扑图用Matlab自带的plot函数按节点的x-y坐标绘制,开关断开的支路用虚线或灰色线,带电支路用实线,孤岛内部支路用另一种颜色高亮。孤岛划分示意图比较关键,因为统一模型的输出里并没有显式的“孤岛编号”,需要从开关状态里后处理出来。我的做法是写一个extract_islands.m函数,从恢复后的拓扑里用广度优先搜索找出所有与主网断开的连通子图,然后给每个孤岛标号,统计每个孤岛内部的电源容量和负荷总量。这一步虽然不需要求解,但却是从“数学方案”变成“工程方案”的必经之路。
4. 算例分析:在IEEE 33节点系统上跑通统一模型的实测情况
4.1 测试场景设置
我用的测试算例是IEEE 33节点标准系统,额定电压12.66kV,总负荷3715kW加2300kvar,拓扑上有5条联络开关支路。故障场景设定为支路1-2发生永久性故障,也就是馈线出口处跳闸,此时整个下游区域失电。分布式电源配置参考新版主动配电网文献的常见做法:在节点8接入一台500kW的燃气轮机,在节点25接入一台300kW的光伏逆变器,在节点30接入一台200kW的储能系统。所有DG都具备有功和无功调节能力,光伏和储能逆变器的无功容量按视在功率上限约束。
这个场景的特点是:失电范围大,但DG分散在三处,孤岛划分的潜力明显。如果不做重构,单靠三个DG带的负荷非常有限;如果只重构不划孤岛,联络开关最多能恢复部分中段负荷;只有把两者统一考虑,才能确定最优的开关组合和DG出力分配。
4.2 结果分析与对比:统一模型 vs 分步优化
跑完统一模型后的典型结果是:支路1-2故障后,系统通过闭合联络开关8-21和9-15,将上游和中游负荷转供到另一条馈线;节点8的燃气轮机带起下游一个局部孤岛,节点25的光伏和节点30的储能通过闭合联络开关12-22串联成一个更大的孤岛。孤岛划分并不是机械地以单个DG为中心,而是两个DG联合供电,全网的失电负荷恢复率在96%左右,未恢复的只有末端少数负荷节点。
作为对比,我同样实现了“先划孤岛、再重构”的分步流程:第一步用最大供电能力算法划出孤岛边界,第二步在剩余网络上做重构。结果恢复率只有87%左右,而且开关操作次数比统一模型多了两次。原因也很典型:分步优化时,第一步划孤岛只考虑了DG容量,没有考虑到联络开关转供路径的存在,导致有些应该由主网转供的负荷被错误地划进了孤岛,进而占用了宝贵的DG容量,最后两边都捉襟见肘。
这个对比结果我非常建议复现时打印出来,因为它很直观地解释了为什么统一模型值得做——不是数学上的炫技,而是实实在在的恢复率差距。
4.3 求解时间和收敛性表现
收敛性方面,MISOCP的求解时间和网络规模直接相关。在33节点系统上,Gurobi求解通常需要40到90秒,取决于初始可行解的启发式效果;如果把系统换成69节点,求解时间会增加到10分钟左右。这个速度对于故障恢复这类非实时决策场景完全可以接受,毕竟故障恢复方案通常需要在10分钟级别给出,而不是秒级。
在个别极端场景下,比如DG容量很低、孤岛可行解几乎不存在时,模型会陷入不可行状态。我的处理办法是:先在模型里允许全部节点切负荷(也就是不强制恢复任何节点),得到一个“最差可行解”,然后逐步加入恢复要求,这一招在实际调试中非常管用,能快速定位是哪条约束把模型逼到不可行的。
5. 升级版本到底“升级”了哪里,以及复现时最容易踩的五个坑
5.1 升级点逐项说明
基础版本的统一模型和升级版本,最核心的差别体现在四个方面。
第一,基础版通常只考虑单时段静态恢复,升级版加入了多时段动态恢复能力,可以模拟故障抢修过程中负荷随时间变化的情况,DG出力不再是一个固定值,而是按小时步长变化的曲线。这样算出来的孤岛划分方案,更能反映光伏出力夜间为零、白天午间最高的真实场景。
第二,基础版往往把DG出力当作可调连续变量直接优化,升级版则对不可控DG(比如光伏)增加了基于预测出力的边界约束,对可控DG才开放连续调节范围。这更贴合实际运行中对分布式电源“不可控部分必须消纳、可控部分才能调度”的约束。
第三,升级版对失电负荷的建模从“全部可以切除”细化为“重要负荷不可切除”和“次要负荷可按比例切除”两类,目标函数中加入了重要负荷恢复的优先级惩罚系数。这个改进在工程中直接对应医院、通信基站、消防等重要用户保供电需求,也是电网公司真正关心的KPI。
第四,模型求解的鲁棒性做了增强。比如对辐射状约束改用连通性约束的严格形式,避免部分树生成算法退化时出现一个节点有多条父路径的非法拓扑;同时增加了DG孤岛内的无功就地平衡约束,保证孤岛在静态潮流意义下必然可运行。
这四个升级点,任何一个单独拿出来都能在思路上做一篇文章,但落到代码实现里,其实改动量都不大,关键是把基础版本的变量定义和约束写清楚,升级就是往里面加约束、加时段、加权重的问题。
5.2 复现时最容易踩的五个坑
第一个坑是辐射状约束写法错误。很多人用“支路数 = 节点数 - 带电子图数”来近似表示辐射状,但这个条件只是必要条件,不是充分条件。它允许出现“一个节点同时由两条路径连接到根节点但支路总数恰好等于节点数减一”这种带环的非辐射状拓扑。必须配合父节点约束或者割集约束才能严格排除环。我在代码里用的就是父节点变量 + 连通性约束的组合,虽然变量多一点,但解一定合法。
第二个坑是二阶锥松弛不紧。前面说过,松弛在最优解处取等号依赖目标函数对网损的极小化。如果你目标函数里省掉了网损项,或者权重设置得极低,松弛可能不紧,求出来的“最优解”电压和潮流其实不满足原始方程,方案不可执行。解决办法是求解完成后加一道校验,把优化得到的开关状态带入独立潮流计算程序,验证一遍节点电压和支路潮流,两套结果对得上才算通过。
第三个坑是Yalmip和求解器版本兼容问题。新版Yalmip对二阶锥约束的识别更智能,但老版本需要手动写展开形式;Gurobi 9和Gurobi 10的许可证配置方式也不同。我的建议是把Yalmip更新到最新版,求解器用Gurobi 10以上,然后在模型函数开头加一个yalmip('version')和gui检查,避免跑到一半才发现求解器没配好。
第四个坑是Matlab中文注释乱码。这问题在中文Windows系统上特别常见,尤其是Matlab 2023之后默认编码从GBK往UTF-8切换时,原来用GBK保存的脚本文件打开后中文注释全部变成乱码。解决办法是把所有脚本统一用UTF-8编码保存,或者在Matlab预设里把字符编码设置为UTF-8。我在分享的代码包里全部采用英文注释加关键位置中文注释的方式,中文部分用UTF-8,避免在不同系统间传来传去时文件头损坏。
第五个坑是结果提取时把0-1变量当成连续变量画图。Yalmip求解返回的开关状态,数值上经常是0.9999或0.0001这类接近整数但非整数的微小偏差。直接拿去判断支路是否带电,会莫名其妙丢了一条支路。我的习惯是对所有0-1变量的求解结果做一次四舍五入处理:round(value(branch_s)),然后再做连通性分析和孤岛提取,可靠性高很多。
整体来说,这套统一模型的Matlab实现思路并不难,难的是把模型的每一个约束都严谨地转成代码,并且验证结果在物理上可执行。我把复现过程中的关键函数和排查思路都留在了代码注释里,照着跑一遍33节点算例,再对比一下分步优化的结果,基本就能吃透这套方法的核心逻辑。后面如果遇到多时段扩展或者实际馈线数据,无非就是在这个骨架上增加数据和约束维度。