一种带有等式约束的状态估计新算法_倪小平

一种带有等式约束的状态估计新算法

倪小平,张步涵

(华中科技大学电气与电子工程学院,武汉430074)

摘要:提出了一种新的带等式约束的状态估计算法,利用等式约束条件来修正加权最小二乘所得的状态量。利用这种算法能有效地利用系统中的一些虚拟零注入量测点数据,并且这种算法整体性不亚于拉格朗日多项式构造的等式约束算法,在增加很少计算量的情况下,达到提高状态估计结果精度的目的,保证了最小二乘状态估计的高效性,有利于状态估计的实时应用。关键词:状态估计;虚拟零注入量测;等式约束;拉格朗日多项式中图分类号:TM 732;TM 744

收稿日期:2001-03-25;修回日期:2001-07-17。

0　引言

随着调度自动化水平的不断提高,人们对状态估计的准确性和可靠性提出了更高的要求。通常情况下,一般采用提高量测系统的冗余度来提高状态估计结果的准确性[1]。这样做,一方面有违于量测系统经济性布置的要求,另一方面随着量测系统冗余度的增加,大大地降低了状态估计的效率(计算量一般随量测个数呈几何级数增加),不利于状态估计的实时性要求。其实,在电力系统中存在许多零注入节点,我们可以在这些节点虚拟一些零注入量测点来达到提高状态估计精度的要求[2,3]。在状态估计中,这种虚拟的零注入量是一种非常精确、可利用的量测类型,并且不必增加量测设备。它的加入可以极大地影响相关节点状态量的拟合趋势,加快算法的收敛速度,有较强的抵御相关量测的残差影响。但是,如何利用这些特殊的虚拟量测在牺牲较小估计效率的前提下来提高估计结果的精度,是一个值得探讨的问题。目前有几种处理办法:①通过提高其权系数来把它作为一种量测加以考虑[4],这种处理在一定程度上损失了虚拟零注入量测的精确性,并且增加了雅可比矩阵的维数,使得计算量增加。此外,由于虚拟量测大的权系数可能会造成系统病态,导致估计结果的不收敛。②通过拉格朗日多项式构造极值函数[2,3],然后导出迭代式。这种方法虽然在一定程度上保证了虚拟零注入量测的精确性,但是在计算过程中增加了过渡性变量λ,使得计算变量增加,增加了较大的计算量。当然,由于它的一些优点[5],现已被用于实时的状态估计中来解决一些等式约束问题。

本文提出了一个新的解决方法,在一定程度上能保证虚拟零注入量测的有效信息,另外也保证了状态估计的计算效率(不会增加太大的计算量)。

1　数学模型

1.1　带等式约束的目标函数

在给定网络接线、支路参数和量测系统的条件下,电力系统状态估计的非线性量测方程可表示为:

z =h (x )+v (1)

式中　z 为m 维量测向量;h (x )为m 维量测函数向

量;x 为n 维状态变量向量;v 为m 维随机量

测误差向量,且有E (v )=0,E (v T v )=R ;R -1

为m ×m 维量测对角权矩阵;m 为量测个数;n 为系统状态量个数。

给定量测向量z 后,带等式约束的状态估计问题可表述为:

min J (x )=[z -h (x )]T

R -1

[z -h (x )](2)s.t.　c (x )=0(3)式中　c (x )为l 维零注入功率等式约束函数向量,

通常l n 。1.2　不考虑等式约束的状态估计模型

当不考虑等式约束时,要使目标函数J (x )最小,则有:

J (x )

x

=0(4)根据式(4)可得估计迭代式为:

Δx (k )

=(H T

R -1

H )-1

H T

R -1

Δz x (k +1)=x (k )+Δx

(k )(5)

式中　H 为m ×n 阶雅可比矩阵,且H =

h (x )

x ;Δz

为m 维计算残差列向量,且有Δz =z -h (x (k )

)。

42

2001年11月10日

N ov.10,2001

1.3　利用等式约束修正估计结果

我们知道,利用式(5)估计出的结果不可能完全满足式(3),因此为了使估计结果更能符合式(3),可以采用式(3)来修正式(5)的结果,具体推导如下。

对于等式约束的极值可通过拉格朗日多项式进行求解,式(2)、式(3)的极值函数可表示为:

L (x )=J (x )+λT

c (x )(6)

式中　L (x )为等式约束下的目标函数;λ为l 维约

束向量。

根据极值求解的步骤,可得极值函数的矩阵表达式[2]

为:

H T R -1H

C T C 0Δx λ=H T R -1

Δz

-Δc

(7)

或

