企业官网建设流程全解析

电网调度控制中心里，PMU（同步相量测量装置）以每秒几十帧甚至上百帧的速度把全网关键节点的功角、电压相量数据推上来，数据量是完全足够的，但直接用这些带噪声的实时量测去判断系统状态，结果会非常不可靠。这时候就需要动态状态估计出马：把系统模型和量测数据结合起来，先预测再校正，输出一条平滑、可信、可预测的状态轨迹。这也是EKF和UKF在电力系统里最重要的实战舞台。

这篇文章就是来拆解基于扩展卡尔曼滤波（EKF）和无迹卡尔曼滤波（UKF）的电力系统动态状态估计的完整实现过程，所有讨论都基于我在Matlab里跑通算例的真实心得。内容包括状态方程和量测方程的建模、两种滤波算法的递推细节、关键参数怎么调、同一算例下EKF和UKF的表现差异，以及我调试过程中踩过的坑。适合正在做电力系统状态估计研究的学生、刚接触卡尔曼滤波想快速落地的工程师，以及任何想搞懂"EKF和UKF到底哪个更适合我的场景"的读者。

1. 为什么说电力系统动态状态估计的关键在非线性滤波

1.1 从静态断面到时变过程：动态估计到底多了一个什么"动态"

先捋一捋静态状态估计和动态状态估计的区别。传统SCADA系统里用的加权最小二乘（WLS）状态估计，本质上是求解一个优化问题：给定某一个时刻的冗余量测，找到一组状态变量（各节点电压幅值和相角）使得量测残差的加权平方和最小。它输出的是一幅"静态断面"，没有利用量测在时间维度的变化规律。换句话说，WLS默认系统是静止的，上一时刻的信息、系统的机电动态模型，全都用不上。

动态状态估计则完全不同。它的基本思路是站在卡尔曼滤波框架里：我不仅知道当前量测值，还知道系统状态的时间演化规律（通常用发电机转子运动方程来描述），那么就能先通过状态方程做一步预测，再用当前量测对预测结果做修正。这样输出的是状态变量在每个采样时刻的条件期望，同时给出估计协方差，可以用来刻画估计的不确定性。这个"预测-校正"结构，在暂态过程中尤其有价值：当系统正在经历扰动、从故障前状态向故障后状态过渡时，SCADA的采样率根本跟不上，而动态状态估计可以凭借模型外推能力，在两拍量测之间给出一个合理的状态过渡估计。

这也是为什么近些年来PMU普及之后，动态状态估计的研究热度明显高了一截：PMU给了高采样率同步量测，等于把"动态估计需要频繁量测校正"这个前提条件补上了。

1.2 经典卡尔曼为什么直接套不上电力系统

如果系统是线性的、噪声是高斯白噪声，标准卡尔曼滤波（KF）是最优的线性无偏估计器，公式简洁、计算量小、理论性质完美。但电力系统的核心量测方程偏偏是非线性的。

你仔细看：PMU量测的电气功率是怎么来的？在经典发电机模型中，发电机电磁功率可以写成：

[ P_e = \frac{E' V_s}{X'} \sin(\delta) ]

功角 (\delta) 是待估计的状态变量，而量测值是功率 (P_e)，量测函数 (h(x) = \frac{E' V_s}{X'} \sin(\delta)) 是正弦函数。再看状态方程，发电机转子运动方程本身也带有 (\sin(\delta)) 项，状态转移也是非线性的。

那么问题来了：标准KF的推导严格依赖"高斯随机变量经过线性变换后仍为高斯分布"这个性质。状态变量是高斯分布，经过状态方程的非线性函数传播之后，就不一定是高斯了，而且你连它的均值、协方差都无法通过简单的矩阵乘法得到。这就是标准KF在非线性系统中失效的根本原因。

EKF和UKF就是在两条不同路径上解决这个问题的：EKF选择把非线性函数在当前估计点做一阶泰勒展开，用线性化的雅可比矩阵代替原函数，强行把问题拉回线性框架里；UKF则选择不线性化，而是采样一组Sigma点，让这些点经过真实非线性函数传播，再用加权统计量近似输出分布的均值和协方差。一个是"把模型变弯为直"，一个是"以点代面做统计近似"，这就是两种算法最本质的思路差异。

