CN103956756B - 一种电力系统低频振荡模态辨识方法 - Google Patents

Abstract

本发明涉及电力系统低频振荡分析技术领域，特别是一种基于改进经验模态分解（EMD）和数位演算法的电力系统低频振荡模态辨识方法。针对EMD分解过程中出现的边界失真现象，提出一种新的端点优化对称延拓法进行改进，利用改进的EMD对电力系统低频振荡信号进行分解，得到固有模态函数（IMF），利用柯尔摩科洛夫-斯米洛夫检验（K-S）法剔除IMF分量中与低频振荡信号相似概率低的伪分量，利用数位演算法对有效IMF分量进行振荡模态瞬时参数的提取，用于电力系统低频振荡模态辨识。该方法能够有效的改善EMD分解的端点效应，准确的剔除伪分量，精确的提取低频振荡信号模态瞬时参数，适用于电力系统等相关部门，用于电力系统低频振荡模态辨识。

Description

一种电力系统低频振荡模态辨识方法

技术领域

本发明涉及电力系统低频振荡分析技术领域，特别是一种基于改进EMD和数位演算法的电力系统低频振荡模态辨识方法。

背景技术

随着电力技术的不断进步，大电网互联已逐步实现，加之快速高放大倍数励磁装置的广泛使用，由此带来的低频振荡问题也越来越成为影响大规模电力系统稳定运行的重要因素，其不仅对电网运行安全上产生影响，同时也成为制约互联电网输送能力的瓶颈之一。如何通过广域测量数据分析电力系统低频振荡的局部动态行为，提取各个振荡模态的瞬时参数，进行广域阻尼控制是近年来的研究热点。

目前针对低频振荡问题所采用的分析方法主要有傅立叶变换、小波变换、Kalman滤波法、矩阵束辨识法、Prony算法和希尔伯特-黄(Hilbert-Huangtransform，HHT)算法等。傅立叶变换无法分析出阻尼特性和局部特性，所以不适于非线性、非平稳信号。小波变换存在频率交叠和自适应基选取问题，只适合于瞬态和非平稳信号。Kalman滤波法能消除噪声的影响，对不同输入信号的适应性较好，但计算精度和收敛速度受初始参数设置的影响很大。矩阵束辨识法能够准确估计系统的振荡模态，并具有较强的抗噪声能力，但若信号存在时变特性，该算法的计算误差较大，无法揭示振荡的动态特性。Prony方法虽然可以提取出振荡信号模式和阻尼等信息，但存在受噪声影响大、计算速度慢和定阶问题不确定等问题。HHT算法是近年发展起来的一种新型的适于非平稳、非线性信号的分析方法，传统的HHT算法受端点效应的影响，虽然能得到振荡模态的瞬时频率、瞬时幅值和衰减因子，难以达到较高的计算精度。但通过对该算法存在的端点效应问题进行改进，利用柯尔摩科洛夫-斯米洛夫检验(K-S)法剔除IMF分量中与低频振荡信号相似概率低的伪分量，可以精确提取出电力系统低频振荡的振荡模式和阻尼特性，有效的实现电力系统低频振荡信号的模态辨识。

发明内容

本发明的目的在于提供一种基于改进EMD和数位演算法的电力系统低频振荡模态辨识方法，该方法能够有效的改善EMD分解的端点效应，准确的剔除伪分量，精确的提取低频振荡信号模态瞬时参数，提高电力系统低频振荡模态辨识能力。

为实现上述目的，本发明的技术方案是：采用端点优化对称延拓法对EMD端点效应进行改进，利用改进之后的EMD对电力系统低频振荡信号进行分解，并利用K-S法对经EMD分解得到的IMF分量进行伪分量的剔除，利用数位演算法对剔除伪IMF分量后的有效IMF分量进行低频振荡模态瞬时参数的提取，得到电力系统低频振荡信号x(t)的瞬时频率、瞬时幅值和衰减因子等信息，其具体步骤如下：

步骤1：对滤波后的电力系统低频振荡信号x(t)进行离散化处理，采样步长取为Δt，得到离散化处理信号：

x_i＝[x₁,x₂,...,x_n]，i＝1,2,...,n

式中，x_i为离散化之后的信号值的集合，x₁,x₂,...,x_n为离散化之后t₁,t₂,...,t_n时刻的信号值；

步骤2：信号x_i两端点值用β、γ替代，构造新的数据序列x'_i：

x'_i＝[x'₀,x'_k,x'_n-1]

式中，x'₀＝β，x'_k＝x_k，x'_n-1＝γ，k＝1,2,...,n-2，i＝0,1,...,n-1；

步骤3：以端点x'₀、x'_n-1为中心对x'_i分别向两端进行对称延拓一个周期得到新的数据序列：

