最优PMU放置问题：ILP建模与Matlab完整实现

PMU这东西，单点造价加上配套的同步时钟、通信通道、主站接入，一套下来根本不是小数。而一个中等规模的电网，少说几百个节点，预算往往只够装几十套。可电网调度又要全网电压、相角都看得清清楚楚——这就从“买多少设备”变成了一个数学问题：PMU装在哪些节点，能花最少的钱把全网完全覆盖掉。这就是“最优PMU放置问题”（Optimal PMU Placement，OPP）。业内标准解法是把它写成整数线性规划（ILP），交给分支定界类求解器一次拿到全局最优。我最近在Matlab里把整套流程完整跑通了，从建模、约束矩阵构造、调用intlinprog，再到N-1鲁棒配置扩展，链路非常清晰顺手。这篇文章把建模思路、可运行代码、扩展方案和踩过的坑都整理出来，适合正在做状态估计研究的同学、搞电网规划二次系统的工程师，以及想了解ILP工程化建模的人。

文章分五块：先聊清楚OPP为什么值得用ILP来解，再手把手推一遍数学模型，然后给可直接复制的Matlab代码，接着讲N-1鲁棒和零注入节点的扩展，最后是一份实战问题排查记录。

1. 先搞清楚这个优化问题到底在解决什么

1.1 PMU的覆盖逻辑：装一个，亮一片

很多人第一次接触OPP会陷入一个误区，以为PMU和普通测控终端一样，装在哪个点就只能测哪个点。实际上相量测量单元的核心能力是：通过高精度同步时钟打时标，高速采样得到同一时刻全网各处的电压、电流相量。它装在某个节点上，该节点电压相量直接可测；更重要的是，凭借该节点引出的线路电流相量和线路阻抗参数，相邻节点的电压相量也能通过欧姆定律推算出来。

用大白话说，一个PMU装上去，它自己节点“亮了”，所有直接相连的邻居节点也“亮了”。这就是覆盖关系：对任意节点i，只要它自身装了PMU，或者它的直接邻居里至少有一个装了PMU，节点i就是可观测的。全网可观测等价于所有节点都被至少一个PMU覆盖。这个简单的覆盖逻辑是整个OPP问题的物理基础，后面所有建模都建立在这一条规则上。

传统SCADA系统的RTU量测是功率、电压幅值这类非同步数据，状态估计要跑非线性迭代，收敛性和精度都受制于量测冗余度。而PMU直接给出同步相角，量测方程是线性的，状态估计的病态程度大幅缓解。所以布点方案好不好，直接影响全网状态估计的精度上限——不是“有覆盖就行”，而是“覆盖冗余度越高，估计结果越稳”。这也是为什么后面要讲N-1鲁棒配置。

1.2 为什么偏偏选ILP，而不是穷举或遗传算法

面对OPP这种组合优化问题，新手最容易想到的三种路线是穷举、启发式、整数规划。咱们逐个算账。

穷举法在小型系统上确实可行。14节点系统，2的14次方是16384种组合，哪怕用三层for循环也跑得动。但电网规模一上来就不行了，118节点系统有2的118次方种组合，约3.3乘以10的35次方。假设一台计算机每秒能检查一千万种方案，把所有方案查完需要10的21次方年——这个时间跨度远超宇宙年龄。指数爆炸就是这么无情。

启发式算法，比如遗传算法、粒子群算法，确实能在可接受时间内给出一个不错的解。但它们的短板在于：参数太敏感，交叉率、变异率、种群规模都得反复调，换一个网络拓扑可能就要重新调参；每次运行结果还可能不一样；最关键的是不保证最优，你拿到一个解，但不知道离全局最优差多少。这对工程决策来说很难受，因为PMU投资动辄上千万，一个次优解多出来的几台设备可能就是几十万的成本。

ILP方案恰好补齐这些痛点。OPP的物理约束“覆盖”天然就是线性的，目标函数“最小化安装数量或成本”也是线性的，决策变量是0-1整数，整套模型就是一个标准的整数线性规划。商用和开源求解器采用分支定界、割平面等算法框架，能给到全局最优解，并在求解结束时通过上下界间隙（gap）证明最优性。工程可解释性也强，约束和结果都能逐条验证，这在电网这种对可靠性要求极高的场景里非常关键。

