CN102074955B - 基于知识发现技术的电力***稳定评估及控制方法 - Google Patents

Abstract

本发明提供了一种电力***暂态稳定评估及实时控制方法，将知识发现相关技术引入电力***稳定评估及控制问题，在构建暂态稳定评估模型的基础上，找出***运行工况和电力***整体稳定性之间的对应关系，掌握电力***内在的运行规律，实现对电力***稳定运行状态的快速评估，用以指导电力***实际的运行管理；对遭受大扰动后失稳的运行状态，能制定及时有效的稳定控制方案，以较小的代价保证***的安全稳定运行。

Description

基于知识发现技术的电力***稳定评估及控制方法

技术领域

本发明属于电力***运行和控制技术领域，具体涉及一种基于知识发现技术的电力***稳定评估及控制方法。

背景技术

在电力***中，稳定问题始终被摆在首位，它决定着电力***能否正常供电。电力***稳定事故中，多数是暂态稳定破坏事故，研究电力***的暂态稳定问题有重要的现实意义。目前主要的研究方法有时域仿真法(或称逐步积分法)、暂态能量函数法(或称直接法)以及这两种方法的混合方法。

利用智能方法进行非模型的电力***暂态稳定判别具有在线计算速度快、容易生成决策的启发规则等优点，可与传统暂态稳定分析方法构成良好的互补。与此同时，计算机技术在电力***中已得到了广泛应用，电力***通过各种监测装置，可以收集到大量运行数据，这为暂态稳定问题直接利用量测数据进行分析提供了充足的技术支撑。当前的突出问题是面对如此丰富的数据，如何***地开发强有力的数据分析工具，准确及时地从原始数据中提取出有用的信息。

将知识发现技术应用于电力***，可更加充分地利用实际运行数据，揭示海量数据背后蕴含的电力***运行特征和规则，找到更加合理的解决问题的方法，为决策的制定提供更加有力的科学依据。

发明内容

本发明提供一种电力***暂态稳定在线评估及控制方法，将知识发现相关技术引入电力***稳定评估及控制问题，在构建暂态稳定评估模型的基础上，找出***运行工况和电力***稳定性之间的对应关系，掌握***内在的运行规律，实现对电力***稳定运行状态的快速识别，用以指导电力***实际的运行管理；对遭受大扰动后失稳的运行状态，提供及时有效的稳定控制方案，以较小的代价保证***的安全稳定运行。

本发明从静态和动态两方面入手，定义稳定评估的输入特征变量，建立稳定评估数学模型。由现场实测电气量或实时仿真数据生成样本。使用数据预处理技术来改善数据质量，降低实际评估计算所用的时间，提高分析的精度和性能。利用支持向量机来构造分类器，进行电力***稳定评估。对失稳状态采用混合型粒子群优化算法搜索可行的紧急控制措施，保证***的稳定运行。

依据本发明的一种电力***稳定评估及控制方法，包括以下步骤：

(1)确定电力***的输入特征变量，构建电力***的暂态稳定评估模型；

(2)获取电力***中电网在线运行数据及其在给定故障集下的短时暂态稳定仿真计算结果，根据步骤(1)确定的输入特征变量，构成训练样本空间，每个样本均由输入特征变量集合x₀和对应的暂稳计算结果y₀两部分构成；

(3)在不损失数据信息的基础上，对初始海量数据样本进行数据预处理，得到新的数据样本集合(x，y)；

(4)根据训练样本集合(x，y)，构造暂态稳定最优分类判别函数：

$f\left(x\right)=\mathrm{sgn}\{(w^{*}\·x)+b^{*}\}=\mathrm{sgn}\{\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\α}}_i^{*}{y}_{i}K(x_i\·x)+b^{*}\}---\left(1\right)$

其中sgn(·)为符号函数，w^＊为一次线性函数的斜率，b^＊为一次线性函数的截距，α_i为与每个样本对应的Lagrange乘子，K(x_i，x)为高斯径向核函数，有

$K(x_i,x)=\exp (-\frac{{|x-x_i|}^2}{{\mathrm{\σ}}^2})$

其中σ为函数的宽度参数，控制了函数的径向作用范围。

(5)对于电力***运行状态未知的运行数据x，将其按步骤(3)经过数据预处理后，代入式(1)，即得到其对应的暂态稳定运行状态f(x)；f(x)的值代表了每个样本到最优分类超平面的距离，在暂态稳定评估问题中，f(x)的值代表了当前运行状态x到暂态稳定分类面的距离，即暂态稳定裕度。f(x)取值为正，表示***稳定，且数值越大***越稳定；反之，取值为负表示***暂态不稳定，且该数绝对值越大表示***越不稳定；

(6)对步骤(5)判定为失稳的电力***运行方式，利用混合型粒子群优化算法，搜索并得到能使电力***恢复暂态稳定运行且切机或切负荷量最小的紧急控制措施。

其中，确定能反映电力***发生故障时暂态特性的输入特征x₀及其对应的暂稳计算结果y₀，构建电力***暂态稳定评估模型；步骤(1)采用的输入特征变量x₀分为静态特征变量和动态特征变量，也即x₀＝(静态特征变量，动态特征变量)，具体包括：

(1)静态特征变量包括

1)整个电力***的发电机的出力和负荷；

2)电力***的最大和最小的发电出力；

3)电力***的最大和最小的负荷；

4)母线电压的最大和最小值；

5)线路最大有功传输功率，及其对应的无功传输功率；

6)关键输电断面和区域联络断面)的传输功率；

7)关键发电机的发电出力；

8)各区域的发电机的出力和负荷；

(2)动态特征变量包括：

以下各式中，δ_i为发电机转子角；ω_i为发电机转速，a_i为发电机转子角加速度，M_i为发电机的惯性常数；P_mi为发电机机械功率；P_ei为发电机电磁功率；V_ki为发电机动能；V_pi为发电机势能；N_G为发电机数且，t₀为故障发生初始时刻，t_cl为故障切除时刻；

1)故障发生初始时刻具有最大功角加速度的发电机的相对转子角度：

$x_1=({\mathrm{\δ}}_{a\max }-{\mathrm{\δ}}_{\mathrm{COI}}){|}_{t_0}$

上式中是相对于***惯性中心的等值转子角，δ_amax是故障发生初始时刻具有最大功角加速度的发电机的转子角；

2)故障切除时刻具有最大动能的发电机的相对转子角度：

$x_2=({\mathrm{\δ}}_{p\max }-{\mathrm{\δ}}_{\mathrm{COI}}){|}_{t_{\mathrm{cl}}}$

上式中δ_pmax是故障切除时刻具有最大动能的发电机的转子角；

3)故障切除时刻领前机与殿后机的功角之差；故障发生初始时刻，具有最大转子角加速度的发电机为领前机；具有最小转子角加速度的发电机为殿后机；

$x_3=\max \left({\mathrm{\δ}}_{a\max }\right)-\min \left({\mathrm{\δ}}_{a\min }\right){|}_{t_{\mathrm{cl}}}$

上式中δ_amax、δ_amin分别是领前机和殿后机的转子角；

4)故障切除时刻最大的发电机转角差：

$x_4=\max \left({\mathrm{\δ}}_i\right)-\min \left({\mathrm{\δ}}_i\right){|}_{t_{\mathrm{cl}}}$

5)故障切除时刻与惯性中心COI相差最大的功角：

$x_5=\max ({\mathrm{\δ}}_i-{\mathrm{\δ}}_{\mathrm{COI}}){|}_{t_{\mathrm{cl}}}$

6)故障切除时刻领前机群和余下机群中的发电机的惯量加权平均角度之差：

$x_6=\underset{i\∈s}{\mathrm{\Σ}}{M}_i{\mathrm{\δ}}_i/\underset{i\∈s}{\mathrm{\Σ}}{M}_i-\underset{j\∈A}{\mathrm{\Σ}}{M}_j{\mathrm{\δ}}_j/\underset{j\∈A}{\mathrm{\Σ}}{M}_j{|}_{t_{\mathrm{cl}}}$

其中，S表示领前机群集合，A表示余下机群集合。对指定时刻的各发电机转角进行降序排序，形如：

$({\mathrm{\δ}}_1,{\mathrm{\δ}}_2,{\mathrm{\δ}}_3,...,{\mathrm{\δ}}_{N_G-2},{\mathrm{\δ}}_{N_G-1},{\mathrm{\δ}}_{N_G}),$ 且 ${\mathrm{\δ}}_1>{\mathrm{\δ}}_2>{\mathrm{\δ}}_3>...>{\mathrm{\δ}}_{N_G-2},{\mathrm{\δ}}_{N_G-1}>{\mathrm{\δ}}_{N_G}$

每两个相邻转角之差构成一个角度间隙，形如：

$({\mathrm{\δ}}_{\mathrm{1,2}},{\mathrm{\δ}}_{\mathrm{2,3}},...,{\mathrm{\δ}}_{N_G-2,N_G-1},{\mathrm{\δ}}_{N_G-1,N_G}),$ 且

${\mathrm{\δ}}_{\mathrm{1,2}}={\mathrm{\δ}}_1-{\mathrm{\δ}}_2,{\mathrm{\δ}}_{\mathrm{2,3}}={\mathrm{\δ}}_2-{\mathrm{\δ}}_3,...,{\mathrm{\δ}}_{N_G-2,N_G-1}={\mathrm{\δ}}_{N_G-2}-{\mathrm{\δ}}_{N_G-1},{\mathrm{\δ}}_{N_G-1,N_G}={\mathrm{\δ}}_{N_G-1}-{\mathrm{\δ}}_{N_G}$

