二阶锥规划求解主动配电网最优潮流：从模型到代码实战

简介：基于二阶锥规划的IEEE33节点主动配电网最优潮流求解程序，面向电力系统优化方向的科研人员与电气工程研究生，可用于多源协同运行、配电网经济运行等场景的算例复现。程序基于MATLAB平台，调用YALMIP和CPLEX求解器实现二阶锥规划，针对IEEE33节点系统开展24小时多时段优化，模型中综合了风电、CB、SVG、OLTC、ESS等常见配电网可控资源。压缩包内含2个文件，以.m主程序为核心，另附说明txt文档，整体约4KB，代码小巧但结构清晰，注释被作者打磨得较为详尽，具备骨灰级注释特点，适合初学者对照算法流程逐步理解。目前已有37人浏览学习，可作为课题研究或课程设计的直接参考。资源中额外提供了参考文献线索，可帮助读者从方法论层面把握优化建模与求解思路。

1. 为什么说二阶锥规划是主动配电网最优潮流的解药

主动配电网里高比例分布式电源接入后，光伏反送、电压越限、潮流倒送经常同时出现。调度要回答的不再是“算个潮流”，而是“怎么调储能和逆变器，让整个网在满足安全约束的前提下跑得最经济”，这就是主动配电网最优潮流问题。问题的核心障碍在于传统交流最优潮流是非凸优化，直接求解容易陷在局部最优里，算出来的策略现场不敢用。二阶锥规划通过把非凸的潮流方程做凸松弛，把最优潮流转成一个可全局求解的凸问题，是目前配电网方向可落地性最强的一条路径。这篇笔记适合做配电网运行策略、规划研究和调度算法的工程师，从模型推导讲到代码实现，把能跑的代码和踩过的坑一起讲清楚。

2. 从非凸到凸：二阶锥规划求解最优潮流的模型演变

2.1 交流最优潮流为什么让传统算法翻车

在配电网里做优化，最直接的想法是拿现成的交流最优潮流模型去套。套过就会知道，那个模型是真的难解。节点注入功率方程长这样：

P_i = V_i Σ V_j (G_ij cosθ_ij + B_ij sinθ_ij)

电压幅值、相角、电导电纳全部耦合在一起，还有三角函数项。这个方程本身不是凸函数，由它构成的可行域也不是凸集。凸规划里“局部最优等于全局最优”的性质在这里不成立，牛顿法、内点法从不同初值出发可能收敛到完全不同的结果，初值给不好甚至直接发散。对于配电网这种动辄几百上千个节点的网络，靠多初值试探来找全局解，计算量完全失控。

有人会退一步做直流潮流最优化，把电压幅值全设成1，忽略无功功率和网损。这在输电网里勉强能用，因为输电网支路阻抗以电抗为主，电压幅值接近1，有功和无功近似解耦。配电网完全不同，支路电阻和电抗在同一数量级，R/X常常在0.5到2之间，甚至更高。光伏和负荷波动带来的无功潮流对节点电压影响非常显著，直流潮流模型算出来的电压分布和网损跟实际差出一大截，拿去调度是危险的。

2.2 DistFlow支路方程：把潮流改写成“可凸”的样子

放射状配电网有一条更顺手的建模路径：DistFlow支路潮流模型，也叫Branch Flow模型。这个模型不直接使用节点电压相量，而是描述每条支路上的有功功率、无功功率、电流平方值和节点电压平方值之间的关系。对一个放射状配网，每条支路(i,j)上有四个核心变量：P_ij代表从节点i流向节点j的有功，Q_ij代表对应的无功，l_ij代表支路电流幅值的平方，v_i代表节点i电压幅值的平方。

DistFlow模型用四条方程描述稳态运行状态。第一条是节点有功平衡：流入节点j的支路有功减去线路损耗，再减去流向子支路的有功，等于节点j的本地注入功率。第二条是无功平衡，形式完全一致。第三条是电压降落方程：末端节点电压平方等于首端电压平方减去支路电阻和电抗上的压降，再加上电流引起的修正项。第四条是支路视在功率的等式定义：P_ij的平方加Q_ij的平方等于首端电压平方乘以电流平方。

