配电网多时间断面优化调度：二阶锥松弛与YALMIP+CPLEX实战

做配电网优化调度这几年，我最大的感受是：潮流计算本身不难，难的是一套框架能同时处理几十个甚至上百个时间断面的联合分析。这次用YALMIP加CPLEX做二阶锥松弛的多时间断面潮流分析，算是我在配电网时序优化这条路上走得比较完整的一次。项目不算新，但把IEEE33和PG69两个经典算例放进同一个框架，统一建模、统一求解、统一收数据的思路，对很多刚入门的同学来说，可以直接“抄作业”。

这篇文章不打算只贴公式，而是把我实际写代码、调求解器、处理数据的过程拆开。从为什么选二阶锥松弛、两个算例的数据怎么准备，到YALMIP里怎么写出能被CPLEX高效求解的模型，再到多时间断面耦合时的坑和排查方案，全都走一遍。

1. 为什么选二阶锥松弛处理配电网潮流

1.1 传统潮流计算的痛点

配电网潮流问题本质上是一个非线性方程组求解问题。经典做法是牛顿拉夫逊法，或者它的各种改进版本，比如前推回代法。这类方法在单一时间断面做潮流计算时很快很成熟，但如果要把潮流放进一个优化问题里反复求解，事情就变了。因为在优化问题里，潮流方程不是“先解出来再带入”，而是作为约束条件被优化器反复调整。你要让优化器去满足一组非线性的潮流等式，这会让整个问题变成非凸、非线性的复杂优化问题，全局最优基本不用想，局部最优解的质量又非常依赖初始点。

在多时间断面的场景下，这个问题被进一步放大。24个时间断面叠加，潮流方程要同时满足，变量数量爆炸式增长，非线性也叠加得更复杂。如果用传统思路做，初始点稍微不合适，求解器就可能在某个断面卡住，要么不收敛，要么收敛到完全不合理的“伪最优解”。

所以在我做这个项目时，第一个关键决策就是把原问题从非凸非线性优化，通过二阶锥松弛转化为一个凸优化问题。这也是目前配电网优化里最主流的数学处理手法之一。

1.2 二阶锥松弛到底做了什么

二阶锥松弛的思想，简单理解就是：把一个难搞的等式约束放宽成一个“相对容易处理”的不等式约束，而且放宽后的可行域是凸的。在配电网潮流问题里，最核心的非线性来源是节点电压和支路电流之间的耦合关系。

考虑一条辐射状支路ij，DistFlow潮流方程有经典形式：

• Pj = Pij - rij * lij - ΣPjk
• Qj = Qij - xij * lij - ΣQjk
• uj = ui - 2(rij Pij + xij Qij) + (rij² + xij²) * lij

其中lij是支路电流的平方，ui是节点电压的平方。看起来这组方程是线性的是不是？问题出在第二组约束：lij和ui之间存在等式关系lij = (Pij² + Qij²) / ui，这才是非线性的根源，也是让整个问题非凸的那个“罪魁祸首”。

二阶锥松弛做的事情特别直白：把这个非凸等式约束放宽为不等式约束，写成向量形式就是：

|| [2Pij; 2Qij; lij-ui] ||₂ ≤ lij + ui

这个约束的意思是，左侧向量的2-范数不超过右侧值，而这种方式定义的可行域正好是一个二阶锥，是凸集。于是这个问题从非凸优化变成了凸的SOCP问题。

我习惯用一个生活化类比来解释这一步：原来的等式约束要求变量必须精确落在一条曲线上，这个曲线上任意两点连成的直线并不都在曲线上，所以它是非凸的。松弛后变成允许落在“曲线包裹的区域”内，这个区域是向外凸的，所以整体是凸集。只要目标函数是凸的，整个问题就成了凸优化，求解器能找到全局最优，而且不需要人为给定一个特别好的初值。

但松弛涉及一个严谨性问题：松弛后得到的解，是否真的满足原始等式？也就是说，最优点会不会恰好落在锥的边界上。对辐射状配电网来说，在目标函数是网损最小化、且电压约束不上下压缩得过紧的前提下，通常能收敛到原问题的精确最优解。我从IEEE33和PG69两个算例的实测结果来看，松弛间隙都足够小，基本上可以认为等价。

2. 算例系统准备与数据预处理

2.1 IEEE33与PG69算例概述

