基于慢同调理论和希尔伯特-黄变换的发电机在线同调识别

宋洪磊，吴俊勇，冀鲁豫

（北京交通大学 电气工程学院，北京 100044）

0 引言

电力系统的快速发展和大区域电网的互联，在提高电力系统传输能力和运行灵活的同时，也给系统的复杂性和安全运行带来了严重的暂态稳定问题，尤其是当电力系统发生大的扰动后，区域电网机群间的失步振荡将会造成大面积停电，并带来极大的经济损失和极坏的社会影响。此时就需要调度人员及时发现网络中的失步机群，从而切断失步机群间的联系，以平息系统间的振荡，所以快速准确地识别系统中的同调机群就变得尤为重要。

同调机群的识别方法有很多，一类是基于特定的扰动，时域仿真得到发电机的功角摇摆曲线识别同调机群，在一段时间尺度下，根据设定阈值区分机组的同调性[1-2]；一类是通过特征值分析，针对系统矩阵的状态矩阵进行识别，例如慢同调方法，文献[3]提出了利用奇异摄动理论和慢同调理论来聚合发电机组，对系统中的各发电机组群进行识别；还有一类是通过预测发电机的摇摆曲线进行小波分析，能够有效地对受扰后的发电机组进行同调性识别[4-5]。但是这些方法大都是离线下对模型进行线性化，对某一特定故障下的发电机进行同调识别，不仅分析计算过程繁琐，也不满足在线分析的要求[6-8]。随着智能电网技术的不断发展，广域测量系统（WAMS）在电网中的应用也日益成熟，通过电网中配置的PMU实时监测发电机的动态特性，并根据WAMS提供的故障后基于实时响应的功角信息，使得发电机群的同调在线识别成为可能。

对于大规模互联电网而言，故障位置、系统结构参数等信息的不明确性，以及失步解列要求的实时快速性，对同调机群的准确识别提出更高要求[9-11]。在电力系统运行中，特别是发生大的扰动后，动态信息具有非线性和非平稳性，且发电机功角轨迹含有很丰富的电力系统动态信息，可以通过提取系统的暂态信息，反映系统的时变特性。由于本文提出利用信号处理工具希尔伯特-黄变换（HHT）方法在线识别发电机群组的同调性，2组或多组具有较强同调性的发电机群组间的弱联系可以理解为系统中的低频振荡模式，所以用HHT方法对发电机的功角曲线进行模式提取，辨识其低频振荡模式，从而实现发电机群组的同调分群。

1 电力系统慢同调理论

在实际运行中，电力系统的动态过程普遍存在着同调现象，所谓同调就是指扰动发生后不同发电机转子具有相同或相近的运动变化趋势，利用聚类的方法将故障后转子摇摆曲线变化趋势相仿的发电机组归为一群。发电机组的同调分群主要表现为具有弱联系的机组间的低频振荡，这是由系统的内在特性和外在影响两部分因素共同决定的。其中内在的影响因素包括网络本身的拓扑结构、发电机惯性和阻尼常数等；而外在的影响因素是指故障发生的类型和位置。慢同调理论是双时间尺度在电力系统中的应用，其由同调机群内部的快动态和非同调机群间的慢动态构成，即局部振荡模式和区域振荡模式，这反映了系统机组的内在固有的振荡模式。从而利用慢同调分群理论对所有发电机进行分群，利用配置的PMU提取主导发电机组的功角信息，为后续HHT的同调识别提供基础和保证。

1.1 同调机群识别

当系统受到大的扰动发生失步振荡时，需要采取紧急解列控制措施保证系统的稳定性，其中同调机群的准确识别是系统安全快速解列的基础和前提。发电机的同调现象，即是故障后各发电机的转子摇摆曲线间的振荡趋势存在相似或一致性。关于同调的定义有多重，最常用的定义如下，在计算时间tє[0，τ]内，如果其中2个发电机i和j的相对转子角的偏差在任意时间都不大于界定的阈值ε（ε≥0），则判定这2个发电机在tє[0， τ]的时间段内同调[12-13]，即：

其中，δi和δj分别为发电机i和j的转子角。若ε=0，则称发电机是完全同调的；若ε＞0，则称发电机是ε相关。如果某一组发电机彼此之间是同调的，则这组发电机是同调机群。

1.2 慢同调分群算法

基于慢同调理论的发电机群同调性的识别，是采用线性化的电力系统模型为基础，提出一种完善的电力系统发电机组慢同调分群算法。忽略凸极效应、励磁绕组暂态及阻尼绕组的作用，同步电机采用经典二阶模型，则发电机转子运动方程式可以表达为[13-15]：