h＝[h_l,h_i,h_r]＝[2x'₀-x'_i,x'_i,2x'_n-1-x'_n-1-i]

式中，i＝0,1,...,n-1；

步骤4：利用三次样条插值法对信号数据序列h进行包络，得到包络线s(x)，优化x'_i的两个端点值β和γ，构建h_i与s_i的偏差评价函数为：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\mathrm{\μ}\∫s^2\left(x\right)dx$

式中，μ为光滑系数，s_i为s(x)对应于h_i的那部分包络线；

步骤5：对偏差评价函数进行离散化，得：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\frac{\mathrm{\μ}}{{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right(h_i\left)\right)}^2$

式中，Δx为采样间隔，▽²为二阶微分算子；

步骤6：将x'₀＝β和x'_n-1＝γ，上式可以变换

$\mathrm{\ϵ}=\frac{1}{2}\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-{x^{\′}}_i)}^2+\frac{\mathrm{\μ}}{2{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right({x^{\′}}_i\left)\right)}^2+{(\mathrm{\β}-x_0)}^2+{(\mathrm{\γ}-x_{n-1})}^2$

当 $\left\{\begin{array}{c}\∂\mathrm{\ϵ}/\∂\mathrm{\β}=0\\ \∂\mathrm{\ϵ}/\∂\mathrm{\γ}=0\end{array}\right.$ 时，ε为最小值，得出β和γ的值；

步骤7：根据数据序列h重新对电力系统低频振荡信号进行对称延拓，分解得到固有模态函数IMF分量c(t)和剩余分量r(t)；

步骤8：如果r(t)为单调函数，则转到步骤9，否则转到步骤2；

步骤9：经验模态分解EMD分解结束，得到固有模态函数IMF分量c(t)；

步骤10：利用柯尔摩科洛夫-斯米洛夫检验K-S法将c(t)与电力系统低频振荡信号x(t)进行对比分析，剔除与x(t)无关的虚假固有模态函数IMF分量，得到余下的有效固有模态函数IMF分量c_i(t)；

步骤11：利用数位演算法对c_i(t)进行各个模态瞬时参数的提取，得到电力系统低频振荡信号x(t)的瞬时频率、瞬时幅值和衰减因子。

相较于现有技术，本发明有以下有益效果：

1、利用端点优化对称延拓法改进之后的EMD进行分解，能使延拓之后的信号包络线端点最大程度逼近原始电力系统低频振荡信号端点，从而有效抑制端点“飞翼”现象；

2、利用K-S法准确的剔除了IMF分量中的伪分量，显著提高IMF分量获取的准确性，同时提高了电力系统低频振荡模态辨识能力；

3、利用数位演算法进行低频振荡模态瞬时参数的提取，能精确得到电力系统低频振荡信号的瞬时频率、瞬时幅值和衰减因子等信息。

附图说明

图1是本发明实施例的工作流程图。

图2是EMD端点效应改进之前EMD分解效果图。

图3是EMD端点效应改进之后EMD分解效果图。

表1是EMD端点效应改进之前电力系统低频振荡信号分析结果。

表2是EMD端点效应改进之后电力系统低频振荡信号分析结果。

具体实施方式

该电力系统低频振荡模态辨识方法结合图1进行说明，并对电力系统低频振荡信号x(t)＝15e^-0.8tcos(3πt)+10e^-0.7tcos(1.6πt)+5e^-0.6tcos(πt)进行EMD改进前后的仿真分析，选取前3个IMF分量进行数据对比分析，得到的效果图示于附图，具体步骤如下：

步骤1：选取经过滤波之后的电力系统低频振荡信号x(t)进行分析，其采样步长取为Δt，将电力系统低频振荡信号进行离散化处理，得到离散化之后的信号如下：

x_i＝[x₁,x₂,...,x_n]，i＝1,2,...,n

式中，x_i为离散化之后的信号值得集合，x₁,x₂,...,x_n为离散化之后t₁,t₂,...,t_n时刻的信号值；

步骤2：将离散化之后的电力系统低频振荡信号x_i两端点值用β、γ替代，构造新的数据序列x'_i：

x'_i＝[x'₀,x'_k,x'_n-1]

其中，x'₀＝β，x'_k＝x_k，x'_n-1＝γ，k＝1,2,...,n-2，i＝0,1,...,n-1；

步骤3：以端点x'₀、x'_n-1为中心对x'_i分别向两端进行对称延拓一个周期得到新的数据序列h：

