CN109787250B - 一种基于多元经验模态分解的电力系统低频振荡模式辨识方法 - Google Patents

Abstract

本发明涉及电力系统领域，尤其涉及一种基于多元经验模态分解的电力系统低频振荡模式辨识方法。包括：利用基于多元经验模态分解的改进算法对电力系统中的低频振荡信号进行预分解处理，在得到对应不同振荡模式的多个本征模函数分量后，利用Teager能量算子的快速响应能力，计算各本征模函数分量的相对能量值，再以能量作为判断依据，筛选出能够反映系统真实振荡情况的主导振荡模式并将虚假噪声部分剔除，最后通过预测误差法计算主导振荡模式所对应的振荡模式参数，即频率及阻尼比，从而完成对电力系统主导振荡模式的辨识。本方法实现了基于广域量测信息的电力系统主导振荡模式的快速、准确、以及高效辨识。

Description

一种基于多元经验模态分解的电力系统低频振荡模式辨识 方法

技术领域

本发明涉及电力系统领域，尤其涉及一种基于多元经验模态分解的电力系统低频振荡模式辨识方法。

背景技术

近年来，区域电网互联规模不断扩大、大容量远距离交直流输电不断增加、可再生能源大规模接入，使得区域间低频振荡已成为限制互联电网输电能力、威胁电网安全稳定运行的重要因素之一。因此，研究在全国联网、可再生能源大规模接入背景下的电力系统低频振荡辨识方法具有十分重要的现实意义和工程实用价值。

目前，电力系统的动态稳定分析方法主要分为基于电力系统数学模型的动态稳定分析方法和基于广域量测数据的动态稳定分析方法。基于电力系统数学模型的动态稳定分析方法，通常在系统稳定平衡点处线性化，得到系统的线性化状态方程和状态矩阵，根据状态矩阵的特征值和特征向量来估计系统的主导振荡模式、振荡模态、参与因子及同调机群。该方法虽可从原理上，全面、系统地评估电力系统的动态稳定性，但在实际应用中存在如下瓶颈：①电力系统元件级和系统级精确建模仍面临巨大挑战；②基于电力系统模型的稳定分析计算量随电网规模扩大而急剧增加；③基于模型的电力系统动态稳定分析方法实时性较差，难以及时跟踪电力系统运行方式变化而实时评估；④稳定分析结果仅对给定运行点有限邻域有效，当系统运行点远离所设定的运行点后，分析结果不能准确表征系统实际的动态稳定性等。因此，基于模型的电力系统动态稳定分析方法适用于电网规划、运行方式安排、安全稳定校核等离线应用场合；目前基于广域量测数据的动态稳定分析方法主要有傅立叶变换法、小波变换算法、Prony算法等。这些方法对于电力系统小干扰稳定性的研究都取得了一些成果，但仍各有欠缺之处：傅立叶变换通过频谱分析得到信号频率，但无法分析信号的阻尼特性；小波算法在傅里叶变换基础上通过小波基函数给信号加窗，能有效反映信号的时频特性，但小波变换算法因频域混叠使其分辨率有限，且小波基函数不易确定；Prony算法是一种多项式线性拟合的辨识方法，该方法首先估计可逼近实测数据的多项式，然后通过求解多项式的根来辨识系统主导振荡模式的振荡频率和阻尼比，由于该方法为线性拟合算法，因而对噪声敏感，抗噪能力差，且该算法只局限于处理电力系统大扰动信号，不适用处理电力系统的小扰动(类噪声)信号。

发明内容

针对上述现有技术中存在的问题，本发明提供了一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，该方法克服了传统基于经验模态分解的主导振荡模式辨识方法对于复杂的多机系统辨识效率低下的不足，实现了基于广域量测信息的电力系统主导振荡模式的快速、准确、高效辨识，提升了电力系统动态稳定在线评估的快速性、准确性和鲁棒性。

为了解决基于传统经验模态分解辨识方法对于多个信道之间振荡模式难以校准，无法很好地反映整个系统内区域间振荡关系的问题，

为了解决背景技术中关于基于随机子空间的电力系统主导振荡辨识耗时较多的不足，本发明实施例提出一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，以实现电力系统主导振荡模式的快速辨识，为电网运行调度人员提供快速、可靠的电网运行状态信息以改善电力系统的动态稳定性。

本发明所采用的技术方案是：

一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述方法包括：

1)通过基于多元经验模态分解的改进算法对电力系统中的初始振荡信号进行与分解处理，得到多个包含系统不同振荡模式的本征模函数分量；

2)利用Teager能量算子的快速响应能力，计算各本征模函数分量在采样点处的相对能量值；