其中，i=1，2，…，n；Mi为第 i台发电机转子的运动惯性时间常数；δi为第i台发电机的转子角；ωi为第i台发电机的转速；ω0为同步转速；Pmi为第i台发电机的机械输入功率，是一个常数；Pei为第i台发电机的电磁输出功率；Di为阻尼常数；Ei、Ej为第 i台和第j台发电机暂态电抗后的电势相量Ei、Ej的幅值；Yij=Gij＋jBij为并入负荷阻抗和发电机暂态电抗后，收缩到发电机内电势节点的导纳矩阵元素；n为发电机台数。

忽略阻尼常数的影响，将公式在各个发电机功率及角度的稳态值Pei0、δi0附近进行线性化，可以得到：

其中，M=diag[M1，M2，…，Mn]为系统中发电机惯性时间常数的对角阵；K为系统状态矩阵。

求取状态矩阵A的全部特征值，每一个特征值λi表示为发电机系统的第i个模式状态，与之对应的时间特性为eλit。特征值的实部表示系统模态的衰减时间常数，虚部代表慢同调的振荡频率。其慢同调分群算法的步骤具体如下。

a.确定最优的慢同调聚合分群数。由于慢同调理论是以双时间尺度为基础，所以根据电力系统的多时间尺度的特征原理，将得到的全部特征值，按照振荡频率绝对值大小排序，即为发电机的台数。 则m=，m的值越小说明系统的时间尺度特征越明显，对于慢同调分区的区域性也越明显。 这样可得到慢模式组 σr=｛λ1，λ2，…，λr｝和对应的右特征向量空间 U=｛u1，u2，…，ur｝，即可确定最优的慢同调聚合分群数。

b.确定慢同调的参考机组。对于特征向量U=｛u1，u2，…，ur｝进行列主元高斯消去法，得到 r组线性无关的向量，其对应的r行选定为进行慢同调聚合的参考机组，并得到慢模式组σr所对应的参考机组向量U1。

c.计算分群矩阵L。对于慢模式对应的特征向量空间U中，参考机组向量为U1，则剩下的n-r行特征向量记为 U2，即特征向量矩阵 U=[U1U2]T，定义分群矩阵L=U2U1-1。

d.形成慢同调分群。根据分群矩阵L每行的最大值，将其他发电机组分配到相应的参考机组中，形成慢同调的分群区域。

对于最优慢同调聚合分群数r，即是系统的慢模式数，每一组内的任何2个状态变量之间都是关于慢模式组σr同调的。在利用慢同调分群算法进行研究时，分群数的选择会影响发电机组的聚合。

对于分群矩阵L=U2U1-1，定义系统发电机的分群数为g，当g＞r时，分群矩阵L是分群模式的一种，不能对发电机组进行准确分群；当g=r时，分群矩阵L即是分群模式的唯一情况，可以根据分群矩阵L对发电机组进行分群；当g＜r时，此时要保证范数‖LU1-U2‖2最小，从而确定发电机组的合理分群数。因此，对于系统慢模式和发电机组的分群之间存在的关系，只有g≤r时，才能保证系统的优化分群。

1.3 同调机组的能观性

通过先前研究中对于慢同调分群算法中分群矩阵L的计算，根据分群矩阵中每行最大值将非参考发电机组分配到相应的参考机组中，形成慢同调的分群区域。由于系统中某些发电机在电气距离和时间惯性常数上的不同，利用慢同调分群算法对系统完成分群时，由于分群矩阵中每行的最大值存在差异，提出根据此差异区分强同调机组和弱同调机组。首先将得到的分群矩阵中所有最大值按照从小到大进行排序，即｛L1-max，L2-max，…，Li-max｝，并根据多时间尺度特性识别强同调机组和弱同调机组，设时间尺度阈值为η，则公式可以表示为：

阈值表示在分群矩阵中每行最大值之间的差异尺度，以此阈值为尺度界定区分强同调机组和弱同调机组。所以得到的分群矩阵中｛L1-max，L2-max，…，Lr-max｝所对应的发电机为弱同调机组，即这些发电机与参考机组是弱相关性；而｛L（r＋1）-max，…，Li-max｝所对应的发电机为强同调机组，即这些发电机与参考机组是强相关性的。

在考虑经济性和能观性的基础上，功角信号的数据提取主要是针对慢同调参考机组和弱相关性的机组，因为参考机组是每个同调机群中的代表机组，对其余相关发电机的慢同调聚合，而弱相关性的机组是考虑到由于故障条件等差异引起的分群不明确性。

