CN111896792A - 基于能量算子和变分模态分解的非稳态电压闪变检测方法 - Google Patents

Abstract

本发明公开了一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，步骤1：利用建立的五点改进能量算子，基于非稳态电压闪变信号的调制模型，获取非稳态电压闪变包络信号的分量v(n)；步骤2：采用变分模态分解方法对v(n)进行分析，基于能量损失系数和能量差确定自适应确定模态个数K值；步骤3：通过自适应变分模态将对v(n)分解得到K个IMF分量，对每个分量进行Hilbert变换，求取每个分量的瞬时频率f_j和瞬时幅值m_j；步骤4：基于K个IMF分量对应的闪变分量确定闪变的起止时刻和幅值突变时刻，实现非稳态电压闪变检测,本发明易于实施，且误差小。

Description

基于能量算子和变分模态分解的非稳态电压闪变检测方法

技术领域

本发明涉及一种基于能量算子和变分模态分解的非稳态电压闪变检测方法。

背景技术

随着电力电子技术的发展，越来越多的非线性和冲击性负荷注入电力系统中，导致大量短时电压波动和闪变等电能质量问题，电压闪变更多的呈现出非稳态、时变的特点，而目前大部分研究都集中在稳态情况下的闪变检测，对非稳态情况下的电压闪变研究较少。

闪变包络的精确提取是分析检测闪变的前提，非稳态电压闪变信号由于其闪变幅值可能突变，其包络信号会多转折点，传统的的能量算子会在转折点出产生很大的误差，而Hilbert变换提取包络的方法会在转折点出产生边缘效应，严重影响闪变参数的检测。

短时傅里叶变换(Short-Time Fourier Transform，STFT)、S变换和小波变换等方法被应用于非稳态时变信号分析，但STFT的时频窗口固定，只适合于分析特征尺度大致相同的过程，对于多尺度过程和突变过程效果较差，而S变换计算量大，运行时间长，小波变换则存在小波基选择困难的问题。

希尔伯特-黄变换(Hilbert-Huang Transform，HHT)能有效处理非线性非平稳信号，可同时在时间和频率达到很高的精度，实现非稳态时变电压闪变信号的测量，而且避免了小波基选取的问题，但HHT中的核心部分经验模态分解(Empirical ModeDecomposition，EMD)方法存在的模式混淆问题，有时无法准确分解出所有闪变包络分量。为了有效解决EMD存在的问题，Dragomiretskiy等人在2014年提出了变分模态分解(Variational Mode Decomposition，VMD)。与EMD相比，该方法通过预设模式数实现多分量信号的准确分离，而且分解精度更高，现主要被应用于机械故障诊断中。

因此，有必要设计一种用于分析非稳态电压闪变的方法。

发明内容

本发明所要解决的技术问题是提供一种基于能量算子和变分模态分解的非稳态电压闪变检测方法易于实施，且检测的精度高，误差小。

发明的技术解决方案如下：

一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，包括以下步骤：

步骤1：利用建立的五点改进能量算子，基于非稳态电压闪变信号的调制模型，获取非稳态电压闪变包络信号的分量v(n)；

步骤2：采用变分模态分解方法对v(n)进行分析，基于能量损失系数和能量差确定自适应确定模态个数K值；

步骤3：通过自适应变分模态将对v(n)分解得到K个IMF分量，对每个分量进行Hilbert变换，求取每个分量的瞬时频率fj和瞬时幅值m_j；

步骤4：基于K个IMF分量对应的闪变分量确定闪变的起止时刻和幅值突变时刻，实现非稳态电压闪变检测；

v(n)与每个分量的瞬时频率f_j、瞬时幅值m_j的关系为：

式中，g_j(n)为取值0或1的离散阶跃函数；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角。

非稳态电压闪变信号的调制模型为：