前三条方程都是线性的，只有第四条是二次等式，非凸性全部集中在这条约束里。

2.3 二阶锥松弛：把等式放宽成锥约束

二阶锥规划的关键操作就在这里：把最后那条等式约束放宽成不等式。原本要求P_ij²加Q_ij²严格等于v_i乘以l_ij，现在允许它小于等于。这样做物理含义很直观：支路视在功率的平方不超过首端电压平方和电流平方的乘积，相当于允许模型暂时低估线路电流。这个放宽后的约束可以直接改写成标准二阶锥形式：

‖ [2P_ij; 2Q_ij; v_i - l_ij] ‖₂ ≤ v_i + l_ij

这是一个旋转二阶锥约束，可以用现成的凸优化求解器高效处理。

但松弛不能白做，要回答一个问题：松弛之后的最优解还能不能映射回物理上真实可行的潮流？答案是：在放射状配电网且目标函数满足一定单调性条件时，最优解会自动落到锥边界上，也就是不等式取到等号，松弛是紧的。实际操作中，大多数以最小化网损、最小化购电成本为目标的最优潮流问题都能满足这个条件。遇到松弛不紧的情况，通常需要检查目标函数或者是否加了一些把解往锥内部推的约束。

下面这张表可以直观看到两个模型的差异：

松弛带来一个额外的好处：求解器能给出全局最优解的下界，这个下界对算法验证极其重要。用潮流计算器去核对SOCP给出的解时，如果偏差大，至少知道是模型本身的问题还是数值问题，而不是在非凸空间里瞎猜。

3. 用Python+CVXPY从零搭一个二阶锥最优潮流

3.1 网络数据准备：3节点放射网的最小案例

理论讲完必须落到代码。我习惯用最小可运行案例验证建模思路，再往真实网络扩展。这里用一个3节点放射状配电网：节点0是变电站出口母线，视为平衡节点，电压固定在1.0标幺；节点1和节点2带负荷，节点1还带了一组分布式光伏。

网络数据全部采用标幺值，基准容量取1 MVA，基准电压取10 kV。支路参数放在字典里，每个支路包含首端节点、末端节点、电阻和电抗。负荷和光伏出力也存成字典，键是节点编号，值是有功和无功的标幺值。注意光伏出力在节点1是0.25，负荷是0.50，净注入是-0.25，意味着这个节点仍然从电网取电，但比原来少了。

数据准备的关键是保持支路方向一致：首端必须是靠近变电站的一端，末端是远离变电站的一端。这样后面构建节点功率平衡时才不会把方向搞反。真实网络的配电网数据通常来自GIS导出的拓扑表，第一步就要做方向规范化，这是个容易踩坑的点。

3.2 变量、约束、目标：CVXPY建模三要素

变量定义分两组。支路变量有三个字典：P、Q、I，键是支路编号，值分别是支路首端有功、首端无功和电流幅值平方。节点变量有一个字典V，键是节点编号，值是电压幅值平方。为什么要用电压平方而不是电压幅值？因为DistFlow方程本身就是关于v的线性方程，用平方形式避免了开方操作，保持了问题的凸性。

约束构造围绕三个层次展开。第一层是变电站节点电压固定，V[0]等于1.0。第二层是支路方程，包括电压降落方程和二阶锥约束。第三层是节点功率平衡，以本地净注入的形式写进等式约束。这里要特别注意：支路变量P和Q不要设置非负约束，分布式电源接入后，配电网支路出现倒送功率是正常工况，把变量限制为非负会把可行解直接截掉。

目标函数选的是最小化网络总损耗。在DistFlow模型里，每条支路的损耗就是电阻乘以电流平方，对l_ij变量而言是线性项。目标函数线性，约束条件是线性和二阶锥，整个问题就是标准的SOCP。

3.3 完整可运行代码：3节点主动配电网最优潮流