确定最大的角度间隙，此间隙之上的发电机构成领前机群S，此间隙之下的发电机构成余下机群A；

7)故障切除时刻与故障发生初始时刻发电机相对转子角度差值的绝对值之和：

$x_7=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}|{\mathrm{\θ}}_{\mathrm{ci}}-{\mathrm{\θ}}_{0i}|$

上式中，θ_i＝δ_i-δ_COI，θ_c为故障切除时刻发电机相对转子角，θ₀为故障发生初始时刻发电机相对转子角。

8)故障切除时刻与故障发生初始时刻发电机相对转子角度差值的绝对值的最大值：

$x_8=\underset{i}{\max }\left(\right|{\mathrm{\θ}}_{\mathrm{ci}}-{\mathrm{\θ}}_{0i}\left|\right)$

9)故障切除时刻与故障发生初始时刻惯性中心坐标COI转速的变化：

$x_9={\mathrm{\ω}}_{\mathrm{COI}}{|}_{t_{\mathrm{cl}}}-{\mathrm{\ω}}_{\mathrm{COI}}{|}_{t_0}$

上式中 ${\mathrm{\ω}}_{\mathrm{COI}}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{M}_i{\mathrm{\ω}}_i/\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{M}_i$

10)故障切除时刻发电机转子角速度的绝对值之和：

$x_{10}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{|{\mathrm{\ω}}_i-{\mathrm{\ω}}_{\mathrm{COI}}|}_{t_{\mathrm{cl}}}$

11)故障切除时刻发电机转子角速度绝对值的最大值：

$x_{11}=\max \left(\right|{\mathrm{\ω}}_i-{\mathrm{\ω}}_{\mathrm{COI}}\left|\right){|}_{t_{\mathrm{cl}}}$

12)所有发电机转子在故障发生初始时刻相对加速度的最大、小值：

$x_{12-1}=(a_{\max }-a_{\mathrm{COI}}){|}_{t_0}$

$x_{12-2}=(a_{\min }-a_{\mathrm{COI}}){|}_{t_0}$

其中是a_COI惯性中心加速度a_max、a_min分别是转子角加速度的最大、小值；

13)所有发电机转子在故障发生初始时刻加速度的均值：

$x_{13}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{a}_i{|}_{t_0}$

14)所有发电机转子在故障发生初始时刻加速度的均方根误差：

$x_{14}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{(a_i-x_{13})}^2}{|}_{t_0}$

15)所有发电机转子在故障发生初始时刻相对加速度的平均值：

$x_{15}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}(a_i-a_{\mathrm{COI}}){|}_{t_0}$

16)所有发电机转子在故障发生初始时刻相对加速度的方差：

$x_{16}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{((a_i-a_{\mathrm{COI}})-x_{15})}^2}{|}_{t_0}$

17)领前机与殿后机在故障瞬间和故障发生初始时刻加速度之差：

$x_{17-1}=\max \left(a_i\right)-\min \left(a_i\right){|}_{t_0}$

$x_{17-2}=\max \left(a_i\right)-\min \left(a_i\right){|}_{t_{\mathrm{cl}}}$

18)故障切除时刻所有转子动能的最大值：

$x_{18}=\max \left(V_{\mathrm{ki}}\right){|}_{t_{\mathrm{cl}}}i=\mathrm{1,2},...,N_G$

19)故障切除时刻具有最大转子角度的发电机p的转子动能：

$x_{19}=\frac{1}{2}{M}_p\×({{\mathrm{\ω}}_p}^2-1){|}_{t_{\mathrm{cl}}}$

20)电力***总的“能量调整”：

$x_{20}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}(P_{\mathrm{ai}}\×{\mathrm{\δ}}_{\mathrm{di}})$

其中，P_ai是第i台发电机在故障发生初始时刻的加速功率；是第i台发电机故障切除时的相对转子角度；

21)故障切除时刻所有转子动能的平均值：

$x_{21}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\frac{1}{2}{M}_i\×({{\mathrm{\ω}}_i}^2-1){|}_{t_{\mathrm{cl}}}$

22)故障前所有发电机机械输入功率的平均值：

$x_{22}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{p}_{\mathrm{mi}}$

23)故障前所有发电机机械输入功率的方差：

$x_{23}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{(P_{\mathrm{mi}}-x_{22})}^2}$

24)发电机在故障发生初始时刻加速功率的均值：

$x_{24}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}(P_{\mathrm{mi}}-P_{\mathrm{ei}}){|}_{t_{\mathrm{cl}}}$

25)发电机在故障发生初始时刻相对加速功率的方差：

$x_{25}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{(\mathrm{\Δ}{P}_i-x_{24})}^2}{|}_{t_0}$

26)发电机在故障发生初始时刻相对加速功率的均值：

$x_{26}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\frac{{\mathrm{\ΔP}}_i}{P_{\mathrm{mi}}}{|}_{t_0}$

27)发电机在故障发生初始时刻相对加速功率的方差：

$x_{27}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{(\frac{{\mathrm{\ΔP}}_i}{P_{\mathrm{mi}}}-x_{26})}^2}{|}_{t_0}$

28)故障切除时刻发电机的总动能：

$x_{28}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\left|V_{\mathrm{KEi}}\right|{|}_{t_{\text{cl}}}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\left|\frac{1}{2}{M}_i{{\mathrm{\ω}}_i}^2\right|{|}_{t_{\mathrm{cl}}}$

29)故障切除时刻发电机的最大动能差：

$x_{29-1}=(\max \left(V_{\mathrm{ki}}\right)-\min \left(V_{\mathrm{ki}}\right)){|}_{t_{\mathrm{cl}}}$

$x_{29-2}=(\max \left({\stackrel{\·}{V}}_{\mathrm{ki}}\right)-\min \left({\stackrel{\·}{V}}_{\mathrm{ki}}\right)){|}_{t_{\mathrm{cl}}}$

30)故障发生初始时刻单台发电机所受的最大有功功率冲击：

$x_{30}=\max (P_{i1}-P_{i0}){|}_{t_0}i=\mathrm{1,2},...,N_G$

式中P_i1，P_i0分别为发电机i在故障瞬间时和故障前的有功输出；

31)故障发生初始时刻单台发电机所受的最小有功功率冲击：

$x_{31}=\min (P_{i1}-P_{i0}){|}_{t_0}i=\mathrm{1,2},...,N_G$

32)故障发生初始时刻单台发电机所受的最大归一化有功功率冲击：

$x_{32}=\underset{i}{\max }\frac{P_{i1}-P_{i0}}{P_{\mathrm{ie}}}{|}_{t_0}i=\mathrm{1,2},...,N_G$

式中P_ie为发电机i的额定有功。

33)故障发生初始时刻单台发电机所受的最小归一化有功功率冲击：

$x_{33}=\min \frac{P_{i1}-P_{i0}}{P_{\mathrm{ie}}}{|}_{t_0}i=\mathrm{1,2},...,N_G$

34)表征故障对***冲击大小的性能特征变量：

$x_{34}=\left(\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{M}_i\right|P_{\mathrm{di}}\left|\right)/\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{M}_i$

上式中是每台发电机相对于***惯性中心的减速功率；P_ci是故障切除瞬间发电机输出功率；P_0i是故障前发电机输出功率；

35)事故清除前一瞬间发电机加速功率之和：

$x_{35}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{P}_{\mathrm{acc},i}^{\mathrm{BF}}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}(P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}})=P_{\mathrm{COI}}^{\mathrm{BF}}$

为事故清除后一瞬间发电机i的加速功率；为事故清除后一瞬间惯性中心的加速功率；分别为事故清除后瞬间发电机的机械功率和电磁功率；

36)事故清除后一瞬间发电机加速功率绝对值的最大值：

$x_{36}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{PF}}\left|\right\}=\underset{i}{\max }\left\{\right|P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}\left|\right\}$

37)事故清除前一瞬间发电机的加速功率与发电机惯性时间常数比值之和：

$x_{37}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{P}_{\mathrm{acc},i}/M_i=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\left(\frac{P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}}}{M_i}\right)$

38)事故清除前一瞬间发电机的加速功率与发电机惯性时间常数比值的最大值：

$x_{38}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{BF}}/M_i\left|\right\}=\underset{i}{\max }\left\{\right|\frac{P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}}}{M_i}\left|\right\}$

39)事故清除后一瞬间发电机的加速功率与发电机惯性时间常数比值之和：

$x_{39}=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}{P}_{\mathrm{acc},i}^{\mathrm{PF}}/M_i=\underset{i=1}{\overset{N_G}{\mathrm{\Σ}}}\left(\frac{P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}}{M_i}\right)$

40)事故清除后一瞬间发电机的加速功率与发电机惯性时间常数比值的最大值：

$x_{40}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{PF}}/M_i\left|\right\}=\underset{i}{\max }\left\{\right|\frac{P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}}{M_i}\left|\right\}.$

其中，所述的海量数据预处理过程包括以下步骤：

(1)数据规范化处理

$x=\frac{x_0-x_{0\max }}{x_{0\max }-x_{0\min }}$

式中x₀、x分别为转换前后的数据样本集合；x_0max为样本最大值；x_0min为样本最小值；

(2)样本规约

1)根据混合F统计值，确定最佳聚类数C；