式中：u(t)表示非稳态电压闪变信号，U₀为基波电压幅值，ω₀、θ₀分别为基波电压角频率和初相角；A(t)表示幅度值变化的信号；v(t)为闪变包络信号；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角；g_j(t)为取值0或1的阶跃函数；

非稳态电压闪变信号的离散调制模型表达式为

式中：u(n)表示离散非稳态电压闪变信号；ω₀＝2πf₀/f_s，ω_j＝2πf_i/f_s，f₀为基波频率，f_s为采样频率，f_j为第j项闪变包络信号的频率，g_j(n)为取值0或1的离散阶跃函数；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角。

上式中u(n)即为非稳态电压闪变信号，v(n)为闪变包络信号的分量，瞬时频率和瞬时幅值分别为f_j和m_j。由于本专利研究的是非稳态电压闪变信号，因此f_j和m_j其可能随时间发生变化，通过瞬时频率曲线和瞬时幅值曲线就可以得到非稳态电压闪变信号的组成，从而达到分析非稳态电压闪变信号的目的。

知晓各个闪变分量的瞬时频率和瞬时幅值曲线信息就能得到闪变包络信号的分量v(n)，从而得到非稳态电压闪变信号u(n)。v(n)是非稳态电压闪变信号u(n)的一部分，v(n)加上1，再乘以载波电压U₀cos(ω₀n+θ₀)就能得到非稳态电压闪变信号。(相位信息在闪变计算闪变值时不需要用到，因此不用考虑相位，第j项闪变分量的相位直接用θj表示就可以)。

基于五点改进能量算子的包络提取表达式为：

其中，v(n)为非稳态电压闪变包络信号的分量；

ψ表示Teager能量算子函数，其表达式为

ψ[x(n)]＝x²(n)-x(n+1)x(n-1) (4)

所述基于能量损失系数和能量差确定自适应确定模态个数K值的具体实现步骤如下：

(1)初始化VMD的模态数M＝1；

(2)令M＝M+1，进行VMD分解；

(3)计算能量损失系数e；

式中u为输入的原始非稳态电压闪变信号，u_k表示第k个模态函数，∑u_k表示重建后的信号；

模态函数也即IMF分量，VMD算法原理是将原始输入信号u分解为M个模态函数，并在各模态函数之和与原始信号相等的要求下，使每个模态函数的估计带宽之和最小，比如当模态数M＝5时，对信号进行VMD分解，会有5个IMF分量，也即5个模态函数，5个模态函数之和就等于原始输入信号u。

||·||表示范数，是数学中的一种基本概念，||·||₂表示2-范数，||·||²表示范数的平方。

(4)重复步骤(2)～(3)，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；μ₁设为0.01；

(5)将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin；

原始信号的能量值E计算式如下

式中：u(i)表示为非稳态电压闪变信号；N为序列长度；不同模态分解层数的能量值之和大小不同，定义各模态分量能量之和与原信号能量的能量差参数计算式如下

式中，E_j表示为第j个模态分量的能量。

步骤3中，通过自适应变分模态分解得到K个IMF分量，对每个分量进行Hilbert变换，以第i个闪变包络分量为例，瞬时幅值求取表达式为：

式中

是u_i(t)的共轭信号，相位求取表达式为：

瞬时频率求取表达式为：

步骤4中，所述得到的非稳态电压闪变信号组成由瞬时频率和瞬时幅值曲线确定】首先根据模态分解个数K确定非稳态电压闪变含有K个闪变分量，也即有K个瞬时频率曲线，通过瞬时频率曲线的开始和结束的时间点确定每个闪变分量的闪变包络信号的起止时刻，然后判断瞬时频率曲线对应的瞬时幅值曲线幅值是否发生变化，若有则结合闪变包络信号的分量v(n)突变点处最大值对应的时刻来确定幅值突变时刻。

本发明认为通过所提的算法对该闪变分量前后检测得到的瞬时幅值不一样即认为幅值发生了变化，如前后的差值大于1％即可认为发生了变化。