IEEE33节点系统是配电网研究入门的“hello world”。33个节点、37条支路，其中32条是分段支路，5条是联络开关，基准电压12.66kV，总的有功负荷3715kW，无功负荷2300kVar。经典参数在网上很容易找到，哪种版本的bus数据、line数据都能直接用。

PG69节点系统实际上是69个节点的配电馈线系统，基准电压同样是12.66kV，但拓扑比IEEE33复杂不少，带了很多支线和分支。总负荷一般取有功约3802.19kW，无功约2694.10kVar。因为节点多、支路多，它在验证算法扩展性上很合适，尤其是多时间断面加进去之后，变量规模能到几千甚至上万，更能看出问题建模和求解器配置的重要性。

两个算例放在同一个框架里，模块化做好之后，切换算例其实就是改一份数据文件的事情。这个设计在后面补充储能、分布式电源时特别方便，新增设备类型只需要扩展对应的结构体就行，不需要去动主循环和约束拼接的代码。

2.2 数据处理时的几个坑

数据预处理这一块，看似简单，实则踩坑最多。

第一是单位问题。原始IEEE33和PG69数据里电阻电抗单位是欧姆，功率单位是kW和kVar，电压基准是12.66kV。在YALMIP里建模时，如果直接用原始单位，数值跨度可能很大，CPLEX内部的数值稳定性会受影响。我的做法是把功率统一到标幺值，基准容量S_base取1MVA，电压标幺化，然后把所有阻抗和功率都换算成标幺值。这个步骤做好了，后面求解基本不会出现“数值病态”的问题。

第二是节点编号的1-index和0-index问题。MATLAB本身的索引从1开始，但很多网上流传的算例数据里节点编号从0开始，支路首末端节点也经常搞混。这个看似小问题，实际最容易导致网络不连通或潮流方向出错。我处理时索性写了一个数据检查函数，加载数据后先验证节点编号范围，再检查支路两端是否都在节点列表内，最后用图论工具确认网络的连通性。

第三是负荷基准。多时间断面分析时，各个节点的负荷不是恒定的，我要把每个节点的负荷按照日负荷曲线缩放。很多算例自带的负荷数据是单点的峰值负荷，所以在做24断面的时候，我会用一个典型日负荷系数曲线去乘每个节点的负荷值。这里有一个容易忽略的地方：IEEE33的节点负荷分布差异很大，有些节点接近0，有些接近400kW，缩放比例要全局统一，不能有的节点用这个曲线、有的用新曲线，否则时间断面的对比就失去意义了。

3. YALMIP + CPLEX 建模实操

3.1 建模前的变量定义

YALMIP最大的优势是建模层抽象得干净，你不用去关心CPLEX内部的数据结构。它把目标函数和约束都转成MATLAB的符号表达式，再交给底层求解器。但这也意味着——你得在建模前把变量想清楚，不然写起来会乱。

对于单时间断面的潮流问题，我需要定义三类主要变量：

• 节点电压平方向量u（n个节点，每个节点1个变量）
• 支路电流平方向量l（m条支路，每条支路1个变量）
• 支路有功和无功向量P、Q（m条支路各1个变量）

在多时间断面版本里，我给每个变量都加上时间维索引，比如P变成mt矩阵，第t列表示第t个断面。YALMIP支持直接定义matrix变量，然后在for循环里按列索引来写约束，这样做可以降低约束拼接的复杂度。

变量定义代码大概是这样的：

nb = 33;   % 节点数
nl = 32;   % 支路数（辐射状考虑）
nt = 24;   % 时间断面数

u = sdpvar(nb, nt, 'full');     % 电压平方
l = sdpvar(nl, nt, 'full');     % 电流平方
P = sdpvar(nl, nt, 'full');     % 支路有功
Q = sdpvar(nl, nt, 'full');     % 支路无功

这里要说明一点：如果只考虑辐射状网络，支路数就是节点数减一，总有功负荷和无功负荷是给定的参数。如果网络带联络开关，结构会更复杂，但为了建模方便，通常把联络开关当作可选的额外支路来处理。

3.2 目标函数与约束的YALMIP表达

目标函数我选择网损最小化，这也是配电网经济调度的最常见目标。全系统网损是所有支路上的电流平方乘以电阻之和。在标幺制里，第t个断面的网损就是：

loss_t = sum(r_all .* l(:, t));
total_loss = sum(loss_t); % 把所有时间断面的网损加起来