3)对各本征模函数分量在所有采样点处的相对能量值进行积分求和，得到整个本征模函数分量的相对能量，并以能量为判断依据，筛选出能够反映系统真实振荡情况的主导振荡模式并将虚假噪声部分剔除；

4)通过预测误差法计算主导振荡模式所对应的振荡模式参数，即频率及阻尼比。

在步骤1)之前，所述方法还包括：

从广域量测系统中获取电力系统的广域量测信息，在广域量测信息中提取出电力系统区域间联络线的有功功率信号作为初始振荡信号。

所述多元经验模态分解对于初始信号的处理过程具体为：

1)将初始信号x(t)在多维空间中沿多个方向进行一次投影，得到多个一次映射；

2)将多个一次映射沿初始信号方向进行二次投影，得到多个二次映射；

3)通过多元样条插值法拟合出各个映射的包络线，并计算所有包络线的均值得到初始输入信号的局部均值m(t)；

4)利用局部均值m(t)不断对初始信号x(t)进行分解，依次得到多个本征模函数分量。

所述一次投影时投影方向向量的选取方式具体为：

等角度采样，即在多维空间中建立单位超球面，分别以等角度差在超球面上建立类似于经线和纬线的曲线，将所有超球面的球心与超球面上曲线交点的连线作为一次投影的单位方向向量。

所述本征模函数的获取方法具体为：

将多通道信号与局部均值作差，得到差值d(t)，即：

x(t)-m(t)＝d(t)

判断差值d(t)是否可以作为一个本征模函数，需要满足以下两个条件：

①整个数据段内，极值点包括极大值和极小值，数目和过零点数目需相等或仅相差一个；

②在任何一点，由极大值确定的上包络与极小值确定的下包络计算出的局部均值为零。

但在实际工况下，绝大多数信号较为复杂，并不能绝对性地满足上述条件，故将上述条件转化为易于实现判别的数量判据，即：

式中N表示信号的采样点个数，P表示分解过程中的迭代次数。根据经验，若SD的取值在0.2至0.3之间时，认为该次迭代的结果可以作为一个本征模函数分量，并将其作为新的输入信号，再次与局部均值作差，重新进行迭代；反之，若SD的取值不在0.2至0.3之间，则将d(t)作为新的输入信号，重新进行迭代，K表示一个变量，k＝0、1、2…N。

重复按照上述过程不断地筛选出本征模函数，直至迭代过程中信号的二次映射接近单调时，信号沿各个方向向量的投影极值点个数小于3个，认为分解过程不会再产生新的本征模函数，多元经验模态分解的分解过程结束。

所述能量算子的计算方法具体为：

对于一个正弦信号，设其离散形式的表达式为

x_n＝A cos(Ω_n+φ)

式中n为信号采样点的总个数，A表示信号的幅值，Ω_n为信号在采样点处的瞬时频率，φ为初相；则该信号在采样点处Teager能量算子的表达式为：

上式中：ψ(x_n)为信号在采样点处的能量值，A表示信号的幅值，Ω为信号在采样点处的瞬时频率，x_n、x_n-1、x_n+1为信号在采样点n、n-1、n+1处的值；

然而在实际工程应用中，需要处理的信号通常不是标准的正弦信号，设一个频率及幅值均随时间变化的信号为

x_n＝a_ncosφ_n

上式中：x_n为信号在采样点n处的值，a_n为信号幅值函数，φ_n为相位函数；

其瞬时频率可表示为

上式中：n是时间刻度，与时间t等价，因为n是离散形式的时间刻度，

就是对相位函数φ_n求时间的导数；

在信号中选择3个相邻的点来构造一个新的对称差分信号

上式中：a_n为信号幅值函数，φ_n为相位函数，Ω为瞬时频率，s_n为n处的能量算子；则s_n的Teager能量算子可以表示为

所述本征模函数的相对能量值计算方法具体为：

利用Teager能量算子计算得到各个本征模函数分量在所有采样点处的瞬时相对能量值后，将所有采样点处的瞬时相对能量值进行积分求和，得到各个本征模函数分量的相对能量值。即

式中n为各本征模函数分量信号的采样点个数，x_n为信号在采样点n处的值，E为IMF 能量总和。

所述主导振荡模式的筛选原则具体为：

低频振荡信号经多元经验模态分解后，若某个或某几个本征模函数分量的相对能量明显高于其他本征模函数分量相对能量，则认为该本征模函数分量中包含着系统的主导振荡模式。其他能量相对较低的本征模函数则被认为包含着系统的虚假振荡模式。

所述低频振荡模式参数的计算方式具体为：

将振荡信号经多元经验模态分解后得到的本征模函数信号分量写成如下形式：

其中，u_t为输入变量，y_t为输出变量，e_t为随机噪声，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数，