基于所述的瞬时频率和瞬时幅值得到非稳态信号的组成成分，只要得到了每个分量的瞬时频率和瞬时幅值，就能通过非稳态电压闪变信号的离散调制模型得到非稳态信号组成成分。

首先，从非稳态电压闪变信号的离散调制模型表达式可知，K个闪变分量之和对应的就是表达式里的v(n)，也即闪变包络信号的分量，v(n)加上1，再乘以载波电压U₀cos(ω₀n+θ₀)就能得到非稳态电压闪变信号。

每个闪变分量都有一个起止时间(起始、停止时间)，也就是K个闪变分量会有相应的K个闪变起止时间，这个只是每个闪变分量对应的起止时间。而整个非稳态电压闪变信号闪变开始时间就是所有闪变分量中开始最早的时间，整个非稳态电压闪变信号停止最晚的时间。只要知晓了每个闪变分量对应的起止时间就能对引起闪变的原因具体化分析，更具研究意义，因此对于整体的闪变发生起止时间本专利不做具体说明。

步骤3通过自适应变分模态将对v(n)分解得到K个IMF分量，对每个分量进行Hilbert变换，求取每个分量的瞬时频率和瞬时幅值，而步骤4是根据步骤3得到每个分量的瞬时频率和瞬时幅值信息，确定每个闪变分量的起止时间以及幅值信息，从而实现非稳态电压闪变的分析。

有益效果：

本发明根据非稳态闪变信号多转折点的特点，提出一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其核心特点在于采用基于五点改进能量算子的闪变包络提取方法，相比传统能量算子计算量小、精度高、抗噪声干扰能力有显著提高，能有效减少非稳态闪变包络信号在转折点处的误差，通过自适应变分模态分解方法，能准确获得各闪变分量的时频信息，闪变的起止时刻和闪变包络幅值突变时刻由闪变包络信号的分量v(n)突变点以及瞬时频率和瞬时幅值曲线的突变点综合得到，根据五点改进能量算子的特点，其起止时刻误差小于3/fs。

另外，基于能量损失系数和能量差值自适应确定VMD的模态个数，通过自适应变分模态分解方法对闪变包络进行分解，求取各分量的瞬时频率和瞬时幅值，并通过瞬时频率曲线、瞬时幅值曲线和闪变包络信号的分量v(n)的突变信息确定闪变的起止时刻，灵活性好，具有自适应的特点。

本发明针对非稳态电压闪变包络信号多转折点的特点，构建一种基于五点改进能量算子的闪变包络提取方法，相比与传统的能量算子提取闪变包络，抗噪声能力更强，在转折点出的误差大大减小，能有效适应非稳态电压闪变包络信号的准确提取。

本发明通过基于所构建自适应变分模态分解方法，结合五点改进能量算子对非稳态电压闪变包络信号进行检测，能准确获得各闪变分量的时频信息，而且可以获得闪变发生的起止时刻信息。相比传统的HHT方法，本发明提出的算法检测精度更高，且避免其存在的模式混淆问题，具有重大的工程应用价值。

附图说明

图1为本发明的实施流程图；

图2为非稳态电压闪变包络信号仿真结果示意图；其中图2(a)为实际闪变包络信号波形图，图2(b)为采用传统能量算子提取闪变包络信号波形图；图2(c)为采用五点改进能量算子提取闪变包络信号波形图；

图3为单一闪变包络分量的测量结果；图3(a)和(b)分别为瞬时幅值曲线和瞬时频率曲线；

图4为闪变包络幅值发生变化时的测量结果示意图；图4(a)和(b)分别为瞬时幅值曲线和瞬时频率曲线；

图5为多频调制的闪变包络的测量结果示意图；图5(a)和(b)分别为瞬时幅值曲线和瞬时频率曲线。

具体实施方式

以下将结合附图和具体实施例对本发明做进一步详细说明：

实施例1：