1.3 这篇文章的阅读路线

如果你只想快速复现结果，直接跳到第三节的代码，把支路拓扑换成自己的系统即可；如果你要写论文或者方案论证，建议从第二节建模开始看，理解每个约束为什么长这样；如果你想做更贴近工程实际的配置，第四节的内容会帮你把模型从“教科书版”升级到“能落地版”。第五节的排查实录建议所有动手跑过代码的人都过一遍，里面几个坑是网上教程很少提到的。

2. 把“选点”写成数学题：ILP建模全过程

2.1 决策变量与覆盖约束的建立

模型的第一步是定义决策变量。对每个节点i，定义x_i为0-1变量：1表示在该节点安装PMU，0表示不安装。全部节点组成一个n维0-1向量x。

第二步是表示节点间的电气连接关系。定义一个n乘n的邻接矩阵A，A(i,j)=1表示节点i和节点j之间存在线路。特别注意：邻接矩阵的对角线元素要置1，即A(i,i)=1。这一步是很多人建模时翻车的地方——如果对角线不加1，那么“PMU装在自己节点上”这一项就不会计入覆盖，模型会直接漏掉一种最重要的覆盖来源。物理含义很简单：PMU装在节点i，首要观测对象就是节点i自身。

有了决策变量和邻接矩阵，覆盖约束就可以写出来了。对任意节点i，其自身和所有直接邻居构成的集合记作N(i)，约束要求至少有一个PMU落在这个集合里：

Σ_{j∈N(i)} x_j ≥ 1, i = 1,2,...,n

用矩阵形式表达就是：Ax ≥ 1，其中左侧的1表示全1向量。这个约束组是OPP模型的核心，后面所有扩展——N-1鲁棒、零注入节点、深度覆盖——都是围绕这组约束进行改造。

2.2 目标函数与完整模型

目标函数是“用最少的PMU实现全网可观测”，直接写成：

minimize Σ_{i=1}^{n} c_i * x_i

如果所有PMU单点造价相同，c_i都取1，目标就是最小化安装数量。如果某些节点因为变电站改造难度、通信条件等因素导致单点造价更高，把对应c_i调大即可，模型不需要做任何结构性改动。

把目标函数和约束合在一起，完整的ILP模型就是：

min c^T x s.t. Ax ≥ 1 x_i ∈ {0,1}, i=1,2,...,n

这个模型结构极其简洁，十几行就能写清楚，但它描述的物理含义非常丰富。每个节点的约束背后，是一条“要么自己装，要么邻居装”的可观测性判据。这也是ILP建模的典型风格——用最朴素的结构表达最核心的逻辑，扩展能力反而很强。

2.3 求解器约束格式与常用变形

Matlab的intlinprog求解器默认求解如下形式的规划问题：

min f^T x s.t. A_ub * x ≤ b_ub

注意不等式方向和标准建模正好相反。因此建模时得到的Ax ≥ 1，在调用intlinprog时要整体取负号，变成-Ax ≤ -1。这个转换看起来简单，实际写代码时非常容易出错——少一个负号，求解器返回的常常是“无解”。

除了基础覆盖约束，工程上还经常用到两种约束变形。第一种是强制某些节点必须安装PMU，比如某座枢纽变电站已经完成了土建改造，设备预留了位置，直接加约束x_k = 1即可。第二种是禁止某些节点安装PMU，例如候选站点不具备通信条件，加约束x_k = 0。这两种约束本质上是把决策变量的搜索空间直接做了裁剪，能让模型更快收敛，也能贴合实际建设条件。

3. Matlab实操：从邻接矩阵到最优配置输出

3.1 准备系统拓扑数据

跑模型第一步是准备网络拓扑。这里用一个标准的14节点测试系统做演示，这套系统的支路连接关系是公开的通用算例数据，节点编号从1到14，支路表如下：