import cvxpy as cp import numpy as np # ---------- 网络数据（标幺值，基准容量1 MVA，基准电压10 kV） ---------- # 支路表: id -> (首端节点, 末端节点, 电阻, 电抗) lines = { 0: (0, 1, 0.010, 0.030), 1: (0, 2, 0.015, 0.040), } # 负荷: 节点 -> (有功, 无功) loads = { 1: (0.50, 0.20), 2: (0.30, 0.10), } # 光伏出力: 节点 -> (有功, 无功) pv = { 1: (0.25, 0.00), 2: (0.10, 0.00), } buses = [0, 1, 2] # ---------- 优化变量 ---------- P = {i: cp.Variable() for i in lines} # 支路首端有功 Q = {i: cp.Variable() for i in lines} # 支路首端无功 I = {i: cp.Variable(nonneg=True) for i in lines} # 支路电流平方 V = {b: cp.Variable() for b in buses} # 节点电压平方 cons = [] # 变电站出口电压固定为1.0 cons.append(V[0] == 1.0) # ---------- 支路方程 ---------- for i, (f, t, r, x) in lines.items(): # 电压降落方程 cons.append(V[t] == V[f] - 2 * (r * P[i] + x * Q[i]) + (r * r + x * x) * I[i]) # 二阶锥松弛: || [2P, 2Q, Vf-I] || <= Vf + I cons.append(cp.norm(cp.hstack([2 * P[i], 2 * Q[i], V[f] - I[i]])) <= V[f] + I[i]) # ---------- 节点功率平衡 ---------- parent_of = {t: i for i, (f, t, r, x) in lines.items()} children_of = {b: [] for b in buses} for i, (f, t, r, x) in lines.items(): children_of[f].append(i) for b in buses: if b == 0: continue # 变电站节点不建平衡方程，功率由目标函数决定 p_net = 0.0 q_net = 0.0 # 从父支路流入本节点的功率（要扣除线路损耗） if b in parent_of: i = parent_of[b] f, t, r, x = lines[i] p_net += P[i] - r * I[i] q_net += Q[i] - x * I[i] # 流出到子支路的功率 for i in children_of[b]: f, t, r, x = lines[i] p_net -= P[i] q_net -= Q[i] # 本地净注入 = 光伏 - 负荷 p_gen = pv.get(b, (0, 0))[0] - loads.get(b, (0, 0))[0] q_gen = pv.get(b, (0, 0))[1] - loads.get(b, (0, 0))[1] cons.append(p_net == p_gen) cons.append(q_net == q_gen) # ---------- 节点电压运行约束（0.95 ~ 1.05 pu） ---------- for b in buses: if b == 0: continue cons.append(V[b] >= 0.95 ** 2) cons.append(V[b] <= 1.05 ** 2) # ---------- 目标：最小化网络损耗 ---------- loss = sum(r * I[i] for i, (f, t, r, x) in lines.items()) prob = cp.Problem(cp.Minimize(loss), cons) prob.solve(solver=cp.CLARABEL, verbose=False) # ---------- 输出 ---------- print(f"求解状态: {prob.status}") print(f"最优网络损耗: {loss.value:.6f} pu ({loss.value * 1.0:.3f} MW)") print(f"变电站出口有功: {sum(P[i].value for i in lines):.4f} pu") for b in buses: print(f"节点 {b} 电压: {np.sqrt(V[b].value):.4f} pu") for i, (f, t, r, x) in lines.items(): print(f"支路 {i}: P={P[i].value:.4f} pu, Q={Q[i].value:.4f} pu, " f"I={I[i].value:.4f} pu")

这段代码有三个关键点要说明。第一，V[f]在二阶锥约束里是支路首端电压平方，不是末端，方向写反会得到完全错误的结果。第二，节点功率平衡里流入项扣掉了损耗r乘I，这是DistFlow的标准写法，漏掉损耗会让功率不守恒。第三，光伏和负荷的符号约定必须统一，代码里采用“注入为正”，所有输出结果也按这个方向解读。

运行这个3节点例子，最优网损通常在0.002 pu量级，也就是几十kW的水平，节点电压会落在0.98到1.02之间。光伏出力大于本地负荷时，支路有功出现负值，代表功率倒送，这是主动配电网的正常工况。

3.4 参数怎么调：从3节点扩展到真实配电网