为对任意的时间t；

对于给定的n_a，n_b＞0，有：

其中

是固定但未知的待辨识参数向量，记为θ，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数；

从而信号也可等价地表示为如下形式：

上式中：y_t为输出变量，e_t为随机噪声，其中

是固定但未知的待辨识参数向量，记为θ，

为对任意的时间t；

其中：

上式中：u_t为输入变量，y_t为输出变量，

为对任意的时间t；

因此，想要求得各本征模函数分量所对应的振荡模式，可以通过计算出参数向量θ中的具体数值，然后建立如下特征方程：

上式中：q为一个常数，

是固定但未知的待辨识参数向量；

计算满足此方程的所有特征根λ_i＝δ_i+jω_i，便可得到每一对共轭复根所对应的系统振荡模式。

其中振荡阻尼比与δ_i相关，具体计算方式为：

振荡频率与ω_i相关，具体计算方式为：

上式中：ω_i特征根的虚部，δ_i特征根的实部，ξ为阻尼比，f为振荡频率；

所述预测误差法估计参数的过程具体为：

预测误差法在确定待辨识系统的模型之后，需要对待辨识的参数θ建立一个损失函数，记作J(θ)，然后运用数值解法对J(θ)进行极小化求解，从而得到待辨识参数θ的具体数值。因此，辨识过程可以抽象成对参数向量θ进行最优估计的过程，即

首先设定初始值θ⁽⁰⁾，参数向量满足如下的一般递归形式：

θ^(k+1)＝θ^(k)+Δθ^(k)＝θ^(k)+γb(J，θ^(k))

式中b(J，θ^(k))是一个修正因子，它可以在每一步迭代中使之前估计值得到轻微改善，使估计值在反复迭代过程中逐渐达到最优。b(J，θ^(k))是一个修正因子，b表示一个函数，关于 J，θ^(k)的函数关系；γ为一个系数；

采用牛顿-拉夫逊法进行迭代，使得预测误差法在迭代过程中的修正因子b(J，θ^(k))可以通过计算J(θ)在当前时刻参数估值θ^(k)的二次逼近来得到，该方法具体迭代形式如下：

θ^(k+1)＝θ^(k)-α_k(V″_n(θ^(k)))^-1V′_n(θ^(k))

其中，α_k代表步长，是一个关于k的递减函数；V′_n(θ^(k))代表损失函数J(θ)在θ^(k)处的梯度；V′_n(θ^(k))是损失函数J(θ)在θ^(k)处的海森矩阵。

在整个迭代过程中，对于预先设定的精度ε₀，当|Δθ^(k)|≤ε₀时，停止计算，即得最终的参数估计值θ^(k+1)；否则用所得的θ^(k+1)来代替θ^(k)，重复计算，直到|Δθ^(k)|满足精度要求为止。

本发明提供的技术方案的有益效果是：

1、本发明实现了基于广域量测信息的电力系统主导振荡模式的快速辨识；可为电网运行调度人员提供更加快速、准确的系统主导振荡模式评估信息，有利于电网运行调度人员及时掌握系统的动态稳定性；

2、本发明解决了基于传统经验模态分解辨识方法对于多个信道之间振荡模式难以校准，无法很好地反映整个系统内区域间振荡关系的问题，实现了电力系统主导振荡模式的快速、准确及鲁棒辨识；

3、本发明可实现从电力系统的广域故障信号和类噪声信号中快速、准确的辨识出系统的主导振荡模式；

4、本发明提升了电力系统动态稳定在线评估的预警响应速度。

附图说明

图1为一种基于多元经验模态分解的电力系统低频振荡模式辨识方法的流程图；

图2为三维坐标系示意图；

图3为采样点集在超球面的分布情况示意图

图4为16机68节点测试系统图；

图5为各区域间联络线有功功率振荡曲线图；

图6为联络线1-2的有功功率振荡曲线分别经过经验模态分解和多元经验模态分解后得到的前5个IMF分量波形图；

图7为各个有功功率信号经多元经验模态分解后IMF1的波形图；

图8为各个有功功率信号经多元经验模态分解后IMF2的波形图；

图9为各个有功功率信号经多元经验模态分解后IMF3的波形图；

图10为各个IMF相对能量图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面对本发明实施方式作进一步地详细描述。

利用基于多元经验模态分解(Multivariate Empirical Mode Decomposition，MEMD)的改进算法对电力系统中的低频振荡信号进行预分解处理，在得到对应不同振荡模式的多个本征模函数Intrinsic Mode Function，IMF分量后，利用Teager能量算子的快速响应能力，计算各信号分量的相对能量，再以能量作为判断依据，筛选出能够反映系统真实振荡情况的主导振荡模式并将虚假噪声部分剔除，最后通过预测误差法Prediction ErrorMethod，PEM计算主导振荡模式所对应的振荡模式参数，即频率及阻尼比，从而完成对电力系统主导振荡模式的辨识。