一种基于五点改进能量算子和自适应变分模态分解的非稳态电压闪变检测方法，其包括：

步骤1，利用建立的五点改进能量算子，针对电压闪变信号的调制模型，提取并得到非稳态电压闪变包络信号的分量v(n)；

步骤2，采用变分模态分解对电压闪变包络信号的分量v(n)进行分析，初始化模态个数M，设为M＝1，令M＝M+1，循环计算能量损失系数e，设置阈值μ₁，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin。

步骤3，通过自适应变分模态分解得到K个IMF分量，对每个分量进行Hilbert变换，求取瞬时频率和瞬时幅值。

步骤4，通过瞬时频率和瞬时幅值曲线的综合得到非稳态信号的组成成分，并通过其中的突变信息结合闪变包络信号的分量v(n)确定闪变的起止时刻。

更优选地，步骤1中所述的非稳态电压闪变信号的调制模型为：

式中：U₀为基波电压幅值，ω₀、θ₀为基波电压角频率和初相角；A(t)表示幅度值变化的信号；v(t)为闪变包络信号；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角；g_j(t)为取值0或1的阶跃函数，以单一闪变频率(幅值、角频率、相角分别为m₁、ω₁、θ₁)调制为例，考虑0～10s的闪变信号，当g₁(t)＝1时，式(1)即为稳态下的电压闪变信号的调制模型；当

表示该闪变信号开始于第2.5s，结束于第5s，持续时间为2.5s。

更优选地，非稳态电压闪变信号的离散调制模型表达式为

式中：ω₀＝2πf₀/f_s；ω_j＝2πf_j/f_s。

更优选地，所述五点改进能量算子其包络提取表达式为：

式中ψ表示Teager能量算子函数，其表达式为

ψ[x(n)]＝x²(n)-x(n+1)x(n-1) (4)

更优选地，所述的VMD旨在将输入信号f分解为K个模态函数u_k(t)，并在各模态函数之和与信号/相等的要求下，使每个模态函数的估计带宽之和最小，上述约束变分问题可表述为

式中，{u_k}为各模态分量的集合，{ω_k}为各模态对应的中心频率集合。

引入Lagrange乘子λ和二次惩罚因子α使式(5)转化成非约束性的变分问题，其表达式为

VMD采用交替方向乘子法求解式(6)，即通过

和

交替更新寻求上式的拉格朗日鞍点。其中

取值可表述为

利用Parseval/Plancherel傅里叶等距变换将式(7)转变到频域，可得二次优化结果为

同样将中心频率的取值问题转换到频域，可得中心频率的计算式为

更优选地，VMD的具体实现过程如下：

(1)初始化

和

且令n＝0；

(2)令n＝n+1，执行整个循环；

(3)取k＝1∶K，根据式(4)和式(5)更新u_k和ω_k；

(4)根据下式更新λ，有

(5)给定判别精度e＞0，当收敛满足式(11)结束整个循环输出结果，否则返回步骤(2)。

由VMD算法原理可知，其处理信号时需要预先设定模态分解个数K，受实际情况的限制，K值通常难以准确设定，当取值较小会导致多个分量分解在同一模态中，导致信息提取不完全，反之，K取值较大会造成模态分量的中心频率重叠，产生虚假分量。

更优选地，所述的自适应变分模态分解特征在于自适应确定K值，

具体实现步骤如下：

(1)初始化VMD的模态数M＝1；

(2)令M＝M+1，进行VMD分解；

(3)计算能量损失系数e；

式中f为输入的原始非稳态电压闪变信号，u_k表示第k个模态函数，∑u_k表示重建后的信号。

(4)重复步骤(2)～(3)，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；

(5)将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin。

更优选地，通过自适应变分模态分解得到K个IMF分量，对每个分量进行Hilbert变换，以第i个闪变包络分量为例，瞬时幅值求取表达式为：

式中

是u_i(t)的共轭信号，相位求取表达式为：

瞬时频率求取表达式为：

