CN104766142A - 基于eeac和轨迹灵敏度的暂态稳定约束最优潮流计算方法 - Google Patents

Abstract

本发明公开了基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，属于电力系统暂态分析的技术领域。所述方法：以常规最优潮流解为初始运行点计算暂态初值并设置预想故障集；对预想故障集中的各故障进行暂态时域仿真、EEAC等值操作，对故障进行分类，将满足故障稳定裕度要求的故障从预想故障集中剔除；计算不满足稳定裕度要求的故障的轨迹灵敏度，给出不同类型下各故障的裕度灵敏度，并根据各故障的裕度灵敏度构建暂态稳定性约束；将所述暂态稳定性约束嵌入常规最优潮流模型后得到含有简单暂态稳定约束的最优潮流模型，避免了只采用功角作为暂态稳定判据时无法量化暂态稳定性，减少了计算复杂性并解决了收敛困难的问题。

Description

基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法

技术领域

本发明公开了基于EEAC(Extended Equal Area Criterion，扩展等面积准则)和轨迹灵敏度的暂态稳定约束最优潮流计算方法，属于电力系统暂态分析的技术领域。

背景技术

随着国民经济的快速发展和人民生活水平的不断提高，电力工业迅猛发展，电力新技术、新措施逐步投入应用，给社会各方面带来了显著的效益。但是，由于电力体制改革的不断深入，电网运营商为了追求更高的经济效益，使电力系统的运行点不断逼近安全稳定边界。各区域电网间互联加强以及新技术的应用，导致电网规模和复杂性日益增大，使电力系统的安全稳定问题越来越突出。电力系统中任任意元件发生故障，都有可能引起事故扩大。如果电网结构不合理，或者安全自动装置不够健全，都有可能使系统陷入稳定危机，造成稳定破坏或大面积停电，给国民经济带来重大损失。现代社会一旦发生大停电，其经济损失可能远远超过市场竞争可能带来的经济收益甚至危及国家和社会安全。

含暂态稳定约束的最优潮流(Transient Stability Constrained Optimal powerFlow,TSCOPF)可以兼顾系统的暂态稳定性和经济性，实现了电力系统暂态稳定离线预防控制，成为解决这一问题的有效手段。

TSCOPF抽象为数学模型后，是一个复杂的大规模非线性微分代数混合约束优化问题。目前计算的基本解法分为两类。

第一类方法为联立求解法，一般是将暂态微分方程差分化变为代数方程后加入最优潮流中联立求解。如文献一《暂态稳定约束下的最优潮流》(中国电机工程学报，2005年第25卷第12期第12页)所述。这类方法原理清晰，但是问题规模随着故障数目和系统规模的增大，呈几何级数增长，易出现维数灾的问题。由于差分化加入最优潮流中后，最优化问题庞大，容易出现收敛困难的问题。另外，这类方法的暂态稳定判据一般为功角阈值，文献二《暂态稳定预防控制和优化新进展》(电力系统自动化，2004年第28卷第10期)所述，该判据存在两个理论上的问题，一是角度稳定是否就意味着系统稳定，二是如何选取最大相对摇摆角。

第二类方法是分解迭代求解法，通常将TSCOPF问题分解为暂态稳定评估和最优潮流两个子问题，将暂态稳定分析得到的约束式嵌入最优潮流问题中迭代求解。这类方法可有效处理多故障问题、减少计算维数。如文献三《基于单机无穷大母线等值和轨迹灵敏度的暂态稳定约束最优潮流》(中国电机工程学报，2011年第31卷第13期)所述，这类解法可以定义系统暂态稳定裕度，可得到系统暂态稳定量化信息，便于优化。这类方法通常采用功角门槛值作为暂态稳定约束，没有量化稳定裕度。文献三定义了一般不稳定和稳定情况下的暂态稳定裕度，构造稳定性约束，极度失稳情况下仍采用功角门槛值构造暂态稳定性约束，采用二次规划算法进行迭代计算。最后的算例结果仅使多个故障达到稳定，并未满足一定的稳定性要求。由于各种随机因素的存在，必须留有一定的暂态稳定裕度。

发明内容

本发明所要解决的技术问题是针对上述背景技术的不足，提供了基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，将暂态稳定约束最优潮流分解为暂态稳定评估、灵敏地分析、最有潮流三个子问题进行迭代求解，解决了现有暂态稳定约束最优潮流计算复杂、难以收敛，以功角阈值作为稳定判据不准确以及选取最大相对摇摆角难的技术问题。

本发明为实现上述发明目的采用如下技术方案：

基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，包括如下步骤：