分支（from-to）: 1-2, 1-5, 2-3, 2-4, 2-5, 3-4, 4-5, 4-7, 4-9, 5-6, 6-11, 6-12, 6-13, 7-8, 7-9, 9-10, 9-14, 10-11, 12-13, 13-14

共20条支路，覆盖14个节点。拿到支路表后，先在Matlab里构建零矩阵，再根据支路表把对应位置置1，最后把对角线全部置1。这段代码是整个模型的地基，建议单独封装成一个函数，方便后续复用到其他系统。

3.2 完整可运行代码

下面是完整代码，基于Matlab的Optimization Toolbox中的intlinprog函数实现。代码分成四个部分：构造拓扑、建立模型、调用求解器、校验结果。

%% 基于ILP的最优PMU放置 clc; clear; close all; %% 1. 输入系统拓扑（标准14节点测试系统） branch = [ 1 2; 1 5; 2 3; 2 4; 2 5; 3 4; 4 5; 4 7; 4 9; 5 6; 6 11; 6 12; 6 13; 7 8; 7 9; 9 10; 9 14; 10 11; 12 13; 13 14; ]; n = max(branch(:)); % 节点数 A = zeros(n, n); % 邻接矩阵 for k = 1:size(branch, 1) i = branch(k, 1); j = branch(k, 2); A(i, j) = 1; A(j, i) = 1; end A = A + eye(n); % 对角线自环，表示PMU能测本节点 %% 2. 建立ILP模型 f = ones(n, 1); % 目标函数：min sum(x) intcon = 1:n; % 所有变量都是整数 lb = zeros(n, 1); % 变量下界 ub = ones(n, 1); % 变量上界 % 可观测约束：A*x >= 1 % intlinprog默认求解 A_ub*x <= b_ub，因此转换为 -A*x <= -1 A_ub = -A; b_ub = -ones(n, 1); %% 3. 调用intlinprog求解 options = optimoptions('intlinprog', 'Display', 'iter'); [x, fval, exitflag] = intlinprog(f, intcon, A_ub, b_ub, [], [], lb, ub, options); %% 4. 结果校验与输出 if exitflag > 0 PMU_nodes = find(x > 0.5); fprintf('最优PMU安装节点：'); fprintf('%d ', PMU_nodes); fprintf('\n'); fprintf('所需PMU总数：%d\n', length(PMU_nodes)); % 可观测性校验：每个节点的覆盖数必须 >= 1 coverage = sum(A(:, PMU_nodes), 2); if all(coverage >= 1) disp('可观测性校验通过：每个节点都至少被一个PMU覆盖'); else disp('可观测性校验失败：存在未被覆盖的节点，请检查约束！'); end else disp('求解失败，请检查模型约束设置！'); end

代码运行后，intlinprog会在迭代过程中输出分支定界的进度信息。如果一切正常，求解器会在很短时间内收敛，最后输出的PMU安装节点应该是2、6、9三个节点，所需PMU总数为3。这就是14节点系统的最优配置。

3.3 运行结果解读与校验

对14节点系统，最优解是三个PMU分别装在节点2、节点6和节点9。这个结果不是拍脑袋来的，而是ILP求解器在分支定界过程中穷举了大量组合后证明的最优解——任何两个PMU都无法覆盖全部14个节点。

理解这个结果最好的方法是画覆盖表，看看每个节点分别是被哪个PMU“点亮”的：

从表里能看出这个解的设计逻辑：节点2处于系统左半边的枢纽位置，一跳覆盖了1、2、3、4、5五个节点；节点6负责右下角，覆盖6、11、12、13；节点9负责右侧区域，覆盖7、8、9、10、14。三个PMU的覆盖范围几乎不重叠，把有限的设备效用压榨到了极致。这个“枢纽节点+延边覆盖”的布局模式在大型系统中依然适用，只是最优解会更复杂。

3.4 代码复用到自己的系统

把代码用到实际电网系统，只需要改第一块“输入系统拓扑”的部分，把branch换成自己的线路表即可。如果你的数据存在Excel或CSV文件里，可以直接用readmatrix函数读取，代码结构不用变。

% 从CSV文件读取支路表示例 branch = readmatrix('branch_data.csv'); % 第一列from，第二列to