得到的非稳态电压闪变信号组成由瞬时频率和瞬时幅值曲线确定，首先根据模态分解个数K确定非稳态电压闪变含有K个闪变分量，也即有K个瞬时频率曲线，通过瞬时频率曲线的开始和结束的时间点确定闪变包络信号的起止时刻，然后判断瞬时频率曲线对应的瞬时幅值曲线幅值是否发生变化，若有则结合闪变包络信号的分量v(n)突变点处最大值对应的时刻来确定幅值突变时刻，根据五点改进能量算子的特点，其误差小于3/fs。

本发明针对非稳态电压闪变信号的准确检测开展研究，首先建立的五点改进能量算子，针对电压闪变信号的调制模型，提取并得到非稳态电压闪变包络信号的分量v(n)；构建基于能量损失系数和能量差的方法自适应确定模态个数K值，通过对电压闪变包络信号的分量v(n)进行自适应变分模态分解(AVMD)得到K个IMF分量，对每个分量进行Hilbert变换，求取瞬时频率和瞬时幅值，并通过其中的突变信息结合闪变包络信号的分量v(n)确定闪变的起止时刻和幅值突变时刻。

在执行本发明之前，首先需要做如下三个工作：

步骤s1，建立非稳态电压闪变信号下的调制模型；

建立的非稳态电压闪变信号的调制模型为：

式中：u(t)表示非稳态电压闪变信号，U₀为基波电压幅值，ω₀、θ₀为基波电压角频率和初相角；A(t)表示幅度值变化的信号；v(t)为闪变包络信号；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角；g_j(t)为取值0或1的阶跃函数，以单一闪变频率(幅值、角频率、相角分别为m₁、ω₁、θ₁)调制为例，考虑0～10s的闪变信号，当g₁(t)＝1时，式(1)即为稳态下的电压闪变信号的调制模型；当

表示该闪变信号开始于第2.5s，结束于第5s，持续时间为2.5s。

非稳态电压闪变信号的离散调制模型表达式为

式中：u(n)表示离散非稳态电压闪变信号；ω₀＝2πf₀/f_s，ω_j＝2πf_j/f_s，f₀为基波频率，f_s为采样频率，f_j为第j项闪变包络信号的频率，g_j(n)为取值0或1的离散阶跃函数。

步骤s2，建立五点改进能量算子提取非稳态电压闪变包络信号方法。

由于非稳态电压闪变包络有很多转折点，采用Hilbert变换提取闪变包络会有很大的边缘效应，而单纯采用现有Teager能量算子提取电压闪变包络信号会在转折点出产生很大误差，为此，本发明人提出建立五点改进能量算子提取非稳态电压闪变包络信号，以便降低转折处的误差，提升抗噪声性能，进而准确追踪电压闪变包络信号的幅值和频率变化。

设连续时间信号为x(t)，其能量算子可表示为

ψ[x(t)]＝[x′(t)]²-x(t)x″(t) (3)

式中：x′(t)、x″(t)分别表示x(t)的一阶和二阶导数。

设其离散信号为x(n)，则能量算子运算结果为

ψ[x(n)]＝x²(n)-x(n+1)x(n-1) (4)

式中：n为任意离散采样点，n-1表示采样点时刻n的前一个采样点值，n+1表示后一个采样点值。

对离散号x(n)与y(n)的乘积进行能量算子变换有

ψ[x(n)y(n)]＝x²(n)ψ[y(n)]+y²(n)ψ[x(n)]-ψ[x(n)]ψ[y(n)] (5)

因此，对电压闪变信号的离散表达式中的u(n)进行能量算子运算，可得

由于实际电网电压闪变包络信号的幅值U_i一般小于工频电压幅值U₀的0.1倍，即调幅波系数m_i通常小于0.1，且闪变频率f_i的变化范围为0.05～35Hz，电网工频f₀等于50Hz，而闪变分量h值一般比较小，因此，忽略误差因子E(n)影响，闪变包络v(n)可近似为