为了提高同调机群识别的准确性和快速性，也为了保证PMU配置的合理性和有效性，本文提出根据参与因子和相关因子提取具有特征信息的机组进行监测。设ui和vi分别是特征值λi的右特征向量和左特征向量，则相应的模式状态i的参与因子Pi定义为[17]：

参与因子Pi表示状态模式中变量的相对参与程度。参与因子的值越大，则对模式状态的观测性和控制的灵敏度也越大。其中，右特征向量反映了系统的振荡模式，表征了状态各变量的相对运动方向，反映了对特征值λi的可观性，左特征向量反映了对特征值λi的可控性。定于相关因子Cki是表示第k个状态量与第i个特征根的相关程度，即：

其中，vki、uki分别为 vi、ui的第 k 个元素。 相关因子值的大小反映了状态量对于相应特征值的可观性和可控性。

因此，相关因子可反映每台发电机与振荡模式的强相关性程度。综合利用相关因子和慢同调分群对PMU合理配置，提取具有代表性的发电机的振荡信息量，可以保证识别的快速性和准确性。

2 基于HHT的同调识别

HHT是近年来发展起来的一种针对信号自适应的时频分析方法，该方法不仅在时域和频域上具有较高的分辨率，还能够自适应地将信号分解，得到不同特征时间尺度上的固有模态函数。HHT方法主要是利用经验模式分解（EMD）对非线性、非平稳的数据进行线性化和平稳化处理，获取发电机功角信号的本征模态函数（IMF），再通过选取IMF并对其进行Hilbert变换，得到IMF分量的瞬时频率和瞬时幅值，利用合成的IMF分量的Hilbert谱分布，实现对信号数据的时频分析，即时间-频率-振幅的三维谱分析，可以很直观地反映低频振荡过程的动态特性和系统的非线性时变特性。

当系统失步运行时，通过配置的PMU提取主导参考机组和弱相关性机组的功角响应信号 x（t）=A（t）×cos φ（t），经过 EMD得到一组线性组合的IMF分量和一个余项的和[16]，即：

然后利用Hilbert变换求出固有模态信号的解析信号，分别得到分量信号的瞬时幅值A（t）和瞬时频率f（t），进而得到信号的 Hilbert谱和 Hilbert边界谱。对每个IMF进行Hilbert变换得：

则Hilbert幅值谱即是上式的展开式，可以记作：

还可以得到边际谱为：

以上的EMD和相应的Hilbert谱信号的分析称为HHT，瞬时频率及Hilbert谱中包含着系统丰富的功角失稳信息，例如在故障发生时刻和系统发生失步时都会在时频图中表现为突变，从而实现对发电机暂态失稳信息的提取，实现对低频振荡的辨识，从而对发电机的在线同调识别提供更加丰富和有价值的信息。并且由于瞬时幅值反映振荡的阻尼特性，瞬时频率则反映各振荡模式时变特性，可以非常精确地描述信号的能量随着频率的分布状况，对各发电机间的振荡模式相互作用刻画得更加准确。

3 同调识别算法分析

本文提出的方法是慢同调分群算法和HHT振荡分析相结合对扰动后的发电机组同调识别，综合考虑了发电机的固有振荡模式和外在扰动影响。首先在离线情况下，对发电机的强相关同调机组和弱相关机组分群合理配置PMU以提取具有代表意义和重要的发电机组信息，其中根据振荡模式状态的参与因子确定慢同调分群得出的同调群组中具有可观测性的发电机组，并且提取预分群配置PMU的实时功角信息；再利用信号处理工具HHT对得到的动态信息进行分析，实现发电机群组的快速准确识别。方法的具体步骤为：

a.在离线的情况下，利用慢同调分群算法对系统的发电机群组分群，并且根据最优的慢同调聚合分群数确定系统的发电机固有振动模式；

b.在考虑同调机组能观性的基础上，根据慢同调分群矩阵的相关系数和相关振荡模式的参与因子，确定预分群组中发电机的PMU信息提取；

c.当系统发生扰动后，提取广域量测系统提供的实时功角数据信息，利用HHT进行边际谱振荡模式能量谱分析，得到失步机群并进行同调识别。

4 算例分析

4.1 WSCC-3机9节点系统

本算例的系统拓扑结构如图1所示，本系统是由3台发电机和3个负荷组成的简单系统。

图1 WSCC系统的拓扑结构图Fig.1 Topology of WSCC system