这里r_all是支路电阻向量，l为电流平方，总网损就是各支路i的电流平方乘以电阻，再对所有支路求和，然后把24个时间断面加起来。因为电流平方乘以电阻就是功率形式的损耗，所以这个表达式是线性的——没错，目标函数反而是整个模型里最“轻松”的部分。

约束条件分三块。第一是DistFlow潮流方程的节点功率平衡，第二是根节点电压固定（通常设为1.0pu），第三是二阶锥松弛约束本身。

节点功率平衡约束写成YALMIP格式就是：

Constraints = [];
for t = 1:nt
    for k = 1:nl
        i = from(k); j = to(k);
        % 支路k的功率注入：从i到j
        Constraints = [Constraints, ...
            P(k,t) - r(k)*l(k,t) - sum(P(find(from==j),t)) == ...
            load_P(j,t)];
        % 这里load_P(j,t)是节点j在t时段的负荷有功
        % 实际需要把支路k发出的功率等于节点j的负荷加上下游支路功率
    end
end

注意，上面的写法只是示意。实际操作中，为了减少for循环嵌套带来的计算负担，我通常把DistFlow约束写成矩阵形式，或者用YALMIP的repmat和索引技巧批量生成约束。因为节点数和时间断面数相乘之后，约束规模会显著上升，for循环写得不好会让建模阶段本身就慢到怀疑人生。

第二块电压约束，包括根节点电压固定，以及所有节点电压的上下限约束：

Constraints = [Constraints, u(1, :) == 1.0];
Constraints = [Constraints, 0.95^2 <= u <= 1.05^2];

上下限取0.95到1.05（标幺值），对应电压偏差±5%。这里用平方是因为u是电压幅值的平方，所以电压上限1.05对应u上限1.1025。

第三块二阶锥约束是核心，也是YALMIP里最容易写错的地方。我的写法是循环每条支路、每个时间断面，用cone函数构建二阶锥约束：

for k = 1:nl
    for t = 1:nt
        ui_kt = u(from(k), t);
        uj_kt = u(to(k), t);
        lij_kt = l(k, t);
        Pij_kt = P(k, t);
        Qij_kt = Q(k, t);
        % 二阶锥约束：|| [2P; 2Q; lij-ui] ||_2 <= lij + uj
        Constraints = [Constraints, ...
            cone([2*Pij_kt; 2*Qij_kt; lij_kt - uj_kt], lij_kt + uj_kt)];
    end
end

这里要注意的是，cone函数第一个参数是锥内向量，第二个参数是标量上限。如果写成cone(v, t)，含义是norm(v,2) <= t。YALMIP会自动把这个约束识别为二阶锥约束，并传递给CPLEX的SOCP接口。

3.3 求解器调用与参数设置

模型拼好后，调用CPLEX求解的方式异常简单：

ops = sdpsettings('solver', 'cplex', 'verbose', 2);
ops.cplex.mip.tolerances.integrality = 1e-6;  % 如果是混合整数模型
ops.cplex.simplex.tolerances.optimality = 1e-6;
ops.cplex.simplex.tolerances.feasibility = 1e-6;

sol = optimize(Constraints, total_loss, ops);

CPLEX处理SOCP默认用的是barrier算法，对于数千变量的凸问题，求解速度很快。我的实测中，IEEE33单断面SOCP基本在1秒内收敛，24断面也在10秒量级。PG69的24断面可能要到几十秒到几分钟，这个要看具体的负荷曲线和约束紧度。

关于CPLEX参数，我重点调了三个：

• optimality tolerance （对偶可行性容差）和 feasibility tolerance （原始可行性容差）：默认值一般是1e-6，如果遇到数值问题，可以适当放宽到1e-5或1e-4试试。
• barrier crossover ：CPLEX解完barrier之后默认会做crossover，把最优解转成基本可行解，这个步骤在某些病态模型上会额外耗时。如果只是跑潮流分析，不需要对偶信息，可以直接设置 ops.cplex.solutiontype = 2 ，强制使用barrier解而不做crossover，可以省不少时间。
• display frequency ：控制输出日志密度，前端不关心时可以把verbose调低。

4. 多时间断面：从单断面走向时序优化

4.1 时间耦合约束怎么写

单断面的SOCP潮流好写，但多时间断面的核心难点不在“每个断面都要满足潮流”，而在“断面之间是有联系的”。最典型的联系就是储能设备：t时段的充电行为，会直接影响t+1时段的荷电状态。也就是说，这是个时序耦合问题，不是简单的24个独立子问题拼凑。