2)计算每个特征变量对稳定结果的信息增益，作为聚类分析中的特征权重W_i；对规范化的样本数据进行聚类分析，将每个样本划分到与其距离最小的聚类中心对应的聚类子集中，形成C个数据样本子集；

(3)维度规约

对样本数据进行对数中心化变换和行向量中心化，通过非线性K-L变换，将样本集的维度(即特征变量数)，由最初的m降到p(p＜m)；

1)对原始数据X＝(x_ij)，i＝1，...，n，j＝1，...，m(n为样本数，m为输入特征变量数)进行对数变换，变换后的数据集记为W＝{w_ij}，有

w_ij＝ln x_ij

2)行向量中心化，令

$y_{\mathrm{ij}}=w_{\mathrm{ij}}-\frac{1}{m}\underset{j=1}{\overset{m}{\mathrm{\Σ}}}{w}_{\mathrm{ij}}$

3)计算矩阵Y＝(y_ij)_n×m的协方差阵及其特征根和与特征值对应的标准正交特征向量

4)对特征值进行降序排序，有λ₁＞λ₂＞…＞λ_m；

5)根据累计贡献率选择累计贡献率＞累计贡献率阀值的前p(p＜m)个特征值，一般的，累计贡献率阀值取0.85～0.95；

6)K-L变换，形成新的数据样本集合

F_i＝Γ_iWi＝1，2，...，p

Γ_i是与上述协方差阵第i个特征值λ_i对应的标准正交特征向量。

其中，还包括下列步骤：

(1)构造电力***紧急控制目标函数：

Max(f(稳定裕度，MVA_G，MVA_L))                  (2)

当***不稳定时：

f＝-abs(稳定裕度)                             (3)

***稳定时：

MVA_Gi＝P_gi×r_gi i＝1，...，N_G

MVA_Lj＝P_lj×r_lj j＝1，...，N_L

上式中weight是为了平衡稳定裕度取值和控制量取值的权重指数，一般取值为100；MVA_G是发电机的切机量，MVA_L是切负荷量，P_gi是被切发电机的有功出力，P_lj是被切负荷的有功负荷。r_gi是发电机切出力比例，r_lj是切负荷比例，满足：

r_gi＝0或1                                     (5)

0≤r_lj≤1                                     (6)

(2)定义待切除的发电机和负荷，包括故障线路送端侧的N_G台发电机和受端侧的N_L个负荷：

(3)初始化N个输入粒子每个粒子是一个(N_G+N_L)维矢量，前N_G维数据代表N_G台待切发电机的状态r_g，在实际的紧急控制中，发电机只有完全切和完全不切两种状态，这个N_G维数据的取值是非0即1；后N_L维数据代表N_L个待切除负荷的切除比例r_l，每一维的取值是一个[0，1]之间的实数；在允许的范围内，随机设置各粒子的初始位置和速度

(4)一个粒子对应一种控制措施，利用电力***分析软件(如PSASP)进行暂态稳定计算，判断加入这个粒子对应的控制措施后，***是否稳定；根据判稳结果，按式(3)或(4)计算每个粒子对应的目标函数；

(5)利用混合粒子群法搜寻最优的控制方案。对第k次迭代，根据式(7)更新每个粒子的速度，对每个粒子的的第1～N_G发电机维，根据式(8)更新对应维度位置；对每个粒子的的第(N_G+1)～(N_G+N_L)负荷维，根据式(9)更新对应维度位置取值

$v_i^{(k+1)}=w\·v_i^{\left(k\right)}+a_1\×\mathrm{rand}\left(\mathrm{0,1}\right)\×({\mathrm{pbest}}_i^{\left(k\right)}-x_i^{\left(k\right)})+a_2\×\mathrm{Rand}\left(\mathrm{0,1}\right)\×({\mathrm{gbest}}^{\left(k\right)}-x_i^{\left(k\right)})---\left(7\right)$

$\mathrm{If}(\mathrm{rand}\left(\right)<S\left(v_i^{(k+1)}\right))\mathrm{then}{x}_i^{(k+1)}=1$

$\mathrm{else}{x}_i^{(k+1)}=0---\left(8\right)$

$x_i^{(k+1)}=x_i^{\left(k\right)}+v_i^{(k+1)}---\left(9\right)$

式中w是惯性权重，a₁和a₂为两个学习因子，rand(0，1)和Rand(0，1)是两个均匀分布在(0，1)之间的随机数，i＝1，2，...，m。w定义方式如下：

$w=w_{\max }-\frac{w_{\max }-w_{\min }}{{\mathrm{iter}}_{\max }}\&times;\mathrm{iter}---\left(10\right)$

其中w_max＝0.9，w_min＝0.4，iter是当前的迭代次数，最大迭代数设为100。

本发明的有益效果是：

利用知识发现理论和技术及有效的数据预处理技术，能够自动从海量的电力运行数据中自动提取运行工况和电力***整体稳定性之间的对应关系，发现电力***内在的运行规律，用以指导电力***实际的运行管理和调度决策。该方法计算速度快，并能给出***暂态稳定性的量化指标，对电力***在线动态安全评估的在线实现具有重要的实用意义。

其中数据预处理技术的应用，在有效降低个体分类预测数据空间的同时，能较好地保持原始数据的分类效果，这在数据量庞大的情况下表现得尤其突出。

本发明提出的混合型粒子群优化算法可以实现快速可行的紧急控制措施搜索，可以在保证***稳定的同时，使付出的控制代价最小，得到的计算结果符合实际情况，整个计算过程自动完成。

附图说明

下面结合附图对本发明进一步说明。

图1是依据本发明的最优分类面。

图2示出了电力***稳定评估及控制流程。

图3示出了暂稳控制混合型粒子群优化算法计算粒子定义实例。

图4示出了确定稳定控制量的混合型粒子群优化算法流程。

图5示出了跨区电网局部区域示意。

图6示出了跨区电网稳控问题的混合型粒子群优化算法目标函数收敛曲线。

图7示出了跨区电网稳控问题的元件切除量收敛曲线。

图8示出了跨区电网无控制措施的功角曲线。

图9示出了跨区电网有控制措施的功角曲线。

具体实施方式

定义

暂态稳定评估电力***暂态稳定评估属于大规模非线性动力学***的动态分析问题。暂态稳定指互联电力***中的同步发电机在遭受严重扰动(如输电线路短路故障)后保持同步的能力。目前主要的研究方法有时域仿真法(或称逐步积分法)、暂态能量函数法(或称直接法)以及这两种方法的混合方法。

知识发现是从数据中识别出有效的、新颖的、潜在有用的，以及最终可理解的模式的高级过程。知识发现过程可分为：问题定义、数据准备和预处理、知识发现，以及结果的解释和评估。

数据预处理数据预处理就是在不损失数据信息的基础上，对数据进行优化处理，提高数据质量，满足知识发现算法对数据的特定要求，提高其后知识发现过程的精度和性能。

聚类分析聚类是把物理或抽象的数据对象按照相似性分成若干类别的过程，即“物以类聚”，它的目的是使属于同一类别的个体之间的距离尽可能的小而属于不同类别的个体间的距离尽可能的大。

K-L(Karhunen-Loeve)变换K-L变换是一种把原来多个指标转化为少数几个互不相关的综合指标的多元统计方法，可以达到数据化简、揭示变量间的关系和进行统计解释的目的，为进一步分析数据的性质和统计特征提供重要信息。

支持向量机支持向量机方法(Support Vector Machines，简称SVMs)是建立在统计学***面，使得超平面与不同类样本集之间的距离最大，即具有最大的泛化能力，从而作为分类未知样本的标准。具体可用图1的两维情况说明。

在图1中，方形点和圆形点分别代表两类样本，H为分类线，H₁、H₂分别为过各类中离分类线最近的样本且平行于分类线的直线，它们之间的距离r叫做分类间隔。所谓最优分类线就是要求分类线不但能将两类数据正确分开，而且使分类间隔最大。分类线方程可写为x·w+b＝0，对它进行归一化，使得对线性可分的样本集(x_i，y_i)，x_i(x_i∈Rⁿ))是输入特征向量，y_i(y_i∈R)是分类目标属性，i＝1，...，N，N是样本个数，满足：

y_i[(w·x_i)+b]-1≥0，i＝1，…，N                  (1)

上式中w为一次线性函数的斜率，b为一次线性函数的截距。此时分类间隔r＝2/||w||。使分类间隔r最大等价于使||w||²最小。满足条件(2-3)且使最小的分类面叫做最优分类面，H₁、H₂上的训练样本点被称作支持向量。

求解最优分类面的问题可以表示成如下的约束优化问题：

Min(Φ(w))＝Min(||w||²/2)＝Min((w·w)/2)         (2)

s.t.y_i(w·x_i+b)-1≥0，i＝1，...，N

定义Lagrange函数：

$L(w,b,\mathrm{\&alpha;})=(w\&CenterDot;w)/2-\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i[y_i(w\&CenterDot;x_i+b)-1]---\left(3\right)$

其中α_i＞0为Lagrange系数。利用Lagrange优化方法可以把上述最优分类面问题转化为其对偶问题(4)，即求解α₁...α_N，使

$\mathrm{Max}\left(Q\left(\mathrm{\&alpha;}\right)\right)=\mathrm{Max}(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i-\frac{1}{2}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i{\mathrm{\&alpha;}}_j{y}_i{y}_j(x_i\&CenterDot;x_j))---\left(4\right)$

$s.t.\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i{y}_i=0,{\mathrm{\&alpha;}}_i\&GreaterEqual;0,i=1,...,N$