实施例1：

本发明一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，参见图1，该方法包括以下步骤：

101：通过基于多元经验模态分解的改进算法对电力系统中的初始振荡信号进行与分解处理，得到多个包含系统不同振荡模式的本征模函数分量；

102：利用Teager能量算子的快速响应能力，计算各本征模函数分量在采样点处的相对能量值；

103：对各本征模函数分量在所有采样点处的相对能量值进行积分求和，得到整个本征模函数分量的相对能量，并以能量为判断依据，筛选出能够反映系统真实振荡情况的主导振荡模式并将虚假噪声部分剔除；

104：通过预测误差法计算主导振荡模式所对应的振荡模式参数，即频率及阻尼比。

在步骤101之前，该方法还包括：

从广域量测系统中获取电力系统的广域量测信息，在广域量测信息中提取出电力系统区域间联络线的有功功率信号作为初始振荡信号。

其中，步骤101中的多元经验模态分解对于初始信号的处理过程具体为：

1)将初始信号x(t)在多维空间中沿多个方向进行一次投影，得到多个一次映射；

2)将多个一次映射沿初始信号方向进行二次投影，得到多个二次映射；

3)通过多元样条插值法拟合出各个映射的包络线，并计算所有包络线的均值得到初始输入信号的局部均值m(t)；

4)利用局部均值m(t)不断对初始信号x(t)进行分解，依次得到多个本征模函数分量。

进一步地，一次投影时投影方向向量的选取方式具体为：

等角度采样，即在多维空间中建立单位超球面，分别以等角度差在超球面上建立类似于经线和纬线的曲线，将所有超球面的球心与超球面上曲线交点的连线作为一次投影的单位方向向量。

进一步地，本征模函数的相对能量值计算方法具体为：

利用Teager能量算子计算得到各个本征模函数分量在所有采样点处的瞬时相对能量值后，将所有采样点处的瞬时相对能量值进行积分求和，得到各个本征模函数分量的相对能量值。

进一步地，主导振荡模式的筛选原则具体为：

低频振荡信号经多元经验模态分解后，若某个或某几个本征模函数分量的相对能量明显高于其他本征模函数分量相对能量，则认为该本征模函数分量中包含着系统的主导振荡模式。其他能量相对较低的本征模函数则被认为包含着系统的虚假振荡模式。

综上所述，本发明实施例通过上述步骤101-步骤104实现了基于广域量测信息的电力系统主导振荡模式的快速辨识；可为电网运行调度人员提供更加快速、准确的系统主导振荡模式评估信息，有利于电网运行调度人员及时掌握系统的动态稳定性。

实施例2：

下面结合具体的计算公式、实例对实施例1中的方案进行进一步地介绍，详见下文描述：

201：从广域量测系统，本领域技术人员公知的技术术语，在此不做赘述，中获取电力系统的广域量测信息，将多机系统发生低频振荡后的各区域间联络线上的有功功率振荡信号作为待辨识的多通道输入信号，对其进行多元经验模态分解，包括：

1)将多通道信号在多维空间中沿多个方向进行一次投影，其投影方向向量集的选取方式为等角度采样，由于高维空间具有不可视性，本节以三维空间为例，将多通道信号的投影方向向量在三维空间中的分布等价地表示其在单位二维超球面上点的分布，并对二维超球面上采样点的选取过程进行详细说明。

如图2所示，在三维坐标系中，多通道信号的方向为z轴正向，单位二维超球面的球心位于坐标系原点o，点q是单位二维超球面与xoy平面交线上的一点，设线段oq与x轴的夹角为φ。保持夹角φ的大小不变，将线段oq绕原点旋转，方向由z轴正向旋转至z轴负向，可以得到单位超球面上的一条经线，设旋转过程中oq与z轴的夹角为θ，在夹角θ变化的同时，以等角度差进行采样，便可得到该条经线上的均匀采样点集，再将这条经线绕z轴转动一周，在夹角φ的变化的同时，以等角度差进行采样，即可得到超球面的均匀采样点集。采样点集在超球面的分布情况示意图参见图3。

多通道信号沿各个投影方向向量的一次映射记作

e^qθx(t)(e^qθ)^* (1)

上式中：x(t)为多通道输入信号，t为时间，θ＝kπ/K，k＝1，…K，其中K是φ在变化范围内均匀采样的个数，e是变量的名称。

2)将一次映射统一沿z轴正向进行二次投影，设z轴正方向单位向量为i，得到的二次映射记作：