在YALMIP里，这种时间耦合实现起来非常直白。以储能系统为例，SOC变量是一个向量，第t+1时段的SOC由第t时段的SOC、充电功率和放电功率共同决定：

SOC = sdpvar(nt+1, 1, 'full');
Pch = sdpvar(nt, 1, 'full');   % 充电功率
Pdis = sdpvar(nt, 1, 'full');  % 放电功率
eff_ch = 0.95;                 % 充电效率
eff_dis = 0.95;                % 放电效率
dt = 1;                        % 时段间隔，小时

Constraints = [Constraints, SOC(1) == 0.5]; % 初始SOC为50%
for t = 1:nt
    Constraints = [Constraints, ...
        SOC(t+1) == SOC(t) + (eff_ch*Pch(t) - Pdis(t)/eff_dis) * dt];
    Constraints = [Constraints, 0 <= Pch(t) <= 0.5]; % 最大充放电功率0.5MW
    Constraints = [Constraints, 0 <= Pdis(t) <= 0.5];
    Constraints = [Constraints, 0 <= SOC(t+1) <= 1]; % SOC在0~1之间
end

同时，储能在节点的注入功率可以简单建模为Pdis减Pch的差值，再加上该节点的原始负荷。这样潮流方程里的负荷项变成了净负荷项，就完成了储能与潮流模型的耦合。

除了储能，风机、光伏等分布式电源的出力也可以按时间序列给定，这时DG出力是已知参数，不需要额外优化变量。如果DG出力可调（比如微型燃气轮机），那就要增加对应时段的出力变量，并加上爬坡约束：

Constraints = [Constraints, Pgen(t+1) - Pgen(t) <= ramp_up];
Constraints = [Constraints, Pgen(t) - Pgen(t+1) <= ramp_down];

这些都是非常简单实用的时间耦合约束，能把“多断面潮流”从静态快照升级成一个真正的时序优化问题。

4.2 储能与分布式电源建模注意事项

在储能和DG接入多断面潮流时，有两个特别容易被忽略的细节。

第一个是“节点注入方向”的符号约定。我在最初写存储建模时，把充电功率写成正，放电功率写成负，结果潮流方程里负荷项怎么都不对——因为充电时储能是从电网吸收功率的，净负荷应该增大；放电时是向电网注入功率的，净负荷应该减小。这个问题最好在建模前就用文档把符号约定固定下来，并且做一次单节点测试来验证，不然排查起来非常痛苦。

第二个是并联DG与储能在同一节点的处理。如果同一个节点上既有光伏出力，又有储能充放电，又有一个峰值负荷，净负荷的计算必须把三者统一加起来。我在项目里用了一个load_net矩阵来存储所有节点在所有时间断面的净负荷，这样在拼DistFlow方程时，只需要引用load_net(j,t)，而不是分别处理每个设备类型。这个方法强烈推荐，它能极大简化约束代码。

4.3 计算结果怎么解读

求出优化结果后，并不是画几张图就完了。我一般会做三件事验证结果合理性。

第一件，检查二阶锥松弛的紧度。计算每一支路每一断面的潮流结果，验证不等式约束是否取等。如果松弛间隙过大，说明原始非凸问题没有被等价松弛，可能是电压约束太紧导致锥约束没有“贴到”边界上。这时候要看是否真的需要原问题的最优解，还是可以接受有界的近似解。

第二件，验证潮流方程是否真的满足。用求得的最优解回代DistFlow方程，检查等式平衡是否在容差范围内。CPLEX的SOCP默认解的最优性和可行性有严格保证，但数值误差还是可能存在，特别是大规模问题时。

第三件，可视化展示电压曲线和网损曲线。IEEE33和PG69的结果差异很明显——优化前PG69在高峰时段的末端电压跌落比IEEE33严重得多，接入储能和DG后的改善幅度也更显著。这类表格或曲线图，比任何分析文字都有说服力。

5. 常见问题与排查技巧实录

5.1 求解器报错速查表

实际跑模型时，我遇到过几类非常典型的问题。这里整理成一个速查表，基本覆盖YALMIP加CPLEX日常能见到的报错。

5.2 数值稳定性与初始点问题

CPLEX作为商业求解器，数值稳定性已经做得相当好，但面对上千个变量、几千条约束的SOCP模型，缩放不好照样会出问题。