求解上述问题后得到的最优分类函数是

$w^{*}=\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i^{*}{y}_i{x}_i---\left(5\right)$

b^＊＝y_k-w^T x_kx_k                                  (6)

$f\left(x\right)=\mathrm{sgn}\{(w^{*}\&CenterDot;x)+b^{*}\}=\mathrm{sgn}\{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i^{*}{y}_i(x_i\&CenterDot;x)+b^{*}\}---\left(7\right)$

sgn(·)为符号函数。b^*是分类阈值，可以用任一个支持向量(满足(1)中的等号)求得，或通过两类中任意一对支持向量取中值求得。

上面的方法在保证训练样本全部被正确分类，通过最大化分类间隔来获得最好的推广性能。考虑到存在某些数据不能满足式(1)的条件，即线性不可分，此时需在约束条件中加入一个松弛量ξi，问题变成：

Min(1/2w^Tw+CΣξi)                               (8)

s.t.yi(w^Tx_i+b)≥1-ξi，i＝1，...，N

ξi≥0，i＝1，...，N

C是某个指定的常数，起控制对错分样本惩罚的程度的作用。此处的优化问题的解与可分情况下的解几乎完全相同，只是条件变为

0≤α_i≤C，i＝1，...，N                       (9)

对非线性问题，可以通过满足Mercer条件的非线性变换K(x，x′)将输入空间变换到一个高维空间，然后在这个新空间中求取最优分类面。这种非线性变换是通过定义适当的内积函数实现的。此时式(4)的优化函数变为：

$\mathrm{Max}\left(Q\left(\mathrm{\&alpha;}\right)\right)=\mathrm{Max}(\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i-\frac{1}{2}\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i{\mathrm{\&alpha;}}_j{y}_i{y}_{j}K(x_i,x_j))---\left(10\right)$

相应的判别函数式(7)也变为

$f\left(x\right)=\mathrm{sgn}\{(w^{*}\&CenterDot;x)+b^{*}\}=\mathrm{sgn}\{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i^{*}{y}_{i}K(x_i\&CenterDot;x)+b^{*}\}---\left(11\right)$

式中内积核函数K(x_i，x)采用高斯径向核函数

$K(x_i,x)=\exp (-\frac{{|x-x_i|}^2}{{\mathrm{\&sigma;}}^2})---\left(12\right)$

式中σ为函数的宽度参数，控制了函数的径向作用范围。

暂态稳定控制在电力***遭受大扰动后，及时进行有效的暂态稳定控制，往往能以较小的代价保证***的安全稳定运行。常用的控制措施主要包括切机、切负荷、快关汽门、电气制动、强制补偿等。

粒子群优化粒子群优化方法(Particle Swarm Optimization，简称PSO)本质是一种多代理算法，研究由简单个体组成的群落与环境以及个体之间的互动行为。粒子群的迁移过程是有方向性的，搜索过程中用到反馈原理并可采用并行计算技术，具有较高的搜索效率，能以较大概率找到问题的全局最优解，计算效率比传统随机方法高、收敛速度快。

在PSO算法中，用粒子的位置表示待优化问题的解，每个粒子性能的优劣程度取决于待优化问题目标函数确定的适应值，每个粒子由一个速度决定其飞行方向和速率大小。设在一个d维的目标搜索空间中，有N个粒子组成一个群体，在第k次迭代时刻粒子i的位置表示为相应的飞行速度表示为开始执行PSO算法时，首先需随机初始化N个粒子的位置和速度，然后通过迭代寻找最优解。将粒子代入问题定义的目标函数，粒子的优劣是根据目标函数适应值的大小来衡量的。在每一次迭代中，每个粒子通过跟踪两个极值来更新自己在解空间中的位置与速度：一个极值是粒子本身在迭代过程中迄今搜索到的自身最优解，称为个体极值，其位置表示为pbest_i＝(pbest_i1，pbest_i2，...，pbest_id)，另一个极值是整个粒子群到目前为止找到的最优解，称为全局极值，其位置表示为gbest＝(gbest₁，gbest₂，...，gbest_d)。

在第k+1次迭代计算时，粒子x_i根据下列规则来更新自己的速度和位置：

(1)若x取值连续，

$v_i^{(k+1)}=w\&CenterDot;v_i^{\left(k\right)}+a_1\&times;\mathrm{rand}\left(\mathrm{0,1}\right)\&times;({\mathrm{pbest}}_i^{\left(k\right)}-x_i^{\left(k\right)})+a_2\&times;\mathrm{Rand}\left(\mathrm{0,1}\right)\&times;({\mathrm{gbest}}^{\left(k\right)}-x_i^{\left(k\right)})---\left(13\right)$

$x_i^{(k+1)}=x_i^{\left(k\right)}+v_I^{(k+1)}---\left(14\right)$

式中w是惯性权重，a₁和a₂为两个学习因子，rand(0，1)和Rand(0，1)是两个均匀分布在(0，1)之间的随机数，i＝1，2，...，m。在(14)中，速度v_i被看作是对位置的修正量，更为合适。另外，粒子速度v_i在每一维都被一个最大速度v_max所限制。如果当前粒子的加速度导致它在某一维的速度超过该维上的最大速度，则该维的速度被限制为最大速度。

(2)若x取值离散，粒子位置更新公式(14)被下面的式(15)所代替，

$\mathrm{If}(\mathrm{rand}\left(\right)<S\left(v_i^{(k+1)}\right))\mathrm{then}{x}_i^{(k+1)}=1$

$\mathrm{else}{x}_i^{(k+1)}=0---\left(15\right)$

在二进制PSO算法中，粒子的每一维x_id和个体极值的每一维pbest_id被限制为1或0，而对速度v_id不作这种限制。用速度来更新位置时，如果v_id高一些，粒子的位置x_id更有可能选1；反之，若v_id低一点，则x_id选0。式(15)中阈值在[0，1]之间的函数S(v)是sigmoid函数，S(v)＝1/(1+e^-v)。

图2是本发明稳定评估及控制流程图。本发明的稳定评估及控制方法包括如下步骤：

步骤一定义输入特征变量，建立电力***暂态稳定评估模型。电力***暂态稳定评估模型特征分为稳态特征和动态特征，动态特征分为功角特征、角速度特征、角加速度特征和功率、能量特征。本发明定义的特征变量如下：

步骤二数据采集，样本生成。本发明中，利用时域仿真对实际电网***模型设置各种预想事故，进行大量仿真计算，利用初始稳态潮流数据和故障发生后1～2秒有限时间内的动态信息，产生仿真样本，然后从中提取需要的特征量数据。

步骤三数据预处理。采集到的电网运行特征数据数据量纲不一致，各特征数据取值大小差别大。计算过程中取值小的特征易被取值大的数据淹没。另外，实测数据随时间累积，数据量大。通过数据预处理技术，可有效提高数据质量，降低后续计算所需时间。本发明主要采用如下数据预处理技术：

1.数据规范化处理将特征数据按比例缩放，使之落入一个小的特定区间(如0.0到1.0)。规范化处理可以有效防止具有较大初始值域的特征与具有较小初始值域的特征相比出现的权重过大问题，可以避免各特征变量的量纲和数量级对计算结果的影响。对原始样本数据，本发明按如下公式进行规范化处理：

$x=\frac{x_0-x_{0\max }}{x_{0\max }-x_{0\min }}---\left(16\right)$

式中x₀、x分别为转换前后的数据样本集合；x_0max为样本最大值；x_0min为样本最小值

2.数据归约对数据集进行压缩，得到小得多的数据表示，但能够产生同样的分析结果。一般地，对于以二维表形式存在的数据，分别从横向样本维度和纵向特征维度先后进行规约。

(1)样本规约本发明利用聚类分析技术对运行数据进行样本规约。可在保证分析结果不变的情况下，将原始大数据集分成若干较小的数据集，新生成的小数据集，各自包含的样本彼此相似，各小数据集之间彼此相异。各小数据集可独立并行处理，以提高计算效率。聚类问题可表示为一个迭代寻优问题。

设由n个样本、m个特征变量构成的规范化后的样本空间X＝{x_ij}，x_ij是样本j在特征变量i上的取值，i＝1，2，...，m，j＝1，2，...，n。设模糊聚类中心矩阵为S＝{s_hi}_c×m，s_hi是第h类样本中第i个变量的聚类中心，h＝1，2，...，c。

设将样本空间X聚成C类。最佳的聚类值C根据混合F统计值确定。

$\mathrm{Mixed}\_F=\frac{m}{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}1/F\left(i\right)}---\left(17\right)$

$F\left(i\right)=\frac{\underset{h=1}{\overset{c}{\mathrm{\&Sigma;}}}n{(s_{\mathrm{hi}}-{\stackrel{\&OverBar;}{s}}_i)}^2(n-c)}{\underset{h=1}{\overset{c}{\mathrm{\&Sigma;}}}\underset{j=1}{\overset{n_h}{\mathrm{\&Sigma;}}}{(x_{\mathrm{ijh}}-s_{\mathrm{hi}})}^2(c-1)}(i=\mathrm{1,2},...,m)---\left(18\right)$

式中n_h是第h类样本数；是第i个变量聚类中心的平均值；x_ijh是第h聚类中第j个样本的第i个变量值。Mixed_F取值越大，聚类效果越好。

定义模糊聚类矩阵U＝{u_hj}，u_hj为样本j归属于类别h的相对隶属度，满足

0≤u_hj≤1