离散求导一般通过差分方式进行，从泰勒展开式出发，若函数f(x₀)在含x₀的展开区间内有n+1阶导数，则在该区间内可泰勒展开为各级项与拉格朗日余弦之后，如下式

因此，

联立公式(9)～(12)有

根据式(13)构造差分方程，

对y(n)进行能量算子运算，有

因此，

当kΩ＜π/4，即f/F_s＜1/(8k)时，可近似有sin(kΩ)＝kΩ，

本发明提出的改进能量算子是针对式(6)中的提取因子D(n)和误差因子E(n)做进一步的计算，忽略其中的幅值较小项，有

因此，忽略误差因子E(n)后的改进能量算子提取的闪变包络v(n)的表达式为

将式(16)代入到式(19)，有

因此，闪变包络的提取公式为

本发明提出的五点改进能量算子残差余弦为O(Δn⁴)，相比Teager能量算子的残差余弦O(Δn)有更高的精度，使用五点求导能减少两点差分求导时的偶然误差，抗噪性能更强，相比与式(7)中传统的能量算子提取闪变包络公式，无需平方根运算，计算速度更快，能有效减少非稳态闪变包络信号转折点出的误差，提高检测的准确性。

步骤s3，基于能量损失系数和能量差构建自适应确定模态个数K值方法。

由VMD算法原理可知，其处理信号时需要预先设定模态分解个数K，受实际情况的限制，K值通常难以准确设定，当取值较小会导致多个分量分解在同一模态中，导致信息提取不完全，反之，K取值较大会造成模态分量的中心频率重叠，产生虚假分量。

基于能量损失系数和能量差构建自适应确定模态个数K值的具体实现步骤如下：

(1)初始化VMD的模态数M＝1；

(2)令M＝M+1，进行VMD分解；

(3)计算能量损失系数e；

式中f为输入的原始非稳态电压闪变信号，u_k表示第k个模态函数，∑u_k表示重建后的信号。

(4)重复步骤(2)～(3)，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；

(5)将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin。

做好上述准备工作后，进入本发明的实施过程，实施流程如图1所示，包括如下步骤：

步骤1，利用步骤s2建立五点改进能量算子，针对步骤s1中建立的非稳态电压闪变信号下的调制模型，提取非稳态电压闪变包络信号的分量，其表达式如式(21)所示。

步骤2，采用步骤S3中建立的基于能量损失系数和能量差构建自适应确定模态个数K值方法。

(1)初始化VMD的模态数M＝1；

(2)令M＝M+1，进行VMD分解；

(3)计算能量损失系数e；

(4)重复步骤(2)～(3)，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；

(5)将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin。

步骤3，通过对电压闪变包络信号的分量v(n)进行自适应变分模态分解(AVMD)得到K个IMF分量，对每个分量进行Hilbert变换，求取瞬时频率和瞬时幅值。

通过自适应变分模态分解得到K个IMF分量，对每个分量进行Hilbert变换，以第i个闪变包络分量为例，瞬时幅值求取表达式为：

式中

是u_i(t)的共轭信号，相位求取表达式为：

瞬时频率求取表达式为：

步骤4，通过瞬时频率和瞬时幅值曲线的综合得到非稳态信号的组成成分，并通过其中的突变信息结合闪变包络信号的分量v(n)确定闪变的起止时刻。

首先根据模态分解个数K确定非稳态电压闪变含有K个闪变分量，也即有K个瞬时频率曲线，通过瞬时频率曲线的开始和结束的时间点确定闪变包络信号的起止时刻，然后判断瞬时频率曲线对应的瞬时幅值曲线幅值是否发生变化，若有则结合闪变包络信号的分量v(n)突变点处最大值对应的时刻来确定幅值突变时刻，根据五点改进能量算子的特点，其误差小于3/fs。

本发明的效果通过如下仿真实验可以获得，具体如下：