1.3 EKF和UKF各自的实际门槛在哪里

搞清楚思路差异后，实际应用的门槛差异就浮出来了。EKF的实现门槛在雅可比矩阵的推导：状态转移矩阵 (F_k) 和量测矩阵 (H_k) 都要算导数，手推容易出错，尤其状态变量一多、方程一复杂，解析雅可比很容易漏项。如果改用数值差分算雅可比，又会引入步长选择和数值误差的问题。

UKF的门槛不在推导上，而在Sigma点参数的选择与数值稳定性上：权重怎么算、Sigma点怎么生成、协方差矩阵怎么保证半正定。一旦这些细节处理不好，滤波发散的速度比EKF还快。

我在实际对比测试中还有一个体会：EKF在系统运行点附近偏离不大的情况下，精度完全够用，而且计算量明显小；但如果系统受到较大扰动、功角变化剧烈，EKF在每一个时刻的局部线性化误差会持续累积，这时候UKF的优势就很明显了。所以选型时，先回答一个问题：你要估计的是稳态小幅波动场景，还是暂态大幅摆动场景？这两个场景的最佳答案往往是不同的。

2. 状态方程与量测方程：两种滤波算法共用的基础

2.1 发电机经典二阶模型与离散化处理

动态状态估计的第一步，一定是把受控系统的连续时间模型写清楚。在电力系统动态状态估计研究里，最常用的基础模型是经典二阶发电机模型（摆动方程），它对单机无穷大系统和多机系统的理论研究都非常常见。连续时间形式如下：

[ \begin{aligned} \frac{d\delta}{dt} &= \omega_0 (\omega - 1) \ \frac{d\omega}{dt} &= \frac{1}{2H} \left( P_m - P_e - D(\omega - 1) \right) \end{aligned} ]

其中 (\delta) 是发电机功角，(\omega) 是发电机转速标幺值，(\omega_0) 是同步转速，(H) 是惯性时间常数，(D) 是阻尼系数，(P_m) 是机械功率，(P_e) 是电磁功率。状态变量取 (x = [\delta, \omega]^T)。

这个连续模型在计算机里无法直接递推，需要离散化。最常用的方法是欧拉法或四阶龙格-库塔法（RK4）。我一般建议用RK4做状态外推，因为EKF的精读很大程度上取决于状态预测的准确性，如果预测本身就有较大的离散化误差，后面的量测修正很难完全补偿回来。RK4虽然每步要多算三次函数评估，但在Matlab的数据规模下，这点计算开销几乎可以忽略。欧拉法只有在步长非常小的情况下才够用，电力系统动态仿真的步长如果取5毫秒以内，欧拉法还能接受；一旦步长到10毫秒以上，误差就开始明显影响了。

注意状态方程中的非线性来源：(P_e = \frac{E' V_s}{X'} \sin(\delta)) 里有 (\sin(\delta))。这意味着即使量测是线性的（直接量测功角），状态转移本身也不是线性的，标准KF依然无法直接套用。这一点经常被初学者忽略，以为"只要量测方程是线性的就能用KF"，实际上状态方程的非线性同样需要处理。

2.2 量测方程：为什么PMU量测直接带来非线性

量测方程的选择直接影响EKF要求导、UKF要凑Sigma点传播的复杂度。常见做法有两种，我分开说。

第一种，直接把PMU输出的功角和转速当量测，量测方程是：

[ z_k = \begin{bmatrix} \delta_k \ \omega_k \end{bmatrix} + v_k ]

这种写法量测矩阵 (H) 是常数矩阵，非常方便，但问题在于：你实际上已经假设量测没有任何非线性变换，这个假设在仿真中可行，在真实工程中却很少见。真实PMU输出的电气量通常是电压相量、电流相量、有功功率等，需要经过相量计算才能还原出功角和转速，稍微处理不当就会引入额外的转换误差。

第二种更贴近物理实际的做法，把发电机电磁功率作为量测：

[ z_k = P_e^{(m)} = \frac{E' V_s}{X'} \sin(\delta_k) + v_k ]