$\underset{h=1}{\overset{c}{\mathrm{\&Sigma;}}}{u}_{\mathrm{hj}}=1---\left(19\right)$

$\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{u}_{\mathrm{hj}}>0$

设模糊聚类中心矩阵为S＝{s_ih}，s_ih是聚类中心h在特征变量i上的规范化取值，i＝1，2，...，m，h＝1，2，...，c。考虑到各特征变量对聚类作用的差异，设权重向量W＝(w₁，...，w_m)^T，且满足

$\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{w}_i=1---\left(20\right)$

模糊聚类问题可表示为一个迭代寻优问题：

$\min \{F\left(u_{\mathrm{hj}}\right)=\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}\underset{h=1}{\overset{c}{\mathrm{\&Sigma;}}}{\left[u_{\mathrm{hj}}\right|\left|w_i(x_{\mathrm{ij}}-s_{\mathrm{ih}})\right|\left|\right]}^2\}---\left(21\right)$

且满足约束条件(19)和(20)。上述问题的解，可通过构造拉格朗日函数，得到U阵和S阵的迭代表达式：

$u_{\mathrm{hj}}=\frac{1}{\underset{k=1}{\overset{c}{\mathrm{\&Sigma;}}}\frac{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\left[w_i(x_{\mathrm{ij}}-s_{\mathrm{ih}})\right]}^2}{\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\left[w_i(x_{\mathrm{ij}}-s_{\mathrm{ik}})\right]}^2}}---\left(22\right)$

$s_{\mathrm{ih}}=\frac{\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{u}_{\mathrm{hj}}^2{w}_i^2{x}_{\mathrm{ij}}}{\underset{j=1}{\overset{n}{\mathrm{\&Sigma;}}}{u}_{\mathrm{hj}}^2{w}_i^2}---\left(23\right)$

对于暂态稳定评估而言，是一个两(或多)模式分类问题，各特征变量对聚类分析的作用程度具有内在规律性。在信息论中，具有最高信息增益或最大熵压缩的特征能使分类所需的信息量最小，对分类问题的贡献最大，本发明以基于熵的信息增益作为特征权重的度量，可以避免权重W是人为给定造成的计算结果受主观因素影响的问题。

设X是包含n个数据样本的集合，在暂稳评估问题中，全体样本可分为稳定与不稳定两类，记为C_i(i＝1，2)。n_i是X中属于类C_i的样本数，特征变量A有v个不同的取值{a₁，a₂，...，a_v}，可以用属性A将X划分为v个子集{X₁，X₂，...，X_v}，其中，X_j包含X中这样的样本，它们在属性A上取值为a_j。特征变量A的信息增益表示为：

$W_A=\mathrm{Gain}\left(A\right)$

$=-\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}{p}_i{\log }_2\left(p_i\right)-\underset{j=1}{\overset{v}{\mathrm{\&Sigma;}}}\frac{n_{1,j}+n_{2,j}}{n}(-\underset{i=1}{\overset{2}{\mathrm{\&Sigma;}}}{p}_{\mathrm{ij}}{\log }_2\left(p_{\mathrm{ij}}\right))---\left(24\right)$

式中p_i是任意样本属于C_i的概率，p_ij是子集X_j中的样本属于类C_i的概率， |·|表示集合的样本数，n_ij是子集X_j中属于类C_i的样本数。

求解上述优化问题，得到模糊聚类中心矩阵S。对于每个运行样本x_i(i＝1，...，n)，根据其与各聚类中心矢量的距离，将其划分到与之距离最小的聚类中心对应的聚类子集中。

(2)维度规约删除或减少不相关的特征属性，达到提高计算效率的目的。本发明利用改进的非线性K-L变换法进行暂态稳定评估的维度规约。K-L变换是一种把原来多个指标转化为少数几个互不相关的综合指标的多元统计方法，可以达到数据化简、揭示变量间的关系和进行统计解释的目的。

在实际问题中，特征变量之间以及综合指标与原始数据之间常呈非线性关系，因此有必要对传统K-L变换中的“线性化”进行改进。对于由n个样本，m个特征构成的数据集X＝(x_ij)_n×m，本发明采用“对数中心化”的非线性K-L变换法，主要做法如下：

(1)对原始数据x_ij进行对数变换，变换后的数据集记为W＝{w_ij}，有

w_ij＝lnx_ij                            (25)

(2)行向量中心化，令

$y_{\mathrm{ij}}=w_{\mathrm{ij}}-\frac{1}{m}\underset{j=1}{\overset{m}{\mathrm{\&Sigma;}}}{w}_{\mathrm{ij}}---\left(26\right)$

(3)计算矩阵Y＝(y_ij)_n×m的协方差矩阵V，求解V的特征值和对应的标准正交特征向量，根据贡献率，确定满足累计贡献率的前p(p＜m)个综合指标；

(4)建立综合指标方程，计算综合指标值，形成新的数据集

F_i＝Γ_iW i＝1，2，...，p            (27)

上式中，Γ_i是与上述协方差阵第i个特征值λ_i对应的标准正交特征向量。

利用K-L变换选取的输入特征变量，可将数据集的维度由m降到p(p＜m)，有助于提高暂态稳定评估的计算效率。

步骤四暂态稳定评估本发明采用支持向量机技术进行面向暂态稳定评估的知识发现。暂态稳定评估的目的是根据收集到的***运行数据，确定对应的***运行状态是稳定的还是不稳定的。对收集到的历史运行数据，首先利用支持向量机算法，确定区分样本空间中稳和不稳运行点的最优分类线(或面)数学表达式(11)。

对于运行状态未知的运行数据x，将其代入式(11)，即可得到其运行状态。支持向量机的输出代表了测试数据到最有超平面的距离，输出越大表明到超平面距离越远。令支持向量机为正时表示***稳定，则该数越大***越稳定，反之亦然。

步骤五暂态稳定控制从数学的角度看，暂态稳定控制属于离散控制范畴；从最优化学科的角度看，可归结为一个点极值问题，即对提高电力***暂态稳定性这一控制目标，如何给出最佳的控制量。

(1)本发明所解决的稳控问题的数学模型可描述为：

1)目标函数

Max(f(稳定裕度，MVA_G，MVA_L))                         (28)

目标函数f是***稳定裕度、切机量MVA_G和切负荷量MVA_L的函数，是根据稳控量和实施稳定控制措施后的***状态确定的。根据以上暂态稳定评估程序获得***运行状态，在***不稳定时，目标函数侧重于稳定裕度的变化；当***状态恢复到稳定后，目标是保证***稳定的同时，找到切除量最小的控制方式，故在这时的目标函数中计及切机切负荷量的因素。目标函数定义如式(29)和式(30)：

***不稳定时：

f＝-abs(稳定裕度)                                    (29)

***稳定时：

MVA_Gi＝P_gi×r_gi，i＝1，...，N_G

MVA_Lj＝P_lj×r_lj，j＝1，...，N_L

上式中weight是为了平衡稳定裕度取值和控制量取值的权重指数，一般取值为100；MVA_G是发电机的切机量，MVA_L是切负荷量，P_gi是被切发电机的有功出力，P_lj是被切负荷的有功负荷。r_gi是发电机减出力比例，r_lj是切负荷比例。

2)约束条件

r_gi＝0或1                                            (31)

0≤r_lj≤1                                            (32)

在寻优之前，需定义待求问题的解，即PSO算法中的粒子。对于暂稳控制量确定问题，每个粒子被表示成为一个(NG+NL)维的向量。第1～N_G维是发电机维，第N_G+1～N_G+N_L维是负荷维。因为在实际的紧急控制中，发电机只有完全切和完全不切两种状态，没有切多少的问题，所以这N_G维数据的取值是非0即1。后N_L维数据代表N_L个待切除负荷的切除比例，每一维的取值是一个[0，1]之间的实数。

(2)暂态稳定控制过程

1)定义待切除的发电机和负荷。原则上，***中的所有发电机和负荷都可作为被控对象，但这样做，一方面不符合实际情况，另一方面也会造成巨大的计算负担。为此必须针对故障点的位置，明确切除元件的范围。原则上是切除故障线路的送端侧发电机和受端侧负荷。

2)初始化N个输入粒子每个粒子是一个(N_G+N_L)维矢量。在允许的范围内，随机设置各粒子的初始位置和速度

3)一个粒子对应一种控制措施，对应暂稳计算程序中的扰动输入文件(如PSASP软件中的ST.S12)。根据每个粒子的内容，写扰动输入文件，运行暂稳程序，判断加入这个粒子对应的控制措施后，***是否稳定。

4)根据判稳结果，按式(29)或(30)计算每个粒子对应的目标函数。

5)当k＝1时，即第一次迭代，每个粒子的坐标设置为其当前位置，计算出其相应的个体极值(即个体极值点的目标函数值)；当k＞1时，则计算粒子的目标函数值，如果好于该粒子当前的个体极值，则将pbest设置为该粒子的位置，即并更新当前的个体极值；对每个粒子，当k＝1时，即第一次迭代，全局极值(即全局极值点的目标函数值)是所有个体极值中最好的，gbest设置为该最好粒子的当前位置；当k＞1时，如果所有粒子的个体极值中最好的好于当前的全局极值，则将gbest设置为该粒子的位置，并更新当前的全局极值；

6)根据式(13)更新每个粒子的速度，根据式(14)和(15)更新每个粒子的位置；如果 则如果则

7)判读问题的循环次数是否超过了最大的次数限制，如果超过，则停止迭代，否则转向步骤3)；