h＝[h_l,h_i,h_r]＝[2x'₀-x'_i,x'_i,2x'_n-1-x'_n-1-i]

其中，i＝0,1,...,n-1；

步骤4：利用三次样条插值法对经过对称延拓之后的电力系统低频振荡信号数据序列h进行包络，得到包络线s(x)，优化x'_i的两个端点值x'₀(β)和x'_n-1(γ)，构建h_i与s_i的偏差评价函数：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\mathrm{\μ}\∫s^2\left(x\right)dx$

式中，μ为光滑系数，s_i为s(x)对应于h_i的那部分包络线；

步骤5：对上式进行离散化，得：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\frac{\mathrm{\μ}}{{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right(h_i\left)\right)}^2$

式中，Δx为采样间隔，▽²为二阶微分算子；

步骤6：由于经过延拓之后的数据序列h与包络线s(x)取决于端点值x'₀＝β和x'_n-1＝γ，故偏差评价函数同样取决于β和γ，上式可以变换为：

$\mathrm{\ϵ}=\frac{1}{2}\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-{x^{\′}}_i)}^2+\frac{\mathrm{\μ}}{2{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right({x^{\′}}_i\left)\right)}^2+{(\mathrm{\β}-x_0)}^2+{(\mathrm{\γ}-x_{n-1})}^2$

当 $\left\{\begin{array}{c}\∂\mathrm{\ϵ}/\∂\mathrm{\β}=0\\ \∂\mathrm{\ϵ}/\∂\mathrm{\γ}=0\end{array}\right.$ 时，ε为最小值，求出β和γ的值；

步骤7：根据数据序列h重新对电力系统低频振荡信号进行对称延拓，分解得到IMF分量c(t)和剩余分量r(t)；

步骤8：如果r(t)为单调函数，则转到步骤9，否则转到步骤2；

步骤9：得到IMF分量c(t)；

步骤10：利用K-S法将IMF与电力系统低频振荡信号x(t)进行对比分析，取离散化之后的电力系统低频振荡信号和任意一阶IMF分量N₁、N₂为x_i、y_i的时间信号序列点数，令x_i、y_i的累积分布函数为f(x)、r(x)，则：

f(x)＝N₁(i)/N₁

r(x)＝N₂(i)/N₂

式中，N₁(i)、N₂(i)分别为x_i、y_i中比x_i、y_i小的数据点数目。作累积分布函数f(x)、r(x)差绝对值的最大值D：

$D=\underset{-\mathrm{\∞}\≤x\≤\mathrm{\∞}}{max}|f\left(x\right)-r\left(x\right)|$

将差绝对值的最大值D转换为相似概率ρ：

$\mathrm{\ρ}\left(D\right)=Q_{ks}\left(\right(\sqrt{N_e}+0.12+\frac{0.11}{\sqrt{N_e}}\left)D\right)$

式中， $N_e=\frac{N_1{N}_2}{N_1+N_2},$ 令 $\mathrm{\λ}=(\sqrt{N_e}+0.12+\frac{0.11}{\sqrt{N_e}})D,$ Q_ks(λ)是K-S分布函数，则：

$\mathrm{\ρ}\left(D\right)=Q_{ks}\left(\mathrm{\λ}\right)=2\underset{j=1}{\overset{\mathrm{\∞}}{\Σ}}{(-1)}^{j-1}{e}^{-2j^2{\mathrm{\λ}}^2}$

当λ→0时，Q_ks(λ)→1；λ→∞时，Q_ks(λ)→0。由此得知，当累积分布函数在统计上趋于相近时，得到的相似程度越高，概率趋于1，反之趋近于0。据此原则剔除与电力系统低频振荡信号x(t)无关的虚假IMF(伪IMF)分量；

步骤11：利用数位演算法对c_i(t)进行各个模态瞬时参数的提取，对余下的有效IMF分量c_i(t)进行瞬时参数的提取，令其单一的离散低频振荡信号x(k)的表达式为：

x(k)＝Xe^-αkΔtcos(2πfkΔt+φ)

式中，k为离散信号采样点数，X为信号初始幅值，α为信号衰减因子，f为信号频率，Δt为采样间隔，φ为信号初始相位；

根据欧拉公式将x(k)展开得：

$x\left(k\right)=\frac{X}{2}{e}^{j\mathrm{\φ}}{e}^{-\mathrm{\α}k\mathrm{\Δ}t}{e}^{j2\mathrm{\π}fk\mathrm{\Δ}t}+\frac{X}{2}{e}^{-j\mathrm{\φ}}{e}^{-\mathrm{\α}k\mathrm{\Δ}t}{e}^{j2\mathrm{\π}fk\mathrm{\Δ}t}$

令a＝e^j2πfΔt,b＝e^-αΔt,z＝Re[a]，则：

x(k)＝Ab^ka^k+A^*b^ka^-k

$b=\sqrt{\frac{x^2(k+2)-x(k+1)x(k+3)}{x^2(k+1)-x\left(k\right)x(k+2)}}$

$z=\frac{bx\left(k\right)+b^{-1}x(k+2)}{2x(k+1)}$

其中，A^*是A的共轭函数。当信号长度大于4并且满足单一的离散低频振荡信号x(k)的表达式时，根据上述式子可以计算出瞬时频率f、瞬时幅值S和衰减因子α如下：

$f=\frac{\arccos z}{2\mathrm{\π}\mathrm{\Δ}t}$

S＝Xe^-αkΔt＝2|A|e^-αkΔt

$\mathrm{\α}=\frac{\ln b}{\mathrm{\Δ}t}$

为了让本领域人员更容易理解本发明，对电力系统低频振荡信号x(t)＝15e^-0.8tcos(3πt)+10e^-0.7tcos(1.6πt)+5e^-0.6tcos(πt)进行EMD改进前后的仿真分析，得到效果图图2和图3；如图2所示，为EMD端点效应改进之前EMD分解效果图；如图3所示，为EMD端点效应改进之后EMD分解效果图；并且选取前3个IMF分量进行数据对比分析，得到表1和表2；如表1所示，为EMD端点效应改进之前电力系统低频振荡信号分析结果；如表2所示，为EMD端点效应改进之后电力系统低频振荡信号分析结果。

以上是本发明的较佳实施例，凡依本发明技术方案所作的改变，所产生的功能作用未超出本发明技术方案的范围时，均属于本发明的保护范围。

表1

表2

Claims (1)

1.一种电力系统低频振荡模态辨识方法，其特征在于：

步骤1：对滤波后的电力系统低频振荡信号x(t)进行离散化处理，采样步长取为Δt，得到离散化处理信号：

x_i＝[x₁,x₂,...,x_n]，i＝1,2,...,n

式中，x_i为离散化之后的信号值的集合，x₁,x₂,...,x_n为离散化之后t₁,t₂,...,t_n时刻的信号值；

步骤2：信号x_i两端点值用β、γ替代，构造新的数据序列x'_i：

x'_i＝[x'₀,x'_k,x'_n-1]