这时量测方程就是非线性的，EKF里需要计算 (H_k = \partial h/\partial x = [\frac{E' V_s}{X'} \cos(\delta_k), 0])。我在仿真里更倾向于用第二种做法，因为它的非线性特征明显，最能体现EKF和UKF在非线性处理能力上的差异，也更容易暴露出算法的数值问题。如果一上来就用常数 (H) 矩阵，你会发现EKF和UKF的结果几乎一样，算法的差异反而不容易看出来。

顺便提一句噪声假设：量测噪声 (v_k) 通常假设为均值为零、协方差为 (R) 的高斯白噪声。这在仿真里用randn直接生成即可，但在真实系统里，PMU量测噪声并不总是严格白色高斯的，可能存在时序相关。这时候EKF和UKF都还勉强能用，但估计结果会偏乐观，协方差会低估。如果项目对不确定性估计的准确性要求高，可以考虑在噪声模型里加入相关性描述，不过这已经超出这篇文章的讨论范围了。

2.3 噪声参数与初始协方差的取值策略

动态状态估计里，过程噪声协方差 (Q) 和量测噪声协方差 (R) 的设置非常关键，而且很多教科书对这个问题的讲法过于理想化。实际调参经验如下。

(R) 的取值相对容易：如果量测是仿真生成的，直接在真实电压电流上叠加标准差为 (\sigma_v) 的高斯噪声，那么 (R = \sigma_v^2) 就是最优选择。PMU的幅值测量误差典型范围在0.1%~0.2%左右，相角测量误差在0.1度左右，据此可以反推功率量测的噪声方差数量级。

(Q) 就麻烦多了。过程噪声代表的是模型不准的部分：发电机模型中 (P_m) 的实际波动、参数误差、离散化误差等，这些很难精确度量。我的经验是，(Q) 的值宁可设得偏大一点，也不要偏小。(Q) 偏小会让滤波器过信模型，量测来了也不愿意修正，最终导致估计滞后于真实状态；(Q) 偏大则会让滤波器更信任量测，虽然噪声多一些，但至少能跟上真值的变化。实际操作时，我一般是先给一个较大初值，观察滤波轨迹是否跟得上真值，然后逐步减小，直到RMSE不再明显下降为止。

初始协方差 (P_0) 同样重要。它代表我们对初始状态估计不确定性的认知。如果初值设得与真值偏差很大，但 (P_0) 又设得很小，滤波器在起初几步会"固执己见"，估计值很难快速收敛到真值附近。我通常把 (P_0) 的对角元素设成初始不确定性的平方，比如功角初始标准差估计5度，那 (P_0(1,1)) 就取25（角度制下）对应值，或者直接取状态典型幅值平方的若干倍。

3. EKF的Matlab实现：雅可比矩阵是核心也是坑

3.1 雅可比矩阵用解析推导还是数值差分

EKF区别于KF的地方，就是在预测和更新两个环节里各插入了一次线性化。预测时需要对状态方程求状态转移矩阵 (F_k)，更新时需要对量测方程求量测矩阵 (H_k)。这一步的雅可比矩阵推导，是EKF编程里最容易出错的环节。

以量测方程 (h(x) = \frac{E' V_s}{X'} \sin(\delta)) 为例，解析雅可比很简单：

[ H_k = \left[ \frac{\partial h}{\partial \delta}, \frac{\partial h}{\partial \omega} \right] = \left[ \frac{E' V_s}{X'} \cos(\delta_k), 0 \right] ]

但如果你的系统状态变量是10个节点的功角和转速，量测又包含节点注入功率、支路潮流，那么每个量测对每个状态的偏导数展开来会有几十项，手推极易出错，而且错了还很难检查出来。

这时候我建议采用数值差分的方式验证解析雅可比：先用中心差分公式（如 (\frac{f(x+\epsilon) - f(x-\epsilon)}{2\epsilon})）计算一个数值雅可比，再和你的解析结果逐元素对比。步长 (\epsilon) 一般取 (\sqrt{eps} \approx 1.49 \times 10^{-8}) 量级与状态幅值之积。如果两者差异在 (10^{-6}) 量级，说明解析推导正确。这个方法我几乎每次写新系统都会用一遍，哪怕已经很有经验。