8)输出全局极优值，这个值就是各机和各负荷的应切除量。

以下是本发明的一个实施例，以某跨区电网进行仿真实验做实施例，进一步说明如下：

(1)样本生成

采用某跨区电网年度运行方式计算用实际数据，电网规模为4936个节点，660个发电机母线，1646个负荷母线，4290条交流线路和3条直流线路。在某一运行方式下，计算各条线路的传输功率，取传输功率最大的前10条线路作为关键线路；同样地，取出力最大的前8台发电机作为关键发电机。合计定义了86个静态特征变量和42个暂稳特征变量。现场提供了15种潮流运行方式，给定的故障列表中包含了39组三相短路故障，通过切除故障线路切除故障。随机取14种潮流运行方式，计算其对应的14×39＝546个暂稳作业的80秒全过程仿真，提取故障前后的静/动态特征变量和对应的暂稳计算结果(稳定状态)，形成训练数据。利用剩下的一组运行方式数据，计算其对应的39个暂稳作业的80秒全过程仿真，获取相应的特征向量值，构成测试数据。值得一提的是，测试数据不参与SVMs的训练过程，属于新的未知的运行数据。

(2)数据预处理

1)按照式(16)，对收集到的样本数据进行规范化处理；

2)样本规约设聚类数C从2变化至10，计算Mixed_F，结果见表1。

表1不同聚类数对应的Mixed_F值

从表1中可以看出，当聚类数为3时，Mixed_F对应的值最大，所以可将训练集样本聚成3类。聚类后，每一类中的样本数分别为147、187和212，样本空间分别是原数据集的73.08％，65.75％和6117％。对于测试集中的每个运行样本，根据其与各聚类中心矢量的距离，划分到与之距离最小的聚类中心对应的聚类子集中，每一组测试集中的样本数分别为7、16和16，这样就形成了3组训练-测试子集。

3)维度规约对经过样本规约的数据接着进行维度规约。在非线性法中，各数据子集中的特征数目由128降至25，31和28个。

(3)暂态稳定评估

对3组训练子集(分别包含147、187和212个训练样本)进行训练，各自构造暂稳评估预测模型。3组训练子集模型分别12、100、14个稳定支持向量和6、11、7个不稳定支持向量组成。支持向量是在预测时真正起作用的数据，它的数量远小于训练样本的数量。对三组子集的训练-测试可以并行执行，取其中最耗时一组的计算时间作为最终计算时间。

较之时域仿真法，基于知识发现技术的暂态稳定评估，最大的优势在于其计算速度快。利用两种方法，先后对国调数据进行仿真计算，结果如表2和表3所示：

表2进行数据预处理后的稳定预测结果

表3知识发现法和时域仿真法计算时间比较

从表2和表3可以看出，得到***运行是否稳定的结论，时域仿真计算需仿真40秒动态过程，计算耗时116.64秒；若应用知识发现方法，仅需13.812秒。预测准确性见表4。

表4基于SVMs模型的暂态稳定预测结果

SVMs值为正表示当前暂稳作业稳定，为负表示不稳定，且SVMs值的绝对值越大，表示对应的暂稳作业距SVMs分类超平面的距离越大；稳定类别1表示稳定，2表示不稳定。在暂稳评估问题中，误分指将稳定样本判定为不稳定，漏分指把不稳定样本判定为稳定。如表4所示：有1个作业被漏分。漏分作业的SVMs值等于0.5008，接近超平面，即y＝0.0。对那些处于临界稳定边界的样本，可用挑选出来，利用时域仿真法进行详细的暂态稳定分析，以确保全体数据最终分析结果的可靠性。

(4)暂态稳定控制

图5是跨区电网局部区域示意图，图中的圆圈表示各发电区域或用电区域，有向线段表示区域间的联络线及其潮流方向。其中区域1和2之间通过直流线3600连接。设在区域1和10之间的联络线54和58上，靠近区域10处，发生三相短路接地永久性故障，故障发生后0.2秒，断开线路：

区域10和区域11之间的联络线53和57；

区域1和区域7之间的联络线56；

区域1内的线路59。

仿真计算表明，开断上述线路并不能使***恢复稳定。现采用混合PSO算法搜索能保证***稳定的控制措施。

1)控制区域确定

对于实际运行大***，若将***中的所有发电机和负荷都作为切除对象，即PSO算法中粒子的每一维数据对应***中的每一台发电机和负荷，将会降低算法的执行效率。在进行PSO搜索之前，必须明确搜索空间的范围，并且这个范围要尽可能的小。在本发明中，搜索空间由被切除的元件(造成功率不平衡的送端发电机和受端负荷)组成。在本算例中，最终确定功率送端区域(11，33，34，37)的发电机和功率受端区域2的负荷母线构成PSO算法的初始搜索空间。

2)PSO算法数据初始化

初始考虑不切负荷，如果仅靠切除发电机还不能使***恢复稳定，则再考虑切除负荷。因此，计算粒子中，发电机对应位初始置1，负荷对应位初始置0

混合PSO算法中的控制参数选择为：粒子数N＝20，学习因子a₁＝a₂＝2.0，速度限值v_max＝1.0，v_min＝-1.0，迭代次数＝100。惯性权重w采用式(33)的定义方式：

$w=w_{\max }-\frac{w_{\max }-w_{\min }}{{\mathrm{iter}}_{\max }}\&times;\mathrm{iter}---\left(33\right)$

式中w_max＝0.9，w_min＝0.4，iter是当前的迭代次数，最大迭代数设为100。

3)稳定控制策略确定

执行所确定步骤，得到使电网恢复稳定的最终切机和切负荷量方案。设定所有的稳控措施在故障发生0.22秒之后动作。表5给出了运行PSO算法得到的最终稳控结果。

表5使跨区电网***恢复稳定的切机和切负荷量

图6表示目标函数迭代收敛曲线，图7表示元件切除量收敛曲线。从图6和图7中可见，利用混合PSO算法能很快地确定紧急控制调节方案，在本实施例中，对跨区电网***迭代10次就能收敛到稳定值。图8和图9表示加入控制措施前后，发电机功角的时域响应结果。从图9可见，加入稳控措施后，所在***能在2秒内恢复稳定，不会发生如图8所示的功角失稳现象。在本发明的PSO算法中，每一次迭代中的每一个粒子代表一种稳控措施，它的适应度函数值，也就是稳控措施对***稳定性的影响，都要通过暂态稳定计算得到。对于一个象跨区电网这样的大***来说，得到最终最优解，需要多次调用暂态稳定计算程序。如果每次都进行40秒全过程时域仿真，那么整个计算所需的时间将是非常可观的。为了提高计算效率，在对跨区电网的稳控仿真中，本发明采用SVMs法取代时域仿真，进行运行状态的快速判断。将得到的SVMs值作为稳定裕度代入式(29)或(30)，计算适应度目标函数，作为进一步迭代的依据。每个粒子，每次迭代，若用时域仿真法判稳，需用34秒；若用SVMs法判稳，仅需不到2秒的时间。设粒子群大小为15，进行10次迭代，获得最终解，运行时域仿真耗时34秒×15×10＝5100秒＝85分钟，而运行SVMs法需要2秒×15×10＝300秒＝5分钟。

Claims (4)

1.一种电力***稳定评估及控制方法，其特征在于包括以下步骤：

(1)确定电力***的输入特征变量，构建电力***的暂态稳定评估模型；

(2)获取电力***中电网在线运行数据及其在给定故障集下的短时暂态稳定仿真计算结果，根据步骤(1)确定的输入特征变量，构成训练样本空间，每个样本均由输入特征变量集合x₀和对应的暂稳计算结果y₀两部分构成；

(3)在不损失数据信息的基础上，对初始海量数据样本进行数据预处理，得到新的数据样本集合(x,y)；

(4)根据训练样本集合(x,y)，构造暂态稳定最优分类判别函数：

$f\left(x\right)=\mathrm{sgn}\{(w^{*}\&CenterDot;x)+b^{*}\}=\mathrm{sgn}\{\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{\mathrm{\&alpha;}}_i^{*}{y}_{i}K(x_i\&CenterDot;x)+b^{*}\}---\left(1\right)$

其中sgn(·)为符号函数，w^*为一次线性函数的斜率，b^*为一次线性函数的截距，为与每个样本对应的Lagrange乘子，y_i是暂稳计算结果变量，K(x_i,x)为高斯径向核函数，有

$K(x_i,x)=\exp (-\frac{{|x-x_i|}^2}{{\mathrm{\&sigma;}}^2})$

其中σ为函数的宽度参数,控制了函数的径向作用范围；

(5)对于电力***运行状态未知的运行数据x，将其按步骤(3)经过数据预处理后，代入式(1)，即得到其对应的暂态稳定运行状态f(x)；f(x)的值代表了每个样本到最优分类超平面的距离，在暂态稳定评估问题中，f(x)的值代表了当前运行状态x到暂态稳定分类面的距离，即暂态稳定裕度，f(x)取值为正，表示***稳定，且数值越大***越稳定；反之，取值为负表示***暂态不稳定，且该数绝对值越大表示***越不稳定；

(6)对步骤(5)判定为失稳的电力***运行方式，利用混合型粒子群优化算法，搜索并得到能使电力***恢复暂态稳定运行且切机或切负荷量最小的紧急控制措施。