1、五点改进能量算子和传统能量算子提取非稳态电压闪变包络仿真

为了验证本文提出算法的有效性和准确性，基于MATLAB平台完成本文算法的仿真实验与分析。设载波电压的幅值U₀为1pu，电网工频f₀为50Hz，闪变包络频率fi为17Hz，闪变包络幅值为

由图2可知，对同一非稳态电压闪变信号，五点改进能量算子提取的闪变包络相比与传统能量算子提取的包络在转折点出的误差显著减小，其整体波动性更小，与实际非稳态电压闪变包络曲线更加接近，由图2(c)可知，幅值突变时刻发生在2.5s和5s处，误差小于3/fs。

2、针对单一频率的非稳态电压闪变包络信号的仿真

设电网工频f₀为50Hz，幅值U₀为1p.u.(p.u.为标称值)，考虑到对照度波动作用最为显著的是8.8Hz闪变包络的正弦波波动电压，本文先以闪变包络信号的频率fi为8.8Hz闪变包络信号进行仿真，设闪变包络系数mi为2.5×10^-3pu，闪变开始于第2s，结束于第8s，仿真检测结果如图3所示。

图3为不同闪变起止时刻的仿真实验结果，其中(a)为闪变包络信号的瞬时幅值曲线，(b)为闪变包络信号的瞬时频率曲线。

由图3(a)可以看出闪变包络信号且在2s开始和8s结束，幅值为2.5×10^-3pu，图3(b)可以看出闪变包络频率为8.8Hz，说明图3中信号在2s时发生闪变和8s时闪变结束。本文提出的方法不仅能准确测量闪变参数，而且可准确测量闪变发生的起止时刻。

3、针对幅值发生变化时的非稳态电压闪变包络信号的仿真

对式(7)所示的电压闪变信号，设载波电压的幅值U₀为1pu，电网工频f₀为50Hz，闪变包络频率fi为17Hz，闪变包络幅值设置同式(26)仿真实验结果如图4所示，图4(a)为闪变包络的幅值变化曲线，从图4(a)可以看出闪变包络幅值在2.5s、5.0s时发生变化，由图4(b)可知闪变包络的频率为17Hz。可见，本文提出的基于自适应变分模态分解的电压闪变包络参数检测方法对闪变包络幅值发生变化的电压闪变信号也能实现准确测量。

4、针对多频率包络调制时的非稳态电压闪变包络信号的仿真

实际电压波动分量中可能含有多种频率成分的正弦波，为验证本文算法在多频闪变包络调制情况下的有效性，设置含有多频率成分调制的电压闪变信号参数如表2所示，其中，25Hz闪变包络分量发生的起止时间分别设定为2.0s和6.0s。

表1多频率闪变包络分量的测量误差

5、与传统检测方法比较

将本发明提出的基于五点改进能量算子和自适应变分模态分解方法、HHT方法与基于4项3阶Nuttall窗双谱线插值FFT方法的电压包络参数检测方法的仿真实验结果进行对比，3种仿真实验的输入参数一致，闪变起止时刻随机产生，测量结果如表2所示。

表2 3种测量方法的仿真结果对比

由表3误差分析结果可知，自适应变分模态分解方法的幅值相对误差整体小于HHT方法和4项3阶Nuttall窗双谱线插值FFT方法。可见，基于自适应变分模态分解的电压闪变包络参数检测方法具有较高的测量精度。此外，基于自适应变分模态分解的包络参数检测方法可准确判断闪变发生的起止时刻，优于不具备时域特性的FFT方法

虽然本发明已以较佳实施例公开如上，但实施例并不限定本发明。在不脱离本发明之精神和范围内，所做的任何等效变化或润饰，同样属于本发明之保护范围。因此本发明的保护范围应当以本申请的权利要求所界定的内容为标准。

Claims (6)

1.一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，包括以下步骤：