真实网络至少几百个节点，数据结构要改成按列表组织，三个基础参数要确认：电压基准值、容量基准值和支路方向。我一般把电压基准设为配电网额定电压，容量基准设为10 MVA或者1 MVA，具体看量级。支路方向处理上，配电网拓扑通常从变电站开始做广度优先搜索，确定每个节点的父节点，再统一生成lines表。

求解器的选择也有讲究。CVXPY默认带的CLARABEL对SOCP问题非常稳，小规模秒出结果。如果网络规模大，ECOS和MOSEK都是好选择，其中MOSEK对数值病态问题的鲁棒性最好，工程上大量案例用MOSEK跑数万节点的配电网SOCP不会出问题。预算有限就选CLARABEL，绝大多数场景足够了。

4. 主动配电网的设备建模：光伏、储能和可调无功怎么塞进SOCP

4.1 光伏逆变器：有功出力约束和无功补偿

主动配电网和传统配电网最大的区别是分布式电源具备主动调节能力。光伏逆变器的有功出力不是固定值，而是可以在一定范围内调节，通常允许弃光。建模时，节点注入的有功需要改写成一个变量而不是常量：

# 光伏有功出力变量: 0 到 1.2 倍额定装机之间可调 p_pv_min = 0.0 p_pv_max = 1.2 # 单位为pu，假设装机1.2MVA P_pv = cp.Variable(nonneg=True) cons.append(P_pv <= p_pv_max) cons.append(P_pv >= p_pv_min) # 无功出力受逆变器视在容量约束 Q_pv = cp.Variable() cons.append(cp.norm(cp.hstack([P_pv, Q_pv])) <= 1.2)

无功容量约束实际上是一个圆盘约束，用cp.norm写成二阶锥形式，还能融入SOCP框架。这里要注意逆变器容量约束通常写成P_pv²加Q_pv²不超过S²，这是凸约束可以直接加，不需要额外近似。

光伏出力的上限还要乘一个气象系数或者预测曲线，真实场景里p_pv_max是一个时段序列，不是常量。在静态单时段问题里，先按当前时刻光照强度折算一个上限值，就能满足大部分分析需要。

4.2 储能系统：跨时段耦合的SOC约束

储能是主动配电网最有价值的调节资源，但建模会引入跨时段耦合，单时段最优潮流直接解不了。原因是储能电池的荷电状态（SOC）满足递推关系：下一时刻的SOC等于当前SOC加上充电功率减去放电功率再乘效率系数。这个递推关系让每个时段的优化不能独立求解，必须把多个时段放在同一个问题里联立求解，这就成了动态最优潮流。

在SOCP框架里，储能模型也不复杂，关键是把充放电功率拆成充电和放电两个非负变量，并且加上两者不同时为正的约束。这个约束虽然是非线性的，但因为目标函数通常会让储能尽量有效率地运行，实践中可以省略互补约束，只靠目标函数里的充放电成本项来避免同时充放电的荒谬结果。

储能约束的SOCP写法：

# 储能变量 SOC = cp.Variable(T + 1, nonneg=True) # T个时段，多一个是初始时刻 P_ch = cp.Variable(T, nonneg=True) # 充电功率 P_dis = cp.Variable(T, nonneg=True) # 放电功率 Q_ess = cp.Variable(T) # 无功出力 # SOC递推 cons.append(SOC[0] == 0.5) # 初始SOC 50% for t in range(T): cons.append(SOC[t+1] == SOC[t] + P_ch[t] * eta_ch * dt / E_cap - P_dis[t] / eta_dis * dt / E_cap) cons.append(SOC[t+1] >= 0.1) cons.append(SOC[t+1] <= 0.9) cons.append(P_ch[t] <= P_ess_max) cons.append(P_dis[t] <= P_ess_max) cons.append(cp.norm(cp.hstack([P_dis[t] - P_ch[t], Q_ess[t]])) <= S_ess_max)

SOC递推公式里的dt和E_cap要换算成统一的标幺值体系，不然数值会直接崩掉。举个例子，假设储能容量2 MWh，基准功率1 MW，时间步长1小时，dt除以E_cap就是1除以2等于0.5，SOC方程的系数就是这个量级。如果基准值选得不合适，SOC的数值在0到1之间变化，但充放电功率的数值可能是几十，两个数量级差太远，求解器数值稳定性会很差。