2.如权利要求1所述的方法，其特征是确定能反映电力***发生故障时暂态特性的输入特征x₀及其对应的暂稳计算结果y₀，构建电力***暂态稳定评估模型；步骤(1)采用的输入特征变量x₀分为静态特征变量和动态特征变量，也即x₀＝(静态特征变量，动态特征变量)，具体包括：

(1)静态特征变量包括

1)整个电力***的发电机的出力和负荷；

2)电力***的最大和最小的发电出力；

3)电力***的最大和最小的负荷；

4)母线电压的最大和最小值；

5)线路最大有功传输功率，及其对应的无功传输功率；

6)关键输电断面和区域联络断面的传输功率；

7)关键发电机的发电出力；

8)各区域的发电机的出力和负荷；

(2)动态特征变量包括：

以下各式中，δ_i为发电机转子角；ω_i为发电机转速，a_i为发电机转子角加速度，M_i为发电机的惯性常数；P_mi为发电机机械功率；P_ei为发电机电磁功率；V_ki为发电机动能；V_pi为发电机势能；N_G为发电机数目，t₀为故障发生初始时刻，t_cl为故障切除时刻；

1)故障发生初始时刻具有最大功角加速度的发电机的相对转子角度：

$x_1=({\mathrm{\&delta;}}_{a\max }-{\mathrm{\&delta;}}_{\mathrm{COI}}){|}_{t_0}$

上式中是相对于***惯性中心的等值转子角，δ_amax是故障发生初始时刻具有最大功角加速度的发电机的转子角；

2)故障切除时刻具有最大动能的发电机的相对转子角度：

$x_2=({\mathrm{\&delta;}}_{p\max }-{\mathrm{\&delta;}}_{\mathrm{COI}}){|}_{t_{\mathrm{cl}}}$

上式中δ_pmax是故障切除时刻具有最大动能的发电机的转子角；

3)故障切除时刻领前机与殿后机的功角之差；故障发生初始时刻，具有最大转子角加速度的发电机为领前机；具有最小转子角加速度的发电机为殿后机；

$x_3=\max \left({\mathrm{\&delta;}}_{a\max }\right)-\min \left({\mathrm{\&delta;}}_{a\min }\right){|}_{t_{\mathrm{cl}}}$

上式中δ_amax、δ_amin分别是领前机和殿后机的转子角；

4)故障切除时刻最大的发电机转角差：

$x_4=\max \left({\mathrm{\&delta;}}_i\right)-\min \left({\mathrm{\&delta;}}_i\right){|}_{t_{\mathrm{cl}}}$

5)故障切除时刻与惯性中心COI相差最大的功角：

$x_5=\max ({\mathrm{\&delta;}}_i-{\mathrm{\&delta;}}_{\mathrm{COI}}){|}_{t_{\mathrm{cl}}}$

6)故障切除时刻领前机群和余下机群中的发电机的惯量加权平均角度之差：

其中，S表示领前机群集合，A表示余下机群集合；对指定时刻的各发电机转角进行降序排序，形如：

且 ${\mathrm{\&delta;}}_1>{\mathrm{\&delta;}}_2>{\mathrm{\&delta;}}_3>\&CenterDot;\&CenterDot;\&CenterDot;>{\mathrm{\&delta;}}_{N_G-2},{\mathrm{\&delta;}}_{N_G-1}>{\mathrm{\&delta;}}_{N_G}$

每两个相邻转角之差构成一个角度间隙，形如：

且

${\mathrm{\&delta;}}_{\mathrm{1,2}}={\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_2,{\mathrm{\&delta;}}_{\mathrm{2,3}}={\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_3,\&CenterDot;\&CenterDot;\&CenterDot;,{\mathrm{\&delta;}}_{N_G-2,N_G-1}={\mathrm{\&delta;}}_{N_G-2}-{\mathrm{\&delta;}}_{N_G-1},{\mathrm{\&delta;}}_{N_G-1,N_G}={\mathrm{\&delta;}}_{N_G-1}-{\mathrm{\&delta;}}_{N_G}$

确定最大的角度间隙，此间隙之上的发电机构成领前机群S，此间隙之下的发电机构成余下机群A；

7)故障切除时刻与故障发生初始时刻发电机相对转子角度差值的绝对值之和：

$x_7=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}|{\mathrm{\&theta;}}_{\mathrm{ci}}-{\mathrm{\&theta;}}_{0i}|$

上式中，θ_i＝δ_i-δ_COI，θ_ci为故障切除时刻发电机相对转子角，θ_0i为故障发生初始时刻发电机相对转子角；

8)故障切除时刻与故障发生初始时刻发电机相对转子角度差值的绝对值的最大值：

$x_8=\underset{i}{\max }\left(\right|{\mathrm{\&theta;}}_{\mathrm{ci}}-{\mathrm{\&theta;}}_{0i}\left|\right)$

9)故障切除时刻与故障发生初始时刻惯性中心坐标COI转速的变化：

$x_9={\mathrm{\&omega;}}_{\mathrm{COI}}{|}_{t_{\mathrm{cl}}}-{\mathrm{\&omega;}}_{\mathrm{COI}}{|}_{t_0}$

上式中 ${\mathrm{\&omega;}}_{\mathrm{COI}}=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{M}_i{\mathrm{\&omega;}}_i/\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{M}_i$

10)故障切除时刻发电机转子角速度的绝对值之和：

$x_{10}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{|{\mathrm{\&omega;}}_i-{\mathrm{\&omega;}}_{\mathrm{COI}}|}_{t_{\mathrm{cl}}}$

11)故障切除时刻发电机转子角速度绝对值的最大值：

$x_{11}=\max \left(\right|{\mathrm{\&omega;}}_i-{\mathrm{\&omega;}}_{\mathrm{COI}}\left|\right){|}_{t_{\mathrm{cl}}}$

12)所有发电机转子在故障发生初始时刻相对加速度的最大、小值：

$x_{12-1}=(a_{\max }-a_{\mathrm{COI}}){|}_{t_0}$

$x_{12-2}=(a_{\min }-a_{\mathrm{COI}}){|}_{t_0}$

其中是a_COI惯性中心加速度a_max、a_min分别是转子角加速度的最大、小值；

13)所有发电机转子在故障发生初始时刻加速度的均值：

$x_{13}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{a}_i{|}_{t_0}$

14)所有发电机转子在故障发生初始时刻加速度的均方根误差：

$x_{14}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{(a_i-x_{13})}^2}{|}_{t_0}$

15)所有发电机转子在故障发生初始时刻相对加速度的平均值：

$x_{15}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}(a_i-a_{\mathrm{COI}}){|}_{t_0}$

16)所有发电机转子在故障发生初始时刻相对加速度的方差：

$x_{16}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{((a_i-a_{\mathrm{COI}})-x_{15})}^2}{|}_{t_0}$

17)领前机与殿后机在故障瞬间和故障发生初始时刻加速度之差：

$x_{17-1}=\max \left(a_i\right)-\min \left(a_i\right){|}_{t_0}$

$x_{17-2}=\max \left(a_i\right)-\min \left(a_i\right){|}_{t_{\mathrm{cl}}}$

18)故障切除时刻所有转子动能的最大值：

$x_{18}=\max \left(V_{\mathrm{ki}}\right){|}_{t_{\mathrm{cl}}}i=\mathrm{1,2},...,N_G$

19)故障切除时刻具有最大转子角度的发电机p的转子动能：

$x_{19}=\frac{1}{2}{M}_p\&times;({{\mathrm{\&omega;}}_p}^2-1){|}_{t_{\mathrm{cl}}}$

20)电力***总的“能量调整”：

$x_{20}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}(P_{\mathrm{ai}}\&times;{\mathrm{\&delta;}}_{\mathrm{di}})$

其中，P_ai是第i台发电机在故障发生初始时刻的加速功率；是第i台发电机故障切除时的相对转子角度；

21)故障切除时刻所有转子动能的平均值：

$x_{21}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\frac{1}{2}{M}_i\&times;({{\mathrm{\&omega;}}_i}^2-1){|}_{t_{\mathrm{cl}}}$

22)故障前所有发电机机械输入功率的平均值：

$x_{22}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{p}_{\mathrm{mi}}$

23)故障前所有发电机机械输入功率的方差：

$x_{23}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{(P_{\mathrm{mi}}-x_{22})}^2}$

24)发电机在故障发生初始时刻相对加速功率的均值：

$x_{24}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}(P_{\mathrm{mi}}-P_{\mathrm{ei}}){|}_{t_{\mathrm{cl}}}$

25)发电机在故障发生初始时刻相对加速功率的方差：

$x_{25}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{(\mathrm{\&Delta;}{P}_i-x_{24})}^2}{|}_{t_0}$

26)发电机在故障发生初始时刻相对加速功率的均值：

$x_{26}=\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\frac{\mathrm{\&Delta;}{P}_i}{P_{\mathrm{mi}}}{|}_{t_0}$

27)发电机在故障发生初始时刻相对加速功率的方差：

$x_{27}=\sqrt{\frac{1}{N_G}\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{(\frac{\mathrm{\&Delta;}{P}_i}{P_{\mathrm{mi}}}-x_{26})}^2}{|}_{t_0}$

28)故障切除时刻发电机的总动能：

$x_{28}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\left|V_{\mathrm{ki}}\right|{|}_{t_{\mathrm{cl}}}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\left|\frac{1}{2}{M}_i{{\mathrm{\&omega;}}_i}^2\right|{|}_{t_{\mathrm{cl}}}$

29)故障切除时刻发电机的最大动能差：

