多元变分模态分解的电力系统强迫振荡源定位方法及装置

1.本发明涉及电力系统领域，尤其涉及一种多元变分模态分解(multivariate variational mode decomposition,mvmd)的电力系统强迫振荡源定位方法及装置。


背景技术：

2.电力系统的强迫振荡是一种由持续周期性扰动引起的功率振荡现象
1.。当电网强迫振荡频率与系统固有振荡频率接近时，将导致电网功率大幅振荡，带来严重的安全稳定事故
[2]-[4]
。已有研究表明：电网中可诱发强迫振荡的扰动源种类繁多，除传统的励磁系统故障、原动机油门异动等，风电、光伏等具有强随机性的新能源接入也可诱发强迫振荡
[5]
。由于强迫振荡具有起振快、失去振荡源后快速平息的特点，实际工程常用抑制电网强迫振荡的措施是快速定位振荡源并对其进行隔离。因此，快速、准确定位强迫振荡源是抑制电网强迫振荡的首要任务
[6]
。
[0003]
目前，随着同步相量测量单元(synchrophasor measurement unit，pmu)在电力系统中的广泛配置，使得对电网的实时动态监测变成了可能。在此背景下，耗散能量流(dissipation energy flow，def)方法作为能量法的一种，凭借其具有计算精度高，耗时短、完全数据驱动性的特点，已广泛应用于电力系统强迫振荡源定位。文献
[7]
首次提出了def的计算方法，该方法仅需获取发电机pmu的量测数据，即可精准定位振荡源。
[0004]
然而，在应用def方法时需对广域量测信息进行处理，以避免无关信息影响振荡源定位的精度。因此，如何挖掘出广域量测信息中的强迫振荡信息量是现阶段振荡源定位亟需解决的重点和难点。文献
[8]
采用带通滤波(band-pass filter,bf)方法提取出强迫振荡时域分量计算耗散能量。文献
[9]
基于同步压缩变换(synchrosqueezing wavelet transform,swt)实现了强迫振荡信号的提取-重构，提高了强迫振荡源定位的精度。文献
[10]
采用连续小波变换(continuous wavelet transform,cwt)处理广域量测信息，利用提取出的表征强迫振荡模式的小波系数定位强迫振荡源。文献
[11]
采用短时傅里叶变换(short-time fourier transform，stft)处理信号并计算发电机耗散能量谱，将耗散能量概念从时域转换到频域，提高了计算效率。文献
[12]
采用多元经验模态分解同步提取出多通道广域量测信息中的强迫振荡本征模函数(intrinsic mode functions，imf)分量，有效提升了计算精度。
[0005]
上述信号处理方法与def结合后均取得了较好的定位效果。然而，随着电网互联规模的不断扩大，系统日趋复杂，在实际采集信号时，广域量测系统(wide area measurement system，wams)总站采集的信号是各pmu子站量测信息的汇总，各通道采集的信号包含的振荡信息量不尽相同，具有很强的多通道特性。如何在兼顾准确性与快速性的同时，从众多wams数据中提取出隐含的强迫振荡模式分量是现阶段强迫振荡源定位亟需解决的问题。
[0006]
参考文献
[0007]
[1]follum j，pierre j w.detection of periodic forced oscillations in power systems[j].ieee transactions on power systems，2016，31(3)：1-11.
[0008]
[2]chen lei，min yong，chen yiping，et al.evaluation of generator damping using oscillation energy dissipation and the connection with modal analysis[j].ieee transactions on power systems,2014，29(3)：1393-1402.
[0009]
[3]kosterev d n，taylor c w，mittelstadt w a.model validation for the august 10，1996wscc system outage[j].ieee transactions on power systems，1999，14(3)：967-979.
[0010]
[4]ghorbaniparvar m.survey on forced oscillations in power system[j].journal of modern power systems&clean energy，2017，5(5)：671-682.
[0011]
[5]栾某德.电力系统强迫振荡扰动源定位研究[d].浙江：浙江大学，2020.
[0012]
[6]曹娜，赵旭，于群.基于暂态能量流的双馈风电机组强迫振荡源定位[j].电力系统自动化，2020，44(10)：103-110.
[0013]
[7]陈磊，闵勇，胡伟.基于振荡能量的低频振荡分析与振荡源定位(一)理论基础与能量流计算[j].电力系统自动化，2012，36(3)：22-27.
[0014]
[8]chen lei，min yong，chen yiping，et al.online monitoring of generator damping using dissipation energy flow computed from ambient data[j].iet generation transmission&distribution，2017，11(18)：4430-4435.
[0015]
[9]estevez p g，marchi p，galarza c，et al.non-stationary power system forced oscillation analysis using synchrosqueezing transform[j].ieee transactions on power systems，2020，36(2)：1583-1593.
[0016]
[10]姜涛，李孟豪，李雪，等.电力系统强迫振荡源的时频域定位方法[j].电力系统自动化，2021，45(9)：98-106.
[0017]
[11]姜涛，高浛，李雪，等.电力系统强迫振荡源定位的耗散能量谱方法[j/ol].中国电机工程学报：1-12[2022-03-01].http://kns.cnki.net/kcms/detail/11.2107.tm.20220224.1015.003.html.
[0018]
[12]姜涛，刘博涵，李雪，等.基于多元经验模态分解的电力系统强迫振荡源定位[j/ol].中国电机工程学报：1-11[2022-01-11].http://kns.cnki.net/kcms/detail/11.2107.tm.20211214.1555.006.html.