在最终交付的代码里，我自己的习惯是：能解析推导的场合就用解析式，因为数值差分在每次递推都要额外调用多次非线性函数，在大规模系统里会拖慢计算速度；但解析式旁边一定要留一个数值差分对比函数用来做单元测试。

3.2 递推主循环的完整流程

EKF的递推主循环代码结构并不复杂，核心就五步。下面是参考实现，注意函数签名和具体模型强相关，这里展示的是结构骨架。

function [x_est_arr, P_arr] = ekf_estimation(z_arr, dt, params) % z_arr: 量测序列，每一列是一个量测向量 % dt: 采样时间 % params: 系统参数结构体 n = length(params.x0); x_est = params.x0; % 初始状态估计 P = params.P0; % 初始协方差 Q = params.Q; R = params.R; x_est_arr = zeros(n, length(z_arr)); P_arr = zeros(n, n, length(z_arr)); for k = 1:length(z_arr) % 1. 状态外推（用RK4离散化） x_pred = rk4_state_equation(x_est, dt, params); % 2. 状态转移矩阵线性化（在x_est处求雅可比） Fk = compute_F(x_est, dt, params); % 3. 协方差外推 P_pred = Fk * P * Fk' + Q; % 4. 量测线性化（在x_pred处求雅可比） Hk = compute_H(x_pred, params); % 5. 滤波更新 Kk = P_pred * Hk' / (Hk * P_pred * Hk' + R); innovation = z_arr(:, k) - h_function(x_pred, params); x_est = x_pred + Kk * innovation; P = (eye(n) - Kk * Hk) * P_pred; P = (P + P') / 2; % 强制对称化 x_est_arr(:, k) = x_est; P_arr(:, :, k) = P; end end

这里有三个值得展开的细节。第一个是compute_F的线性化时机：我是在 (x_k) 处而不是 (x_{k+1}) 处求雅可比，严格来说这是标准EKF的写法，属于一阶精度；如果对精度要求更高，可以考虑迭代EKF或在中间点线性化。

第二个是 (\sin) 项导致的步长敏感问题：状态方程里带 (\sin(\delta))，如果 (dt) 较大，RK4离散化后的状态转移与连续模型的偏差就会变大，这些偏差没有体现在 (Q) 里的话，滤波器会误以为预测很准，从而降低对量测的信赖。所以在设置 (Q) 时，一定要把离散化误差算进过程噪声里。