$x_{29-1}=(\max \left(V_{\mathrm{ki}}\right)-\min \left(V_{\mathrm{ki}}\right)){|}_{t_{\mathrm{cl}}}$

$x_{29-2}=(\max \left({\stackrel{\&CenterDot;}{V}}_{\mathrm{ki}}\right)-\min \left({\stackrel{\&CenterDot;}{V}}_{\mathrm{ki}}\right)){|}_{t_{\mathrm{cl}}}$

30)故障发生初始时刻单台发电机所受的最大有功功率冲击：

$x_{30}=\max (P_{i1}-P_{i0}){|}_{t_0},i=\mathrm{1,2},...,N_G$

式中P_i1，P_i0分别为发电机i在故障瞬间时和故障前的有功输出；

31)故障发生初始时刻单台发电机所受的最小有功功率冲击：

$x_{31}=\min (P_{i1}-P_{i0}){|}_{t_0},i=\mathrm{1,2},...,N_G$

32)故障发生初始时刻单台发电机所受的最大归一化有功功率冲击：

$x_{32}=\underset{i}{\max }\frac{P_{i1}-P_{i0}}{P_{\mathrm{ie}}}{|}_{t_0},i=\mathrm{1,2},...,N_G$

式中P_ie为发电机i的额定有功；

33)故障发生初始时刻单台发电机所受的最小归一化有功功率冲击：

$x_{33}=\min \frac{P_{i1}-P_{i0}}{P_{\mathrm{ie}}}{|}_{t_0},i=\mathrm{1,2},...,N_G$

34)表征故障对***冲击大小的性能特征变量：

$x_{34}=\left(\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{M}_i\right|P_{\mathrm{di}}\left|\right)/\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{M}_i$

上式中 $P_{\mathrm{di}}=(P_{\mathrm{ci}}-P_{0i})/M_i-\left[\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}(P_{\mathrm{ci}}-P_{0i})\right]/\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{M}_i,$ 是每台发电机相对于***惯性中心的减速功率；P_ci是故障切除瞬间发电机输出功率；P_0i是故障前发电机输出功率；

35)事故清除前一瞬间发电机加速功率之和：

$x_{35}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{P}_{\mathrm{acc},i}^{\mathrm{BF}}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}(P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}})=P_{\mathrm{COI}}^{\mathrm{BF}}$

为事故清除前一瞬间发电机i的加速功率；为事故清除前一瞬间惯性中心的加速功率；分别为事故清除前一瞬间发电机的机械功率和电磁功率；

36)事故清除后一瞬间发电机加速功率绝对值的最大值：

$x_{36}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{PF}}\left|\right\}=\underset{i}{\max }\left\{\right|P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}\left|\right\}$

式中，分别为事故清除后一瞬间发电机的机械功率和电磁功率；

37)事故清除前一瞬间发电机的加速功率与发电机惯性时间常数比值之和：

$x_{37}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{P}_{\mathrm{acc},i}^{\mathrm{BF}}/M_i=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\left(\frac{P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}}}{M_i}\right)$

38)事故清除前一瞬间发电机的加速功率与发电机惯性时间常数比值的最大值：

$x_{38}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{BF}}/M_i\left|\right\}=\underset{i}{\max }\left\{\right|\frac{P_{\mathrm{mi}}^{\mathrm{BF}}-P_{\mathrm{ei}}^{\mathrm{BF}}}{M_i}\left|\right\}$

39)事故清除后一瞬间发电机的加速功率与发电机惯性时间常数比值之和：

$x_{39}=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}{P}_{\mathrm{acc},i}^{\mathrm{PF}}/M_i=\underset{i=1}{\overset{N_G}{\mathrm{\&Sigma;}}}\left(\frac{P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}}{M_i}\right)$

40)事故清除后一瞬间发电机的加速功率与发电机惯性时间常数比值的最大

值：

$x_{40}=\underset{i}{\max }\left\{\right|P_{\mathrm{acc},i}^{\mathrm{PF}}/M_i\left|\right\}=\underset{i}{\max }\left\{\right|\frac{P_{\mathrm{mi}}^{\mathrm{PF}}-P_{\mathrm{ei}}^{\mathrm{PF}}}{M_i}\left|\right\}.$

3.如权利要求1所述的方法，其特征在于其中所述的海量数据预处理过程包括以下步骤：

(1)数据规范化处理

$x=\frac{x_0-x_{0\max }}{x_{0\max }-x_{0\min }}$

式中x₀、x分别为转换前后的数据样本集合；x_0max为样本最大值；x_0min为样本最小值；

(2)样本规约

1)根据混合F统计值，确定最佳聚类数C；

2)计算每个特征变量对稳定结果的信息增益，作为聚类分析中的特征权重W_i；对规范化的样本数据进行聚类分析，将每个样本划分到与其距离最小的聚类中心对应的聚类子集中，形成C个数据样本子集；

(3)维度规约

对样本数据进行对数中心化变换和行向量中心化，通过非线性K-L变换，将样本集的维度、即特征变量数，由最初的m降到，p<m；

1)对原始数据X＝(x_ij),i＝1,...,n,j＝1,...,m，n为样本数，m为输入特征变量数，进行对数变换，变换后的数据集记为W＝{w_ij}，有

w_ij＝lnx_ij

2)行向量中心化，令

$y_{\mathrm{ij}}=w_{\mathrm{ij}}-\frac{1}{m}\underset{j=1}{\overset{m}{\mathrm{\&Sigma;}}}{w}_{\mathrm{ij}}$

3)计算矩阵Y＝(y_ij)_n×m的协方差阵及其特征根和与特征值对应的标准正交特征向量

4)对特征值进行降序排序，有λ₁＞λ₂＞…＞λ_m；

5)根据累计贡献率选择累计贡献率>累计贡献率阀值的前p(p<m)个特征值，一般的，累计贡献率阀值取0.85～0.95

6)K-L变换，形成新的数据样本集合

$F_i={\mathrm{\&Gamma;}}_{i}W,i=\mathrm{1,2},...,p$

Γ_i是与上述协方差阵第i个特征值λ_i对应的标准正交特征向量。

4.如权利要求1所述的方法，其特征在于还包括下列步骤：

(1)构造电力***紧急控制目标函数：

Max(f(稳定裕度,MVA_G,MVA_L))         (2)

当***不稳定时：

f＝－abs(稳定裕度)             (3)

***稳定时：

MVA_Gi＝P_gi×r_gi i＝1,...,N_G

MVA_Lj＝P_lj×r_lj j＝1,...,N_L

上式中weight是为了平衡稳定裕度取值和控制量取值的权重指数，一般取值为100；MVA_G是发电机的切机量，MVA_L是切负荷量，P_gi是被切发电机的有功出力，P_lj是被切负荷的有功负荷；r_gi是发电机切出力比例，r_lj是切负荷比例，满足：

r_gi＝0或1           (5)

0≤r_lj≤1           (6)

(2)定义待切除的发电机和负荷，包括故障线路送端侧的N_G台发电机和受端侧的N_L个负荷；

(3)初始化N个输入粒子每个粒子是一个(N_G+N_L)维矢量，前N_G维数据代表N_G台待切发电机的状态r_g，在实际的紧急控制中，发电机只有完全切和完全不切两种状态，这个N_G维数据的取值是非0即1；后N_L维数据代表N_L个待切除负荷的切除比例r_l，每一维的取值是一个[0,1]之间的实数；在允许的范围内，随机设置各粒子的初始位置和速度

(4)一个粒子对应一种控制措施，利用电力***分析软件进行暂态稳定计算，判断加入这个粒子对应的控制措施后，***是否稳定；根据判稳结果，按式(3)或(4)计算每个粒子对应的目标函数；

(5)利用混合粒子群法搜寻最优的控制方案；对第k次迭代，根据式(7)更新每个粒子的速度，对每个粒子的的第1～N_G发电机维，根据式(8)更新对应维度位置；对每个粒子的的第(N_G+1)～(N_G+N_L)负荷维，根据式(9)更新对应维度位置取值

$v_i^{(k+1)}=w\&CenterDot;v_i^{\left(k\right)}+a_1\&times;\mathrm{rand}\left(\mathrm{0,1}\right)\&times;({\mathrm{pbest}}_i^{\left(k\right)}-x_i^{\left(k\right)})+a_2\&times;\mathrm{Rand}\left(\mathrm{0,1}\right)\&times;({\mathrm{gbest}}^{\left(k\right)}-x_i^{\left(k\right)})---\left(7\right)$

$\mathrm{If}(\mathrm{rand}\left(\right)<S\left(v_i^{(k+1)}\right))\mathrm{then}{x}_i^{(k+1)}=1$

$\mathrm{else}{x}_i^{(k+1)}=0---\left(8\right)$

$x_i^{(k+1)}=x_i^{\left(k\right)}+v_i^{(k+1)}---\left(9\right)$

式中，为第k、k+1次迭代时刻粒子i的位置，为第k、k+1次迭代时刻粒子i的飞行速度，为个体极值，gbest^(k)为全局极值，S(v)是sigmoid函数，w是惯性权重，a₁和a₂为两个学习因子，rand(0,1)和Rand(0,1)是两个均匀分布在(0,1)之间的随机数，i＝1,2,…,m，w定义方式如下：

$w=w_{\max }-\frac{w_{\max }-w_{\min }}{{\mathrm{iter}}_{\max }}\&times;\mathrm{iter}---\left(10\right)$

其中w_max＝0.9，w_min＝0.4，iter是当前的迭代次数，最大迭代数设为100。