A.以常规最优潮流解为初始运行点计算暂态初值并设置预想故障集；

B.对预想故障集中的各故障进行暂态时域仿真、EEAC等值操作，对故障进行分类，将满足故障稳定裕度要求的故障从预想故障集中剔除，仅当预想故障集中所有故障都被剔除时直接输出最优解；

C.对于步骤B中不满足稳定裕度要求的故障：计算各故障的轨迹灵敏度，给出不同类型下各故障的裕度灵敏度，并根据各故障的裕度灵敏度构建暂态稳定性约束；

D.将所述暂态稳定性约束嵌入常规最优潮流模型后得到含有简单暂态稳定约束的最优潮流模型，重新计算常规潮流解，返回步骤A。

进一步的，作为所述暂态稳定约束最优潮流计算方法的优化方案，步骤C中所述暂态稳定性约束为：η为各故障稳定裕度构成的向量，为各故障裕度灵敏度构成的向量，u为系统控制变量构成的向量，Δu^T为系统控制变量在各故障下的变化量构成的向量，η_min为各故障稳定裕度最小值构成的向量，

各故障稳定裕度构成的向量η包括：一般失稳故障的稳定裕度η_u、稳定故障的稳定裕度η_s、极度失稳故障的稳定裕度η_eu：

${\mathrm{\η}}_u=-\frac{1}{2}{M}_E{\left({\mathrm{\ω}}_E\left(t_u\right)\right)}^2,$

${\mathrm{\η}}_s=-{\∫}_{{\mathrm{\δ}}_E\left(t_r\right)}^{{\mathrm{\δ}}_E\left(t_u\right)}{P}_a\mathrm{d\δ}\≈\frac{1}{2}{P}_a\left(t_r\right)({\mathrm{\δ}}_E\left(t_r\right)-{\mathrm{\δ}}_E\left(t_u\right)),$

η_eu＝-P_a(t_min)＝-min{P_mE(t)-P_eE(t),t>t_cl}，

其中：t_u为等值后系统达到DSP点的时刻，M_E为发电机的等值惯性时间常数，ω_E(t_u)、δ_E(t_u)分别为t_u时刻发电机的转速、功角，t_r为等值后系统首次达到FEP点的时刻，P_a(t_r)、δ_E(t_r)为分别为t_r时刻的不平衡功率、发电机功角，t_min为不平衡功率达到最小值的时刻，P_a(t_min)为t_min时刻的不平衡功率值，P_mE(t)、P_eE(t)分别为等值后发电机在t时刻的机械功率、电磁功率，t_cl为故障清除时刻。

再进一步的，作为所述暂态稳定约束最优潮流计算方法的优化方案，步骤B中所述一般失稳故障的稳定裕度η_u，在其暂态稳定性约束中增加校正因子λ： ${\mathrm{\λ\η}}_u+\frac{{\∂\mathrm{\η}}_u}{\∂u}{{\mathrm{\Δu}}_{{\mathrm{\η}}_u}}^T>{\mathrm{\λ\η}}_{u\min },$

其中：λ＝0.8～0.9，为系统控制变量在一般失稳故障下变化量构成的向量，η_umin为一般失稳故障的稳定裕度最小值。

进一步的，作为所述暂态稳定约束最优潮流计算方法的优化方案，步骤B中所述EEAC等值操作具体为：将对各故障进行暂态时域仿真得到的功角曲线进行EEAC等值，发电机等值后的模型为：