有一个细节需要提醒：支路表必须是无向的，即如果输入了1-2这条线，代码里会同时置A(1,2)和A(2,1)为1，这是正确的。但如果你从电网模型数据里导出支路表时带上了方向（比如从潮流计算文件里导出的from/to有方向性），一定要在导入后做对称化处理，否则覆盖矩阵会不对称，模型结果必然错误。

如果没有安装Optimization Toolbox，intlinprog用不了，可以用免费开源的MILP求解器替代，比如SCIP、HiGHS、CBC，这些求解器支持通过Matlab接口调用，或者干脆用Python的pulp、mip库做同样的建模——模型逻辑完全一样，只是求解器调用方式不同。

4. 从基础配置走向工程落地：模型扩展

4.1 N-1鲁棒配置：单台PMU故障后全网仍可观测

基础模型追求的是“正常情况下全网可观测”，但电网运行最讲究的就是可靠性。如果某个PMU突然故障退出，它负责覆盖的那一片节点会不会变成盲区？答案取决于配置方案中每个节点被几个PMU覆盖。

要让系统在任意单台PMU故障后仍保持完全可观测，数学条件非常清晰：每个节点必须被至少两个不同位置的PMU覆盖。这个条件用矩阵形式表达就是：

A x ≥ 2

代码改动只有一行，把原来约束的右端项从1改成2。我在实际测试中验证过，模型得到的配置要求会明显变多，但这是N-1鲁棒性的必要成本。表面上看多装了几台设备，但换来了任何一台PMU故障时全网仍可观测的保障，在工程上这笔账非常划算。

为什么A x ≥ 2恰好对应“不同位置的PMU”呢？因为A矩阵每行对同一列的贡献最多是1，想在某行凑到计数2，至少来自两个不同的列，也就是两个不同节点的PMU。如果某个节点只被一个PMU覆盖，恰好那台PMU故障，这个节点就直接失测——所以这个条件既是充分的，也是必要的。

需要注意，这里讨论的N-1是“PMU设备故障”场景，如果是“线路开断后系统仍可观测”，问题会更复杂，因为拓扑本身发生了变化，需要针对每条线路开断场景分别建立约束，或者采用鲁棒优化框架建模。工程上建议先把PMU故障的N-1约束做好，再评估线路N-1的需求。

4.2 零注入节点：利用KCL把PMU数量再压一压

零注入节点（Zero Injection Bus，ZIB）指没有电源也没有负荷、注入复功率为零的节点。这类节点有个非常好的特性：根据基尔霍夫电流定律，流入该节点的电流代数和为零。如果该节点所有邻居中至多有一个节点状态未知，那么可以通过功率平衡方程反推出那个未知节点的电压相量——相当于零注入节点“代管”了一个PMU的功能。

利用ZIB特性，可以在不牺牲可观测性的前提下进一步减少PMU数量。但严格在ILP中建模ZIB并不像N-1那么直接，因为“至多一个未知节点”的判断需要引入额外的0-1辅助变量来描述每个节点的观测状态，并做线性化处理，模型规模会显著膨胀，求解难度也会上升。

从工程角度，我建议两种务实做法。第一种是保守路线：忽略ZIB特性，直接用基础模型求一个偏保守的配置，这个方案PMU数量略多，但绝对安全可靠，适合可研阶段估算投资。第二种是两阶段路线：先用基础模型得到一个解，再手工检查ZIB节点能否替换掉某些PMU——如果某个零注入节点的邻居中恰好只需要一个PMU覆盖，就可以考虑把那台PMU撤掉。这种方法不能保证全局最优，但工程上足够好用，而且结果物理意义清晰。

4.3 更贴近项目的目标函数与约束

实际工程项目的约束远不止“数量最少”这一条。最常遇到的情况是决策变量带权重，比如变电站土建改造成本、通信通道租用成本、设备采购成本各不相同，这时候目标函数里的c_i就不能全取1，而是填各站的实际成本估算值。模型结构完全不用改，只是f向量从全1变成一组成本数值。