4.3 可调无功设备：并联电容器和静止无功补偿器

配电网里传统的无功调节设备主要有并联电容器组和有载调压变压器。电容器组是离散设备，投切状态是整数变量，放进SOCP会让问题变成混合整数二阶锥规划（MISOCP）。这倒不是不能解，只是计算代价明显上升，特别是网络规模大时求解时间会增加一到两个数量级。

工程上的替代做法是：先把电容器组当成连续变量求解，得到最优解后再用启发式规则把连续值映射回最近的离散档位，然后固定档位重新求解一次。这个两阶段法的误差通常在一个档位的容量范围内，对工程调度来说够用。我一般会在代码里加一个后处理函数，把连续的Q_cap值圆整到最近的投切档位。

静止无功补偿器（SVG）是连续调节设备，建模就简单得多，直接作为一个无功注入变量，加上容量上下限约束就好。一个配电网里SVG的调节速度和精度都比电容器组高得多，适合应对光伏波动带来的快速电压变化，在做动态最优潮流时优先用SVG做无功支撑。

4.4 目标函数加惩罚项：怎么调权重不翻车

设备接入后，目标函数往往不能只有网损。比如要最小化购电成本、弃光惩罚、储能老化成本等。多目标加权求和时，权重的标定直接影响解的质量。网损以pu为单位通常数值很小（0.001量级），弃光惩罚如果直接写上百倍的系数，等于强制保弃光为零，这个方向不一定是对的。

我常用的做法是先跑一次纯网损最优得到网损量级，再跑一次纯弃光最优得到弃光量级，然后用两个量级之比设定初始权重，再按现场电价反推修正。这样权重至少在一个合理的数量级内，不会被某个目标单方面支配。权重系数确定后，目标函数仍然是线性的或者凸二次型，问题性质不变，求解器不需要换。

5. 二阶锥规划求解最优潮流的五个高频坑：现象、原因、解法

5.1 求解器报infeasible：标幺值系统不统一

现象：代码看起来没问题，约束都对，但求解器直接返回不可行，甚至报“numerical issues”。这是最常见的翻车现场。

原因：网络数据里有的量是标幺值，有的量是有名值，混在一起用。比如支路电阻给了欧姆值，电压基准却用的10 kV，没有除以基准阻抗。在DistFlow模型里，电压平方、功率、阻抗必须全部在同一个标幺值体系下，混用会导致约束条件数量级混乱，可行域被压成一个空集。

解决：写一个统一的数据转换函数，把所有有名值转成标幺值后再进入优化模型。基准阻抗等于基准电压平方除以基准容量，这个值一算出来，整个网络参数就统一了。调试时打印每条支路的r和x，确认范围在0.0001到0.1之间，超过这个范围基本就是标幺值出问题了。

5.2 二阶锥约束写错方向：不等号方向颠倒

现象：求解器能算出结果，但结果明显不合理，比如电压接近下限时网损反而为零，或者支路功率出现荒谬的数值。

原因：旋转锥约束有两种常见写法，一个是‖[2P; 2Q; v_i-l_ij]‖≤v_i+l_ij，另一个等价写法是P²加Q²小于等于v_i乘l_ij。有人会习惯性写成大于等于，这个方向一错，问题就从凸优化变成了非凸优化，求解器给的所谓最优解毫无物理意义。

解决：写完之后做一次独立校验，取任意一条支路，把P、Q、v_i、l_ij的最优值代回到P²加Q²和v_i乘l_ij里，打印两者的差值。正常松弛紧的情况下两者相差不超过1e-6，如果差得很大或者方向反了，立刻检查约束代码。

5.3 结果和真实潮流对不上：忽略松弛非紧

现象：SOCP最优解算出来电压分布和网损都不错，但拿这个解去跑交流潮流计算器，节点电压和支路功率全对不上，误差超过5%。

原因：目标函数或者约束条件把最优解推到了锥的内部而不是边界。典型场景是最小化购电成本时，如果节点电价设定异常，求解器倾向于让某些支路电流虚高，这时锥松弛不是紧的，SOCP解在物理上是不可实现的。