第三个是协方差对称化。理论上 (P) 经过递推公式仍然保持对称，但由于Matlab的数值舍入误差，长时间运行后 (P) 可能轻微失去对称性。加上P = (P + P') / 2这一行，成本极低，却可以避免很多后续奇异问题。

3.3 保证滤波稳定的细节处理

EKF在实际运行中比教科书上更容易出数值问题，一个常见表现就是协方差矩阵逐渐失去正定性，甚至变成负定，然后卡尔曼增益算出一个极不合理的值，估计轨迹瞬间"飞掉"。

要防止这种情况，我的经验有三条。第一，做奇异值检查：每隔几步对 (P) 做一次特征值分解或用eig(P)检查最小特征值，如果接近零或为负，就在对角上加一个极小的正则项 (\epsilon I)，一般取 (10^{-10}) 量级，原理类似岭回归。第二，量测更新时不要直接用inv(H * P_pred * H' + R)求矩阵逆，而应该用Matlab的\运算符解线性方程组：Kk = P_pred * Hk' / (Hk * P_pred * Hk' + R)。中间矩阵可能接近奇异，显式求逆会放大数值误差。第三，如果量测冗余度很高，多个量测对状态都有强约束，那么更新步的增益矩阵 (K) 可能过度放大某项量测残差，导致振荡。这时候可以引入渐消因子降级处理，不过这是自适应滤波的话题了，这里点到为止。

4. UKF的Matlab实现：用Sigma点绕过求导

4.1 无迹变换的权重设置原理

UKF的核心是无迹变换（Unscented Transform）。它的想法很直接：一个高斯分布虽然通过非线性函数之后不再高斯，但我们不去计算它的解析形式，而是从原始分布里精心挑选一组样本点（Sigma点），让这些点经过非线性函数，再用传播后的点重构高斯分布的均值和协方差。

对于 (n) 维状态变量，需要生成 (2n+1) 个Sigma点。假设当前状态均值是 (\bar{x})，协方差是 (P)，那么每个Sigma点按下式生成：

[ \begin{aligned} \chi^{(0)} &= \bar{x} \ \chi^{(i)} &= \bar{x} + \left(\sqrt{(n+\lambda)P}\right)_i, \quad i=1,...,n \ \chi^{(i+n)} &= \bar{x} - \left(\sqrt{(n+\lambda)P}\right)_i, \quad i=1,...,n \end{aligned} ]

其中 (\lambda = \alpha^2 (n + \kappa) - n) 是一个尺度参数。这里的 ((\sqrt{(n+\lambda)P})_i) 表示协方差矩阵开平方后取第 (i) 列。在Matlab里，可以通过Cholesky分解得到协方差矩阵的下三角根矩阵，然后取列向量。要注意：chol(P)默认返回上三角，所以要么转置，要么直接用(chol(P))'取列，这两种写法容易搞混，是新手经常踩的坑。

对应权重的计算公式为：

[ \begin{aligned} W^{(0)}_m &= \frac{\lambda}{n+\lambda} \ W^{(0)}_c &= \frac{\lambda}{n+\lambda} + (1 - \alpha^2 + \beta) \ W^{(i)}_m &= W^{(i)}_c = \frac{1}{2(n+\lambda)}, \quad i=1,...,2n \end{aligned} ]

注意区分均值权重 (W_m) 和协方差权重 (W_c)，两者只有第一个点的取值不同，后面 (2n) 个点完全一样。这个细节忘了的话，算出来的协方差会明显偏大或偏小，滤波性能直接受影响。

4.2 时间更新与量测更新的代码骨架

UKF的递推循环和EKF很相似，只是把"线性化"换成了"传播Sigma点"。代码骨架如下：

function [x_est_arr, P_arr] = ukf_estimation(z_arr, dt, params) n = length(params.x0); x_est = params.x0; P = params.P0; Q = params.Q; R = params.R; alpha = 1e-3; beta = 2; kappa = 0; lambda = alpha^2 * (n + kappa) - n; [Wm, Wc] = ut_weights(n, alpha, beta, kappa); x_est_arr = zeros(n, length(z_arr)); P_arr = zeros(n, n, length(z_arr)); for k = 1:length(z_arr) % 1. 生成Sigma点 Xsig = generate_sigma_points(x_est, P, lambda); % 2. 时间更新：Sigma点经过状态方程传播 Xsig_pred = zeros(n, 2*n+1); for i = 1:size(Xsig, 2) Xsig_pred(:, i) = rk4_state_equation(Xsig(:, i), dt, params); end x_pred = sum(Wm .* Xsig_pred, 2); P_pred = zeros(n, n); for i = 1:size(Xsig, 2) diff = Xsig_pred(:, i) - x_pred; P_pred = P_pred + Wc(i) * (diff * diff'); end P_pred = P_pred + Q; % 3. 量测更新：Sigma点经过量测方程传播 Zsig = zeros(size(z_arr, 1), 2*n+1); for i = 1:size(Xsig_pred, 2) Zsig(:, i) = h_function(Xsig_pred(:, i), params); end z_pred = sum(Wm .* Zsig, 2); Pzz = zeros(size(Zsig, 1), size(Zsig, 1)); for i = 1:size(Zsig, 2) diff_z = Zsig(:, i) - z_pred; Pzz = Pzz + Wc(i) * (diff_z * diff_z'); end Pzz = Pzz + R; Pxz = zeros(n, size(Zsig, 1)); for i = 1:size(Zsig, 2) diff_x = Xsig_pred(:, i) - x_pred; diff_z = Zsig(:, i) - z_pred; Pxz = Pxz + Wc(i) * (diff_x * diff_z'); end % 4. 更新 Kk = Pxz / Pzz; innovation = z_arr(:, k) - z_pred; x_est = x_pred + Kk * innovation; P = P_pred - Kk * Pzz * Kk'; P = (P + P') / 2; x_est_arr(:, k) = x_est; P_arr(:, :, k) = P; end end

这段代码里最有意思的点在于：UKF完全不需要雅可比矩阵，所以不存在手推导数出错的问题。你只需要保证两个非线性函数rk4_state_equation和h_function的输入输出维度正确，Sigma点传播正确，剩下的都是矩阵运算。这也是我向新手优先推荐UKF的原因：至少你不会踩"雅可比推导错误"这样一个隐蔽的大坑。

4.3 alpha、beta、kappa到底怎么调

UKF的调参问题在几乎所有教科书里都被一笔带过，但它实际影响很大。三个参数的作用和典型取值如下表所示：

我在实践中基本固定 (\beta = 2)，因为假设状态分布接近高斯，这是最优选择。(\kappa) 取0；当 (n) 较大时，(3-n) 的取值会导致权重为负，容易破坏协方差正定性，不建议在电力系统这种维数不高的场景冒险。真正需要细调的是 (\alpha)：(\alpha) 太小时，Sigma点过度集中在均值附近，非线性传播的结构信息丢失，UKF退化得几乎像EKF；(\alpha) 太大时，Sigma点分布过宽，协方差估计偏大，滤波器增益异常。我通常从 (1 \times 10^{-3}) 开始，逐步调到 (0.1) 左右，观察RMSE曲线来定最优值。

另外一个与参数同等重要的细节是：每次生成Sigma点时都需要对 (P) 做Cholesky分解，如果 (P) 不是严格正定的，chol会直接报错。所以在上一步更新之后加上对称化和微小正则化，对UKF来说不是可选项，而是必备操作。

5. 同一算例下EKF与UKF的实测对比结果

5.1 测试条件与评价指标

要比较两种算法，必须放在同一组数据和同样的起哄条件下跑。我在单机无穷大系统模型上做了测试：经典二阶模型，状态变量为功角和转速偏差，量测为发电机电磁功率叠加高斯白噪声。系统参数取 (H = 5) 秒，(D = 2)，(X' = 0.3) 标幺，(E' = 1.0) 标幺，(V_s = 0.995) 标幺。仿真时长5秒，采样步长0.01秒，共500个采样点。真值轨迹通过在系统模型上施加一个机械功率阶跃（扰动场景）和一个小幅随机波动（稳态场景）分别生成，量测噪声标准差设为真值功率幅值的1%。两种算法的 (Q)、(R)、(P_0) 初始值完全一致。

评价指标用均方根误差（RMSE）和平均计算耗时。RMSE对功角和转速分开统计，用标幺值或角度制都行，关键是两种算法用同一单位。

5.2 精度、耗时与收敛性的对比表

我跑了多组随机量测噪声重复实验，取平均后的结果大致如下：

三个结论从数据里很清晰。第一，在稳态场景下，EKF和UKF的精度差异其实没有想象中那么大，RMSE差距在10%以内，工程上可以认为两者都够用。第二，在暂态阶跃扰动场景下，UKF的精度优势明显放大，功角RMSE比EKF低了大约34%。原因是扰动导致的功角摆开，让EKF在每次线性化点的偏差都偏大，误差逐步累积；而UKF用Sigma点传播真实非线性函数，对大幅摆动的跟踪能力更强。第三，UKF每步耗时大约是EKF的1.5到1.8倍，主要开销在 (2n+1) 个Sigma点逐一做状态外推和量测函数计算上。

还有一个关键观察：UKF在初始误差较大的情况下收敛速度更快。我把初始功角误设到真值差10度，UKF在约0.2秒内拉回到真值附近，EKF则需要接近0.5秒。原因还是那个：UKF对非线性函数的统计近似更准确，初期的增益计算更合理。

5.3 从结果反推适用场景

有了这组数据，选型思路就很明确了。如果你的项目是电力系统稳态运行状态下的在线监测，量测质量尚可、运行点偏移不大，EKF完全胜任，而且代码量和计算开销都更小，边际成本最低。如果你的目标是暂态稳定评估、扰动后状态快速跟踪，或者系统模型本身有较强的非线性，UKF带来的精度收益是值得额外计算开销的。另外，如果项目后续要把状态估计与机电暂态仿真闭环结合，模型本身会频繁在非线性区域运行，UKF的鲁棒性优势会被进一步放大。

当然，这组数据是在单机模型下得到的。多机系统的量测和状态变量维度更高，UKF的Sigma点数量随维度线性增长，计算量增幅还能接受，但EKF雅可比矩阵的推导复杂度会快速上升，那时候我会更倾向于UKF，因为它不需要为新增的每台发电机重新手推导数。

6. 我在编写调试过程中踩过的坑与最终建议

6.1 滤波发散：最常见也最难排查的问题

我在最初跑EKF时碰到过非常典型的发散现象：前二百步状态估计看起来一切正常，RMSE小幅波动且整体收敛，但某个时刻起估计轨迹突然跳变，随后一直偏离真值并振荡扩大。当时第一反应是代码写错了，但反复检查公式和矩阵维度都没问题。后来逐项排查，发现问题出在过程噪声 (Q) 设置得太小：真实系统里机械功率 (P_m) 有一个我没有建模的缓慢时变，这个未建模动态让模型预测持续偏小，量测残差却因为 (Q) 小而得不到足够的修正权重，于是偏差慢慢积累，直到某个量测残差被放大器增益一次性放大才爆发出来。

这个案例的教训是：滤波发散不一定意味着代码逻辑错误，更常见的原因是"过程噪声协方差与未建模动态不匹配"。处理办法也很直接：给 (Q) 增加一个对角项来吸收未建模动态。在 (P_m) 缓慢变化场景里，可以在状态方程里把 (P_m) 扩展成状态变量并赋予一个小的随机游走过程噪声，效果立竿见影。这也是我建议用"模型误差分析"来设置 (Q) 而不是拍脑袋的原因。

UKF发散的原因则通常更集中在协方差数值问题上。有一次我让UKF跑长时间的连续仿真，中途协方差矩阵的某个特征值变成负的，chol直接报错，导致整个仿真中断。查了半天发现是在计算权重时，第一个Sigma点的协方差权重因为 (\beta = 2) 修正后变得很小，加上浮点舍入，累积若干步之后矩阵失去了正定性。解决办法是在更新步之后强制做P = (P + P') / 2，并且每隔50步对eig(P)做一次检查，最小特征值低于 (10^{-12}) 就直接加上正则项。

6.2 代码框架建议与后续扩展方向

最后给一个代码组织上的建议：无论EKF还是UKF，建议把"系统模型"和"滤波器本体"完全分离。写两个独立的函数文件：system_model.m里定义状态方程、量测方程、雅可比矩阵；ekf_filter.m和ukf_filter.m里只做纯滤波递推，不包含任何电力系统具体参数。这样换一个IEEE节点系统，只需要改系统模型参数和函数，滤波器代码一行不用动。我在多期项目里靠这个分工节省了大量返工时间。

后续可以扩展的方向也顺手列一下。一是做EKF/UKF的自适应版本，比如用新息序列实时调节 (Q) 和 (R)，可以应对噪声统计特性变化的情况。二是把扩展卡尔曼平滑器（EKS）和相应的无迹平滑器接在滤波后面做离线校正，能够进一步提升状态估计精度，尤其在PMU数据后处理的场景下性价比很高。三是把UKF的Sigma点思想换成正交滤波、容积卡尔曼滤波（CKF）或粒子滤波，对比它们在更强非线性条件下的表现。这些方向我在自己的项目里都试过，每一步都有不少值得展开的坑和心得，回头有机会再单独写一篇。

跑完这一个算例、踩过一轮发散和数值稳定性的坑之后，我的切身体会是：EKF胜在结构简单、计算高效、代码好调试，适合线性化偏差可控的常规场景；UKF胜在对非线性的描述更准、在暂态场景下更稳，代价是计算量和参数调优的复杂度更高。工具本身没有绝对的优劣，先想清楚你要解决的问题处在哪种非线性强度下，再选型，会比先选算法再适配问题高效得多。