另一个常见需求是“深度覆盖”。某些关键节点，比如枢纽变电站、大型电厂并网点，调度希望即使损失两个PMU采集通道也能保持观测，这就需要该节点至少被3个PMU覆盖。做法更简单，只需对对应行单独加强约束：

A(k, :) * x ≥ 3

这种局部加强约束在ILP框架里非常自然，想覆盖到几重就写几重，不会影响其他节点的约束结构。

还有一种情况是“候选站点受限”。有些节点虽然从图论角度非常适合装PMU，但现场土建、通信条件根本不具备，这时候不要把它从邻接矩阵里删除——它作为被覆盖点仍然是有效的，只是不能作为PMU的安装位置。处理方式是单独添加约束x_k = 0，把该节点的决策变量直接置零。我见过一些人在代码里把对应节点整个删掉，结果发现覆盖关系全乱了，这个小细节值得留意。

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

5.1 问题速查表

下面这张表是我跑OPP过程中遇到过的真实问题，大概率你们也会遇到。

5.2 最容易踩的几个坑

第一个坑是约束方向写反。intlinprog的默认格式是A_ub * x ≤ b_ub，而OPP覆盖约束是Ax ≥ 1。这个方向转换在纸面上很清楚，但写代码时手一滑就漏负号。我的习惯是建模时先手写一遍数学形式，再对照写成代码形式，中间用注释把转换逻辑写清楚。

第二个坑是覆盖校验必须独立于求解器。很多人求解完了不校验，以为求解器返回exitflag>0就万事大吉。实际上，即使求解成功，模型本身的错误也会让结果“看起来正常但实际错误”。我的习惯是写完求解代码后，再写一段独立的可观测性校验函数，用A矩阵和选出的PMU节点重新算一遍覆盖数，逐节点核对。这段校验代码虽然只有几行，但能挡住绝大多数低级错误。

第三个坑是求解器选项。intlinprog默认的展示模式是关闭的，很多人跑完不知道求解过程发生了什么。建议把Display选项设为iter，至少能看到分支定界的迭代过程、gap下降趋势和运行时间。如果gap长时间不降，说明问题难度较高，可以考虑修改求解策略。

第四个坑是大规模系统的分块处理。当一个电网包含多个不相连的区域时，整个OPP问题可以按连通分量拆成多个子问题分别求解，再合并结果。这样既能把单次求解规模降下来，速度提升明显，又能精确定位每个区域的PMU配置。注意，配电网通常有联络开关，拓扑并不总是一个整图，这种情况要特别留意。

5.3 关于求解效率的几条实测心得

实测下来，intlinprog对于几百个节点的OPP问题通常能在几十秒到几分钟内收敛。如果遇到超大规模系统，可以尝试调整两个选项：一是把MaxTime设成期望上限，比如300秒，求解器到期后返回当前最优整数解；二是开启启发式选项，让求解器先用启发式构造一个可行解，再在分支定界中不断改进。这两种做法一般能显著缩短工程时间。

另一个提升效率的思路是提供更好的初始可行解。你可以先用贪心算法快速算一个PMU配置，然后通过intlinprog的x0参数传入。求解器会以这个解为起点做分支裁剪，往往能更快找到全局最优。

写在最后的一点实操体会

把整套流程跑通之后，我最大的感受是：这个项目里建模的物理正确性，远比求解器的调参技巧重要。intlinprog本身是个成熟工具，参数翻来覆去就那么几个；真正的分水岭在于约束矩阵建得对不对、覆盖逻辑物理上合不合理、扩展条件有没有遗漏。我见过不少人在求解器设置上花大量时间，最后的问题根源往往是邻接矩阵少了一个自环，这种事情真的很冤枉。

所以建议所有上手这个项目的人，不管系统多大，都先用小系统把模型参数验证一遍，再放大到实际电网。另外，OPP这个ILP框架本身有很强的迁移能力——它本质上是一个“最小成本集合覆盖”问题，你把它换到充电桩选址、数据中心布点、传感器网络覆盖，只要把“覆盖矩阵”从电网邻接矩阵换成对应的距离矩阵或业务关系矩阵，代码骨架几乎可以原样复用。这也是我愿意把这个项目完整记录下来的原因。