配电网的无功优化，听起来像是教科书里才会出现的名词，但做实际工程的人都知道，这几乎是每天都要面对的问题。尤其当分布式光伏、充电桩大规模接入后，台区电压波动、线损超标这些问题越来越频繁，传统“固定补偿+变压器档位”的做法已经快跟不上节奏了。我之前在做一个园区配电网改造项目时，就被电压越限和网损偏高两个问题折磨过很久，后来系统梳理才发现，核心还是落在无功优化上面。
这篇文章就围绕“IEEE33节点配电网的多目标无功优化”这个载体，把整套思路、模型、算法和实操中容易踩的坑完整记录下来。先说明一下，IEEE33节点系统是配电网研究里最经典的标准算例之一，几乎所有做配电网潮流、重构、优化的人都会拿它练手，因为数据开源、结构合理、复现性强，非常适合用来验证算法效果。不管你是做电力系统研究的在校学生，还是刚接触配电网优化的一线工程师，这篇文章都能给你一条清晰的实操路线。
1. 项目概述与可行性拆解
1.1 无功优化到底在解决什么问题
讲无功优化之前，得先聊清楚无功功率本身。电动机、变压器、电抗器这些感性设备在工作时需要建立磁场，这部分能量不做功，却在电网和负荷之间来回流动，就是无功功率。无功不是“没用”，但它的流动会占用线路容量，增大电流，进而带来额外的有功损耗和电压降落。
配电网的电压偏低，绝大多数时候不是有功不够，而是无功缺了。无功不足导致电压支撑能力下降，尤其是线路末端，电压跌落得厉害。反过来看，无功过剩（比如光伏大发但就地消纳不掉）又会让电压抬升甚至越上限。所以无功优化的本质，就是通过调节无功源（电容器、电抗器、静止无功发生器、变压器分接头等）的出力，让电压合格、网损最小、经济性最好。
IEEE33节点系统虽然只是一个标准算例，但它完整复现了实际配电网的长度特征：线路总长约3.7公里，首端到末端跨越了多个分支，末端节点电压在重负荷时下降明显。这正好用来测试无功优化对“长线路末端电压”这个真实痛点的改善效果。
1.2 为什么选择IEEE33节点系统作为研究载体
IEEE33节点系统有33个节点、32条支路，还有一个联络开关，额定电压12.66千伏，基准容量10兆伏安。它最大的优势是数据公开透明，参数早就经过多方验证，不需要自己搭线路模型，省去了大量前置工作。
更重要的是，IEEE33节点在重载工况下天生就有电压偏移问题。比如在基准负载下，如果不加任何补偿，末端节点电压会降到0.903左右（标幺值），低于0.95的合格线。这种“天生有病”的特性，反而成了优化算法最好的试验田：优化前电压越限，优化后电压恢复，效果一目了然。
还有一个选择上的考量是通用性。IEEE33节点不只是学术论文里的常客，很多商业软件（比如某电力系统仿真平台）也把它作为内嵌算例。用这套系统做验证，结论很容易被同行复现和认可，这也是我推荐它作为第一台“练手平台”的原因。
1.3 多目标优化的核心矛盾
所谓多目标，在这类问题里通常指的是两个甚至三个互相打架的目标。项目中列了两个主要目标：节点电压偏移最小、网损最小。这两个目标真的“一致”吗？不一定。降低网损最直接的办法是无功就地平衡，也就是让无功功率在负荷附近得到补偿，减少线路上的无功流动；而电压调整则更看重补偿点位置，有时为了把末端电压撑住，需要在网络末端多投无功，哪怕这个位置的补偿对全网降损贡献不大。
这两个目标之间存在权衡关系：网损最低的解不一定电压最平，电压最平的解不一定网损最小。所以叫“多目标优化”——不是在找唯一的全局最优解，而是在找一组互不支配的Pareto解集，让决策者根据实际偏好选解。
做项目的时候如果只把两个目标用权重系数加成一个单目标，简单是简单，但问题也很明显：权重怎么定？0.7对0.3还是0.5对0.5？量纲不统一时权重更敏感。这个隐藏在起步阶段不太容易发现，等算法跑完，看到结果才会明白——原来换一组权重，补偿方案完全变了。
2. 理论核心与多目标无功优化模型构建
2.1 目标函数的数学刻画与量纲处理
在这个项目里，优化目标写成了两部分叠加的形式。网损目标很直观，就是所有支路的有功损耗之和，公式表达为：
P_loss = Σ G_ij (U_i² + U_j² - 2U_i U_j cosθ_ij)
其中G_ij是支路电导，U_i和U_j是节点电压幅值，θ_ij是节点间的相角差。这个量纲是千瓦或兆瓦。
节点电压偏移目标则更偏向“电能质量”维度。标准的电压偏差指标是各节点电压偏离额定值的程度，通常写成：
ΔU = Σ ((U_i - U_n) / U_n)²
U_n是额定电压。这个量纲是标幺值的平方，无量纲。问题来了：网损是几十千瓦量级，电压偏移是0.001到0.1量级，直接相加的话，电压偏移项会被网损项淹没，优化结果基本只照顾网损，电压质量完全不管。
所以做归一化很关键。一种通用做法是把网损除以其基准值（优化前的网损），把电压偏移除以允许的最大偏移值（比如0.05的平方乘以节点数），让两个目标都落在0到1区间再叠加。我在实操中体会，归一化处理不等于设权重，它解决的是量纲问题，权重则是人为偏好的表达，两者不能混为一谈。
2.2 约束条件与变量类型分析
约束条件是这个项目最容易出错的地方。功率平衡约束要求每个节点的注入功率必须等于负荷功率与支路功率之和，这是潮流计算本身的性质，不是额外加的。节点电压约束要求各节点电压在0.95到1.05标幺值之间。补偿容量约束则限定了电容器投切档位的上限。
决策变量的选取直接决定了模型复杂度。我建议从最简单的场景入手：只考虑在固定的候选节点上安装并联电容器组，变量就是每个补偿节点投入的电容量（以kvar为单位）。电容器分组投切时变量是离散的；如果按照连续值处理，算出来是浮点数，工程上还得圆整到实际档位。很多文章不写这一步，但实际工程里离散圆整可能让结果明显变差，所以在前处理时就要考虑档位粒度。
还有一种做法是把变压器分接头档位也作为决策变量。但这会带来一个问题：分接头调节是非线性的，会改变潮流方程中的变比参数，导致模型复杂度直线上升。我前期的项目就吃过这个亏——变量一多，搜索空间膨胀，算法收敛性立刻变差。建议先把电容器补偿做好，再考虑要不要加分接头。
2.3 为什么潮流计算是非做不可的环节
没有潮流计算的多目标无功优化，就像闭着眼睛开车——目标函数都算不出来，算法怎么可能收敛。潮流计算的意义在于给出系统真实的运行状态：各节点电压幅值和相角、各支路的功率流向和损耗。有了这些数据，才能评估任意一组补偿方案下的网损和电压偏移。
配电网和输电网在电气特征上有显著区别：配电网支路电阻相对大，R/X比值高，PQ分解法容易失效；输电网常用的牛顿-拉夫逊法在配网里也可能因为雅可比矩阵病态而难以收敛。所以配电网潮流计算通常采用前推回代法，这种方法对辐射状网络有天然的适应性和极快的收敛速度，是工程上最稳妥的选择。
前推回代法本质是“猜初始值、逐轮修正”的过程：假设各节点电压为额定值，从末端往首端推支路电流，再从首端往末端回推电压，反复迭代，直到两次迭代的电压差小于允许误差。这是整份代码里幂计算量最大、也最需要细心的地方。
3. 基于IEEE33节点系统的潮流计算与实现要点
3.1 前推回代法的通用计算步骤
前推回代法落地时，标准流程分三步。第一步，初始化所有节点电压为根节点电压。
第二步是前推支路电流。从末端节点开始，根据当前节点电压和节点注入功率计算支路电流。注入功率等于该节点的负荷功率加上所有下游支路的功率之和。这里要特别注意功率方向约定，负荷是流出节点，光伏等分布式电源是注入节点，方向反了潮流结果就会错。
第三步是回代更新电压。从首端往末端走，用上一级的节点电压减去该支路电流在本支路阻抗上的压降，得到下一级节点电压。公式可以写成U_j = U_i - Z_ij * I_ij，这里的Z_ij是复数阻抗，I_ij是刚刚前推算出的支路电流复数。
反复迭代第二步和第三步，直到两次更新的节点电压最大差值小于收敛阈值。阈值通常取0.0001标幺值就行。IEEE33节点系统规模不大，这个收敛速度通常会在几毫秒到几十毫秒之间完成，完全不用担心计算时间，整套系统的计算瓶颈反而在优化算法的多次潮流调用上。
3.2 IEEE33节点系统参数与数据预处理
IEEE33节点系统有三列关键数据：节点到首端的连接关系、支路阻抗（电阻和电抗）、各节点的有功和无功负荷。这些数据可以从公开文献中找到原文表格，我不建议直接手抄，容易抄错，最好找一份带CSV格式的数据文件导入。
数据导入后必须做一项关键预处理：确认基准值。IEEE33节点系统默认基准电压为12.66千伏，基准功率为10兆伏安。所有负荷数据都是有名值（千瓦、千乏），在写潮流程序时要转换成标幺值。除以基准功率就得到标幺值，比如节点18的负荷是90千瓦和40千乏，标幺值就是0.009和0.004。
有一个非常容易踩的坑是联络开关的处理。IEEE33节点系统原始拓扑中有一个常开的联络开关，做无功优化时默认网络处于辐射状运行状态，联络开关是断开的。如果误把联络开关合上，网络变成环网，前推回代法根本不适用，潮流立马出问题。所以在写程序时，支路数据里要明确排除开关支路，只保留32条闭合支路。
3.3 潮流计算结果的校验方法
潮流程序写完以后不要急着接优化算法，先做一次不加补偿的基准潮流计算，把结果和IEEE33节点的公开结果对比。校验重点看三个数据：根节点注入功率是否约为3.98兆瓦和2.3兆乏；末端节点电压是否在0.903左右（标幺值）；全网有功网损是否约为0.202兆瓦，无功网损约为0.135兆乏。
这三个数字几乎出现在所有IEEE33节点的文献中，如果算出的结果差太多，说明代码有bug。在这台“仪器”校准之前，后面所有优化结果都不可信。每次改动潮流代码后也必须回归测试，我因为这个吃过一次亏：后来优化结果怎么调都不对，排查了半天发现是支路数据矩阵录入时漏了一行。
4. 多目标优化算法选型与实现细节
4.1 从单目标加权到多目标Pareto的思路演变
在做多目标优化之前，得先说清楚为什么不能直接做单目标加权。加权法简单粗暴：把网损和电压偏移乘以各自权重再加起来，变成一个标量，然后跑一次常规优化算法。这个方法能用，但问题很多。
最大的问题是两个目标的量纲和数值范围差异巨大，权重很难解释和调整。其次是权重一旦定下来，搜索方向就锁死了，很难发现Pareto前沿上那些“极致点”——比如网损很大但电压分布极平、或者电压一般但网损极小的解。实际工程中，不同阶段对目标的侧重是不同的，比如迎峰度夏时电压更优先，平时则更看重经济运行，如果没有整条Pareto前沿供选择，每次目标变化都要重新优化。
所以我最终选择用NSGA-II（带精英策略的非支配排序遗传算法）来做多目标优化。它不要求把多个目标合成一个，而是直接寻求一组互不支配的Pareto解集，解集中每个解都在两个目标上有不同的权衡表现，决策人可以结合运行需求选解。
4.2 NSGA-II的核心流程拆解
NSGA-II的流程可以拆成五步。第一步，初始化种群。假设种群规模为100，每个个体是一组补偿方案——比如在候选节点上各投入多少千乏的电容器，编码成实数向量或整数向量。每个个体都要调用一次潮流计算，得到自身的网损和电压偏移两个目标值。
第二步，非支配排序。每个个体按照“支配”关系进行层级分类，Pareto前沿上的解属于第一层，被更少数个体支配的属于第二层，依次往后。这一步的目的是把优劣关系理出来，让算法优先保留前沿解。
第三步，计算拥挤度距离。同属一个前沿层的解，还要看它们在目标空间中是否拥挤。落在稀疏区域的解更有价值，因为它能帮助维持解的多样性，避免结果挤在一小段曲线上。
第四步是遗传操作：锦标赛选择、模拟二进制交叉、多项式变异。这一步我建议标准差取0.1左右，变异率不要超过0.1。变异率过大会把解集打成随机搜索，过小则容易早熟。第五步，精英保留策略，将父代和子代合并，整体排序后取前N个个体作为新一代种群。循环迭代直到达到最大代数，通常取100到200代。
4.3 决策变量编码与约束处理
决策变量如何编码，直接影响算法的搜索效率。IEEE33节点系统有32个负荷节点，理论上每个节点都能装电容器，但工程上不可能每个节点都装。常见做法是预选几个候选补偿节点，比如选取电压偏低最严重的末端节点、负荷集中的节点、网络中间节点，通常选3到5个。
编码方式可以采用实数编码或者整数编码。电容器是分组投切的，我在项目中取值区间设置为0到300千乏，步长30千乏（即10档），这样编码是一个正整数，更贴近实际设备。
约束处理主要有两个层面。第一是补偿容量上下限，这个简单，编码时强制在区间内取整。第二是节点电压上下限，这个比较麻烦。如果一个个体经过潮流计算后，某些节点电压超出0.95到1.05范围，应该怎么办？我建议罚函数法：给这个个体的目标函数加上一个很大的惩罚项，让它在非支配排序中自动“沉底”。这样做比直接删除个体更好，因为接近可行的解还能通过遗传变异进化到可行域里。
5. MATLAB实现与关键代码解析
5.1 主程序框架搭建
我当时的开发环境是MATLAB R2021a，但其实逻辑和Python完全互通。主程序的框架结构可以分成数据加载、参数配置、潮流函数封装、优化算法主循环、结果可视化几个模块。潮流计算封装成函数是最关键的一步，函数输入是补偿方案向量，返回两个目标函数值。
下面给出一份精简但完整的潮流计算代码框架，可以在这个基础上扩展优化算法：
function [V, Ploss, Qloss] = powerflow(branch, bus, comp) % branch: 支路数据 [首端节点, 末端节点, 电阻, 电抗] % bus: 负荷数据 [节点编号, 有功, 无功] % comp: 补偿无功数组，对应各节点注入的无功 n = size(bus, 1); V = ones(n, 1); % 电压初始化（标幺值） err = 1; itermax = 50; iter = 0; % 支路功率存储 S_branch = zeros(size(branch, 1), 1); while err > 1e-6 && iter < itermax V_old = V; % 前推：从末端往首端推电流 I = zeros(size(branch, 1), 1); for k = size(branch, 1):-1:1 end_node = branch(k, 2); start_node = branch(k, 1); % 节点注入功率 = 负荷 - 补偿 + 下游支路功率 S_load = bus(end_node, 2) + 1j * (bus(end_node, 3) - comp(end_node)); % 加上下游支路的功率（此处需要累加） S_inj = S_load; % 找所有以 end_node 为首端的支路 downstream = find(branch(:, 1) == end_node); for d = 1:length(downstream) S_inj = S_inj + S_branch(downstream(d)); end I(k) = conj(S_inj) / conj(V(end_node)); S_branch(k) = S_inj; % 保存支路功率 end % 回代：从首端往末端更新电压 for k = 1:size(branch, 1) start_node = branch(k, 1); end_node = branch(k, 2); Z = branch(k, 3) + 1j * branch(k, 4); V(end_node) = V(start_node) - Z * I(k); end err = max(abs(V - V_old)); iter = iter + 1; end % 计算网损 Ploss = 0; Qloss = 0; for k = 1:size(branch, 1) n1 = branch(k, 1); n2 = branch(k, 2); Z = branch(k, 3) + 1j * branch(k, 4); I_k = (V(n1) - V(n2)) / Z; Ploss = Ploss + abs(I_k)^2 * branch(k, 3); Qloss = Qloss + abs(I_k)^2 * branch(k, 4); end % 电压偏移指标 V_dev = sum((V - 1).^2); end上面这段代码的前推回代部分没有提前聚合节点注入功率，严格来说需要先做一次拓扑遍历算出节点分层顺序，才能正确累加下游功率。我在代码里标注了关键点，实际工程实现时建议先用深度优先搜索给节点排个序，而不是直接按支路索引循环。
5.2 多目标优化主循环代码
NSGA-II主循环可以直接用MATLAB全局优化工具箱里的gamultiobj函数，它内置了NSGA-II的核心机制。也可以自己写非支配排序，但工程量不小，前期验证时用内置函数足够。
以下是一个调用gamultiobj求解无功优化问题的示例：
% 目标函数封装 fun = @(x) multi_objective(x, branch, bus); % 变量个数：这里假设有3个补偿节点 nvars = 3; % 变量下界和上界（0到10档，每档30kvar） lb = zeros(1, nvars); ub = 10 * ones(1, nvars); % 整数约束 IntCon = [1, 2, 3]; options = optimoptions('gamultiobj', ... 'PopulationSize', 100, ... 'MaxGenerations', 150, ... 'Display', 'iter', ... 'PlotFcn', @gaplotpareto); % 求解 [x, fval] = gamultiobj(fun, nvars, [], [], [], [], lb, ub, [], IntCon, options);这里multi_objective函数内部需要调用潮流函数，返回两个目标值向量：P_loss（标幺值或有名值）和V_dev。需要特别注意的是，目标函数中必须禁止显示输出，否则每评估一个个体就会打印一大堆中间数据，优化过程会慢到让人崩溃。
5.3 结果可视化与方案比选
优化结束后，得到一组Pareto解集，但很多情况下这组解有几十个，不可能全部拿到工程现场去执行。如何从中挑出最终方案？我的习惯做法是先画Pareto前沿散点图，横轴是网损（kW），纵轴是电压偏移指标。前沿曲线通常是L形，左端网损小但电压偏移大，右端电压偏移小但网损大，中间区域两者相对均衡。
取终点解时我会定一个偏好：如果当前系统电压水平普遍偏低，优先选电压偏移较小但网损稍有增加的解；如果系统电压尚可但线损超标，则选网损更小的解。同时还要检查这个解对应的补偿方案是否合理——比如补偿是否集中在极少数节点上，单点补偿量是否超过了实际设备容量。
我在项目中做了一版以电压质量为优先的选解，最终选出的方案让最小节点电压从0.903提升到0.972，网损从202千瓦降到142千瓦左右，两个目标同时改善。这说明初试方案下电压和网损的矛盾没那么尖锐，因为初始工况太恶劣，随便补一补都能同时改善。等优化到一定程度，两个目标才开始出现明显权衡，这时候才真正体现Pareto前沿的价值。
6. 常见问题与排查技巧实录
6.1 潮流计算不收敛怎么办
潮流计算不收敛是最常见的问题，几乎每个新手都会遇到。典型的症状是迭代次数跑满但误差始终降不下来，或电压出现负值或NaN。
首先是检查支路数据方向。前推回代法要求每条支路的首端节点必须是距离根节点更近的那一侧，如果数据录入时方向反了，电流前推和电压回代的路径就乱了，迭代必然发散。
其次是检查负荷单位。IEEE33节点系统原始数据有的文献用千瓦、千乏，有的文献用兆瓦、兆乏。如果不统一成标幺值，数值大小可能差出三个数量级，导致潮流完全跑飞。我在代码里统一用标幺值计算，每录入一份新数据都会先打印汇总负荷，看总负荷是否为3.98兆瓦左右。
还有一条是看迭代初值。如果根节点电压设定为1.0标幺值，IEEE33节点系统通常是可行的。但有些改进的算例（比如加入分布式电源后）需要把初值设为1.05，否则可能导致迭代发散。建议初始电压直接设为1.0，并限制迭代中电压幅值不低于0.5，如果出现低于0.5的节点说明模型有错误，要停住排查。
6.2 多目标结果普遍很差怎么排查
如果发现所有个体的网损都很高或者电压偏移很大，不要先怀疑算法，而要先检查目标函数是否正确。一个典型错误是补偿变量的正负号：电容器发出无功功率，所以节点注入无功是负荷无功减去补偿量，即Q_inj = Q_load - Q_comp。如果代码里写反了，变成加法，那补偿越多电压越低，结果自然越差。
另一个容易忽视的点是潮流收敛阈值。阈值设得太宽松（比如0.01）会导致目标函数值有较大噪声，同一个补偿方案重复调用潮流得到的目标值略有波动，NSGA-II误以为是不同的解，最终前沿非常杂乱。建议阈值至少取1e-6，这样目标函数可重复性才有保证。
还有一个是检查候选补偿节点选得是否合理。如果选的节点本身就是电源节点附近（比如节点1或2），补电容对整条线路的电压和网损贡献极小，算法再怎么优化也不可能把末端电压拉起来。候选节点应该选在电压最薄弱的末端或负荷密集的节点。
6.3 算法早熟收敛与参数调优技巧
NSGA-II在求解无功优化时偶尔会出现早熟收敛，即所有个体都集中到一小块区域，Pareto前沿非常短，缺少多样性。这种情况多半是交叉率、变异率设置不合理或者种群规模太小。
种群规模建议不要低于100，代数至少100代。如果发现前沿形状不完整，先试着把种群规模加大到200、代数加到200，看前沿是否延展。如果前沿已经延展但很稀疏，可以尝试增大变异率。如果前沿形状可以但目标值次优，问题大概率出在目标函数里。
还有一个实用技巧：对离散变量做约束时，不要直接在优化算法内用取整函数，而是把补偿变量的上下界改成整数档位的索引，解码时才乘上每档容量。这样做的好处是让遗传算子始终在有效整数域内搜索，避免算法浪费大量算力去评估那些无效解。
7. 项目扩展与工程落地思考
以IEEE33节点系统为起点跑通多目标无功优化后，再往实际工程走，会面临几个质的飞越。第一个问题是系统规模。IEEE33节点只有33个节点，但实际配电网动辄几百上千个节点，前推回代法依然能用，但单次潮流计算时间会从毫秒级涨到秒级，优化算法调用数万次潮流就变得非常耗时。这时需要引入基于灵敏度的补偿点筛选或图卷积网络代理模型，先把候选补偿点缩小到十几个，再跑优化。
第二个问题是动态时变性。IEEE33节点是静态负荷，但实际负荷曲线24小时都在跳。白天光伏大发时无功可能过剩，晚上负荷高峰时无功又紧缺。静态优化显然不够，要做多时段动态无功优化，目标函数变成全天网损总和或全天电压合格率，变量数乘24倍，复杂度直线上升。我在工程中见过不少团队把这套静态模型直接套到动态场景里，结果频繁调节设备，严重影响设备寿命。
第三个问题是分布式电源的引入。IEEE33节点原始算例里没有光伏和风电，但有文献扩展版本把分布式电源接在节点18和节点33等处。加分布式电源后，系统出现双向潮流，电压分布更复杂，无功优化的目标里还需要考虑分布式电源的无功支撑能力。这是当前研究最火的方向之一，也是IEEE33节点系统仍然被反复使用的核心原因。
这些扩展方向不是一次性能做完的。我的体会是先把静态多目标无功优化吃透，把Pareto前沿、潮流计算、约束处理这些基本功打牢，再逐步往动态、分布式的方向加复杂度，每一步都回到IEEE33节点上去验证，这样踩坑成本最低。
最后再分享一个小技巧：多目标优化结果出来后，不要只看数字，把电压分布曲线画出来对比优化前后，会非常直观。我在项目里自始至终保留了一张电压分布对比图，汇报、写报告、定位问题都用它，比一堆表格管用得多。