步骤1：利用建立的五点改进能量算子，基于非稳态电压闪变信号的调制模型，获取非稳态电压闪变包络信号的分量v(n)；

步骤2：采用变分模态分解方法对v(n)进行分析，基于能量损失系数和能量差确定自适应确定模态个数K值；

步骤3：通过自适应变分模态将对v(n)分解得到K个IMF分量，对每个分量进行Hilbert变换，求取每个分量的瞬时频率f_j和瞬时幅值m_j；

步骤4：基于K个IMF分量对应的闪变分量确定闪变的起止时刻和幅值突变时刻，实现非稳态电压闪变检测；

v(n)与每个分量的瞬时频率f_j、瞬时幅值m_j的关系为：

式中，g_j(n)为取值0或1的离散阶跃函数；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角。

2.根据权利要求1所述的一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，非稳态电压闪变信号的调制模型为：

式中：u(t)表示非稳态电压闪变信号，U₀为基波电压幅值，ω₀、θ₀分别为基波电压角频率和初相角；A(t)表示幅度值变化的信号；v(t)为闪变包络信号；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角；g_j(t)为取值0或1的阶跃函数；

非稳态电压闪变信号的离散调制模型表达式为

式中：u(n)表示离散非稳态电压闪变信号；ω₀＝2πf₀/f_s，ω_j＝2πf_i/f_s，f₀为基波频率，f_s为采样频率，f_j为第j项闪变包络信号的频率，g_j(n)为取值0或1的离散阶跃函数；h为波动分量的项数；m_j、ω_j、θ_j分别为第j项波动分量的幅值系数、角频率和相角。

3.根据权利要求2所述的一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，基于五点改进能量算子的包络提取表达式为：

其中，v(n)为非稳态电压闪变包络信号的分量；

ψ表示Teager能量算子函数，其表达式为

ψ[x(n)]＝x²(n)-x(n+1)x(n-1) (4)

4.根据权利要求1所述的一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，所述基于能量损失系数和能量差确定自适应确定模态个数K值的具体实现步骤如下：

(1)初始化VMD的模态数M＝1；

(2)令M＝M+1，进行VMD分解；

(3)计算能量损失系数e；

式中u为输入的原始非稳态电压闪变信号，u_k表示第k个模态函数，∑u_k表示重建后的信号；

(4)重复步骤(2)～(3)，通过迭代运算直至损失系数e小于设定的阀值μ₁，得到更新后的M值；μ₁设为0.01；

(5)将最大模态分解个数设为M，先求取原始信号的能量值E，然后计算当模态分解个数设为j时的累加能量E_j(j＝1，2，...，M)，并求取能量差值ρ，选取能量差ρ最小时的j_ρmin值作为最佳模态分解个数K，K＝j_ρmin；

原始信号的能量值E计算式如下

式中：u(i)表示为非稳态电压闪变信号；N为序列长度；不同模态分解层数的能量值之和大小不同，定义各模态分量能量之和与原信号能量的能量差参数计算式如下

式中，E_j表示为第j个模态分量的能量。

5.根据权利要求1所述的一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，步骤3中，通过自适应变分模态分解得到K个IMF分量，对每个分量进行Hilbert变换，以第i个闪变包络分量为例，瞬时幅值求取表达式为：

式中

是u_i(t)的共轭信号，相位求取表达式为：

瞬时频率求取表达式为：

6.根据权利要求1-5任一项所述的一种基于能量算子和变分模态分解的非稳态电压闪变检测方法，其特征在于，步骤4中，首先根据模态分解个数K确定非稳态电压闪变含有K个闪变分量，也即有K个瞬时频率曲线，通过瞬时频率曲线的开始和结束的时间点确定每个闪变分量的闪变包络信号的起止时刻，然后判断瞬时频率曲线对应的瞬时幅值曲线幅值是否发生变化，若有则结合闪变包络信号的分量v(n)突变点处最大值对应的时刻来确定幅值突变时刻。

Images