第26卷　第13期2002年7月10日电　力　系　统　自　动　化A utom ati on of E lectric Pow er System s V o l .26　N o.13July 10,2002基于内点法的含暂态稳定约束的最优潮流计算袁　越1,久保川淳司2,佐佐木博司1,宋永华3(1.广岛大学,日本;2.广岛工业大学,日本;31西安交通大学电力工程系,陕西省西安市710049)摘要:建立了含暂态稳定约束的最优潮流的数学模型,模型中考虑了多个预想事故。

提出了一种基于原—对偶内点法的含暂态稳定约束的最优潮流算法。

通过充分开发修正矩阵的稀疏性,并在求解时采用稀疏技巧,开发出了高性能的计算程序。

在日本60H z 电力网的10机模型系统的优化计算结果表明,所提算法不仅具有强大的处理等式约束和不等式约束的能力,而且具有良好的收敛性,能够有效地解决考虑多个预想事故时的含暂态稳定约束的最优潮流问题。

关键词:最优潮流;原—对偶内点法;暂态稳定分析中图分类号:TM 711;TM 744收稿日期:2002201228。

0　引言自从20世纪60年代法国的Carpen tier 提出最初的最优潮流模型以来,广大学者对最优潮流问题进行了大量的研究。

然而,由于传统的最优潮流没有考虑暂态稳定约束,在其得出的最优运行方式下系统可能会遇到暂态稳定性问题。

特别是在电力市场化运营机制下,系统不可能再在保守的方式下运行,如何能够把系统的安全性和经济性融为一体,就显得更为重要。

为此,近年来,关于包含暂态稳定约束的最优潮流(SCO PF )的研究引起了各国学者广泛的兴趣。

迄今为止,对于包含暂稳约束的最优潮流,提出了两种求解方法:一种是在传统的最优潮流中直接加入暂稳约束条件,然后采用与一般的最优潮流问题相同的算法来求解[1,2];另一种方法是基于Euclidean 空间变换,把一个含有大量约束的SCO PF 优化问题转换为与一般的最优潮流相同规模的优化问题[3]。

在本文中,我们称第1种方法为SCO PF 求解的“直接法”,第2种方法为SCO PF 求解的“间接法”。

“间接法”的长处在于降低了优化问题的规模,而“直接法”则具有可以借鉴和采用各种发展成熟的暂态稳定分析方法的优点。

目前,包含暂稳约束的最优潮流还处于发展阶段,特别是“间接法”才得到小系统的验证。

此外,所有关于SCO PF 的研究还仅局限于考虑一个预想事故。

显然,如果求取整个系统的既安全又经济的运行方式,仅考虑一个预想事故是不够的。

为此,本文建立了考虑多个预想事故时SCO PF 的数学模型,并且提出了一种求解方法。

从数学表达式来看,SCO PF 问题属于非线性规划问题。

实际上,任何求解一般的最优潮流的算法都可以用于求解SCO PF 问题。

不同于文献[1～3]中的算法,本文提出了一种基于原—对偶内点法的SCO PF 求解算法。

作为一种功能强大的优化算法,它已经成功地解决了许多带大量约束的大规模电力系统的优化问题[4,5]。

事实上,本文研究表明,原—对偶内点法在求解考虑多个预想事故的SCO PF 问题上同样可以达到令人满意的性能。

1　含暂态稳定约束的最优潮流模型本文采用多机电力系统的经典数学模型,各发电机用x d ′后的恒定电势E ′来模拟,负荷用恒定阻抗模型。

发电机的转子运动方程为[6]:∆i =Ξi -Ξ0M i Ξi =Ξ0(-D i Ξi +P m i -P e i )(1)　P e i =E i ′2G ii ′+∑n gj =1j ≠iE i ′E j ′B ij ′sin (∆i -∆j )+E i ′E j ′G ij ′co s (∆i -∆j )式中:i ∈S G ;S G 为发电机节点集合;P e i 为发电机的电磁功率;Y ij ′=G ij ′+j B ij ′(i ,j =1,2,…,n g )为发电机内电势节点的自导纳(i =j )和互导纳(i ≠j )。

为便于表示发电机的摇摆特性,取系统的惯性中心(CO I )作为参考。

CO I 的角度定义为[7]:∆CO I =∑n gi =1Mi∆i∑n gi =1Mi(2)1.1　目标函数本文采用发电燃料总费用作为目标函数,机组的燃料特性采用二次函数关系式:41F =∑i ∈S Gfi(P g i )f i (Pg i )=a i +b i P g i +c i2P2g i(3)1.2　等式约束条件不同于一般的最优潮流,SCO PF 的等式约束条件中,除了基本潮流方程式外,还加入了用于暂稳计算的转子运动方程式和初值方程式。

a .潮流方程式本文采用极坐标形式的潮流方程:P g i -P l i -V 2i G ii - Vi∑j ∈IVj(G ij co s Ηij +B ij sin Ηij )=0Q r i -Q l i +V 2i B ii- Vi∑j ∈IV j (G ij sin Ηij -B ij co s Ηij )=0(4)式中:i ∈S N ;S N 为总节点集合;I 为与节点i 相联的节点集合。

b .转子运动方程式如何处理微分方程形式的转子运动方程式是实现SCO PF 的一个关键。

本文继承“直接法”的思想——将微分方程转化成相应的等值差分方程。

在数值积分方法的选择上,文献[8]认为隐式梯形积分法是比较理想的方法。