首先利用本文提出的慢同调分群算法将整个系统的发电机组分群。通过建立的系统状态矩阵得到3 个特征值分别为 0.007 8、-75.503 1、-178.464 6，且根据系统多时间尺度的特征原理，识别出本系统的慢模式有2个，即本系统的最优分群数为2个，得到慢同调参考机组为发电机1和发电机2（与母线1和2直接相连的发电机称为发电机1和发电机2，其他类似），则分群矩阵 L=[0.3089 0.6914]，识别出发电机2和发电机3为慢同调机群，发电机1为单独一群。由于在分群矩阵中得到发电机2和发电机3属于二级相关同调，考虑到本算例的发电机数量较少，可以同时在3个发电机处配置PMU，监测系统的动态信息曲线。

a.系统发生故障，不发生失步状态。

故障形式为线路4-6处在2 s发生三相瞬时短路故障，持续时间为0.2 s，仿真时间为10 s，计算步长为0.01 s，所用的软件为PSASP。仿真得到的发电机功角曲线如图2所示，图中δ21为发电机2与发电机1间功角差，其他类似。

图2 发电机在不失步状态下的功角曲线Fig.2 Power angle trajectory of generator under synchronous state

本文提取仿真后的数据作为PMU的信息输入，对其进行分析。当故障切除后，根据提取的发电机功角曲线可以看出系统没有出现失步现象，但是很难辨识出同调机群。这里提出利用HHT对发电机功角曲线进行分析，由振荡频率识别同调机群。首先对每条功角曲线进行EMD线性化和平稳化处理，分解得到各发电机功角的IMF分量，再进行Hilbert变换，得到发电机的瞬时频率和瞬时振幅，进而可求得信号的Hilbert谱和Hilbert边际谱，实现了对信号的时频分析。发电机功角的Hilbert边际谱的幅值表示某个频率成分的能量大小，能量的值越大，则说明该频率的振动波在整个时间范围轴上出现的可能性就越大，精确反映了信号的振荡频率和能量幅值随时间变化的规律，具有很好的时频分辨率，表征了信号的局部特征。具体的各发电机的HHT边际谱分析如图3所示。

图3 发电机功角HHT边际谱Fig.3 HHT marginal spectrum of generator power angle

由图3可以看出共有3个振荡模式频率，即为0.1875 Hz、0.5625 Hz和 1.3125 Hz。 其中振荡模式频率在0.1875 Hz时，发电机1和发电机2为同调机群，发电机3为单独一群；振荡模式频率为1.3125 Hz时，发电机2和发电机3为同调机群，发电机1为单独一群。其中当振荡频率为1.3125 Hz时HHT的边际谱振荡模式频率能量最大，说明此频率为系统的主导振荡频率，可以由此得到系统的分群模式为：发电机2和发电机3为同调机群，而发电机1为单独一群，这和利用线性化的慢同调分群算法得到的结果是一致的。

b.系统发生故障，发生失步状态。

改变故障线路和延长故障切除的时间，当线路8-9处在2 s发生三相瞬时短路故障，持续时间为0.25 s，仿真时间为10 s，计算步长为0.01 s。仿真得到各发电机功角曲线如图4所示。

图4 发电机在失步状态下的功角曲线Fig.4 Power angle trajectory of generator under out-of-step state

图5为发电机各功角的HHT边际能量谱。可以看出，发电机功角δ21和δ23在0.1875 Hz处的振幅能量非常明显，是当前故障模式下的主导振荡频率，由于本系统的发电机数目较少，振荡模式的相互影响较小，比较容易识别主导振荡频率。如果在较多发电机系统中，各发电机会相互影响，识别出的振荡模式可能会为多模式的合成。通过主导振荡模式能量识别出系统在失步情况下，发电机2为失步机群，发电机1和发电机3为同调机群，这与慢同调分群算法得到的同调分群不一样，这说明同调机群的划分和故障位置及严重程度有关，通过对慢同调分群算法的研究合理配置PMU，提取有意义的信号特征实现同调机群的划分。

图5 发电机功角HHT边际谱Fig.5 HHT marginal spectrum of generator power angle

4.2 新英格兰10机39节点系统

为进一步验证本文方法对发电机同调识别的有效性和可行性，以新英格兰10机39节点系统为例进行仿真验证。新英格兰10机39节点系统包括10台发电机、39条母线和19条负荷母线。其中，39号机的母线节点设定平衡节点，30～38号发电机的母线节点为PQ节点。