$\left\{\begin{array}{c}{M}_E\frac{{\mathrm{d\ω}}_E}{\mathrm{dt}}=P_a=P_{\mathrm{mE}}-P_{\mathrm{eE}}\\ \frac{d{\mathrm{\δ}}_E}{\mathrm{dt}}={\mathrm{\ω}}_E\end{array}\right.,$

其中：δ_E、ω_E、P_mE、P_eE分别为发电机的等值功角、等值转速、等值机械功率、等值电磁功率，P_a为等值后的不平衡功率。

再进一步的，作为所述暂态稳定约束最优潮流计算方法的优化方案，步骤B中所述对故障进行分类，采用如下判据：

一般失稳故障的判据： $\left\{\begin{array}{c}{P}_a\left(t_u\right)=P_{\mathrm{mE}}\left(t_u\right)-P_{\mathrm{eE}}\left(t_u\right)=0\\ {\stackrel{\·}{P}}_a\left(t_u\right)=\frac{{\mathrm{dP}}_a}{\mathrm{dt}}>0\\ {\mathrm{\ω}}_E\left(t\right)>0,\∀t>t_0\end{array}\right.,$

稳定故障的判据： $\left\{\begin{array}{c}{P}_a\left(t_r\right)=P_{\mathrm{mE}}\left(t_r\right)-P_{\mathrm{eE}}\left(t_r\right)<0\\ {\mathrm{\&omega;}}_E\left(t_r\right)=0\\ {\mathrm{\&omega;}}_E\left(t\right)>0,{\&ForAll;t}_0<t<t_r\end{array}\right.,$

极度失稳故障的判据： $P_a\left(t\right)>0,{\&ForAll;t}_0<t<t_{\mathrm{end}},$

其中：P_a(t_u)、P_mE(t_u)、P_eE(t_u)分别为t_u时刻的不平衡功率、发电机机械功率、发电机电磁功率，P_a(t)、ω_E(t)分别为t时刻的不平衡功率、发电机转速，t₀为故障开始时刻，P_mE(t_r)、P_eE(t_r)、ω_E(t_r)分别为t_r时刻发电机的机械功率、电磁功率、转速，t_end为暂态时域仿真结束时刻。

更进一步的，作为所述暂态稳定约束最优潮流计算方法的优化方案，步骤D中采用非线性原对偶内点法求解含有简单暂态稳定约束的最优潮流模型。

本发明采用上述技术方案，具有以下有益效果：

(1)采用EEAC法进行暂态稳定评估，考虑不平衡功率以及发电机相角构建了极度失稳、一般失稳和稳定情况下的暂态稳定判据，进而量化表示各故障的暂态稳定裕度，避免了只采用功角作为暂态稳定判据时无法量化暂态稳定性以及构建的暂态稳定约束不可靠的问题；

(2)针对各故障计算轨迹灵敏度，进而确定各故障的裕度灵敏度，再由轨迹灵敏度以及裕度灵敏度确定二次规划问题中暂态稳定约束的表达式，指导最优控制方向；

(3)本发明构建的将所述暂态稳定性约束嵌入常规最优潮流模型后得到含有简单暂态稳定约束的最优潮流模型，将一般失稳故障调整为稳定故障，将极度失稳故障调整为一般失稳故障，在以后的迭代中再将一般失稳故障调整为稳定故障，以满足给定的稳定性要求，实现系统由轻微失稳向稳定的调整，调整的过程中控制量变化较小，调整精确性较高且过度调整不明显；

(4)将暂态稳定约束最优潮流分解为暂态稳定评估、灵敏地分析、最优潮流三个子问题进行迭代求解，构建的各故障暂态稳定性约束嵌入最常规潮流计算后，其计算量与常规最优潮流计算几乎相同，减少了计算复杂性并解决了收敛困难的问题。

附图说明

图1是EEAC等值后一般失稳故障的等值功角曲线示意图；

图2是EEAC等值后稳定故障的等值功角曲线示意图；

图3是EEAC等值后极度失稳故障的等值功角曲线示意图；

图4是计算流程图。

具体实施方式

下面结合附图对发明的技术方案进行详细说明。本发明的思路是：将暂态稳定约束最优潮流分为三个子问题，包括暂态稳定评估问题、灵敏度分析问题和最优潮流问题，并将三个子问题迭代求解。如图4所示，具体按照以下步骤：

步骤1、基于EEAC对各预想故障进行暂态稳定评估：

步骤101、计算系统常规最优潮流，即计算不含暂态稳定约束的最优潮流，以此作为初始运行点。

建立发电机模型，负荷模型，计算暂态初值，对各故障进行暂态时域仿真。暂态稳定计算的数学模型为：

$G(y\left(t_0^{-}\right),u)=0---\left(1\right)$

$L(x\left(t_0^{-}\right),y\left(t_0^{-}\right),u)=0---\left(2\right)$

$\stackrel{\&CenterDot;}{x}\left(t\right)-F(x\left(t\right),y\left(t\right),u)=0,t\&Element;T---\left(3\right)$

G(x(t),y(t),u)＝0 t∈T  (4)

式(1)至式(4)中，x为系统状态量向量，包括发电机功角和转速等；y为系统代数量向量，包括系统节点电压的幅值和相角等；u为系统控制量向量，包括发电机有功出力和无功出力等；为故障发生前一时刻， 为故障发生时刻，t_cl为故障清除时刻，t_end为暂态仿真结束时刻；式(1)为潮流方程；式(2)为发电机暂态初值方程；式(3)和式(4)描述系统的暂态运行过程，为一组微分代数方程。

然后将仿真后功角曲线进行EEAC等值。如文献四《运动稳定性量化理论—非自治非线性多刚体系统的稳定性分析》(南京：江苏科学技术出版社，1999)所述，系统发电机等值后模型为：

$\left\{\begin{array}{c}{M}_E\frac{{\mathrm{d\&omega;}}_E}{\mathrm{dt}}=P_a=P_{\mathrm{mE}}-P_{\mathrm{eE}}\\ \frac{d{\mathrm{\&delta;}}_E}{\mathrm{dt}}={\mathrm{\&omega;}}_E\end{array}\right.---\left(5\right)$

式(5)中：M_E为等值惯性时间常数；δ_E为等值功角；ω_E为等值转速；P_mE为等值机械功率；P_eE为等值电磁功率；P_a为等值后的不平衡功率。

步骤102、EEAC等值后，将故障分为一般失稳、稳定、极度失稳三种情况。

在仿真时间内，若系统无DSP点，只有FEP点，则系统稳定；若系统达到DSP点，则系统一般失稳；若系统既无DSP点也无FEP点，则系统极度失稳。下面具体给出各稳定情况的判据以及裕度表达式：

(1)一般失稳故障的判据为：

$\left\{\begin{array}{c}{P}_a\left(t_u\right)=P_{\mathrm{mE}}\left(t_u\right)-P_{\mathrm{eE}}\left(t_u\right)=0\\ {\stackrel{\&CenterDot;}{P}}_a\left(t_u\right)=\frac{{\mathrm{dP}}_a}{\mathrm{dt}}>0\\ {\mathrm{\&omega;}}_E\left(t\right)>0,\&ForAll;t>t_0\end{array}\right.---\left(6\right)$

式(6)中：t_u为EEAC等值后系统达到DSP点的时刻，一般失稳故障的等值功角曲线示意图如图1所示，其稳定裕度定义为：

${\mathrm{\&eta;}}_u=-\frac{1}{2}{M}_E{\left({\mathrm{\&omega;}}_E\left(t_u\right)\right)}^2---\left(7\right)$

(2)系统稳定的判据为：

$\left\{\begin{array}{c}{P}_a\left(t_r\right)=P_{\mathrm{mE}}\left(t_r\right)-P_{\mathrm{eE}}\left(t_r\right)<0\\ {\mathrm{\&omega;}}_E\left(t_r\right)=0\\ {\mathrm{\&omega;}}_E\left(t\right)>0,{\&ForAll;t}_0<t<t_r\end{array}\right.---\left(8\right)$

式(8)中：t_r为系统EEAC等值后第一次达到FEP点的时刻；δ_E(t_r)为系统EEAC等值后第一次达到FEP点时的等值功角，稳定故障的等值功角曲线示意图如图2所示，其稳定裕度可近似表达为：

${\mathrm{\&eta;}}_s=-{\&Integral;}_{{\mathrm{\&delta;}}_E\left(t_r\right)}^{{\mathrm{\&delta;}}_E\left(t_u\right)}{P}_a\mathrm{d\&delta;}\&ap;\frac{1}{2}{P}_a\left(t_r\right)({\mathrm{\&delta;}}_E\left(t_r\right)-{\mathrm{\&delta;}}_E\left(t_u\right))---\left(9\right)$

式(9)中：δ_E(t_u)是系统临界失稳时在t_u时刻的等值功角。

(3)系统极度失稳的判据为：

$P_a\left(t\right)>0,{\&ForAll;t}_0<t<t_{\mathrm{end}}---\left(10\right)$

由于极度失稳故障不存在DSP点和FEP点，极度失稳故障的等值功角曲线示意图如图3所示，所以将最小不平衡功率定义为其稳定裕度：

η_eu＝-P_a(t_min)＝-min{P_mE(t)-P_eE(t),t>t_cl}  (11)

式(6)至式(11)中：P_a(t_u)、P_mE(t_u)、P_eE(t_u)分别为t_u时刻的不平衡功率、发电机机械功率、发电机电磁功率，P_a(t)、ω_E(t)分别为t时刻的不平衡功率、发电机转速，t₀为故障开始时刻，P_mE(t_r)、P_eE(t_r)、ω_E(t_r)分别为t_r时刻发电机的机械功率、电磁功率、转速，t_min代表故障后不平衡功率达到最小值的时刻，P_a(t_min)代表故障后不平衡功率的最小值。

步骤2、计算轨迹灵敏度以及稳定裕度灵敏度

步骤201、针对上一步中不满足稳定裕度要求的故障，基于上述暂态时域仿真的数据，计算轨迹灵敏度初值。

当系统控制量发生变化时，会造成潮流分布的改变，所以首先对式(1)求偏导，得到静态潮流结果对于控制量的灵敏度然后对式(2)求偏导，得到暂态初值对于控制量的灵敏度具体公式如下：

$\left\{\begin{array}{c}{G}_y{y}_u\left(t_0^{-}\right)+G_u=0\\ L_x{x}_u\left(t_0^{-}\right)+L_y{y}_u\left(t_0^{-}\right)+L_u=0\end{array}\right.---\left(12\right)$

由于状态变量不会突变，所以其灵敏度也不会突变，即 $x_u\left(t_0^{-}\right)=x_u\left(t_0^{+}\right),$ 而代数变量会发生突变，所以 $y_u\left(t_0^{+}\right)\&NotEqual;y_u\left(t_0^{-}\right),$ 可由式(4)计算得到。所以进行整理后，轨迹灵敏度的初值计算方程为：

$\left\{\begin{array}{c}{x}_u\left(t_0^{+}\right)=L_x^{-1}\left(t_0^{-}\right)L_u\left(t_0^{-}\right)+L_x^{-1}\left(t_0^{-}\right)L_y\left(t_0^{-}\right)G_y^{-1}\left(t_0^{-}\right)G_u\left(t_0^{-}\right)\\ y_u\left(t_0^{+}\right)=-G_y^{-1}\left(t_0^{+}\right)(G_x\left(t_0^{+}\right)x_u\left(t_0^{+}\right)+G_u\left(t_0^{+}\right))\end{array}\right.---\left(13\right)$

步骤202、计算轨迹灵敏度和裕度灵敏度。

轨迹灵敏度方程可由式(3)和式(4)对控制量u求偏导得到：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_u=F_x\left(t\right)x_u+F_y\left(t\right)y_u+F_u\left(t\right)\\ 0=G_x\left(t\right)x_u+G_y\left(t\right)y_u+G_u\left(t\right)\end{array}\right.---\left(14\right)$

式(14)中：F_x(t)、F_y(t)、F_u(t)、G_x(t)、G_y(t)和G_u(t)为随系统轨迹变化的时变偏导矩阵。

为了刻画控制变量对系统稳定裕度的影响，基于上述轨迹灵敏度，引入裕度灵敏度。对于不同的故障类型，有不同的裕度灵敏度表达式。

(1)对于极度失稳故障，其裕度灵敏度为：

$\frac{{\&PartialD;\mathrm{\&eta;}}_{\mathrm{eu}}}{\&PartialD;u}=-\frac{{\&PartialD;P}_{a\min }}{\&PartialD;u}=-(\frac{{\&PartialD;P}_{\mathrm{mE}}\left(t_{\min }\right)}{\&PartialD;u}-\frac{{\&PartialD;P}_{\mathrm{eE}}\left(t_{\min }\right)}{\&PartialD;u})---\left(15\right)$

(2)对于一般失稳故障，其裕度灵敏度为：

$\frac{{\&PartialD;\mathrm{\&eta;}}_u}{\&PartialD;u}-M_E{\mathrm{\&omega;}}_E\left(t_u\right)\frac{{\&PartialD;\mathrm{\&omega;}}_E\left(t_u\right)}{\&PartialD;u}---\left(16\right)$

(3)对于稳定故障，当其稳定裕度较小而不满足稳定性要求时，仍需调整控制变量，所以其裕度灵敏度可表达为：

$\frac{{\&PartialD;\mathrm{\&eta;}}_s}{\&PartialD;u}=\frac{1}{2}\frac{{\&PartialD;P}_a\left(t_r\right)}{\&PartialD;u}({\mathrm{\&delta;}}_E\left(t_r\right)-{\mathrm{\&delta;}}_E\left(t_u\right))+\frac{1}{2}{P}_a\left(t_r\right)(\frac{{\&PartialD;\mathrm{\&delta;}}_E\left(t_r\right)}{\&PartialD;u}-\frac{{\&PartialD;\mathrm{\&delta;}}_E\left(t_u\right)}{\&PartialD;u})---\left(17\right)$

步骤203、基于上述各故障的裕度灵敏度，可构造暂态稳定性约束，可简写为如下形式：

$\mathrm{\&eta;}+\frac{\&PartialD;\mathrm{\&eta;}}{\&PartialD;u}{\mathrm{\&Delta;u}}^T>{\mathrm{\&eta;}}_{\min }---\left(18\right)$

式(18)中：η为各个故障稳定裕度构成的向量，为各故障裕度灵敏度构成的向量，u为系统控制变量构成的向量，Δu^T为系统控制变量在各故障下变化量构成的向量，η_min为各故障稳定裕度最小值构成的向量，各故障稳定裕度构成的向量η包括：式(7)所示一般失稳故障的稳定裕度η_u、式(9)所示稳定故障的稳定裕度η_s、式(11)所示极度失稳故障的稳定裕度η_eu，ω_E(t_u)、δ_E(t_u)分别为t_u时刻发电机的转速、功角，t_r为等值后系统首次达到FEP点的时刻，P_a(t_r)、δ_E(t_r)为分别为t_r时刻的不平衡功率、发电机功角，t_min为不平衡功率达到最小值的时刻，P_a(t_min)为t_min时刻的不平衡功率值，P_mE(t)、P_eE(t)分别为等值后发电机在t时刻的机械功率、电磁功率，t_cl为故障清除时刻。

一般失稳故障的稳定性约束目标是将一般失稳故障调整为稳定。由于极度失稳故障的裕度是根据故障后最小不平衡功率定义的，所以其稳定性约束的目标是将极度失稳调整为一般失稳故障，以后的迭代中再将一般失稳故障调整为稳定。在极度失稳故障和一般失稳故障的稳定性约束中，一般取η_min＝0，在稳定故障中，根据实际电网需求，给定η_min，使发生故障时系统稳定裕度不会过小，符合预防控制的要求。

步骤3、计算含简单暂态稳定约束的最优潮流

步骤301、由上述稳定性约束形成含简单暂态稳定约束的最优潮流模型：

Min f(u)  (19)

$\begin{array}{cc}s.t.& G(y\left(t_0^{-}\right),u)=0\end{array}---\left(20\right)$

$H(y\left(t_0^{-}\right),u)\&le;0---\left(21\right)$

${\widetilde{\mathrm{\&eta;}}}^{\left(k\right)}(x\left(t\right),x_u\left(t\right),y\left(t\right),y_u\left(t\right),u)\&le;0,k=1,...,n_F---\left(22\right)$

式(19)为目标函数，如发电费用函数；式(21)为故障前稳态时需满足的运行约束方程，式(22)为稳定性约束方程，其中k＝1,…,n_F代表第k个故障，n_F为故障的总数目。

为提高计算准确性，对于较为严重的一般失稳故障，迭代时可在稳定性约束中增加一个校正因子λ＝0.8～0.9，构成校正的稳定性约束为：

${\mathrm{\&lambda;\&eta;}}_u+\frac{{\&PartialD;\mathrm{\&eta;}}_u}{\&PartialD;u}{{\mathrm{\&Delta;u}}_{{\mathrm{\&eta;}}_u}}^T>{\mathrm{\&lambda;\&eta;}}_{u\min }---\left(23\right)$

该约束使迭代后系统为轻微失稳。而其余迭代中，校正因子λ＝1。采用轨迹灵敏度计算控制量时，会出现过度调整的现象。由于轨迹灵敏度是基于原轨迹的线性灵敏度，当控制量的变化较小时，应用灵敏度调整稳定裕度误差较小，精确性较高。通过这一处理方法，可使在预防控制过程中，系统从失稳调为稳定时，系统是由轻微失稳调整为稳定，控制量变化较小，调整精确性较高，过度调整不明显。

效果验证：

为验证本发明算法的有效性，针对IEEE10机39节点标准算例进行了测试，分别计算含单一故障和多故障的TSCOPF，并进行比较分析。发电机采用经典二阶模型，负荷采用恒阻抗模型，仿真时间为4s(考虑多摆失稳)，仿真步长为0.01s。flag是标识符，用于记录系统发生故障后的暂态稳定情况，0表示稳定，1表示一般失稳，2表示极度失稳。令符号L(i,j)表示线路i-j在母线i处0s时刻发生三相短路故障，t_cl时刻切除故障线路。下面选取故障A₁和故障A₂，故障A₁为L(29,28)，切除故障时间为0.32s，为极度失稳故障；故障A₂为L(2,1)，切除故障时间为0.24s，为一般失稳故障。故障A₁中38号机为临界群，其余为剩余群；故障A₂中39号机为剩余群，其余为临界群，两者为不同失稳模式的故障。表1和表2分别为单独考虑故障A₁和故障A₂时计算TSCOPF的过程，表3为同时考虑故障A₁和A₂计算TSCOPF的过程，表中CCT表示临界切除时间。本算例分析中将发电机有功出力作为控制变量，取η_min＝0.1。当一般失稳故障裕度小于-1时，取校正因子λ＝0.8。

表1 故障A₁的TSCOPF迭代过程

表1中，由原对偶内点法计算常规最优潮流，得到初始运行点。对故障A₁仿真，判断为极度失稳故障，经过第1次迭代优化后，转为一般失稳故障，是严重的一般失稳故障。为了避免出现过度调整的情况，将校正的稳定性约束嵌入最优潮流计算，得到第2次迭代结果，此时仍然失稳，但裕度较大，经第3次迭代后，系统达到稳定。此时，系统仍不满足稳定性要求，经第4次迭代，达到稳定性要求。

表2 故障A₂的TSCOPF迭代过程

发电机	初始OPF	第1次迭代	第2次迭代	第3次迭代
					30(MW)	243.96	242.94	241.45	250.89
31(MW)	569.62	555.87	548.57	548.2
					32(MW)	646.06	648.39	647.35	647.03
33(MW)	636.67	631.85	629.25	628.25
					34(MW)	513.11	509.66	507.5	506.94
35(MW)	655.63	650.68	648.69	647.38
					36(MW)	563.31	558.3	556.29	555.09
37(MW)	539.11	531.92	530.51	527.27
					38(MW)	837.3	812.29	813.47	809.78
39(MW)	981.08	1041.86	1060.39	1062.3
					flag	1	1	0	0
裕度(p.u.)	-1.4182	-0.2537	0.0355	0.1321
					费用($/h)	61656.05	61689.91	61711.22	61717.5
CCT(s)			0.242	0.244

表3 故障A₁+A₂的TSCOPF迭代过程

故障A₂为多摆失稳故障，表2给出故障A₂的TSCOPF计算过程，与故障A₁类似，经3次迭代可达到稳定性要求。表3中，同时处理含故障A₁和故障A₂这两个不同失稳模式故障的暂态稳定最优潮流问题，经4次迭代可达到稳定性要求。对比表1、2和3可见，考虑多故障的发电费用略高于考虑单故障的发电费用。

由此可见，本文提出的方法是合理有效的。

综上所述，本发明具有以下有益效果：

(1)采用EEAC法进行暂态稳定评估，考虑不平衡功率以及发电机相角构建了极度失稳、一般失稳和稳定情况下的暂态稳定判据，进而量化表示各故障的暂态稳定裕度，避免了只采用功角作为暂态稳定判据时无法量化暂态稳定性以及构建的暂态稳定约束不可靠的问题；

(2)针对各故障计算轨迹灵敏度，进而确定各故障的裕度灵敏度，再由轨迹灵敏度以及裕度灵敏度确定二次规划问题中暂态稳定约束的表达式，指导最优控制方向；

(3)本发明构建的将所述暂态稳定性约束嵌入常规最优潮流模型后得到含有简单暂态稳定约束的最优潮流模型，将一般失稳故障调整为稳定故障，将极度失稳故障调整为一般失稳故障，在以后的迭代中再将一般失稳故障调整为稳定故障，以满足给定的稳定性要求，实现系统由轻微失稳向稳定的调整，调整的过程中控制量变化较小，调整精确性较高且过度调整不明显；

(4)将暂态稳定约束最优潮流分解为暂态稳定评估、灵敏地分析、最优潮流三个子问题进行迭代求解，构建的各故障暂态稳定性约束嵌入最常规潮流计算后，其计算量与常规最优潮流计算几乎相同，减少了计算复杂性并解决了收敛困难的问题。

本发明所述方案不局限于具体实施例中涉及的常规潮流模型，本领域的技术人员在本发明技术方案的指导下，能够将本发明所述方案推广应用到其它常规潮流模型中，其它应用于任意常规潮流模型中的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法都落入本发明的保护范围。

Claims (6)

1.基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，包括如下步骤：

A.以常规最优潮流解为初始运行点计算暂态初值并设置预想故障集；

B.对预想故障集中的各故障进行暂态时域仿真、EEAC等值操作，对故障进行分类，将满足故障稳定裕度要求的故障从预想故障集中剔除，仅当预想故障集中所有故障都被剔除时直接输出最优解；

C.对于步骤B中不满足稳定裕度要求的故障：计算各故障的轨迹灵敏度，给出不同类型下各故障的裕度灵敏度，并根据各故障的裕度灵敏度构建暂态稳定性约束；

D.将所述暂态稳定性约束嵌入常规最优潮流模型后得到含有简单暂态稳定约束的最优潮流模型，重新计算常规潮流解，返回步骤A。

2.根据权利要求1所述的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，步骤C中所述暂态稳定性约束为：η为各故障稳定裕度构成的向量，为各故障裕度灵敏度构成的向量，u为系统控制变量构成的向量，Δu^T为系统控制变量在各故障下的变化量构成的向量，η_min为各故障稳定裕度最小值构成的向量，

各故障稳定裕度构成的向量η包括：一般失稳故障的稳定裕度η_u、稳定故障的稳定裕度η_s、极度失稳故障的稳定裕度η_eu：

${\mathrm{\&eta;}}_u=-\frac{1}{2}{M}_E{\left({\mathrm{\&omega;}}_E\left(t_u\right)\right)}^2,$

${\mathrm{\&eta;}}_s=-{\&Integral;}_{{\mathrm{\&delta;}}_E\left(t_r\right)}^{{\mathrm{\&delta;}}_E\left(t_u\right)}{P}_a\mathrm{d\&delta;}\&ap;\frac{1}{2}{P}_a\left(t_r\right)({\mathrm{\&delta;}}_E\left(t_r\right)-{\mathrm{\&delta;}}_E\left(t_u\right)),$

η_eu＝-P_a(t_min)＝-min{P_mE(t)-P_eE(t),t>t_cl}，

其中：t_u为等值后系统达到DSP点的时刻，M_E为发电机的等值惯性时间常数，ω_E(t_u)、δ_E(t_u)分别为t_u时刻发电机的转速、功角，t_r为等值后系统首次达到FEP点的时刻，P_a(t_r)、δ_E(t_r)为分别为t_r时刻的不平衡功率、发电机功角，t_min为不平衡功率达到最小值的时刻，P_a(t_min)为t_min时刻的不平衡功率值，P_mE(t)、P_eE(t)分别为等值后发电机在t时刻的机械功率、电磁功率，t_cl为故障清除时刻。

3.根据权利要求2所述的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，步骤B中所述一般失稳故障的稳定裕度η_u，在其暂态稳定性约束中增加校正因子λ：

其中：λ＝0.8～0.9，为系统控制变量在一般失稳故障下变化量构成的向量，η_umin为一般失稳故障的稳定裕度最小值。

4.根据权利要求2所述的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，步骤B中所述EEAC等值操作具体为：将对各故障进行暂态时域仿真得到的功角曲线进行EEAC等值，发电机等值后的模型为：

$\left\{\begin{array}{c}{M}_E\frac{d{\mathrm{\&omega;}}_E}{\mathrm{dt}}=P_a=P_{\mathrm{mE}}-P_{\mathrm{eE}}\\ \frac{d{\mathrm{\&delta;}}_E}{\mathrm{dt}}={\mathrm{\&omega;}}_E\end{array}\right.,$

其中：δ_E、ω_E、P_mE、P_eE分别为发电机的等值功角、等值转速、等值机械功率、等值电磁功率，P_a为等值后的不平衡功率。

5.根据权利要求4所述的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，步骤B中所述对故障进行分类，采用如下判据：

一般失稳故障的判据： $\left\{\begin{array}{c}{P}_a\left(t_u\right)=P_{\mathrm{mE}}\left(t_u\right)-P_{\mathrm{eE}}\left(t_u\right)=0\\ {\stackrel{\&CenterDot;}{P}}_a\left(t_u\right)=\frac{d{P}_a}{\mathrm{dt}}>0\\ {\mathrm{\&omega;}}_E\left(t\right)>0,\&ForAll;t>t_0\end{array}\right.,$

稳定故障的判据： $\left\{\begin{array}{c}{P}_a\left(t_r\right)=P_{\mathrm{mE}}\left(t_r\right)-P_{\mathrm{eE}}\left(t_r\right)<0\\ {\mathrm{\&omega;}}_E\left(t_r\right)=0\\ {\mathrm{\&omega;}}_E\left(t\right)>0,\&ForAll;t_0<t<t_r\end{array}\right.,$

极度失稳故障的判据： $P_a\left(t\right)>0,\&ForAll;t_0<t<t_{\mathrm{end}},$

其中：P_a(t_u)、P_mE(t_u)、P_eE(t_u)分别为t_u时刻的不平衡功率、发电机机械功率、发电机电磁功率，P_a(t)、ω_E(t)分别为t时刻的不平衡功率、发电机转速，t₀为故障开始时刻，P_mE(t_r)、P_eE(t_r)、ω_E(t_r)分别为t_r时刻发电机的机械功率、电磁功率、转速，t_end为暂态时域仿真结束时刻。

6.根据权利要求1至5中任意一项所述的基于EEAC和轨迹灵敏度的暂态稳定约束最优潮流计算方法，其特征在于，步骤D中采用非线性原对偶内点法求解含有简单暂态稳定约束的最优潮流模型。