上式中，φ_n＝nπ/N，n＝1，…Nθ_k＝kπ/K，k＝1，…K，，其中K和N分别是θ和φ在变化范围内，均匀采样的个数，x(t)为多通道输入信号，t为时间，e是变量的名称，P是变量的名称。

3)找到所有二次映射的全部极值点

并通过多元样条插值法拟合出各个映射的包络线，记作

计算所有包络线的均值，得到多通道信号局部均值m(t)，即：

上式中，φ_n＝nπ/N，n＝1，…Nθ_k＝kπ/K，k＝1，…K，K和N分别是θ和φ在变化范围内，均匀采样的个数，局部均值m(t)，e是变量的名称。

4)将多通道信号与局部均值作差，得到差值d(t)，即：

x(t)-m(t)＝d(t) (4)

上式中，局部均值m(t)，初始信号x(t)。

判断差值d(t)是否可以作为一个本征模函数，需要满足以下两个条件：

①整个数据段内，极值点包括极大值和极小值，数目和过零点数目需相等或仅相差一个；

②在任何一点，由极大值确定的上包络与极小值确定的下包络计算出的局部均值为零。

但在实际工况下，绝大多数信号较为复杂，并不能绝对性地满足上述条件，故将上述条件转化为易于实现判别的数量判据，即：

式中N表示信号的采样点个数，P表示分解过程中的迭代次数。根据经验，若SD的取值在0.2至0.3之间时，认为该次迭代的结果可以作为一个本征模函数分量，并将其作为新的输入信号，再次与局部均值作差，重新进行迭代；反之，若SD的取值不在0.2至0.3之间，则将d(t)作为新的输入信号，重新进行迭代，K表示一个变量，k＝0、1、2...N。

重复按照上述过程不断地筛选出本征模函数，直至迭代过程中信号的二次映射接近单调时(信号沿各个方向向量的投影极值点个数小于3个)，认为分解过程不会再产生新的本征模函数，多元经验模态分解的分解过程结束。

202：利用Teager能量算子的快速响应能力，计算各本征模函数分量在采样点处的相对能量值：

对于一个正弦信号，设其离散形式的表达式为

x_n＝A cos(Ω_n+φ) (6)

式中n为信号采样点的总个数，A表示信号的幅值，Ω_n为信号在采样点处的瞬时频率，φ为初相；则该信号在采样点处Teager能量算子的表达式为：

上式中：ψ(x_n)为信号在采样点处的能量值，A表示信号的幅值，Ω为信号在采样点处的瞬时频率，x_n、x_n-1、x_n+1为信号在采样点n、n-1、n+1处的值。

然而在实际工程应用中，需要处理的信号通常不是标准的正弦信号，设一个频率及幅值均随时间变化的信号为

x_n＝a_ncosφ_n (8)

上式中：x_n为信号在采样点n处的值，a_n为信号幅值函数，φ_n为相位函数。

其瞬时频率可表示为

上式中：n是时间刻度，与时间t等价，因为n是离散形式的时间刻度，

就是对相位函数φ_n求时间的导数。

在信号中选择3个相邻的点来构造一个新的对称差分信号

上式中：a_n为信号幅值函数，φ_n为相位函数，Ω为瞬时频率，s_n为n处的能量算子。

则s_n的Teager能量算子可以表示为

203：计算本征模函数相对能量值并筛选出主导振荡模式：

1)利用Teager能量算子计算得到各个本征模函数分量在所有采样点处的瞬时相对能量值后，将所有采样点处的瞬时相对能量值进行积分求和，得到各个本征模函数分量的相对能量值。即

式中n为各本征模函数分量信号的采样点个数，x_n为信号在采样点n处的值，E为IMF 能量总和。

2)所述主导振荡模式的筛选原则具体为：

低频振荡信号经多元经验模态分解后，若某个或某几个本征模函数分量的相对能量明显高于其他本征模函数分量相对能量，则认为该本征模函数分量中包含着系统的主导振荡模式。其他能量相对较低的本征模函数则被认为包含着系统的虚假振荡模式。

204：基于预测误差法的振荡模式参数计算：

1)将振荡信号经多元经验模态分解后得到的本征模函数信号分量写成如下形式：

其中，u_t为输入变量，y_t为输出变量，e_t为随机噪声，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数，

为对任意的时间t。

对于给定的n_a，n_b＞0，有：

其中

是固定但未知的待辨识参数向量，记为θ，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数；

从而信号也可等价地表示为如下形式：

上式中：y_t为输出变量，e_t为随机噪声，其中

是固定但未知的待辨识参数向量，记为θ。

为对任意的时间t。

其中：

上式中：u_t为输入变量，y_t为输出变量，

为对任意的时间t。