H T R -1H Δx +C T λ=H T R -1Δz

C Δx =-Δc

(8)

式中　C 为l ×n 维矩阵,且有C = c (x )

;-Δc 为虚

拟零注入点的残差矢量。由于中间变量λ的引入,使得式(7)的矩阵维数增加,增加了较大的计算量,随着虚拟零注入测点个数的增多,尤为明显。为了减小这一不利因素,对式(8)进一步约化。

当网络可观测时,H T R -1

H 矩阵的逆矩阵存在,另外由于C 阵为节点注入型量测的雅可比矩阵元素,其行向量必然线性无关,即C 矩阵的秩等于虚

拟零注入量测数,因此C [H T R -1H ]-1C T

的逆矩阵必然存在。根据上述说明,可对式(8)消去中间变量λ

,有:Δx =[H T R -1H ]-1H T R -1Δz -[H T R -1H ]-1C T ·

[C [H T R -1H ]-1C T ]-1

·

(C [H T R -1H ]-1H T R -1

Δz +Δc )(9)对于式(9),定义

Δx =Δx 1-Δx 2

(10)且有

Δx 1=[H T R -1H ]-1H T R -1

Δz (11)Δx 2=[H T

R -1

H ]-1

C T

[C [H T

R -1

H ]-1

C T ]-1·

(C [H T R -1H ]-1H T R -1Δz +Δc )(12)假如我们已计算出了Δx 1,且有x 1=x 0+Δx 1;x 0表示每次迭代的初值。利用泰勒公式,对Δc 进行展开:Δc (x 0)=Δc (x 0+Δx 1-Δx 1)=Δc (x 1)-C Δx 1(13)根据式(11)、式(13),则式(12)变为:

Δx 2=[H T R -1H ]-1C T [C [H T R -1H ]-1C T ]-1Δc (x 1)

(14)

从式(14)可看出,由于H T R -1

H 的因子表在对

式(11)的求解中已得到,对于C [H T R -1H ]-1C T

的

求解,可利用该因子表。

设H T R -1H 的因子表为LL T ,L 为下三角矩阵,

代入C [H T R -1H ]-1C T

,得:

C [LL T ]-1C T =C [L -1]T L -1C T =

[L -1C T ]T L -1C T =B T B

式中　B =L -1C T 。

根据消去过程可求得B 矩阵,然后可进一步对B T

B 进行因子化。

因此式(14)的求解可表述为:a .先求出矩阵B ,然后令:

B T

B Δt =Δc (x 1)

b .求出Δt ,然后可求出Δx 2:

H T R -1H Δx 2=C T Δt 由于B T

B 仅是一个l ×l 阶的对称正定矩阵,因此对Δt 的求解不会增加太大的计算量和存储内存。

然后利用式(14)的计算结果去修正式(10),得:

x 1′

=x 1-Δx 2(15)从整个推导过程来看,这样处理并不会损失计算精度,但却减少了计算量(相当于把一个高维矩阵解耦成两个低维矩阵进行计算),当矩阵阶数较高并且零注入节点较多时,这种分解会较大地提高计算效率。另外,特别需要注意的是,虽然C 阵与H 阵有相同的形式,但最大的区别是,C 阵只与网络的结构有关,而与量测设备的状态无关,只要网络的结构不发生变化,就不必重新形成C 阵,这一点对状态估计的实时应用有较大的好处。1.4　带等式约束的状态估计模型

根据 1.2节、1.3节的推导,最后得到带等式约束的状态估计模型:

Δx (k )=[H T R -1H ]-1H T R -1Δz (x (k )

)

(16)x (k +1)=x (k )+Δx (k )(17)

　Δx ′(k )=[H T R -1H ]-1C T [C [H T R -1H ]-1C T ]-1

·

Δc (x (k +1))(18)

x ′(k +1)=x (k +1)-Δx ′(k )(19)

假设X x 为状态量的收敛判据,X r 为残差r 的门槛值。其具体实现为:

a .置初值x (0)

。

b .形成矩阵H T R -1H 因子表、C [H T R -1H ]-1·C T

的分解矩阵B 。

c .置迭代次数k =0,并置最大误差为0

。d .计算残差Δz (x (k ))。

e .利用式(16)、式(17)计算Δx (k ),x (k +1)

。f .最大误差为max |Δx (k )|。

g .计算Δc (x (k +1)),并判断max |Δc (x (k +1)

)|是否小于X r ,若是,转至下一步;否则,利用式(18)、式

(19)计算Δx ′(k ),x ′(k +1)

,并计算最大误差。

43

·学术论文·　倪小平等　一种带有等式约束的状态估计新算法