首先利用本文提出的慢同调分群算法将整个系统的发电机组分群。通过建立的系统状态矩阵特征值识别出慢模式为3个，因为时间尺度关系为最小值，即本系统的最优分群数为3个，可以得到发电机30、37和38称为同调机群1，发电机31、32和39称为同调机群2，发电机33、34、35和36称为同调机群3。其中在同调机群1中发电机38为弱相关同调机组，同调机群3中发电机33为弱相关同调机组。根据在参考机组对应的振荡模式为相应特征值的可观性进行判断，在每一个慢同调分群基础上得到各发电机组的参与因子。本算例假定提取同调机群1中发电机30和38的功角信息、同调机群2中发电机31和39的功角信息、同调机群3中发电机33和34的功角信息。

本算例中假设的故障形式为线路16-19处在5 s发生三相瞬时短路故障，持续时间为0.2 s，仿真时间为15 s，故障切除后系统出现失稳。仿真得到的30～38号发电机功角曲线如图6所示（其中39号发电机设定为参考发电机）。

图6 发电机失步状态下的功角曲线Fig.6 Power angle trajectory of generator under out-of-step state

通过功角曲线得到发电机33和34出现失稳现象。为了验证慢同调分群算法和HHT同调识别算法的适用性，提取相关发电机的功角信息，进行边际谱振荡模式能量分析。

图7—11为发电机功角的边际谱能量谱分布，δi_39代表i号发电机和参考发电机间的功角差。

图7 δ30_39的 HHT 边际谱Fig.7 HHT marginal spectrum of δ30_39

图8 δ31_39的 HHT 边际谱Fig.8 HHT marginal spectrum of δ31_39

图9 δ33_39的 HHT 边际谱Fig.9 HHT marginal spectrum of δ33_39

图10 δ34_39的 HHT 边际谱Fig.10 HHT marginal spectrum of δ34_39

图11 δ38_39的 HHT 边际谱Fig.11 HHT marginal spectrum of δ38_39

可看出，系统共存在 0.1875 Hz、0.3125 Hz、1.1875 Hz和1.3125 Hz这4种振荡模式。其中，发电机δ33_39和δ34_39在频率0.1875 Hz处出现最大的边际谱振荡模式能量振幅，且频率的组成较单一，模式能量幅值约为55000～60000，而其他发电机功角的振荡频率都是多种频率的相互影响。各发电机功角的振荡频率和能量幅值整理如表1所示，模式能量幅值约为20～90，与失步机群的能量幅值差别明显。HHT同调识别方法中的边际谱振荡模式能量可以作为一个明显的识别指标，实现对发电机组的同调性辨识。

表1 发电机功角的振荡频率和能量幅值Tab.1 Oscillation frequency and energy amplitude of generator power angle

从表1中可以看出，在非失稳的发电机中都是由2种能量幅值相近的振荡频率组成。为了解释此问题，在相同故障条件下，提取所有的发电机功角进行HHT分析，可以判断出发电机33和34出现较大的低频振荡边际谱能量，是系统的失步机群，这与仿真结果是对应的。在网络拓扑结构上，距离失步发电机较近的发电机31和32都受到主导失稳振荡频率0.1875 Hz的影响，虽然该频率不是这些发电机的主导频率，但是由于系统间存在的非线性作用，导致邻近发电机产生相互混合的振荡模式。通过HHT的分析可以很明显地识别出发电机的失步机群，从而实现对同调机群的合理分群；并且HHT可以直观地反映低频振荡过程的动态特性和系统的非线性时变特性，得到的边际谱振荡模式能量差别明显。

以上结果表明：以慢同调分群算法的离线分群结果为基础，结合参与因子表现的能观性，可以通过有效配置PMU提取必要发电机的功角信息；根据HHT的同调识别可以实现在线的发电机同调分群，这种分群算法是充分考虑了系统网络本身的拓扑结构、发电机惯性以及故障发生类型和位置等外在的影响因素。

5 结论

本文提出了一种新的同调机群识别方法，将慢同调分群算法和新兴的信号处理工具HHT相结合，有效利用最优配置的PMU提取故障后发电机的功角信息，实现发电机群的同调分群，此方法为在线快速的同调机群识别和失步解列控制策略提供了基础和依据，并且可以实现WAMS的实时功角信息的充分利用。本文在WSCC系统和新英格兰39节点系统上的应用结果表明，提出的同调分群识别新方法将网络系统固有的拓扑特征的慢同调分群算法和在线的HHT故障信息特征提取方法相结合，更能全面反映系统在受扰后（包括不失步和失步2种形式）的动态特性和系统的非线性时变特性。本文下一步研究的重点在于利用HHT的振荡模式能量谱的提取，分析在系统发生失步解列时的最佳解列时间，利用信号处理分析实现在线的主动解列控制。