最近在做一个电力系统状态估计相关的仿真任务，标题很直白：用WLS（加权最小二乘法）配合PMU（相量测量单元）来估计系统的电压幅值和相角，然后把估计结果和Newton-Raphson潮流计算出来的状态摆在一起对比。整套代码用Matlab实现。
这个题目我第一眼看过去，就觉得它是个特别适合练手的组合：没有花里胡哨的深度学习，也没有动辄上千页的电力系统手册，就是把电力系统分析里两个最核心的计算工具拉到同一个平台下跑一遍。做完这个项目，你对状态估计、潮流计算、量测建模、权重设计这些概念的理解，会从"课本上见过"变成"真的知道它们在代码里是怎么串起来的"。这篇文章我就把整个实现过程的思路、算法细节、代码设计逻辑和踩过的坑都写清楚，给后面做类似方向的朋友一个可直接参考的路线。
1. 为什么要做这个对比：状态估计在电网里到底解决什么问题
先说清楚一个基本问题：状态估计在电力系统里到底干什么。电网调度中心（EMS）每天要做的事，本质上就是"我手里有一大堆量测数据，但我不完全相信它们，我要通过这些数据推断出电网当前最可能的运行状态"。这里的"运行状态"指的就是每个节点的电压幅值和相角。有了准确的状态，下游的潮流分析、安全评估、经济调度才有意义。
但传统SCADA量测有个天生的短板：它量的是节点注入功率和支路功率，这些量本身就是电压的非线性函数，而且SCADA的刷新率低、精度有限，再加上天气、设备老化、通信干扰等因素，数据里常有噪声甚至坏数据。PMU的出现改变了这个局面：它直接测量电压和电流的相量，也就是能直接给出电压幅值和相角，精度高、带时标、刷新速度快。所以把PMU量测融合到状态估计里，是让估计结果更准、更快的关键方向。
那为什么又要和Newton-Raphson潮流对比？这两个东西其实是"数据来源不同但本质相通"的一对：潮流计算是已知节点注入功率，求解电压分布；状态估计是已知带噪声的量测，反推系统状态。如果量测没有噪声，而且量测配置是完全的，那么WLS状态估计的结果就应该和潮流计算的解非常接近。换句话说，NR潮流解出的电压幅值和相角，可以当作一个"理想参考基准"，用来验证WLS+PMU的估计效果到底准不准。这就是这个题目对比逻辑的核心。
从我的实际体验看，这个项目非常适合三类人：一是电力系统方向的研究生，想快速上手状态估计仿真；二是刚进电网相关岗位的工程师，想搞清楚PMU接入对状态估计的实际影响；三是做综合能源、微电网方向的人，需要搭建一个基础的状态估计验证平台。这个项目不大不小，刚好能把理论、代码、结果分析三块都练到。
2. 核心算法拆解：WLS、NR和PMU量测模型
2.1 WLS状态估计的数学骨架
WLS之所以是状态估计的经典算法，核心在于它的目标函数设计得非常直观：让"量测值"和"量测方程计算值"之间的加权残差平方和最小。
[ \min J(x) = [z - h(x)]^T R^{-1} [z - h(x)] ]
其中，z是量测向量，h(x)是量测函数（把状态变量映射到量测值的函数），R是量测误差协方差矩阵。权重矩阵取R^{-1}，意思很明确：误差方差越小的量测，可信度越高，在目标函数里权重就越大。
这个问题的求解需要用Gauss-Newton迭代。迭代的核心公式是：
[ \Delta x_k = (H^T R^{-1} H)^{-1} H^T R^{-1} [z - h(x_k)] ]
[ x_{k+1} = x_k + \Delta x_k ]
这里的H是量测函数h(x)对状态变量x的雅可比矩阵。G = H^T R^{-1} H叫做增益矩阵，它是整个迭代里最关键也最容易出问题的地方：如果量测配置不好，G可能奇异或病态，迭代就没法收敛。
状态变量x在这个问题里就是所有节点的电压幅值V和相角θ（平衡节点相角通常固定为0）。量测函数h(x)则取决于你有哪些量测：如果是节点注入有功P_i，那h_i(x)就是潮流方程里的P_i表达式；如果是PMU电压相量，那h_i(x)就是V_i和θ_i本身，这组量测方程是线性的，处理起来非常舒服。
2.2 PMU量测模型长什么样
PMU接入状态估计最大的优势在于，它的量测方程是线性的。假设节点i装了PMU，那它直接给出两个量：
[ z_{PMU,i} = \begin{bmatrix} V_i \ \theta_i \end{bmatrix} + e ]
也就是"电压幅值 = 状态变量真值 + 噪声"。写成矩阵形式，就是z = H x + e，其中H矩阵里对应的行是稀疏的单位行。这种线性结构意味着一件事：这些量测对状态变量的"指示"非常直接，不需要经过非线性函数间接推算，所以对收敛性和精度都有正面作用。
实际操作中，通常可以把PMU量测和SCADA功率量测同时放进同一个WLS框架，只是它们的雅可比矩阵行不同：功率量测对应的行是潮流方程对V和θ的偏导数，PMU量测对应的行是常数。权重矩阵R^{-1}里，PMU量测对应的权重因为精度高，会明显大于功率量测。
这里有一个很实用的经验：PMU的权重如果设置得过高，可能出现数值上"压倒"其他量测的情况，导致坏数据检测能力下降；如果设置得过低，又发挥不出PMU的优势。我在项目里一般把PMU电压幅值的标准差设为0.001标幺（常规功率量测设0.01标幺），相角标准差设0.001弧度，算下来权重比大约在100倍量级，效果比较平衡。
2.3 Newton-Raphson潮流和状态估计的关系
Newton-Raphson在电力系统的应用主要指向潮流计算。它的核心思路是用牛顿法求解一组非线性功率平衡方程：
[ P_i = V_i \sum_{j \in i} V_j (G_{ij}\cos\theta_{ij} + B_{ij}\sin\theta_{ij}) ]
[ Q_i = V_i \sum_{j \in i} V_j (G_{ij}\sin\theta_{ij} - B_{ij}\cos\theta_{ij}) ]
给定除平衡节点外的有功注入P和PV节点的电压幅值（以及无功注入Q在PQ节点上），通过迭代求解，得到所有节点的V和θ。NR潮流与WLS状态估计的迭代框架非常像，都需要雅可比矩阵，都是反复求解线性方程组。差别只在于：NR里方程右侧是"给定的注入功率"，而WLS里"给定的"是带噪声的量测向量。
在实际项目里，我通常先把NR潮流跑一遍，把它收敛后的V和θ作为"真值/参考值"来生成模拟量测：在真值基础上叠加高斯噪声，模拟PMU和SCADA的实测数据。然后再拿WLS去做状态估计，最后把WLS的估计结果和NR的真值结果放一起算误差。整个过程环环相扣，非常流畅。
3. Matlab实现全过程：从数据生成到迭代收敛
3.1 仿真系统与量测数据准备
我建议先别直接上大系统，拿一个规模适中的标准节点系统练手就够了。IEEE 14节点系统是公认的"最小区间里能体现电力系统几乎所有典型特征"的模型：有平衡节点、有PV节点、有PQ节点，线路和变压器都有，无环网也有拓扑复杂性。我这次就用14节点系统作为示例。
数据准备分三步：
- 确定系统参数：节点导纳矩阵Y（由线路参数和变压器变比计算）、节点类型、基准功率。
- 跑一次NR潮流，得到各节点电压幅值V_true和相角θ_true。
- 基于V_true和θ_true，计算各注入功率和支路潮流的真值，然后加高斯噪声生成量测。
量测生成时，噪声的标准差设置成多少挺关键。太小了体现不出状态估计的滤波平滑作用，太大了结果会明显偏离真值。我常用的设置是：
- 节点注入有功/无功量测：标准差 = 0.01标幺功率
- 支路有功/无功潮流：标准差 = 0.01标幺功率
- PMU电压幅值量测：标准差 = 0.001标幺
- PMU相角量测：标准差 = 0.001弧度
这个设置比较接近工程实际里SCADA和PMU的精度差别。
3.2 量测配置与权重矩阵设计
量测配置决定了雅可比矩阵H的结构，也决定了问题能不能解。理论上讲，量测数量至少要不小于状态变量数量，而且量测要覆盖到系统的各个区域，否则增益矩阵奇异。
我在14节点系统上的典型配置是：
- 部分节点安装PMU，比如在节点4、7、9上装PMU，直接量测V和θ
- 其余节点和线路配置SCADA功率量测：所有PQ节点的注入有功无功、若干条关键支路的有功无功潮流
这样既模拟了PMU尚未全网覆盖的过渡阶段，也能对比"有PMU"和"无PMU"的估计效果差异。
权重矩阵W = R^{-1}的构建，在Matlab里就是个对角矩阵。每个对角线元素等于对应量测误差方差的倒数。写代码时可以直接用一个diag函数搞定：
% 假设量测误差标准差向量为sigma，量测向量为z R = diag(sigma.^2); W = inv(R); % 等效于 diag(1 ./ sigma.^2)注意不要真的对满秩大矩阵求逆，直接用倒数对角元素就行，避免引入不必要的数值开销。
3.3 WLS迭代核心代码与收敛判据
WLS迭代的核心就是反复计算雅可比矩阵H、增益矩阵G和增量Δx。Matlab伪代码如下：
x = x0; % 初始状态，通常用平启动，V=1.0，theta=0 maxIter = 20; tol = 1e-6; for iter = 1:maxIter [h, H] = measurementFunction(x); % 计算量测函数值和雅可比矩阵 G = H' * W * H; dx = G \ (H' * W * (z - h)); % 解线性方程 x = x + dx; if norm(dx, inf) < tol break; end end其中measurementFunction你要自己实现，它根据当前的状态估计x，计算每个量测的估计值h(x)，以及对应的雅可比矩阵行。功率量测的雅可比推导比较繁琐，建议对照潮流方程逐项求偏导，或者用数值差分校验一遍自己的解析表达式。
关于初始值，我用的是平启动：所有节点V=1.0，θ=0。对于状态估计来说，只要量测配置合理，WLS对初值不像NR潮流那么敏感，但初值太离谱也会影响收敛。
收敛判据我习惯用无穷范数norm(dx, inf) < 1e-6，比二范数更严格更直接。也可以配合检查目标函数J的下降情况，如果J不再明显下降但dx还很大，就要怀疑是否有坏数据或者H算错了。
3.4 一个隐藏的效率优化点
在小系统上不太明显，但如果你后面打算换到118节点甚至更大系统，有一个优化点值得现在做：雅可比矩阵H和增益矩阵G都是稀疏矩阵，Matlab里应该用稀疏矩阵格式存储和计算，而不是满阵。
% 建议用稀疏方式组装 H_sparse = sparse(H); G = H_sparse' * W * H_sparse;这个习惯从一开始就养成，后面扩展到大规模系统会省很多事。
4. 结果对比：误差、收敛性和计算效率
4.1 电压幅值和角度的估计效果对比
完成WLS迭代后，把估计得到的V_est和θ_est与之前NR潮流得到的V_true和θ_true做对比。这里我建议画两类图：一是节点电压幅值对比曲线，二是节点相角对比曲线，每条曲线同时画真值、估计值和（如果有的话）不含PMU时的估计值。
不过光看图不够，还要有定量的误差指标。通常我会计算最大绝对误差、平均绝对误差和均方根误差（RMSE）。以下是我在一次典型仿真中得到的对比结果（数值经过整理，体现趋势）：
| 对比项 | 无PMU时平均误差 | 含PMU时平均误差 | 最大误差（含PMU） |
|---|---|---|---|
| 电压幅值（标幺） | 3.2e-3 | 0.8e-3 | 1.6e-3 |
| 相角（弧度） | 4.5e-4 | 1.2e-4 | 2.8e-4 |
可以看到，PMU加入后电压幅值和相角的估计误差都显著下降。这里面的逻辑很直接：PMU直接量测V和θ，相当于在状态空间里"钉"住了几个关键节点的准确状态，WLS在迭代时不会被其他间接量测的误差带到偏离方向。
另外我还发现一个现象：靠近PMU安装节点的区域，估计误差改善最明显；离得远的节点也有改善，但幅度略小。这说明PMU的布置位置对全网估计精度分布有直接影响，这也是后面可以单独做优化研究的方向。
4.2 收敛行为对比
NR潮流在相同系统上的迭代次数通常在3到5次，WLS状态估计也差不多，一般4到6次收敛到1e-6精度。两者的收敛行为很像，这和它们都使用牛顿类迭代有关。
有一个细节值得注意：如果量测集不含PMU，WLS的收敛路径有时候会走一个"中间偏差点"，导致前面一两次迭代的Δx偏大；加入PMU之后，线性量测项的约束作用会让迭代路径更平稳，几乎不会出现中途大幅偏离的情况。我在代码里打印每次迭代的Δx范数时能明显看到这个差异。
计算效率上，14节点系统在Matlab里跑一次WLS迭代耗时在毫秒级，肉眼基本无感。但这个项目的重点不在绝对耗时，而在于理解量测配置和权重对结果的影响，所以效率倒不是主要限制。
4.3 一个必须做的极端测试
我非常推荐你做一个测试来验证代码逻辑是否正确：把量测噪声全部设为0，也就是"完美量测"情况。此时WLS状态估计的结果应该与NR潮流结果高度一致，误差在数值精度范围内（通常在1e-10到1e-12数量级）。如果这个测试下误差仍然很大，说明代码里有bug，最常见的是雅可比矩阵求错，或者量测函数表达式与真实量测不一致。
这个测试的价值在于它把问题"理想化"了，排除了噪声干扰后，整套算法的正确性可以直接和NR潮流互相印证。我在调试过程中就是靠这个测试定位过两个错误，一个是支路量测函数里忽略了变压器变比，另一个是相角的单位搞混了（度数/弧度）。
5. 常见问题与排查技巧
5.1 不收敛或迭代发散怎么办
这是最常见的坑。排查顺序按优先级来：
- 检查量测数量是否足够，是否存在某个节点既没有功率量测也没有PMU量测，导致增益矩阵奇异。
- 检查雅可比矩阵H的表达式，尤其是PV节点的无功边界处理。如果H某一行和潮流方程不匹配，迭代方向就是错的。
- 检查状态变量与量测函数的单位是否统一。比如相角用弧度，但量测生成时用了度数，就会出现反复震荡无法收敛。
- 尝试改用阻尼牛顿法，在Δx上乘一个小于1的步长因子，看能否稳定下来。
我遇到过一次比较隐蔽的情况：量测里同时包含节点注入功率和该节点支路潮流，二者在数学上是线性相关的（注入功率等于该节点所有支路潮流之和），导致H行之间存在近似线性相关性，增益矩阵条件数变得很大。解决办法是删掉部分冗余量测，或者给相关量测设置略低的权重。
5.2 增益矩阵奇异或条件数过大
增益矩阵G = H^T W H，当量测不足或量测冗余性极差时，G会奇异。判断方法很直接：在Matlab里用cond(G)查看条件数，如果大于1e12就要警惕。
解决办法有两个方向：一是增加量测，特别是增加PMU点位的覆盖；二是调整权重，避免某个量测权重过大导致"伪奇异"——比如把某个PMU的权重设置成其他量测的1万倍，数值上增益矩阵就会偏向那几个对角线元素，其他方向的约束几乎被忽略。这个问题用前面说的"PMU权重比约100倍"的经验值就能避开。
5.3 对比结果时相位基准不一致
这个坑特别容易忽略。NR潮流里平衡节点的相角通常固定为0，但WLS估计时，如果量测里没有直接给出平衡节点的相角约束，估计出来的相角可能整体偏移一个常数。所以做对比之前，必须把估计结果的相角基准对齐。
我的做法是：取同一个节点（一般是平衡节点或某个可信的PMU节点）作为参考，把所有估计相角减去该节点的估计相角，再和真值对比。换句话说，我在WLS迭代时也会固定平衡节点相角为0，或者在迭代完成后用参考节点做平移对齐。否则你画出来的相角对比曲线会整个错开，误差指标也会虚高。
5.4 PMU位置选取的初步策略
如果你要在这个项目基础上进一步扩展PMU布置优化，我建议先从"关键节点优先"开始：优先在带有较大发电机接入的节点、枢纽节点和网络拓扑中心区域装PMU。这个策略虽然不一定最优，但实操中效果通常不错，也容易在结果里看出来哪些节点对估计精度的贡献最大。
我个人的调试心得是：每改一次PMU位置，都重新跑一遍"完美量测"测试和"含噪量测"测试，看两组误差指标的变化。前者反映算法正确性，后者反映工程效果，两个都正常，整个改动才算靠谱。
6. 我在实操中的一点体会
做这个项目最大的收获，不是背下了WLS的迭代公式，而是把"量测-状态-算法"三者之间的关系想通了。以前看课本上讲状态估计，总感觉是一个抽象的最优化问题，真到了代码层面才发现，每一行雅可比矩阵、每一个权重设置、每一次收敛判断，背后都有实实在在的物理含义。
特别是把PMU和WLS放在一起之后，我对"数据质量如何影响算法上限"这件事有了更具体的感知：算法本身再精巧，量测不到位，估计结果就是有瓶颈；反过来，量测信息足、权重设计合理的时候，WLS能把精度推到接近理论极限。
如果你也想做这个方向，我的建议是：先跑通14节点系统的完整流程，把每一段代码和公式对应上，然后大胆改条件——加PMU、加坏数据、换权重、换网络拓扑，让结果发生变化，再去理解为什么。这个过程比单纯把代码跑出图有价值得多。后续如果你有兴趣，还可以在PMU优化配置、坏数据检测与辨识、抗差状态估计这几个方向上继续深挖，这套代码就是很扎实的起步基础。
