💥1 概述
电力系统状态估计正面临不同类型的 异常。这些可能包括由严重测量误差或 通信系统故障。负载或发电的突然变化可能是 根据实现的状态估计方法被视为异常。 此外,将电网视为一个信息物理系统,状态 估计容易受到错误数据注入攻击。现有的 异常分类的方法无法准确分类(区分 之间)上述三种类型的异常,尤其是当它来的时候 区分突然的负载变化和错误的数据注入攻击。 该文提出一种检测异常存在的新算法,对 异常类型和识别异常的来源,即测量 在数据错误的情况下包含严重错误,或与负载相关的总线 经历突然变化,或错误数据针对的状态变量 注入攻击。该算法结合了分析和机器学习 (ML) 方法。第一阶段利用分析方法来检测异常 通过组合 $\chi^2$-测试和异常检测指数来呈现。第二个 阶段利用 ML 对异常类型进行分类并识别其 原产地,特别是对突然负载变化的区分 和虚假数据注入攻击。所提出的基于ML的方法被训练成 独立于网络配置,无需重新训练 网络拓扑更改后的算法。通过实施获得的结果 所提算法在IEEE 14总线测试系统上验证了算法的准确性和 所提算法的有效性。
详细文章见第三部分。
📚2 运行结果
部分代码:
% Used to calculate H matrix, the state variables should be the forcasted function H = H_matrix(num,V,del) %num = 14; ybus = ybusppg(num); % Get YBus.. zdata = zdatas(num); % Get Measurement data.. bpq = bbusppg(num); % Get B data.. nbus = max(max(zdata(:,4)),max(zdata(:,5))); % Get number of buses.. type = zdata(:,2); % Type of measurement, Vi - 1, Pi - 2, Qi - 3, Pij - 4, Qij - 5, Iij - 6.. % z = zdata(:,3); % Measuement values.. fbus = zdata(:,4); % From bus.. tbus = zdata(:,5); % To bus.. % Ri = diag(sig); % Measurement Error.. % V = ones(nbus,1); % Initialize the bus voltages.. % del = zeros(nbus,1); % Initialize the bus angles.. % E = [del(2:end); V]; % State Vector.. G = real(ybus); B = imag(ybus); vi = find(type == 1); % Index of measurements.. ppi = find(type == 2); qi = find(type == 3); pf = find(type == 4); qf = find(type == 5); nvi = length(vi); % Number of Voltage measurements.. npi = length(ppi); % Number of Real Power Injection measurements.. nqi = length(qi); % Number of Reactive Power Injection measurements.. npf = length(pf); % Number of Real Power Flow measurements.. nqf = length(qf); % Number of Reactive Power Flow measurements.. % iter = 1; % tol = 5; % Jacobian.. % H11 - Derivative of V with respect to angles.. All Zeros H11 = zeros(nvi,nbus-1); % H12 - Derivative of V with respect to V.. H12 = zeros(nvi,nbus); for k = 1:nvi for n = 1:nbus if n == k H12(k,n) = 1; end end end % H21 - Derivative of Real Power Injections with Angles.. H21 = zeros(npi,nbus-1); for i = 1:npi m = fbus(ppi(i)); for k = 1:(nbus-1) if k+1 == m for n = 1:nbus H21(i,k) = H21(i,k) + V(m)* V(n)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n))); end H21(i,k) = H21(i,k) - V(m)^2*B(m,m); else H21(i,k) = V(m)* V(k+1)*(G(m,k+1)*sin(del(m)-del(k+1)) - B(m,k+1)*cos(del(m)-del(k+1))); end end end % H22 - Derivative of Real Power Injections with V.. H22 = zeros(npi,nbus); for i = 1:npi m = fbus(ppi(i)); for k = 1:(nbus) if k == m for n = 1:nbus H22(i,k) = H22(i,k) + V(n)*(G(m,n)*cos(del(m)-del(n)) + B(m,n)*sin(del(m)-del(n))); end H22(i,k) = H22(i,k) + V(m)*G(m,m); else H22(i,k) = V(m)*(G(m,k)*cos(del(m)-del(k)) + B(m,k)*sin(del(m)-del(k))); end end end % H31 - Derivative of Reactive Power Injections with Angles.. H31 = zeros(nqi,nbus-1); for i = 1:nqi m = fbus(qi(i)); for k = 1:(nbus-1) if k+1 == m for n = 1:nbus H31(i,k) = H31(i,k) + V(m)* V(n)*(G(m,n)*cos(del(m)-del(n)) + B(m,n)*sin(del(m)-del(n))); end H31(i,k) = H31(i,k) - V(m)^2*G(m,m); else H31(i,k) = V(m)* V(k+1)*(-G(m,k+1)*cos(del(m)-del(k+1)) - B(m,k+1)*sin(del(m)-del(k+1))); end end end % H32 - Derivative of Reactive Power Injections with V.. H32 = zeros(nqi,nbus); for i = 1:nqi m = fbus(qi(i)); for k = 1:(nbus) if k == m for n = 1:nbus H32(i,k) = H32(i,k) + V(n)*(G(m,n)*sin(del(m)-del(n)) - B(m,n)*cos(del(m)-del(n))); end H32(i,k) = H32(i,k) - V(m)*B(m,m); else H32(i,k) = V(m)*(G(m,k)*sin(del(m)-del(k)) - B(m,k)*cos(del(m)-del(k))); end end end % H41 - Derivative of Real Power Flows with Angles.. H41 = zeros(npf,nbus-1); for i = 1:npf m = fbus(pf(i)); n = tbus(pf(i)); for k = 1:(nbus-1) if k+1 == m H41(i,k) = V(m)* V(n)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n))); else if k+1 == n H41(i,k) = -V(m)* V(n)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n))); else H41(i,k) = 0; end end end end % H42 - Derivative of Real Power Flows with V.. H42 = zeros(npf,nbus); for i = 1:npf m = fbus(pf(i)); n = tbus(pf(i)); for k = 1:nbus if k == m H42(i,k) = -V(n)*(-G(m,n)*cos(del(m)-del(n)) - B(m,n)*sin(del(m)-del(n))) - 2*G(m,n)*V(m); else if k == n H42(i,k) = -V(m)*(-G(m,n)*cos(del(m)-del(n)) - B(m,n)*sin(del(m)-del(n))); else H42(i,k) = 0; end end end end % H51 - Derivative of Reactive Power Flows with Angles.. H51 = zeros(nqf,nbus-1); for i = 1:nqf m = fbus(qf(i)); n = tbus(qf(i)); for k = 1:(nbus-1) if k+1 == m H51(i,k) = -V(m)* V(n)*(-G(m,n)*cos(del(m)-del(n)) - B(m,n)*sin(del(m)-del(n))); else if k+1 == n H51(i,k) = V(m)* V(n)*(-G(m,n)*cos(del(m)-del(n)) - B(m,n)*sin(del(m)-del(n))); else H51(i,k) = 0; end end end end % H52 - Derivative of Reactive Power Flows with V.. H52 = zeros(nqf,nbus); for i = 1:nqf m = fbus(qf(i)); n = tbus(qf(i)); for k = 1:nbus if k == m H52(i,k) = -V(n)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n))) - 2*V(m)*(-B(m,n)+ bpq(m,n)); else if k == n H52(i,k) = -V(m)*(-G(m,n)*sin(del(m)-del(n)) + B(m,n)*cos(del(m)-del(n))); else H52(i,k) = 0; end end end end % Measurement Jacobian, H.. H = [H11 H12; H21 H22; H31 H32; H41 H42; H51 H52];
🌈3 Matlab代码及详细文章
🎉4 参考文献
部分理论来源于网络,如有侵权请联系删除。