因此，想要求得各本征模函数分量所对应的振荡模式，可以通过计算出参数向量θ中的具体数值，然后建立如下特征方程：

上式中：q为一个常数，

是固定但未知的待辨识参数向量。

计算满足此方程的所有特征根λ_i＝δ_i+jω_i，便可得到每一对共轭复根所对应的系统振荡模式。

其中振荡阻尼比与δ_i相关，具体计算方式为：

振荡频率与ω_i相关，具体计算方式为：

上式中：ω_i特征根的虚部，δ_i特征根的实部，ξ为阻尼比，f为振荡频率。

2)预测误差法估计参数的过程具体为：

预测误差法在确定待辨识系统的模型之后，需要对待辨识的参数θ建立一个损失函数，记作J(θ)，然后运用数值解法对J(θ)进行极小化求解，从而得到待辨识参数θ的具体数值。因此，辨识过程可以抽象成对参数向量θ进行最优估计的过程，即

首先设定初始值θ⁽⁰⁾，参数向量满足如下的一般递归形式：

θ^(k+1)＝θ^(k)+Δθ^(k)＝θ^(k)+γb(J，θ^(k)) (22)

式中b(J，θ^(k))是一个修正因子，它可以在每一步迭代中使之前估计值得到轻微改善，使估计值在反复迭代过程中逐渐达到最优。b(J，θ^(k))是一个修正因子，b表示一个函数，关于 J，θ^(k)的函数关系；γ为一个系数。

本发明方法采用牛顿-拉夫逊法进行迭代，使得预测误差法在迭代过程中的修正因子 b(J，θ^(k))可以通过计算J(θ)在当前时刻参数估值θ^(k)的二次逼近来得到，该方法具体迭代形式如下：

θ^(k+1)＝θ^(k)-α_k(V″_n(θ^(k)))^-1V′_n(θ^(k)) (23)

其中，α_k代表步长，是一个关于k的递减函数；V′_n(θ^(k))代表损失函数J(θ)在θ^(k)处的梯度；V′_n(θ^(k))是损失函数J(θ)在θ^(k)处的海森矩阵。

在整个迭代过程中，对于预先设定的精度ε₀，当|Δθ^(k)|≤ε₀时，停止计算，即得最终的参数估计值θ^(k+1)；否则用所得的θ^(k+1)来代替θ^(k)，重复计算，直到|Δθ^(k)|满足精度要求为止。

实施例3：

下面结合具体的实例、图4-图10、以及表1-表2，对实施例1和2中的方案进行可行性验证，详见下文描述：

本实例是以16机68节点测试系统的主导振荡模式辨识为例，验证本发明实施例1和2 的有效性，16机68节点测试系统如图4所示。

本实例设置的故障为节点1与节点2之间发生区域间的三相短路故障，故障发生时间为仿真开始之后的0.1s，故障持续时间0.1s，0.1s后节点1侧断路器跳开，0.12s后节点2侧断路器跳开。本实例选取68节点测试系统中各个区域间的联络线上的有功功率作为辨识对象，采样频率为100Hz，仿真总时长60s。

图5是发生三相短路后，1-2，8-9，46-49，1-47，52-42，42-41这6个区域间联络线有功功率振荡曲线。以图5中5s-25s的区域间联络线上的有功功率振荡信号作为本实例的输入，经过归一化处理之后，对图5中6个联络线有功功率振荡信号进行多元经验模态分解。

本实例以联络线1-2为例，在图6中依次给出了其有功功率振荡曲线通过传统的经验模态分解(a)和多元经验模态分解(b)后得到的前5个IMF分量波形图：

对比图6中的两列波形可以看出：相较于传统的经验模态分解，多元经验模态分解对振荡信号的模式尺度分离效果更好，可以更加有效地抑制模式混叠现象的产生。

图7-图9还分别对各个有功功率信号经多元经验模态分解后IMF1，IMF2，IMF3的波形进行了对比：

利用Teager能量算子计算各个IMF分量在采样点处的瞬时能量，进而通过对瞬时能量进行积分的方式得到各IMF的相对能量值。同样以联络线1-2的有功功率为例，表1中给出了该线路的有功功率信号经多元经验模态分解后，前5个IMF分量的相对能量值。

为了能够直观地体现出有功功率信号分解得到的各个IMF之间的相对能量差异，图10 绘制了各个IMF相对能量图。

由图10可知：该线路的有功功率信号经多元经验模态分解得到的前5个IMF中，IMF1 的相对能量明显高于其他IMF，故该IMF被认为含有系统的主导振荡模式，剩余的IMF由于相对能量较低，被认为对应着系统的虚假振荡模式。因此在辨识节点联络线1-2的振荡模式时，只对IMF1进行计算。其他联络线的主导振荡模式的确定方式与该联络线相同。