技术实现要素：

[0019]
本发明提供了一种多元变分模态分解的电力系统强迫振荡源定位方法及装置，本发明无需系统详细模型，仅需发电机广域量测信息即可实现强迫振荡源定位；本方法作为一种新型的多通道完全非递归的模态变分方法，分解稳定性高，可实现多通道量测信息的同步处理，有效提升强迫振荡源定位的精度和效率，详见下文描述：
[0020]
第一方面，一种多元变分模态分解的电力系统强迫振荡源定位方法，所述方法包括：
[0021]
设置mvmd分解模态数量k及惩罚因子a；构建多通道广域量测信息矩阵xi的mvmd变分约束模型；
[0022]
构建增广拉格朗日函数并求解变分约束模型，得到广域量测信息矩阵xi的频域模态矩阵，进而由傅里叶逆变换得到时域imf矩阵；
[0023]
分别计算各imf分量的teager能量算子及能量权重，得到能量权重矩阵ξi；筛选出能量权重最高的imf分量作为强迫振荡分量予以保留；
[0024]
利用筛选得到的强迫振荡分量imf，计算各发电机的耗散能量流，并依据强迫振荡源判据定位系统的强迫振荡源。
[0025]
其中，所述构建多通道广域量测信息矩阵xi的mvmd变分约束模型为：
[0026][0027]
其中，x
k,n
(t)为第n个通道中第k个imf分量，x
k.δp
为有功功率变化量第k个imf分量，x
δp
(t)为有功功率变化量的量测信息，x
k.δq
为无功功率变化量第k个imf分量，x
δq
(t)为无功功率变化量的量测信息，x
k.δlnu
为电压幅值对数变化量第k个imf分量，x
δlnu
(t)为电压幅值对数变化量的量测信息，x
k.δf
(t)为频率变化量第k个imf分量，x
δf
(t)为频率变化量的量测信息。
[0028]
其中，所述构建增广拉格朗日函数并求解变分约束模型为：
[0029]
构建的增广拉格朗日函数如下：
[0030][0031]
式中：λn(t)为拉格朗日乘子；a为惩罚因子；《
·
，
·
》为求取内积，xn(t)为第n个通道的量测信息；
[0032]
采用admm求解变分模型，得估计的频域内模态，即：
[0033][0034]
其中，为第n个通道量测信息的频域分量，为拉格朗日乘子，ω为角频率；获取分解的第k个多通道imf分量的中心频率为：
[0035][0036]
式中，为时域模态的傅里叶变换；
[0037]
获取发电机i的多通道电压幅值对数、频率、有功和无功变化量xi经mvmd分解后的频域模态矩阵，对进行傅里叶逆变换并取实部得到对应的imf分量，对各频域模态进行傅里叶逆变换并取实部，得各时域imf分量及残差。
[0038]
其中，所述分别计算各imf分量的teager能量算子及能量权重，得到能量权重矩阵ξi；
[0039]
筛选出能量权重最高的imf分量作为强迫振荡分量予以保留具体为：
[0040]
能量矩阵ξi为：
[0041][0042]
式中，ξ
k,δp
、ξ
k,δq
、ξ
k,δlnu
、ξ
kδf
分别为有功功率变化量、无功功率变化量、电压幅值对数变化量和频率变化量的第k个imf分量对应的能量权重；
[0043]
若矩阵ξi中各电气量的第k行的能量权重均为最高，则强迫振荡模式位于第k个imf分量中，将xi各电气量第k个imf分量作为强迫振荡分量提取出来，即：
[0044][0045]
其中，imf
ik,δp
(τ)为发电机i有功功率变化量第k个imf分量的第τ个元素，imf
ik,δq
(τ)为发电机i无功功率变化量第k个imf分量的第τ个元素，imf
ik,δlnu
(τ)为发电机i电压幅值对数变化量第k个imf分量的第τ个元素，imf
ik,δf
(τ)为发电机i频率变化量第k个imf分量的第τ个元素。
[0046]
进一步地，所述利用筛选得到的强迫振荡分量imf，计算各发电机的耗散能量流为：
[0047]
耗散能量流w
dp,i
为：
[0048]wdp
＝w
dp,p,f
w
dp,q,u
[0049]
＝∫(2πimf
ik,δp
imf
ik,δf
dt imf
ik,δq
d(imf
ik,δlnu
))
[0050]
其中，imf
ik,δp
为发电机i有功功率变化量第k个imf分量，imf
ik,δf
为发电机i频率变化量第k个imf分量，imf
ik,δq
为发电机i无功功率变化量第k个imf分量，imf
ik,δlnu
为电压幅值对数变化量第k个imf分量。
[0051]
第二方面、一种多元变分模态分解的电力系统强迫振荡源定位装置，所述装置包括：处理器和存储器，所述存储器中存储有程序指令，所述处理器调用存储器中存储的程序指令以使装置执行第一方面中的任一项所述的方法步骤。
[0052]
第三方面、一种计算机可读存储介质，其特征在于，所述计算机可读存储介质存储有计算机程序，所述计算机程序包括程序指令，所述程序指令被处理器执行时使所述处理器执行第一方面中的任一项所述的方法步骤。
[0053]
本发明提供的技术方案的有益效果是：
[0054]
1、本发明定位强迫振荡源时无需系统详细模型，仅根据发电机广域量测信息即可实现强迫振荡源的准确、有效定位，具有较强的工程实用性；
[0055]
2、本发明所提基于mvmd的强迫振荡源定位方法通过同步提取多通道广域量测信息中的强迫振荡imf分量进而计算耗散能量流，可准确定位电力系统的强迫振荡源，有效提高定位精度和效率；
[0056]
3、本发明能够有效兼顾各量测通道信息的相关性，实现各通道间强迫振荡分量的同步提取，提高了现有方法的辨识精度和并行能力；且所提方法对多通道量测信息分解出
的imf分量带宽跨度更小，在强迫振荡频率处的模态更为集中，有效避免了模态混叠现象；4、本方法在较强噪声水平下仍可有效定位出强迫振荡源，具有良好的鲁棒性。
附图说明
[0057]
图1为基于mvmd的强迫振荡源定位方法的流程图；
[0058]
图2为wecc(western electricity coordinating council，西部电力协调委员会)-179节点系统的拓扑示意图；
[0059]
图3为各发电机的机端电气参数图；
[0060]
图4为发电机有功功率和频率imf的分量图；
[0061]
图5为各发电机的imf分量的能量权重图；
[0062]
图6为各发电机的强迫振荡模式imf1分量图；
[0063]
图7为基于mvmd的发电机耗散能量曲线图；
[0064]
图8为基于vmd的发电机耗散能量曲线图；
[0065]
图9为基于memd的发电机耗散能量曲线图。
具体实施方式
[0066]
为使本发明的目的、技术方案和优点更加清楚，下面对本发明实施方式作进一步地详细描述。
[0067]
强迫振荡是威胁电力系统安全稳定运行的重要因素之一，准确、高效定位强迫振荡源是抑制强迫振荡的关键，但现有基于广域量测信息的强迫振荡源定位方法实现较为复杂且效率和精度有待提升。为此，本发明实施例提出了一种基于多元变分模态分解的电力系统强迫振荡源定位方法及装置。
[0068]
该方法首先对电力系统多通道广域量测信息构建变分约束模型，基于变分约束模型进而构建增广拉格朗日函数求解变分约束模型，将多通道量测信息分解出一系列含不同振荡模式的imf分量；然后利用teager能量法筛选出含强迫振荡模式的imf分量，进而利用基于mvmd的耗散能量流计算式，计算各发电机耗散能量流，最后实现强迫振荡源定位。
[0069]
实施例1
[0070]
本发明实施例提供了一种基于多元变分模态分解的电力系统强迫振荡源定位方法，参见图1，下面对本发明实施方式作进一步详细描述，其中101-107为计算方法的详细步骤，详见下文描述：
[0071]
101：从pmu中获取各发电机电压幅值对数、频率、有功功率和无功功率量测信息，以发电机为单位，构建多通道广域量测信息矩阵xi；
[0072]
102：设置mvmd分解模态数量k及惩罚因子a；
[0073]
103：构建多通道广域量测信息矩阵xi的mvmd变分约束模型；
[0074]
104：构建增广拉格朗日函数并求解变分约束模型，得到广域量测信息矩阵xi的频域模态矩阵，进而由傅里叶逆变换得到时域imf矩阵；
[0075]
105：分别计算各imf分量的teager能量算子(teager energy operator，teo)及能量权重，得到能量权重矩阵ξi；
[0076]
106：筛选出能量权重最高的imf分量作为强迫振荡分量予以保留，忽略其余分量；
[0077]
107：利用筛选得到的强迫振荡分量imf，计算各发电机的耗散能量流，并依据强迫振荡源判据定位系统的强迫振荡源。
[0078]
综上所述，以上步骤101-107为基于多元变分模态分解的电力系统强迫振荡源定位方法的详细步骤。通过上述步骤能够准确甄别出广域量测信息中的强迫振荡模式分量，并提升振荡源的定位效率，进而实现电力系统强迫振荡源的快速、准确定位。
[0079]
实施例2
[0080]
下面结合具体的计算公式、实例对实例1中的方案进行进一步地介绍，详见下文描述：
[0081]
201：选取电网中发电机节点pmu采集发电机的母线频率信号、机端电压对数信号、有功功率信号和无功功率信号(实际中可选取受强迫振荡影响较大的机组)，将其作为本发明实施例的输入数据，分别构建多通道量测信息矩阵xi；
[0082]
202：采用mvmd分解多通道量测信息矩阵xi；
[0083]
其中，采用mvmd分解电力系统多通道量测信息计算过程如下：
[0084]
采用mvmd从电力系统的多通道广域量测信息中分解出表征系统强迫振荡模式的imf分量时，首先需以发电机为单元从pmu中获取各发电机有功功率、无功功率、电压幅值和频率的广域量测信息，并计算电压幅值的对数值，进一步计算上述各电气量相对稳态值的变化量，进而构建含发电机有功功率变化量、无功功率变化量、电压幅值对数变化量和频率变化量的多通道广域量测信息矩阵xi：
[0085]
xi＝[x
i,δp x
i,δq x
i,δlnu x
i,δf
]
ꢀꢀꢀ
(1)
[0086]
式中，x
i,δp
、x
i,δq
、x
i,δlnu
、x
i,δf
分别为发电机i的有功功率、无功功率、电压幅值对数和频率的变化量。
[0087]
为了有效提取出多通道量测信息矩阵x的多通道振荡模式imf分量xk(t)，并使得提取后的imf分量能够完整重构出原始信号且各imf分量的带宽之和最小。
[0088]
为此，xk(t)的带宽可通过梯度函数l2平方范数估计得到，即：
[0089][0090]
式中，表示对时间求偏导；ωk为mvmd提取出第k个模态的中心频率。为xk(t)的解析调制信号。
[0091]
为避免单一频率分量ωk应用于导致混频，因此需在多通道振荡模式imf分量xk(t)中找到各通道间的公共频率分量。可通过将的各通道单边频谱以中心频率ωk平移，并采用frobenius范数进一步估计的带宽，则有：
[0092][0093]
式中，为通道n中第k个振荡模态的解析调制信号，k为振荡模态，n为量测信息通道。
[0094]
设定待分解的imf数量为k，进而构建发电机i的多通道广域量测信息矩阵xi的mvmd变分约束模型：
[0095][0096]
其中，x
k,n
(t)为第n个通道中第k个imf分量，x
k.δp
为有功功率变化量第k个imf分量，x
δp
(t)为有功功率变化量的量测信息，x
k.δq
为无功功率变化量第k个imf分量，x
δq
(t)为无功功率变化量的量测信息，x
k.δlnu
为电压幅值对数变化量第k个imf分量，x
δlnu
(t)为电压幅值对数变化量的量测信息，x
k.δf
(t)为频率变化量第k个imf分量，x
δf
(t)为频率变化量的量测信息。
[0097]
对于式(4)所示的多通道量测信息的变分问题，构建的增广拉格朗日函数如下：
[0098][0099]
式中：λn(t)为拉格朗日乘子；a为惩罚因子；《
·
，
·
》为求取内积，xn(t)为第n个通道的量测信息。
[0100]
进一步采用admm求解变分模型，可得估计的频域内模态，即：
[0101][0102]
其中，为第n个通道量测信息的频域分量，为拉格朗日乘子，ω为角频率。
[0103]
对应的分解的第k个多通道imf分量的中心频率为：
[0104][0105]
式中，为时域模态的傅里叶变换。由式(6)、(7)可求得多通道量测信息经mvmd分解得到的频域模态及其对应的中心频率。
[0106]
则发电机i的多通道电压幅值对数、频率、有功和无功变化量xi经mvmd分解后的频域模态矩阵为：
[0107][0108]
式中，x
ik,δp
(ωk)、x
ik,δq
(ωk)、x
ik,δlnu
(ωk)、x
ik,δf
(ωk)分别为发电机i的有功功率变化量、无功功率变化量、电压幅值对数变化量和频率变化量的第k个频域模态；r
i,δp
(ω)、r
i,δq
(ω)、r
i,δlnu
(ω)、r
i,δf
(ω)分别为发电机i有功功率变化量、无功功率变化量、电
压幅值对数变化量和频率变化量的残差。
[0109]
进一步，可对进行傅里叶逆变换并取实部得到对应的imf分量，即：
[0110][0111]
分别对式(8)中各频域模态进行傅里叶逆变换并取实部，得各时域imf分量及残差，即：
[0112][0113]
式中，imf
ik,δp
(τ)、imf
ik,δq
(τ)、imf
ik,δlnu
(τ)、imf
ik,δf
(τ)分别为发电机i有功功率、无功功率、电压幅值对数和频率变化量的第k个imf分量的第τ个元素，r
i,δp
(τ)、r
ij,δq
(τ)、r
i,δlnu
(τ)、r
i,δf
(τ)分别为发电机i有功功率、无功功率、电压幅值对数和频率变化量的残差。
[0114]
203：利用teo能量筛选出多通道量测信息的强迫振荡imf分量；
[0115]
mvmd可将多通道电力系统广域量测信息分解为一系列多通道的imf分量集合，如何在众多imf分量中提取出含强迫振荡模式的imf分量是精准定位强迫振荡源的关键。针对含强迫振荡模式的imf分量的振荡能量通常最大，在众多imf分量中所占能量权重最高这一特点，通过计算各imf分量的teo能量，进而计算各imf分量的能量权重，筛选出多通道量测信息的强迫振荡imf分量。
[0116]
对于单通道imf分量x
k,n
(t)，其teo能量可表示为：
[0117][0118]
式中，e
k,n
为电压幅值信号中第n个通道内第k个模态的teo能量。
[0119]
进一步，计算得各imf分量相对能量权重ξ
k,n
：
[0120][0121]
进一步，由式(12)计算出各imf的teo能量权重ξ
k,n
，得到能量矩阵ξi为：
[0122][0123]
式中，ξ
k,δp
、ξ
k,δq
、ξ
k,δlnu
、ξ
kδf
分别为有功功率变化量、无功功率变化量、电压幅值
对数变化量和频率变化量的第k个imf分量对应的能量权重。
[0124]
进一步，若矩阵ξi中各电气量的第k行的能量权重均为最高，则可认为系统的强迫振荡模式位于第k个imf分量中，可将式(10)中xi各电气量第k个imf分量作为强迫振荡分量提取出来，即：
[0125][0126]
其中，imf
ik,δp
(τ)为发电机i有功功率变化量第k个imf分量的第τ个元素，imf
ik,δq
(τ)为发电机i无功功率变化量第k个imf分量的第τ个元素，imf
ik,δlnu
(τ)为发电机i电压幅值对数变化量第k个imf分量的第τ个元素，imf
ik,δf
(τ)为发电机i频率变化量第k个imf分量的第τ个元素。
[0127]
204：基于mvmd的耗散能量流计算；
[0128]
定义节点i流入节点j的耗散能量为：
[0129][0130]
式中，δp
ij
和δq
ij
分别为从节点i流入节点j的有功功率和无功功率的变化量；δθi为节点i电压相角的变化量；δlnui为节点i的电压幅值对数的变化量；δfi为节点i频率的变化量。
[0131]
进一步根据式(14)提取出的各电气量的变化量的强迫振荡imf分量，可得式(15)中有功功率变换量与频率变化量的乘积为：
[0132]wdp,p,f
＝∫2πimf
ik,δp
imf
ik,δf
dt
ꢀꢀꢀ
(16)
[0133]
同理，式(15)中无功功率变换量与电压幅值对数变化量的乘积为：
[0134]wdp,q,u
＝∫imf
ik,δq
d(imf
ik,δlnu
)
ꢀꢀꢀ
(17)
[0135]
根据式(16)和(17)可得基于mvmd分解的发电机i在系统强迫振荡期间的耗散能量流w
dp,i
为：
[0136][0137]
其中，imf
ik,δp
为发电机i有功功率变化量第k个imf分量，imf
ik,δf
为发电机i频率变化量第k个imf分量，imf
ik,δq
为发电机i无功功率变化量第k个imf分量，imf
ik,δlnu
为电压幅值对数变化量第k个imf分量。
[0138]
式(18)即为基于mvmd分解的耗散能量流计算公式，采用基于mvmd分解的强迫振荡imf分量计算各发电机的耗散能量流有效避免了多通道量测信息中冗余信息对强迫振荡源定位精度的干扰。
[0139]
205：基于mvmd的强迫振荡源定位判据。
[0140]
借鉴传统def定位强迫振荡的判据可知基于mvmd分解的强迫振荡源定位判据：若发电机的w
dp,i
《0且在整个振荡期间内变化趋势持续向下，则可认为该发电机持续向电网中注入振荡能量，为强迫振荡源；若发电机的w
dp,i
》0且在整个振荡期间内变化趋势持续向上或保持不变，则可认为该发电机持续从电网中吸收振荡能量，为非强迫振荡源。
[0141]
综上所述，本发明实例通过上述步骤201-205实现了基于mvmd的强迫振荡源定位方法，通过同步提取多通道广域量测信息中的强迫振荡imf分量进而计算耗散能量流，可准确定位电力系统强迫振荡源，与传统方法def方法相比定位精度和效率更高。
[0142]
实施例3
[0143]
下面结合具体实例，wecc-179节点测试系统拓扑结构如图2所示，该系统共有29台发电机，所有发电机阻尼系数均设置为4，系统各振荡模式均具有良好的阻尼。本节对该测试系统设置了如下扰动：向发电机g79的励磁系统持续注入0.46hz的正弦波激发出电网的强迫振荡模式，仿真时长为40s。图3为强迫振荡期间各发电机有功功率、无功功率、电压幅值对数和频率的仿真数据。
[0144]
根据本发明实施例所提的mvmd分解方法，首先计算图3中各发电机0～40s的有功功率、无功功率、电压幅值对数和频率量测数据相对稳态值的变化量，并以发电机为单元，构建各发电机的多通道广域量测信息矩阵xi＝[x
i,δp x
i,δq x
i,δlnu x
i,δf
]，共得到29个量测信息矩阵。前期大量仿真表明，模式分解数量k取值6、惩罚因子取值1000时mvmd具有较好的分解效果，进一步，构建各发电机多通道广域量测信息矩阵xi的mvmd变分约束模型，通过求解变分约束模型得到各电气量的时域imf分量，各发电机的多通道量测信息矩阵xi均被分解为1个6
×
1200
×
4的3维量测信息矩阵，受篇幅限制，图4给出了发电机g79的有功功率变化量和频率变化量原始量测信息及经mvmd分解后的imf分量。
[0145]
由图4可知：发电机g79的多通道量测信息矩阵中δp和δf经mvmd均分解出6个imf分量，各imf分量均包含不同的振荡模式。为从分解的众多imf分量中提取出表征系统强迫振荡模式的imf分量，进一步根据式(12)分别计算各imf分量的能量权重，并构建各发电机能量权重矩阵ξi。为便于分析，图5进一步给出了经mvmd分解后计算得到的各imf分量的teo能量权重。由图5可知，在经mvmd分解的各发电机有功功率变化量、无功功率变化量、电压幅值对数变化量和频率变化量的imf分量中，imf1能量权重均为最高，其余imf对应能量权重均远小于imf1。能量权重最高的imf分量表征系统的主导强迫振荡模式，故可认为imf1为系统的强迫振荡模式imf分量。因此可将各发电机有功功率变化量、无功功率变化量、电压幅值对数变化量及频率变化量的imf1分量作为强迫振荡分量提取出来。
[0146]
图6给出了提取出的各电气量的强迫振荡imf分量，进一步，将图6中的imf1分量代入式(18)计算各发电机基于mvmd的耗散能量流，结果如图7所示。
[0147]
由图7可知，发电机g79的耗散能量流曲线在0-40s内均呈明显下降趋势，由强迫振荡源判据可知：0-40s时段内g79持续向电网中注入振荡能量，为强迫振荡源；除g79外，其余发电机耗散能量流均无明显下降趋势，表明此时段内其余发电机持续吸收电网中振荡能量，为非振荡源。将振荡源定位结果与扰动设置对比可知：所提方法定位结果与注入扰动的发电机完全一致，有效验证了所提方法定位强迫振荡源的准确性。
[0148]
图8进一步给出了采用vmd的强迫振荡源定位结果，由图8可知：采用vmd计算的耗散能量流中除发电机g79呈明显下降趋势外，发电机g36及g144也具有下降趋势，表明此时g79、g36和g144均向电网中注入振荡能量，显然此辨识结果与注入扰动不符，若采用该振荡源定位结果将导致运行人员出现误判进而错误切除强迫振荡源。对比vmd与mvmd的原理可知：所提mvmd定位强迫振荡源方法能够有效避免单通道定位振荡源方法的不足，综合考虑各通道间量测信息的相关性，实现对多通道量测信息中强迫振荡模式的整体提取。
[0149]
图9给出了基于memd的耗散能量曲线。由图9可知：g79的耗散能量呈振荡下行趋势，但在40s量测时间内g79耗散能量仅减少了0.0116，下降趋势并不明显从而易造成误判，而图7中g79耗散能量减少了0.0653。综合对比mvmd和memd的耗散能量下降趋势可知，所提mvmd方法相对memd方法更能完整地提取出了多通道量测信息中的强迫振荡模式分量，在计算耗散能量时有效避免了冗余信息的干扰，具有更高的强迫振荡源定位精度。
[0150]
实施例4
[0151]
基于同一发明构思，本发明实施例还提供了一种多元变分模态分解的电力系统强迫振荡源定位装置，该装置包括：处理器和存储器，存储器中存储有程序指令，处理器调用存储器中存储的程序指令以使装置执行实施例中的以下方法步骤：
[0152]
设置mvmd分解模态数量k及惩罚因子a；构建多通道广域量测信息矩阵xi的mvmd变分约束模型；
[0153]
构建增广拉格朗日函数并求解变分约束模型，得到广域量测信息矩阵xi的频域模态矩阵，进而由傅里叶逆变换得到时域imf矩阵；
[0154]
分别计算各imf分量的teager能量算子及能量权重，得到能量权重矩阵ξi；筛选出能量权重最高的imf分量作为强迫振荡分量予以保留；
[0155]
利用筛选得到的强迫振荡分量imf，计算各发电机的耗散能量流，并依据强迫振荡源判据定位系统的强迫振荡源。
[0156]
其中，构建多通道广域量测信息矩阵xi的mvmd变分约束模型为：
[0157][0158]
其中，x
k,n
(t)为第n个通道中第k个imf分量，x
k.δp
为有功功率变化量第k个imf分量，x
δp
(t)为有功功率变化量的量测信息，x
k.δq
为无功功率变化量第k个imf分量，x
δq
(t)为无功功率变化量的量测信息，x
k.δlnu
为电压幅值对数变化量第k个imf分量，x
δlnu
(t)为电压幅值对数变化量的量测信息，x
k.δf
(t)为频率变化量第k个imf分量，x
δf
(t)为频率变化量的量测信息。
[0159]
其中，构建增广拉格朗日函数并求解变分约束模型为：
[0160]
构建的增广拉格朗日函数如下：
[0161][0162]
式中：λn(t)为拉格朗日乘子；a为惩罚因子；《
·
，
·
》为求取内积，xn(t)为第n个通道的量测信息；
[0163]
采用admm求解变分模型，得估计的频域内模态，即：
[0164]
[0165]
其中，为第n个通道量测信息的频域分量，为拉格朗日乘子，ω为角频率；获取分解的第k个多通道imf分量的中心频率为：
[0166][0167]
式中，为时域模态的傅里叶变换；
[0168]
获取发电机i的多通道电压幅值对数、频率、有功和无功变化量xi经mvmd分解后的频域模态矩阵，对进行傅里叶逆变换并取实部得到对应的imf分量，对各频域模态进行傅里叶逆变换并取实部，得各时域imf分量及残差。
[0169]
其中，分别计算各imf分量的teager能量算子及能量权重，得到能量权重矩阵ξi；筛选出能量权重最高的imf分量作为强迫振荡分量予以保留具体为：
[0170]
能量矩阵ξi为：
[0171][0172]
式中，ξ
k,δp
、ξ
k,δq
、ξ
k,δlnu
、ξ
kδf
分别为有功功率变化量、无功功率变化量、电压幅值对数变化量和频率变化量的第k个imf分量对应的能量权重；
[0173]
若矩阵ξi中各电气量的第k行的能量权重均为最高，则强迫振荡模式位于第k个imf分量中，将xi各电气量第k个imf分量作为强迫振荡分量提取出来，即：
[0174][0175]
其中，imf
ik,δp
(τ)为发电机i有功功率变化量第k个imf分量的第τ个元素，imf
ik,δq
(τ)为发电机i无功功率变化量第k个imf分量的第τ个元素，imf
ik,δlnu
(τ)为发电机i电压幅值对数变化量第k个imf分量的第τ个元素，imf
ik,δf
(τ)为发电机i频率变化量第k个imf分量的第τ个元素。
[0176]
进一步地，利用筛选得到的强迫振荡分量imf，计算各发电机的耗散能量流为：
[0177]
耗散能量流w
dp,i
为：
[0178]wdp
＝w
dp,p,f
w
dp,q,u
[0179]
＝∫(2πimf
ik,δp
imf
ik,δf
dt imf
ik,δq
d(imf
ik,δlnu
))
[0180]
其中，imf
ik,δp
为发电机i有功功率变化量第k个imf分量，imf
ik,δf
为发电机i频率变化量第k个imf分量，imf
ik,δq
为发电机i无功功率变化量第k个imf分量，imf
ik,δlnu
为电压幅值对数变化量第k个imf分量。
[0181]
基于同一发明构思，本发明实施例还提供了一种计算机可读存储介质，存储介质包括存储的程序，在程序运行时控制存储介质所在的设备执行上述实施例中的方法步骤。
[0182]
该计算机可读存储介质包括但不限于快闪存储器、硬盘、固态硬盘等。
[0183]
这里需要指出的是，以上实施例中的可读存储介质描述是与实施例中的方法描述相对应的，本发明实施例在此不做赘述。
[0184]
在上述实施例中，可以全部或部分地通过软件、硬件、固件或者其任意组合来实现。当使用软件实现时，可以全部或部分地以计算机程序产品的形式实现。计算机程序产品包括一个或多个计算机指令。在计算机上加载和执行所述计算机程序指令时，全部或部分地产生按照本发明实施例的流程或功能。
[0185]
计算机可以是通用计算机、专用计算机、计算机网络、或者其他可编程装置。计算机指令可以存储在计算机可读存储介质中，或者通过计算机可读存储介质进行传输。计算机可读存储介质可以是计算机能够存取的任何可用介质或者是包含一个或多个可用介质集成的服务器、数据中心等数据存储设备。可用介质可以是磁性介质或者半导体介质等。
[0186]
本发明实施例对各器件的型号除做特殊说明的以外，其他器件的型号不做限制，只要能完成上述功能的器件均可。
[0187]
本领域技术人员可以理解附图只是一个优选实施例的示意图，上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。
[0188]
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。