式中，x'₀＝β，x'_k＝x_k，x'_n-1＝γ，k＝1,2,...,n-2，i＝0,1,...,n-1；

步骤3：以端点x'₀、x'_n-1为中心对x'_i分别向两端进行对称延拓一个周期得到新的数据序列：

h＝[h_l,h_i,h_r]＝[2x'₀-x'_i,x'_i,2x'_n-1-x'_n-1-i]

式中，i＝0,1,...,n-1；

步骤4：利用三次样条插值法对信号数据序列h进行包络，得到包络线s(x)，优化x'_i的两个端点值β和γ，构建h_i与s_i的偏差评价函数为：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\mathrm{\μ}\∫s^2\left(x\right)dx$

式中，μ为光滑系数，s_i为s(x)对应于h_i的那部分包络线；

步骤5：对偏差评价函数进行离散化，得：

$\mathrm{\ϵ}=\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-h_i)}^2+\frac{\mathrm{\μ}}{{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right(h_i\left)\right)}^2$

式中，Δx为采样间隔，▽²为二阶微分算子；

步骤6：将x'₀＝β和x'_n-1＝γ，上式可以变换

$\mathrm{\ϵ}=\frac{1}{2}\underset{i=0}{\overset{n-1}{\Σ}}{(s_i-{x^{\′}}_i)}^2+\frac{\mathrm{\μ}}{2{\left(\mathrm{\Δ}x\right)}^3}\underset{i=0}{\overset{n-1}{\Σ}}{\left({\▿}^{2}s\right({x^{\′}}_i\left)\right)}^2+{(\mathrm{\β}-x_0)}^2+{(\mathrm{\γ}-x_{n-1})}^2$

当 $\left\{\begin{array}{c}\∂\mathrm{\ϵ}/\∂\mathrm{\β}=0\\ \∂\mathrm{\ϵ}/\∂\mathrm{\γ}=0\end{array}\right.$ 时，ε为最小值，得出β和γ的值；

步骤7：根据数据序列h重新对电力系统低频振荡信号进行对称延拓，分解得到固有模态函数IMF分量c(t)和剩余分量r(t)；

步骤8：如果r(t)为单调函数，则转到步骤9，否则转到步骤2；

步骤9：经验模态分解EMD分解结束，得到固有模态函数IMF分量c(t)；

步骤10：利用柯尔摩科洛夫-斯米洛夫检验K-S法将c(t)与电力系统低频振荡信号x(t)进行对比分析，剔除与x(t)

无关的虚假固有模态函数IMF分量，得到余下的有效固有模态函数IMF分量c_i(t)；

步骤11：利用数位演算法对c_i(t)进行各个模态瞬时参数的提取，得到电力系统低频振荡信号x(t)的瞬时频率、瞬时幅值和衰减因子。