解决：计算锥间隙。对每条支路计算(P_ij²加Q_ij²)除以(v_i乘l_ij)再减去1，取所有支路的最大值。如果这个值大于1e-4，说明松弛损失了精度。针对这种情况，常用的修法是在目标函数里加一个很小的电流惩罚项，比如0.001乘以所有支路的l_ij之和，这样会鼓励求解器把电流压低，让解回到锥边界上。

5.4 储能SOC和功率数量级差距过大导致求解器震荡

现象：加了储能模型后，求解器迭代次数暴涨，有时几百次不收敛，有时收敛了SOC曲线出现锯齿状跳变。

原因：储能容量基准和功率基准不一致。比如基准功率是1 MW，储能容量是5 MWh，时间步长是15分钟，那么dt除以E_cap等于0.25除以5等于0.05，而充放电功率的标幺值可能在0.8左右，SOC数值在0到1之间，三个量级完全不同，导致线性方程组的条件数变得很差。

解决：统一容量和功率的基准。把储能容量也折算成标幺值，或者直接把SOC递推方程的系数整体放大100倍，让约束矩阵的条件数维持在合理范围。另一个更简单的做法是用秒作为时间单位，确保dt和E_cap的量级差不超过100倍。

5.5 电容器组离散档位处理不当

现象：把电容器组当成连续变量求解，结果出来了，但圆整到离散档位后重新跑潮流，电压越限。

原因：圆整操作改变了无功注入量，直接影响节点电压。如果原来的最优解里，某个电容器的无功出力刚好在档位边界附近，圆整到上一档或下一档都可能让电压越过上下限。

解决：不要直接圆整到最近档位，而是优先选择能让电压留在安全区间的档位。具体做法是：先固定其他设备，把电容器组的每个离散档位依次代入原问题求解一次，选目标函数最优且所有约束满足的档位。网络规模大时可以用灵敏度分析缩小候选档位集合，一般不超过5个档位。这个方法是增加计算量的，一般固定档位后再求解一次SOCP即可。

6. 从单时段到动态最优潮流：验证方法和一个进阶技巧

单时段SOCP最优潮流跑通后，下一步自然就是动态最优潮流，把一天24个或者96个时段放进同一个问题里联立求解。储能SOC的递推约束把各时段耦合在一起，模型不再是一个SOCP在独立运算，而是一整个时间块的SOCP，变量数量直接乘以时段数。求解器的性能在这里拉开差距，同一批数据用CLARABEL可能几分钟，换MOSEK只要几十秒，体感差别很明显。

我验证一个SOCP最优潮流结果的习惯动作有三步。第一步看锥间隙，逐条支路检查松弛紧度，数值大于1e-4就要加惩罚项修正。第二步做交流潮流回代，把SOCP给出的节点有功无功注入值固定住，用传统潮流计算器跑一次完整ACPF，对比电压幅值和支路功率误差，2%以内算正常。第三步做时序一致性检查，动态场景里看储能SOC曲线的连续性，SOC跳变意味着约束写错了或者数值有问题。

一个值得尝试的进阶技巧是：把SOCP当成“热启动器”，用它的解作为非线性交流最优潮流的初值。具体做法是先解SOCP得到各DG的出力和储能策略，然后把这个解作为内点法的初始点去解精确的ACOPF。这样ACOPF的迭代次数会大幅减少，而且大概率能收敛到和SOCP接近的局部最优解。等于用凸问题画了一个靠谱的可行域入口，再用精确模型做收尾。这个方法在多个配电网测试系统上都验证过，比直接跑ACOPF稳定得多，也比只用SOCP结果可信得多。

做主动配电网调度这一年多下来，最大的体会是SOCP最优潮流的价值不在理论研究，而在于它给了一套可以反复使用的计算底座。模型从3节点扩展到实际馈线，坑来来去去就那么几个，把标幺值、松弛紧度、求解器选型这三件事管住，大部分问题都能顺下来。希望这篇笔记能帮你少走一段弯路，也希望你第一次跑通代码时，能体会到“看电压曲线从越限被拉回来”的那种踏实感。

本文还有配套的精品资源，点击获取