做配电网电压无功优化VVO的人一开始往往只盯两个数网损和电压偏移。我也是这么过来的。直到有次仿真跑完方案确实是全网网损最低、电压全部合格但一天96个时段里电压主变有载分接开关OLTC动作了二十多次并联电容器组投切超过三十次。按设备手册里的机械寿命估算这种动作频率撑不了几年就要进大修。从那之后我在任何VVO模型里都会把“设备动作损耗”写进目标函数并且把求解方式从集中式改成了分布式。这篇内容就围绕“考虑设备动作损耗的配电网分布式电压无功优化”来展开梳理建模逻辑、ADMM求解框架以及一套能在Matlab里落地的代码套路。适合看这篇内容的人很明确在做配电网无功优化课题的学生、搞配网自动化工程建设与运维的工程师、以及在Matlab里写过或准备写优化算法的开发者。文中不会堆一堆脱离工程的理论重点是我实际调试这类算法时反复踩过的坑和验证过好用的处理方式。1. 为什么配电网电压无功优化必须把设备动作当成成本来算1.1 频繁动作不是小事OLTC与电容器组的机械寿命账有载调压变压器最金贵的部件就是有载分接开关切换一次分接头机构内部要完成一系列机械动作储能弹簧压缩、释放、触头切换、过渡电阻通断。这一套动作一次触头就磨损一次绝缘油也在电弧作用下慢慢劣化。并联电容器组用的接触器或专用开关同样有机械寿命和电寿命限制。工程上通常用“动作次数上限”来管理这类设备。某地区负荷峰谷差大又赶上分布式光伏春夏季节出力波动快一天内电压从上限附近掉到下限附近非常常见。如果优化目标里没有对动作次数的惩罚求解器会怎么省事怎么来——反正在多个时段里改变OLTC档位能降低网损那就变。结果就是前面说的一天二三十次动作。集中式无功优化算法在纯学术算例上可能看不出问题拿到真实调度环境里调度员首先问你的就是设备还能这样折腾几年1.2 传统VVO目标函数缺了哪一块传统VVO的优化目标通常写成[ \min \quad \sum_{t} P_{loss,t} \alpha \sum_{t} \sum_{i} (V_{i,t} - V_{ref})^2 ]第一项是网损第二项是电压偏移惩罚。这个目标本身没错但它是“无记忆”的t时段的决策只取决于t时段的运行状态完全不关心上一个时段设备做了什么。加了设备动作损耗之后目标函数要改写成[ \min \quad \sum_{t} P_{loss,t} \alpha \sum_{t} \sum_{i} (V_{i,t} - V_{ref})^2 \sum_{k} c_k \sum_{t} |u_{k,t} - u_{k,t-1}| ]其中 (u_{k,t}) 是第k组离散设备OLTC档位或电容器组数在t时段的状态(|u_{k,t} - u_{k,t-1}|) 就是动作次数(c_k) 是折算到目标函数里的单次动作等效成本。加了这一项求解器在让设备动作之前必须先掂量一下这个动作带来的网损下降到底值不值一次机械磨损。1.3 设备动作成本怎么定价才不会失真这是建模里最容易被忽略的一步。动作成本 (c_k) 必须有量纲而且要和网损项的功率/能量单位对齐。我的做法是折算成“等效网损成本”。思路很简单一次OLTC动作给设备寿命造成的损失大约相当于设备在电网上多带几千瓦负荷运行一小时造成的经济损失。比如一台主变的OLTC单次动作折合运维和寿命成本约20元那就把它看成一个“每次动作消耗20元的虚拟功率损耗”电容器组单次动作折合5元那就当成5元。这样目标函数里所有项的单位统一成元网损项用实时电价乘电量动作项用等效成本乘次数加起来就可以直接比较了。如果你不给动作成本赋一个明确数值等于在目标函数里放了一个没有校准的权重结果要么动作惩罚太弱设备乱动要么太强设备一个档位都不肯变电压越限。我调试时会把网损成本和动作成本打印出来逐时段对比确保两者处于同一个数量级优化的行为才合理。2. 分布式优化不是炫技集中式方案在运行场景中的真实痛点2.1 集中式求解的真实工程限制集中式VVO把全网节点电压、支路潮流、负荷数据、DG出力数据全部汇总到一个控制中心然后解一个大规模优化问题。算法上这不是不能做但工程现场有几个现实约束通信代价和时延。配电网节点动辄上千个每个节点实时上送量测控制中心每15分钟做一次滚动优化通信链路压力很大。多利益主体。配电网现在往往有多个分布式电源投资运营主体各自的出力数据属于商业敏感信息。你让所有DG运营商把逆变器无功能力曲线全部报到同一个中心人家不一定愿意。单点故障风险。控制中心挂了全网无功优化就瘫了。分布式架构至少能在局部区域继续维持“次优但合格”的电压控制。2.2 考虑设备动作损耗后集中式计算量更下不来设备动作损耗这个“记忆项”把相邻时段耦合在一起问题从“各时段独立优化”变成“跨时段联合优化”还带着OLTC档位、电容器组数这类整数变量。这是一个混合整数非线性规划MINLP问题工程上常用二阶锥松弛整数分支定界去求解即混合整数二阶锥规划MISOCP。单时段问题规模小还好说一旦时间维度拉到96个时段整数变量数目是设备数乘以时段数再乘以分支定界的复杂度求解时间会非常难看。我的经验是一个中等规模馈线两三百个节点的日内滚动VVO集中式MISOCP在普通办公用机器上可能要跑到分钟级甚至小时级。对调度周期15分钟的应用来说这个速度不可接受。分布式优化的思路就是把这个大问题“拆了”按网络拓扑切成几个区域每个区域只解自己的子问题区域之间只交换边界节点的电压、功率信息。单个子问题的规模小一个数量级整数变量也分散到各个子区域里求解速度自然就上来了。2.3 为什么选ADMM作为分布式框架常见分布式求解方法有三种交替方向乘子法ADMM、目标级联分析法ATC、辅助问题原理法APP。我做过的对比体验如下方法收敛条件实现难度适用场景ADMM线性约束下收敛性好参数敏感性中等低约束写成等式后用增广拉格朗日即可以边界一致性为耦合的大规模凸/近似凸问题ATC收敛性强适合层级结构中高需要上下层目标函数交互设计多个层级、多个利益主体APP需要预处理保证收敛高参数设计复杂理论分析强工程上手门槛高ADMM对配电网VVO最友好区域间的耦合天然是“边界节点电压相等、边界支路功率守恒”这类等式约束正好是ADMM的标准形式。而且子问题可以用现成的优化工具箱去解主循环代码量不大调试也直观。2.4 分布式求解会不会损失最优性这是所有人都会问的问题也确实值得问。对一个凸问题ADMM能收敛到全局最优对含整数变量的问题凸松弛之后本身就做了一个近似ADMM再收敛也只能得到“工程上的满意解”。但VVO这类问题有一个特点网损目标在海运成本函数上比较平缓在一个合理电压解附近很多近似最优解的目标值差距很小。分布式解和集中式解之间的差距通常在百分之几以内而设备动作次数可能差别很大——分布式解由于各区域自扫门前雪动作次数反而更少。我的验证办法很土但很有效同一个算例集中式和分布式各跑一遍对比总成本和电压曲线。差距在可接受范围内就认为分布式实现是可信的。3. 设备动作损耗的建模细节成本系数与离散变量的坑3.1 离散设备状态与“前后时段变位”的线性化OLTC档位和电容器组数本身是整数变量。动作次数是当前时刻状态和上一时刻状态差的绝对值。绝对值不能直接进MISOCP的目标函数需要引入辅助变量 (s_{k,t} \ge 0)并写成两个不等式[ s_{k,t} \ge u_{k,t} - u_{k,t-1}, \quad s_{k,t} \ge u_{k,t-1} - u_{k,t} ]然后把目标函数里的动作代价项替换成 (\sum_k c_k \sum_t s_{k,t})。这样绝对值就被线性化了模型保持混合整数凸结构求解器不用处理任何非光滑函数。这是个很小但很关键的建模动作漏了这一步后面所有求解器都会哭。3.2 电压无功与潮流模型的凸化LinDistFlow配电网潮流方程本身是非线性非凸的。直接上完整交流潮流方程分布式算法基本跑不动。工程标准做法是采用DistFlow方程再忽略高次小项得到LinDistFlow[ p_{ij,t} - \sum_{h:(j,h)} p_{jh,t} P_{j,t}^{load} - P_{j,t}^{DG} ][ q_{ij,t} - \sum_{h:(j,h)} q_{jh,t} Q_{j,t}^{load} - Q_{j,t}^{DG} - q_{j,t}^{shunt} ][ v_{j,t} v_{i,t} - 2(r_{ij} p_{ij,t} x_{ij} q_{ij,t}) (r_{ij}^2 x_{ij}^2) l_{ij,t} ]其中 (v_{i,t}V_{i,t}^2)(l_{ij,t}I_{ij,t}^2)。这里定义的 (l_{ij,t}) 是虚构的平方量还需要下面的二阶锥约束保证物理一致[ l_{ij,t} \ge \frac{p_{ij,t}^2 q_{ij,t}^2}{v_{i,t}} ]这个约束看起来丑但它是凸的二阶锥现代优化求解器处理起来非常高效。把潮流方程从“完整非线性”降到“线性方程二阶锥约束”是整个工程可实现性的基石。注意在实际计算里我通常再用一个下界约束比如 (l_{ij,t} \ge \epsilon) 来避免除零问题(\epsilon) 取1e-6就够。3.3 完整模型长什么样综合上面的讨论一个配电网VVO模型可以整理成目标网损成本 电压偏移惩罚 设备动作损耗成本约束LinDistFlow潮流方程、二阶锥电流约束、节点电压上下限、OLTC分接头可调范围与每时段最大动作次数、电容器组投切组数上下限、DG无功出力范围变量连续变量支路功率、节点电压平方、DG无功、整数变量OLTC档位、电容器组数。这是我说的“基础模型”。在此基础上你可以加三相不平衡、储能或需求响应但核心就这一套。我自己在实际代码里会把电压偏移惩罚设成软约束只有在电压超出上下限时才计入惩罚而不是罚偏离参考值太多。这样模型不会为了把电压往1.0上凑而白白增加设备动作。3.4 一个容易被忽略的边界条件首端电压与OLTC的关系OLTC不是直接告诉你节点电压而是改变主变二次侧电压。建模时我会把OLTC所在支路当作一个“可变电压源”分接头档位对应二次侧电压平方的一个离散增量[ v_{0,t} (V_{base} \Delta V \cdot tap_t)^2 ]每次档位变化 (\Delta V) 一般是额定电压的0.0125到0.02之间。我这里提醒一句如果直接让OLTC档位进潮流方程它是一个整数变量乘以常数再平方的结构这会让模型变成非凸。工程上更稳的做法是预先把每个档位对应的二次侧电压平方值列成一个有限集合用整数变量去选择这个集合里的值。本质还是一个MISOCP可以处理的离散选择问题。4. ADMM框架落地配网区域划分、耦合变量与子问题拆分4.1 区域划分不是随便切两刀区域划分直接决定算法收敛速度和通信复杂度。我的经验是三个原则按馈线分支切、让区域之间只共享少量边界节点、尽量减少跨区无功流动。以一个自建的33节点辐射状配电网为例主干一条三条分支。我会把主变压器和主干上游划为区域A两条较长分支各划为区域B和区域C。区域A和区域B之间只共享一个边界节点该节点下游的支路功率是耦合变量。区域C同理。这样一来每个区域内部三五十个节点变量规模小一个量级实际迭代中通信量也很小——每个区域每轮只需给相邻区域传边界节点的电压平方和边界支路功率外加一个对偶变量标量。4.2 耦合变量与一致性约束在每个区域m内部优化问题可以写成[ \min_{x_m \in \mathcal{X}_m} f_m(x_m) ]其中 (x_m) 包括区域m内部所有设备的变量。区域m和相邻区域n之间的边界条件要满足[ A_m x_m z_m, \quad A_n x_n z_n, \quad z_m z_n ]这里 (z_m, z_n) 就是边界变量。实际里比较方便的是直接取“边界下游支路功率”和“边界节点电压平方”两个量。约束“区域m看到的边界节点电压 区域n看到的边界节点电压”“区域m送出的功率 区域n接收的功率”是等式ADMM可以自然处理。注意一个细节边界变量不要在区域双方都自己“确定”而要有一个公共参考量 (z)让两个区域分别去靠近这个参考量。这一点决定了ADMM会不会振荡。我在初版代码里省掉了公共z直接在两个区域之间互相逼近结果残差曲线一直忽高忽低后来改成公共z之后才稳定。4.3 带动作损耗的子问题怎么进ADMMADMM要求在每轮迭代里依次对每个区域求解如下形式的子问题[ x_m^{(k1)} \arg\min_{x_m \in \mathcal{X}_m} \left( f_m(x_m) \frac{\rho}{2}|A_m x_m - z^{(k)} \lambda_m^{(k)}|_2^2 \right) ]注意到 (f_m(x_m)) 里已经包含了区域m内部的OLTC与电容器动作成本项。动作成本项的跨时段耦合只在区域内部所以在ADMM迭代里它是“局部可解”的不会给区域间通信增加负担。这正好体现了分布式优化的核心优势本地约束本地消化边界约束才需要跨区协商。公共变量z的更新是把所有区域边界变量的平均作为新的共识值[ z^{(k1)} \frac{1}{M}\sum_{m1}^{M} A_m x_m^{(k1)} ]对偶变量更新[ \lambda_m^{(k1)} \lambda_m^{(k)} A_m x_m^{(k1)} - z^{(k1)} ]残差计算[ r_{prim} \sum_m |A_m x_m - z|2, \quad r{dual} \rho |z^{(k1)} - z^{(k)}|_2 ]4.4 收敛判据怎么定收敛判据我一般取原始残差和对偶残差都小于1e-3或者直接看总目标函数连续三轮相对变化小于0.1%。后者在工程上更直观目标函数都稳定了残差还在缓慢波动也不影响决策。还有一个重要经验ADMM的迭代次数和罚参数 (\rho) 强相关不要执着于调一个“万金油” (\rho)动态调整更实际。我会在第6章专门讲。5. Matlab实现结构与关键代码逻辑拆解5.1 程序文件拆成五个模块Matlab代码最忌一个脚本写到头。我这个工程的文件组织方式如下你可以直接抄这个骨架文件名职责输入输出main_vvo_admm.m主循环调度各模块系统参数、负荷/DG时序各区域决策结果、目标函数值build_region_data.m定义各区域拓扑、线路参数、负荷曲线算例基础数据region结构体数组setup_distflow.m生成LinDistFlow方程的系数矩阵区域拓扑与线路参数系数矩阵A、B、Csolve_subproblem.m求解单个区域的MISOCP子问题区域数据、边界共识值、对偶变量、rho区域解、目标值update_boundary.m更新公共边界变量和对偶变量计算残差各区域解、上一个z、rho新z、lambda、残差、自适应rho5.2 子问题求解核心代码以下代码采用Matlab里常用的建模工具箱语法风格函数名我按通用习惯写。sdpvar是连续变量intvar是整数变量optimize是求解入口。你本机装了什么求解器就把solver_name替换成对应的配置function [x, obj] solve_subproblem(reg, z_boundary, lambda_prev, rho) % reg区域数据包含节点数、支路数、OLTC与电容器位置、负荷曲线 % z_boundary上一轮公共边界变量 % lambda_prev上一轮对偶变量 % rhoADMM罚参数 T reg.T; % 时段数 nb reg.nb; % 支路数 nn reg.nn; % 节点数 % 变量声明 p sdpvar(nb, T, full); % 支路有功 q sdpvar(nb, T, full); % 支路无功 l sdpvar(nb, T, full); % 支路电流平方 v sdpvar(nn, T, full); % 节点电压平方 tap intvar(reg.nTap, T); % OLTC档位 cap intvar(reg.nCap, T); % 电容器投入组数 qdg sdpvar(reg.nDg, T); % DG无功出力 % 网损与电压惩罚 loss_cost sum(sum(reg.r .* l)); % reg.r为支路电阻 volt_pen reg.alpha * sum(sum(max(0, v - reg.vmax).^2 ... max(0, reg.vmin - v).^2)); % 设备动作成本用辅助变量线性化 delta_tap intvar(reg.nTap, T); delta_cap intvar(reg.nCap, T); action_cost sum(sum(reg.cTap * delta_tap reg.cCap * delta_cap)); % 边界耦合项ADMM增广拉格朗日 [boundary_val, ~] extract_boundary(reg, p, q, v); coupling 0.5 * rho * norm(boundary_val - z_boundary lambda_prev)^2; % 约束 Constraints []; for t 1:T Constraints [Constraints, distflow_equations(reg, p(:,t), q(:,t), l(:,t), v(:,t))]; Constraints [Constraints, l(:,t) (p(:,t).^2 q(:,t).^2) ./ v(:,t)]; Constraints [Constraints, reg.vmin v(:,t) reg.vmax]; Constraints [Constraints, delta_tap(:,t) tap(:,t) - tap(:,max(1,t-1))]; Constraints [Constraints, delta_tap(:,t) tap(:,max(1,t-1)) - tap(:,t)]; % 电容器同理 Constraints [Constraints, 0 cap(:,t) reg.capMax]; Constraints [Constraints, -reg.qdgMax qdg(:,t) reg.qdgMax]; end Objective loss_cost volt_pen action_cost coupling; ops sdpsettings(verbose, 0, solver, your_solver); sol optimize(Constraints, Objective, ops); x.p value(p); x.q value(q); x.l value(l); x.v value(v); x.tap value(tap); x.cap value(cap); x.qdg value(qdg); x.obj value(Objective); end这是结构示意没有一个一个把DistFlow方程全写出来因为你的区域拓扑不同系数矩阵构建方式会不同。关键是把extract_boundary.m写对它决定区域m暴露给相邻区域的耦合变量是什么。我实际代码里这一步是用一个稀疏选择矩阵完成的效率比循环逐节点拼接高很多。5.3 ADMM主循环代码主循环很短真正干活的是子问题求解器%% main_vvo_admm.m 主循环核心片段 max_iter 200; tol 1e-3; rho 0.1; regions build_region_data(); for m 1:M lambda{m} zeros(boundary_dim(m), 1); end z zeros(total_boundary_dim, 1); for k 1:max_iter % 并行或串行求解各区域子问题 for m 1:M x{m} solve_subproblem(regions{m}, z_m, lambda{m}, rho); Ax{m} extract_boundary_value(regions{m}, x{m}); end % 更新公共边界变量 z z_new average_boundary(Ax); % 更新对偶变量 for m 1:M lambda{m} lambda{m} rho * (Ax{m} - z_new); end % 残差 r_prim 0; for m 1:M r_prim r_prim norm(Ax{m} - z_new); end r_dual rho * norm(z_new - z); % 自适应rhop if r_prim 10 * r_dual rho rho * 2; elseif r_dual 10 * r_prim rho rho / 2; end z z_new; if r_prim tol r_dual tol fprintf(收敛于第%d次迭代\n, k); break; end end注意代码里的z_m是当前区域m对应的那部分公共边界变量不要一次性传整个z给每个区域否则边界信息就“串味”了。真实分布式场景里每个区域只知道自己的边界变量和对偶变量信息交互只发生在相邻区域之间。5.4 单机仿真怎么模拟分布式交互很多读者在Matlab里做分布式算法只是把“分布式”体现在代码结构上本质还是单机循环。这没问题验证算法逻辑足够。但有一个细节不要让每个子问题直接读到其他区域的内部变量。我在代码里刻意用extract_boundary_value把区域内部解“封装”起来主循环只接触边界变量。这培养的是一种好习惯——将来真正部署到多代理系统时把函数调用换成网络接口通信算法主循环不用改。6. 算例测试里的调参心得与常见收敛陷阱6.1 罚参数 (\rho) 到底怎么调ADMM的 (\rho) 影响的是“区域之间协商的刚硬度”。(\rho) 太大边界变量很快达成一致但整目标函数里耦合项权重大会压过网损和设备动作成本项导致最终解偏离真实最优(\rho) 太小边界变量半天谈不拢迭代次数暴涨。我在自建33节点配网算例上的经验是先用一个小 (\rho)比如0.01跑几十轮观察两个残差的比值然后动态调整。经典的自适应准则是若原始残差远大于对偶残差比如10倍以上说明z更新得太慢(\rho) 需要增大若对偶残差远大于原始残差说明z波动太大(\rho) 需要减小。动态调 (\rho) 后我的迭代次数从固定 (\rho) 时的200次以上减少到约80到100次收敛目标函数值基本不变。但要注意每次改 (\rho) 时对偶变量也要相应缩放否则残差计算会突然跳变。我这次偷懒没做缩放结果对偶残差曲线出现了一次毛刺后来才补上。6.2 整数变量导致的抖动问题是第一号麻烦这是ADMM用于MISOCP最大的一道坎。OLTC档位在相邻两次迭代里可能是0.0和1.0之间反复横跳因为网损变化对档位很敏感而ADMM的线性化耦合项又不能完全锁住整数。我的处理办法分两层第一层先做连续松弛把OLTC档位和电容器组数当连续变量跑ADMM等边界变量收敛后把连续值投影到最近的整数档位再固定整数变量做一次潮流校验。这种方式工程上够用但严格来说可能丢最优性。第二层如果一定要在整数状态收敛我会在子问题目标里给整数变量加一个微小的二次邻近项比如 (0.01 \times |tap - tap_{prev}|^2)相当于给整数变量一个“惯性”让它别在相邻整数间来回跳。实际效果上第二层更稳定代价是多调试一个邻近项系数。我把两层都写在代码里用同一个开关控制选哪一层。6.3 初值没给好ADMM也会在原地打转ADMM初值的影响经常被低估。我在第一次实现时把边界变量初值全部给了零结果迭代前期区域之间的功率边界来回震荡很久。后来我把公共电压平方初值设成1.0即电压1.0 p.u.边界功率初值设成0对偶变量初值设成0问题立刻好很多。另外潮流方程里的 (l_{ij}) 初值不要设成0因为我前面提到的二阶锥约束里分母是 (v_i)虽然加了 (\epsilon) 保护但 (l0) 会让那些本不该有大电流的支路出现数值敏感性。给一个0.01的小初值数值稳定性提升明显。6.4 一组典型对比结果下面这组数据来自我手头的一个模拟算例小型33节点辐射状配网一台OLTC主变、两组电容器、两个分布式光伏接入96个时段分布式电源出力按典型日内曲线波动。设备动作成本按1.3节的口径折算OLTC单次20元电容器单次5元电价0.6元/千瓦时。方案网损成本元OLTC动作次数电容器动作次数总等效成本元集中式、不考虑动作损耗213634413021集中式、计入动作损耗2235682395分布式、计入动作损耗2238782418这个结果很能说明问题。不计动作损耗时网损最低但设备每天都在超负荷地“工作”按寿命折算下来总成本最高。计入动作损耗后损失了约100元的网损但省下六百多块设备损耗成本。分布式解和集中式解的总成本差约1%工程上完全可以接受。这个模拟算例同时也验证了我开头说的网损最优不等于总运行成本最优更不等于设备健康。6.5 调试顺序建议别一上来就全网络分布式我踩过的最大一个坑就是“一上来把全网络跑分布式”出现异常时根本不知道是区域子问题建模错、边界提取错、还是ADMM迭代参数错。正确的调试顺序是先在单个区域内做集中式VVO验证潮流方程、动作成本建模和求解器配置都没问题再把网络切成两个区域只保留边界电压一致性暂时去掉动作成本项验证ADMM能否收敛到集中式解加上设备动作成本项检查目标函数里三项权重的比例最后扩展到三个区域及更复杂的DG时序。每一步都留一个集中式求解的“基准解”对照。这样一旦分布式解偏移太多你能立刻定位是哪一层引入的问题而不是对着一个黑盒算法瞎猜。我最后一次做这三区域算例从子问题建模到ADMM收敛调试时间比之前跳着来省了一半不止。