利用预测误差法对每条联络线中筛选出的IMF建立预测误差模型，并通过牛顿-拉夫逊迭代法计算出它们所对应的主导振荡模式的振荡频率及阻尼比，从而完成整个辨识过程。

表2给出特征值分析结果作为参考，同时将本发明所提方法和传统的经验模态分解辨识结果与之进行对比。

通过观察表2中的数据可知，本次故障激发出了模式4和模式14两种低频振荡模式。相较于传统的经验模态分解，利用本专利所提方法得到的振荡模式辨识结果与参考结果之间的误差更小，计算结果更为精确。

本发明实施例对各器件的型号除做特殊说明的以外，其他器件的型号不做限制，只要能完成上述功能的器件均可。

本领域技术人员可以理解附图只是一个优选实施例的示意图，上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

表1各IMF的相对能量

E<sub>1</sub>	E<sub>2</sub>	E<sub>3</sub>	E<sub>4</sub>	E<sub>5</sub>
					0.177	0.0588	0.0317	2.06e-04	2.14e-04

表2不同方法的辨识结果对比

Claims (10)

1.一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述方法包括：

1)通过基于多元经验模态分解的改进算法对电力系统中的初始振荡信号进行分解与处理，得到多个包含系统不同振荡模式的本征模函数分量；

2)利用Teager能量算子的快速响应能力，计算各本征模函数分量在采样点处的相对能量值；

3)对各本征模函数分量在所有采样点处的相对能量值进行积分求和，得到整个本征模函数分量的相对能量，并以能量为判断依据，筛选出能够反映系统真实振荡情况的主导振荡模式并将虚假噪声部分剔除；

4)通过预测误差法计算主导振荡模式所对应的振荡模式参数，即频率及阻尼比。

2.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，在步骤1)之前，所述方法还包括：

从广域量测系统中获取电力系统的广域量测信息，在广域量测信息中提取出电力系统区域间联络线的有功功率信号作为初始振荡信号。

3.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述多元经验模态分解对于初始信号的处理过程具体为：

1)将初始信号x(t)在多维空间中沿多个方向进行一次投影，得到多个一次映射；

2)将多个一次映射沿初始信号方向进行二次投影，得到多个二次映射；

3)通过多元样条插值法拟合出各个映射的包络线，并计算所有包络线的均值得到初始输入信号的局部均值m(t)；

4)利用局部均值m(t)不断对初始信号x(t)进行分解，依次得到多个本征模函数分量。

4.根据权利要求3所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方

法，其特征在于，所述一次投影时投影方向向量的选取方式具体为：

等角度采样，即在多维空间中建立单位超球面，分别以等角度差在超球面上建立类似于经线和纬线的曲线，将所有超球面的球心与超球面上曲线交点的连线作为一次投影的单位方向向量。

5.根据权利要求3所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述本征模函数的获取方法具体为：

将多通道信号与局部均值作差，得到差值d(t)，即：

x(t)-m(t)＝d(t)

上式中：局部均值m(t)，初始信号x(t)，判断差值d(t)；

判断差值d(t)是否可以作为一个本征模函数，需要满足以下两个条件：

①整个数据段内，极值点数目和过零点数目需相等或仅相差一个；

所述极值点包括极大值和极小值；

②在任何一点，由极大值确定的上包络与极小值确定的下包络计算出的局部均值为零；

但在实际工况下，绝大多数信号较为复杂，并不能绝对性地满足上述条件，故将上述条件转化为易于实现判别的数量判据，即：

式中N表示信号的采样点个数，P表示分解过程中的迭代次数；根据经验，若SD的取值在0.2至0.3之间时，认为该次迭代的结果作为一个本征模函数分量，并将其作为新的输入信号，再次与局部均值作差，重新进行迭代；反之，若SD的取值不在0.2至0.3之间，则将d(t)作为新的输入信号，重新进行迭代，K表示一个变量，k＝0、1、2…N；

重复按照上述过程不断地筛选出本征模函数，直至迭代过程中信号的二次映射接近单调时，即，信号沿各个方向向量的投影极值点个数小于3个，认为分解过程不会再产生新的本征模函数，多元经验模态分解的分解过程结束。

6.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述能量算子的计算方法具体为：

对于一个正弦信号，设其离散形式的表达式为

x_n＝Acos(Ω_n+φ)

式中n为信号采样点的总个数，A表示信号的幅值，Ω_n为信号在采样点处的瞬时频率，φ为初相，则该信号在采样点处Teager能量算子的表达式为：

上式中：ψ(x_n)为信号在采样点处的能量值，A表示信号的幅值，Ω为信号在采样点处的瞬时频率，x_n、x_n-1、x_n+1为信号在采样点n、n-1、n+1处的值；