我最常用、也最推荐的做法是把网络参数和功率都转成标幺值。以S_base=1MVA为例，整个系统的潮流结果通常在0.1到10之间，电压平方在0.9到1.1之间，这正好落在CPLEX喜欢的数值范围内。如果你用原始单位（比如功率单位是W），数值会到几万甚至几十万，CPLEX内部计算精度很容易被影响。

另外，在多时间断面模型里，如果时间断面数特别多（比如96个点），可以考虑先用24断面模型跑通，再扩展到96断面。这样排查问题时，问题规模小，出错更容易定位。

5.3 YALMIP建模效率提升技巧

YALMIP建模阶段有一个容易忽略的性能瓶颈：在for循环里不断拼接约束矩阵。如果节点数较多，断面上百，for循环加约束会非常慢。我有一次跑96个断面的PG69，建模阶段就卡了十几分钟，后来发现是约束不断扩充导致矩阵重新分配内存。

解决办法有两个方向。第一个是用YALMIP的约束数组批量生成，把约束放进cell数组，最后一次性合并，避免中间频繁拼接。第二个是尽量利用矩阵向量化思想，把同一类约束在多个时间断面上用repmat或者reshape批量生成，而不是逐断面循环。

实测下来，这两种方式都能把建模时间缩短一个数量级以上。对于大算例，优化求解前的建模效率，往往比求解器求解时间更值得花心思去优化。

6. 工具链选型心得

6.1 为什么选CPLEX而非Gurobi

关于求解器选择，我最初也犹豫过Gurobi还是CPLEX。两个都能解SOCP，性能和稳定性都很好。我最终用CPLEX的原因其实很庸俗——项目环境里已经装了CPLEX，学术版许可也方便，所以直接用。

但如果你要开始一个新项目，我建议Gurobi和CPLEX都测试一下，特别是拿你的模型各跑几个算例比较时间。不同模型对两个求解器的敏感度不一样，尤其SOCP问题，有些场景Gurobi的barrier会更占优势，有些场景CPLEX更快。不要迷信任何一方，拿自己的实际数据测一测最靠谱。

6.2 YALMIP之外的另一个选择

如果不喜欢YALMIP，MATLAB里还有MATPOWER的某些扩展模块，或者用CVX。CVX建模SOCP确实更方便，但CVX对MISOCP的整数变量支持不如YALMIP灵活。如果你的多时间断面模型要加入储能投切、联络开关切换等0-1离散变量，YALMIP的表达能力会明显更强，而且可以直接把CPLEX的MIP参数透传进去调优。

根据我个人的习惯，如果只是简单SOCP，CVX完全可以胜任；一旦涉及混合整数或者更多自定义约束组合，我倾向于YALMIP。

7. 多时间断面扩展的方向

这个项目的模型还可以往几个方向延伸。最直接的，是把单相模型扩展为三相不平衡模型，考虑配电网三相负荷不平衡的情况，这在低压配电网中特别有意义。另外，可以在模型中增加需求响应约束，把可中断负荷或者可平移负荷作为优化变量，让负荷不再是固定的，从而优化整个时间断面组合。

另一个值得提的方向是把机会约束加进来。分布式电源出力的随机性、负荷预测误差，都可以用机会约束来处理，让模型在考虑不确定性的同时还能保持计算效率。这个方向可以把多时间断面潮流分析从“确定性优化”推向“鲁棒优化”，解出来的结果更有工程参考价值。

代码层面，把IEEE33和PG69的模型扩展成PG80甚至更大规模，主要工作量在数据整理，建模框架不用改，这也是当初做好抽象封装的好处。

8. 实操经验收尾

最后分享两个我在这个项目里踩过比较多遍的坑，算是给要动手复现的朋友提个醒。

第一个是别把多时间断面想得太简单。24个断面不是把单断面模型复制24份就行，时间耦合才是灵魂。储能SOC的递推、DG出力的时序曲线、甚至负荷曲线的时间颗粒度，每一个都牵一发动全身。建议先跑通2到3个断面的小规模测试，验证时间耦合逻辑后再扩充。

第二个是数据可视化一定要趁早做。不要把结果导出来之后才开始画图，我习惯在每次求解完成后，立刻画一张电压曲线总览图、网损曲线图，以及松弛紧度热力图。这样调模型、调参数时能快速看到影响，而不是等所有计算都搞完再做后处理。

这个项目的核心路线——DistFlow方程建模、二阶锥松弛、YALMIP建模、CPLEX求解——放到现在依然是配电网时序优化的经典套路。把这条链路走通，再去碰更复杂的鲁棒优化、多目标优化，底子就扎实了。