按照这种方法,对于任一预想事故,有下列一组转子运动方程式:∆t i (k )-∆t -1i (k )-∃t 2Ξt i (k )+Ξt -1i (k )=0Ξt i (k )-Ξt -1i (k )-∃t 21M i-D i Ξt i (k )+ P g i -P te i (k )+1Mi-D i Ξt -1i (k )+　P g i -P t -1e i (k )=0(5)式中:∃t 为积分步长;i ∈S G ;t ∈S T ;k ∈S K ;S T 为积分时刻集合;S K 为预想事故集合。

c .初值方程式为了计算发电机的状态变量初值,引入了如下的初值方程式:E i ′V g i sin (∆0i -Ηg i )-x d i ′P g i =0V 2g i -E i ′V g i co s (∆0i -Ηg i )+x d i ′Q g i =0(6)式中:i ∈S G 。

1.3　不等式约束条件除了与一般的最优潮流相同的运行约束外,SCO PF 中加入了暂态稳定性约束。

a .运行约束运行约束主要包括有功电源出力上下限约束、可调无功电源出力上下限约束、节点电压模值上下限约束、线路通过的最大有功潮流约束:P g i ≤P g i ≤P g i i ∈S GQ r i ≤Q r i ≤Q r i i ∈S R V i ≤V i ≤Vii ∈S NP ij ≤P ij ≤P ij (i ,j )∈S CL(7)式中:S R 为可调无功电源集合;S CL 为线路集合。

b .暂态稳定性约束通过规定各发电机转子角度与系统惯性中心角度∆CO I 之间相对值的上下限,暂态稳定性约束可定义为: ∆≤∆0i -∆0CO I ≤∆∆≤∆t i (k )-∆tCO I (k )≤∆(8)式中:i ∈S G ;t ∈S T ;k ∈S K ;约束限值∆和∆可依据实际系统的运行经验确定。

2　原—对偶内点法最早的内点法可以追溯到20世纪50年代,然而,对于内点法的发展真正具有里程碑意义的是1984年Kar m arkar 提出的具有多项式时间可解性的线性规划内点算法[9]。

随后,无论在理论上还是在实践上,人们都对内点法投入了极大的研究兴趣,并取得了大量的研究成果。

其中,对于求解非线性规划(NL P )问题,文献[10]提出的原—对偶N ew ton 内点法被证明是一种非常有效的方法,能够快速求解大规模非线性规划问题。

在本节中简要说明如何用这种算法求解含暂稳约束的最优潮流问题。

假设SCO PF 中的变量可以表示为一个N 维向量:x ≡x contro l x state T ∈R (N ),于是SCO PF 问题可以表达为一个NL P 问题:m in f (x )(9)s .t .h (x )=0 h (x ):R (N )→R (m)g ≤g (x )≤g g (x ):R (N )→R (r)引入松弛变量向量(l ,u )∈R (r ),式(9)可以转化为: m in f (x )(10) s .t .h (x )=0g (x )-g -l =0g (x )-g +u =0(l ,u )≥0问题(10)的L agrangian 函数可以定义为:L (x ,l ,u ;y ,z ,w ,z ～,w ～)≡f (x )-y Th (x )- z T g (x )-g -l -w T g (x )-g +u - z ～T l -w T u(11)式中:y ∈R (m )和(z ,w ,z ～,w ～)∈R (r )为L agrangian 乘子;z ～=z ;w ～=-w 。

51・优化算法在电力系统中的应用专题・　袁　越等　基于内点法的含暂态稳定约束的最优潮流计算按照摄动KKT (Pertu rbed Karu sh 2Kuhn 2T ucker )一阶最优性条件可得:L x = f (x )- h (x )y - g (x )(z +w )=0L l =L Ze -Λe =0L u =UW e +Λe =0L y =h (x )=0L z =g (x )-g -l =0L w =g (x )-g +u =0(l ,u )≥0,y ≠0,z ≥0,w ≤0(12)式中:(L ,U ,Z ,W )∈R (r ×r )是对角元素分别为l i ,u i ,z i ,w i 的对角阵;e 为单位列向量,e =[1,1, (1)T∈R (r );Λ>0为摄动因子。

对式(12)用N ew ton 法求解,修正方程式为:[ 2h (x )y + 2g (x )(z +w )- 2f (x )]∃x + h (x )∃y + g (x )(∃z +∃w )=L x 0Z ∃l +L ∃z =-L l 0W ∃u +U ∃w =-L u 0( h (x ))T ∃x =-L y 0( g (x ))T ∃x -∃l =-L z 0( g (x ))T ∃x +∃u =-L w 0(13)式中:(L x 0,L l 0,L u 0;L y 0,L z 0,L w 0)为摄动KKT 方程式在展开点的值; 2h (x )和 2g (x )分别为h (x )和g (x )的H essian 矩阵。

从式(13)中消去变量(∃l ,∃u ,∃z ,∃w ),可以得到如下的降阶修正方程式:H ( )J T(x )J (x )0∃x ∃y =-7( ,Λ)　h (x )(14)H ( )= 2h (x )y + 2g (x )(z +w )- 2f (x )+g (x )(U -1W -L -1Z ) g (x )TJ (x )= h (x )T　7( ,Λ)=-L x 0+ g (x )[U -1(WL w 0-L u 0)-L -1(Z L z 0+L l 0)]显然,降阶修正方程式(14)中已经消去了变量不等式约束和函数不等式约束,其阶数取决于变量数和等式约束数,远小于式(13)的阶数。