然而在实际工程应用中，需要处理的信号通常不是标准的正弦信号，设一个频率及幅值均随时间变化的信号为：

x_n＝a_n cosφ_n

上式中：x_n为信号在采样点n处的值，a_n为信号幅值函数，φ_n为相位函数；

其瞬时频率可表示为：

上式中：n是时间刻度，与时间t等价，因为n是离散形式的时间刻度，

就是对相位函数φ_n求时间的导数；

在信号中选择3个相邻的点来构造一个新的对称差分信号：

上式中：a_n为信号幅值函数，φ_n为相位函数，Ω为瞬时频率，s_n为n处的能量算子；

则s_n的Teager能量算子可以表示为：

7.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述本征模函数的相对能量值计算方法具体为：

利用Teager能量算子计算得到各个本征模函数分量在所有采样点处的瞬时相对能量值后，将所有采样点处的瞬时相对能量值进行积分求和，得到各个本征模函数分量的相对能量值；即

式中n为各本征模函数分量信号的采样点个数，x_n为信号在采样点n处的值，E为IMF能量总和。

8.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述主导振荡模式的筛选原则具体为：

低频振荡信号经多元经验模态分解后，若某个或某几个本征模函数分量的相对能量明显高于其他本征模函数分量相对能量，则认为该本征模函数分量中包含着系统的主导振荡模式；其他能量相对较低的本征模函数则被认为包含着系统的虚假振荡模式。

9.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述低频振荡模式参数的计算方式具体为：

将振荡信号经多元经验模态分解后得到的本征模函数信号分量写成如下形式：

其中，u_t为输入变量，y_t为输出变量，e_t为随机噪声，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数，

为对任意的时间t；

对于给定的n_a，n_b＞0，有：

其中

是固定但未知的待辨识参数向量，记为θ，A为输出变量的幅值常数，B为输入变量的幅值常数，q为一个常数；

从而信号也可等价地表示为如下形式：

上式中：y_t为输出变量，e_t为随机噪声，其中

是固定但未知的待辨识参数向量，记为

为对任意的时间t；

其中

上式中：u_t为输入变量，y_t为输出变量，

为对任意的时间t；

因此，想要求得各本征模函数分量所对应的振荡模式，通过计算出参数向量θ中的具体数值，然后建立如下特征方程：

上式中：q为一个常数，

是固定但未知的待辨识参数向量；

计算满足此方程的所有特征根λ_i＝δ_i+jω_i，便可得到每一对共轭复根所对应的系统振荡模式；

其中振荡阻尼比与δ_i相关，具体计算方式为：

振荡频率与ω_i相关，具体计算方式为：

上式中：ω_i特征根的虚部，δ_i特征根的实部，ξ为阻尼比，f为振荡频率。

10.根据权利要求1所述的一种基于多元经验模态分解的电力系统低频振荡模式辨识方法，其特征在于，所述预测误差法估计参数的过程具体为：

预测误差法在确定待辨识系统的模型之后，需要对待辨识的参数θ建立一个损失函数，记作J(θ),然后运用数值解法对J(θ)进行极小化求解，从而得到待辨识参数θ的具体数值；因此，辨识过程可以抽象成对参数向量θ进行最优估计的过程，即

首先设定初始值θ⁽⁰⁾，参数向量满足如下的一般递归形式：

θ^(k+1)＝θ^(k)+Δθ^(k)＝θ^(k)+γb(J,θ^(k))

式中b(J,θ^(k))是一个修正因子，它能够在每一步迭代中使之前估计值得到轻微改善，使估计值在反复迭代过程中逐渐达到最优；b(J,θ^(k))是一个修正因子，b表示一个函数，关于J,θ^(k)的函数关系；γ为一个系数；

采用牛顿—拉夫逊法进行迭代，使得预测误差法在迭代过程中的修正因子b(J,θ^(k))可以通过计算J(θ)在当前时刻参数估值θ^(k)的二次逼近来得到，迭代形式如下：

θ^(k+1)＝θ^(k)-α_k(V″_n(θ^(k)))^-1V′_n(θ^(k))

其中，α_k代表步长，是一个关于k的递减函数；V′_n(θ^(k))代表损失函数J(θ)在θ^(k)处的梯度；V′_n(θ^(k))是损失函数J(θ)在θ^(k)处的海森矩阵；

在整个迭代过程中，对于预先设定的精度ε₀，当|Δθ^(k)|≤ε₀时，停止计算，即得最终的参数估计值θ^(k+1)；否则用所得的θ^(k+1)来代替θ^(k)，重复计算，直到|Δθ^(k)|满足精度要求为止。

Images