快速分解状态估计

2024-07-23

快速分解状态估计（精选4篇）

快速分解状态估计篇1

0 引言

电力系统状态估计是能量管理系统(EMS)的重要组成部分,是高级应用软件执行的基础[1],在该领域无论是理论还是工程实践均有大量的研究成果。随着同步相量测量技术,乃至广域系统在电力系统中的应用,使电压和电流相量直接量测得以实现,该技术对状态估计的影响得到广泛关注[1~16]。

基于广域量测信息所进行的线性估计[7,8],是理想化的,尚不具有实用性。因此,计及SCADA和PMU的混合估计研究更符合实际,也一直是研究热点。由PMU的电压相量测量信息,直接将电压幅值和相角作为量测方程,虽算法简单,易于实现,但对估计贡献效果不明显[9]。通过将电压电相量转换成相应的潮流或电流[9,10]的方式在迭代过程中增加修正方程,扩大了其对信息矩阵的影响,但本质上利用的还是电压相量量测,对PMU间接提供的支路电流相量量测没有充分利用。文献[11]就电流相量对估计影响进行研究,指出其可提高状态估计的精度。在基于电流相量的研究中,采用非线性估计时,将电压、电流相量转换成支路有功、无功潮流形式,或者相邻节点的电压相量;而采用线性量测方程时,SCADA量测转换成等效的电流量测[12],而等效量测的计算需要依赖每次计算得到的状态变量,因此算法退化为迭代算法,且电流相量在转换成直角坐标时会带来了数值稳定的问题[5]。而无论采用哪种方式进行量测转换,均需按误差传递规律计算其相应的间接量测的误差方差分布。文献[13,14]在充分保留传统状态估计研究成果基础上,利用传统估计结果得到的伪量测与PMU的线性量测再次进行估计,但在PMU配置数量不多时,由此带来计算效果并不明显。以上种种研究表明,PMU带来的混合估计使问题求解及模型建立变得更为复杂,计算规模增大。

然而,目前电力系统中基于SCADA的状态估计完全可以实施,PMU配置的目的主要是用于电力系统的快过程监测及安全控制,因此往往配置在关键节点或关键输电断面上,实现系统局部可观,便于对系统运行状况做出快速判断[18]。可见,PMU的功效不在于状态估计,但PMU使系统产生了自动划分,似乎可将复杂系统状态估计自然利用分解协调原理进行,这样无疑可提高状态估计的速度。因此,本文在传统状态估计基础之上,利用PMU对网络进行分割,建立分解协调状态估计模型,并对各子系统参考点及边界点协调处理进行了详细阐述。最后采用IEEE-39及IEEE-118节点验证了模型的有效性。

1 问题的总体描述

在统一参照系下,依据PMU配置所对应的量测,其优越性在于系统实现局部可观,这也是PMU引入电力系统的意义所在。在此背景下,仅基于SCADA系统的状态估计,可以实现局部估计,以及利用PMU特性的分散估计,这就是本文的基本思想。

根据PMU使系统实现局部可观的性质,必然使系统可实施分散的估计,即划分为若干可观测岛以实施分解协调的估计方法。系统划分完全依据PMU的配置地点,以给出划分定义:

PMU可观测岛:有直接电气连接的PMU配置节点集合的割集;

SCADA可观测岛:仅含有依据SCADA系统量测量,且不包含PMU配置节点的系统部分节点构成的割集;

边界节点:PMU可观测岛与SCADA可观测岛交集所对应的节点。

上述的示意图见图1所示。

根据此解耦方式,可将系统节点分为三类:

PMU bus,节点安装有PMU;

Boundary Bus,相邻节点至少有一个安装有PMU;

SCADA Bus,无PMU并且相邻节点也无PMU。

相应节点分别为

将系统量测向量表示为

其中:sz为与SCADA BUS相关的量测,包含支路潮流,节点注入,节点电压幅值;zb I为Boundary Bus的注入量测及零注入伪量测;zt为PMU Bus与Boundary Bus之间联络线上的潮流量测;zp I为PMU Bus的节点注入;zps为PMU Bus相关的相量量测,包括节点电压相量及关联支路电流相量。

依据上述划分,对PMU岛建立量测方程表达式如式(1)。

其中:ae为量测误差,h为量测方程的非线性向量函数。

对式(1)进行状态估计,可得到联络线上潮流zˆt,即相对于SCADA岛边界点的注入,可构造伪量测

伪量测与其它量测一起形成SCADA岛量测方程

根据以上分析,可得到算法流程如下:

(1)预处理:根据PMU节点安装位置,对全网进行分割,得到若干PMU及SCADA观测岛;

(2)对PMU观测岛进行混合量测估计,得到联络线上潮流

(3)对SCADA岛进行估计;

(4)协调各子岛参考节及边界点相角,得到全网统一的状态变量解。

2 几个关键问题的处理

2.1 PMU局部估计模型

PMU岛在整个网络中作为独立部分,可以看做一个广义节点,通过该广义节点的割集电流之和为零。

以典型π型支路为例(图2),若j,k仅有一端j在岛内,则流出岛内的电流

若jk均在岛内,该支路流出岛内的电流为

当该岛内无注入时,流出岛内电流代数和为

为提高估计精度及收敛速度,对无注入的PMU岛可建立带有等式约束的状态估计模型

式中:R-1为量测误差阵。

对式(5)构造拉格朗日函数

式中:λ为拉格朗日乘子。

对式(7)x,λ分别求极值,可得

对该非线性方程组求解,迭代格式如下:

式中:为2×n维矢量,n为岛内待求状态变量个数。

2.2 混合量测匹配问题

在进行状态估计时,需要指定参考节点,令其电压相角为0,其余节点相角为相对于该参考节点的相角差。对于混合量测系统来说,增加了PMU相角量测,而该相角量测参考点为全球定位系统GPS的参考点,两者之间需要进行转换。一个简单的方法是在指定安装PMU节点做为参考节点,岛内其它相角量测以该点为参考进行修正

式中:δi,δi分别为相角量测绝对值及校正后的值。参与状态估计的相角相应的量测方差则由σi2修正为2σi2。

由于参考节点电压相角没有直接参与状态估计,亦即没有量测残差,也就无法采用传统的坏数据识别方法进行识别,而相角误差对状态估计带来的影响相对于其它量测更为明显,对于参考节点相角的坏数据可采用文献[17]的方法进行识别。

2.3 参考点的协调

各观测岛分别有自己的参考节点,并且边界点的状态变量估计值分别由相邻子岛得到多两个不同解,协调层负责将各岛相角统一到同一参考点上并对边界点进行一次校正。

u=[u2u3...un],iu是第i子岛的参考节点相角,子岛1的参考节点做为全网参考节点1u=0。

进而,协调层估计模型可表达为

是伪量测与状态变量的线性表达,eco是误差向量,Hcoθco可详细表述为

θa(-i)pseudo为θa在第i个子岛的估计结果。各伪量测的量测误差σi2为各子岛估计收敛后G-1相应对角元素值。

3 算例分析

本文程序为VC6.0编制,程序运行于cpu1.8GHz,内存512MB的PC机。分别对IEEE39节点及118节点系统进行模拟分析。量测数据为系统真值叠加相应的正态分布随机量测误差形成,量测误差取自文献[16](σvoltage=0.004,σinjection=0.01,σphasor=0.0001)。SCADA量测配置保证统可观且无关键量测,取每条支路两端潮流及节点电压幅值量测。对IEEE-39节点系统,PMU安装位置取自文献[18],选择节点1,3,16,18进行分割,得到4个SCADA估计岛,3个PMU局部估计岛,见图3所示。

算例分别模拟了无坏数据,联络线出现一个坏数据;边界点注入出现坏数据;PMU电压相角出现坏数据的四种情形。以下为39节点系统的详细情况。

case 1:SCADA量测及PMU量测均无坏数据。

对PMU量测采用转换成支路潮流量测的方式。系统共有39个电压幅值,92对潮流量测,16对支路潮流间接量测,总量测数m=255,n=77。表1为集中估计与分布估计结果比较。

分布式估计J(x)要大于集中式估计,这是因为协调层与各子岛之间没有迭代过程产生的,但由于系统未有坏数据,J(x)均小于阈值,在计算时间上,分布式估计要优于集中式估计。

case 2:在联络线上产生一个坏数据

联络线1-2潮流P12真值为-1.18,将其替换为-2.18进行估计,在对节点1形成的局部PMU估计岛进行估计时,J(x)=93.9,使用最大标准残差法进行识别,P12的标准残差Nri=-9.7,坏数据被有效识别;将P12的量测替换为2.18,得到J(x)=93.39,相应的Nri=9.67,坏数据同样被识别。无论坏数据发生在联络线的哪一端,在局部估计岛丰富量测条件下,均可得到校正。

case 3:边界点注入出现坏数据

在节点39处增加一个注入量测,将真值P39=-1.04替换为-2.04进行估计,其中SCADA岛估计的目标函数J(x)=58.19,相应Nri=-7.44,去掉该坏数据,重新估计,目标函数下降为J(x)=1.9。

case 4:PMU节点电压相角出现坏数据。

节点3的相相角真值为θ3=-0.172204 rad,当量测数据增加5°噪声后进行估计。在使用量测变换的方法进行估计时,J(x)=5.4687,ˆθ3=-0.176646 rad,估计误差为0.0044 rad。如仅使用电压相量进行混合估计J(x)=437.07,相角量测标准残差为Nri=-20.9。说明将直接量测转换为多个相关间接量测时,相角误差被分摊到多个关联支路上,降低了对系统的影响。也增加了系统识别坏数据的难度。

在对118节点的模拟中,量测配置仍为支路潮流及节点电压幅值,PMU配置选取联络线38-65两端节点及电机节点24,69,100共5个PMU。系统分割出3个SCADA估计岛及4个PMU估计岛。最大SCADA估计岛为66个节点。从集中与协调估计结果可以看出,在无坏数据的情形下,协调估计目标函数仍然大于集中式估计,但在计算时间上有明显提高。H阵在集中估计时为908×236维,协调估计最大维数482×137,接近集中式估计的一半。

表中协调估计目标函数值略高于集中式估计,而协调估计的整体收敛速度仍然优于集中式估计。

4 结论

利用PMU节点附近量测冗余高及状态量直接可知的特点,基于SCADA系统可观的前提下,充分考虑当前PMU节点安装位置特点,本文提出了一种基于SCADA/PMU混合量测系统的分解协调状态估计方法。经理论研究及算例分析,该算法有如下特点:

(1)雅可比矩阵及信息阵规模下降;

(2)整体计算速度快;

(3)能够处理边界点注入功率;

(4)能够协调多观测岛参考点相角。

充分利用传统估计的成果,不需更改程序,只需增加估计前的分割处理及估计结束之后进行各区域之间相角的协调及校正。

由于算法基于分布式计算模式,计算结果与集中混合估计结果不完全相同,但完全能够满足估计精度的要求,在计算速度上,由于本文算例是在单CPU机上进行串行处理,而各SCADA估计岛是完全可以并行处理的,如在可并行处理多任务多cpu机器上运行,其运算速度还可进一步提高。

快速分解状态估计篇2

为了实现互联电网经济调度和安全稳定运行,有必要进行全网一体化的仿真分析,而获得一体化的潮流断面是开展各类仿真分析的前提和基础。目前,在多区域互联电网中,各区域调度中心只采用本区域内部量测以及外网等值信息进行独立的状态估计[1,2,3,4]。然而,由于外网信息不够准确或者难以获取,影响了区域电网状态估计准确性,导致难以从独立状态估计结果形成全网匹配的潮流断面。分布式状态估计不仅可以获得全网一致收敛的状态估计结果,还能整合各调度中心内现有的异构计算资源,是实现互联电网在线数据整合的有效手段。

目前,分布状态估计相关研究的重点是设计合理的计算模型和高性能算法。文献[5,6,7,8,9]提出了一种双层状态估计算法,首先进行底层各分区状态估计,进而将所得结果和联络线、边界节点信息一同作为量测数据用于上层电网状态估计。文献[10,11]对上述算法作了改进,此类算法适用于由不同电压等级网络构成的互联电网的分布式状态估计。文献[12]提出了一种分布式能量管理系统中的异步分布式状态估计算法,采用笛卡尔积空间的分量解法实现全网状态估计的分解协调计算,其优点在于分区之间只交换边界节点状态,且可以实现异步迭代求解,然而此类算法的收敛条件较为苛刻,在大规模互联电网中应用时可靠性难以保证。文献[13]提出了状态估计的分布式算法,通过节点撕裂法实现区域电网的划分,并将通过无约束优化模型描述的全网状态估计问题转化为一个有约束的优化问题,进而实现分区状态估计问题的解耦和边界节点状态的协调计算。文献[14,15]对该方法进行了多方面的优化,这类算法中采用Guass-Newton法求解全网状态估计问题,协调计算中忽略了非边界节点对其他分区状态估计的影响,从而导致其算法收敛性和收敛速度有所下降。

为了克服上述不足,本文提出了一种新的分布式状态估计算法。首先,基于对无约束优化问题的KKT(Karush-Kuhn-Tucker)条件方程的分解,建立了状态估计分解协调计算模型。进一步,采用JFNG(Jacobian-free Newton-GMRES (generalized minimal residual))算法和逆Broyden拟Newton法求解边界协调方程,从而实现了分区状态估计之间的协调计算。所提出的协调计算方法不仅具有较高的收敛性和收敛速度,而且所需数据接口较为简化,可兼容不同分区状态估计算法。最后,以IEEE标准系统和实际电网系统的测试算例验证了所提出的算法的准确性和高效性。

1 分布式状态估计分解协调计算

1.1 分解协调计算模型

如图1所示,以2分区电力系统为例,采用节点撕裂法,可将其分为3个区域S1,S2,SB。其中,S1和S2为2个区域电网(以下称分区侧,其计算服务简称分区服务);SB为由边界节点构成的上级电网(以下称协调侧,其计算服务简称协调服务);N1,in和N2,in分别为2个分区侧的内部节点;NB为边界节点集合。边界节点一般选择没有量测的节点(若有量测可以通过添加虚拟节点的方法解决[16])。协调侧只负责边界节点的电压幅值和相角的估计。分区将边界节点视为参考节点,不对其状态进行修正,只对分区内部的节点、支路以及与边界节点相邻支路状态进行估计。

以加权最小二乘估计为例(其他状态估计算法也适用[17]),全网集中式的状态估计可表示为:

$\min J (x) = (z - h (x))^{Τ} G (z - h (x)) (1)$

式中:x为待预估的网络状态向量,一般为节点的电压幅值、相角,对于含有n个节点的系统,x的维数为2n-1(平衡节点的相角恒定);h(x)为量测值的理论值,一般为电压幅值、节点注入有功和无功、联络线的有功和无功,其维数是m;z为量测值,与h(x)物理量对应,维数也是m。

针对图1所示系统,根据节点所在分区的不同,可将量测量分为z1和z2,量测值的理论值分为h1和h2,估计量分为x1,x2,xB,权重矩阵分为G11和G22,即

$G = [\begin{matrix} G_{11} \\ G_{22} \end{matrix}] (2)$

则式(1)可以写为:

$\min J (x) = J_{1} (x_{1}, x_{B}) + J_{2} (x_{2}, x_{B}) (3)$

当J(x)取最小值时,对应KKT条件如下:

$\frac{\partial J (x)}{\partial x} = \frac{\partial (J_{1} (x_{1}, x_{B}) + J_{2} (x_{2}, x_{B}))}{\partial (x_{1}, x_{2}, x_{B})} = 0 (4)$

进一步,可对式(4)进行如下分解:

${\begin{cases} L_{1} = \frac{\partial (J_{1} (x_{1}, x_{B}))}{\partial x_{1}} = 0 \\ L_{2} = \frac{\partial (J_{2} (x_{2}, x_{B}))}{\partial x_{2}} = 0 \\ L_{B} = \frac{\partial (J_{1} (x_{1}, x_{B}))}{\partial x_{B}} + \frac{\partial (J_{2} (x_{2}, x_{B}))}{\partial x_{B}} = 0 \end{cases} (5)$

式中:L1和L2分别对应分区1和分区2内部状态估计收敛的条件;LB描述了全网状态估计收敛时边界节点状态应满足的条件。

由式(5)可见,当给定边界节点状态 $\bar{x}_{B}$ ,分区1和分区2可独立求解L1和L2,从而得到x1和x2;由分区获得收敛的状态估计结果 $\bar{x}_{1}$ 和 $\bar{x}_{2}$ 后,协调侧可求解LB,得到边界节点状态xB的调整量。换言之,在L1和L2收敛前提下,通过求解LB可实现全网状态估计分解协调计算。以下分别介绍分区侧和协调侧的计算模型。

1.2 分区侧求解

在上述分解协调计算中,分区侧可将边界节点视为参考节点,分区i需要求解下述方程:

$L_{i} = \frac{\partial (J_{i} (x_{i}, \bar{x}_{B}))}{\partial x_{i}} = 0 (6)$

式中:i为分区编号; $\bar{x}_{B}$ 为协调侧传递给分区侧的边界节点信息,通常为边界节点幅值、相角,其数值在分区侧计算过程中保持恒定; $J_{i} (x_{i}, \bar{x}_{B})$ 为分区的优化目标,可以表示如下:

$J_{i} (x_{i}, \bar{x}_{B}) = (z_{i} - h_{i} (x_{i}, \bar{x}_{B}))^{Τ} G_{i i} (z_{i} - h_{i} (x_{i}, \bar{x}_{B})) (7)$

由式(6)、式(7)可以得到分区i的求解方程,如下:

$Η_{i i}^{Τ} G_{i i} (z_{i} - h_{i} (x_{i}, \bar{x}_{B})) = 0 (8)$

式中: $Η_{i i} = \partial h_{i} (x_{i}, \bar{x}_{B}) / \partial x_{i}, Η_{i i}$ 和 $h_{i} (x_{i}, \bar{x}_{B})$ 的计算均需要分区内部节点状态xi和边界节点状态 $\bar{x}_{B}$ 。

注意到式(6)求解所需的所有量测都是本区域自有的,因此其求解是准确的,并且无需改变原来的状态估计程序。

1.3 协调侧求解

从式(5)可知,协调侧需要满足下述方程:

$L_{B} = \frac{\partial J_{1} (\bar{x}_{1}, x_{B})}{\partial x_{B}} + \frac{\partial J_{2} (\bar{x}_{2}, x_{B})}{\partial x_{B}} = 0 (9)$

式中: $\bar{x}_{1}$ 和 $\bar{x}_{2}$ 为分区侧计算收敛后的节点状态信息。

结合式(7)可以得到其表达式如下:

$\frac{\partial J_{i} (\bar{x}_{i}, x_{B})}{\partial x_{B}} = - 2 (\frac{\partial h_{i} (\bar{x}_{i}, x_{B})}{\partial x_{B}})^{Τ} G_{i i} (z_{i} - h_{i} (\bar{x}_{i}, x_{B})) (10)$

式(10)的求解在分区侧完成内部状态估计后进行,计算量小,所得结果将发送返回协调侧。

式(9)是一个非线性方程组,它的一般形式是:

$F (U) = 0 (11)$

式中:U=xB;F(U)=LB。

对式(11)采用经典Newton法进行迭代求解:

${\begin{cases} A_{k} Δ U_{k} = - F (U_{k}) \\ U_{k + 1} = U_{k} + Δ U_{k} \end{cases} (12)$

式中:Ak=∂F(Uk)/∂Uk为Jacobian矩阵;下标k为迭代求解的次数。

从式(9)～式(12)可以看出,分布式状态估计中边界节点的修正量的求解是算法实现的关键。如可显式地获得Jacobian矩阵Ak,则能采用常规的Newton法求解。然而,在分布式计算环境中,为了简化通信的数据接口,提高算法整合异构计算资源能力,协调侧应避免显式生成该Jacobian矩阵。JFNG算法[18]是一种不精确Newton法,其特点是使用GMRES方法求解Newton迭代中的线性修正方程。文献[18]设计了自适应的预处理方法,可以利用GMRES和Newton迭代中间结果不断更新预处理矩阵使其逼近所求非线性方程组的Jacobian矩阵的逆,进而加速相关迭代收敛。逆Broyden方法[19]是一类拟Newton法,其特点是使用近似Jacobian矩阵(或直接用其逆矩阵)来构建Newton迭代中的线性修正方程,并利用迭代中间结果修正近似Jacobian矩阵(或其逆阵),使其逐渐逼近真实Jacobian矩阵(或其逆阵),从而保证Newton迭代收敛。这2种方法虽然原理不同,但都可避免显式形成Jacobian矩阵,在只能获得方程输入、输出信息的前提下求解非线性方程组。

因此,本文采用具有Jacobian-free特性的JFNG算法[18,20,21]和逆Broyden拟Newton法[19]求解式(11),进而实现分区状态估计的协调计算。

2 算法实现流程

全网状态估计的分解协调计算流程可以用图2表示。

在此流程中,JFNG算法和逆Broyden拟Newton法可动态地形成和更新Jacobian逆矩阵的近似矩阵,并用于式(12)所示修正方程的求解。为了提高协调求解的收敛速度,减少因迭代计算引发的通信,需要改进JFNG算法和逆Broyden拟Newton法使用的Jacobian逆矩阵近似矩阵的初始化方法。可行的方法包括:①由连续在线运行的分布式状态估计过程中得到的有用历史信息,将上一个时间断面协调计算收敛得到的Jacobian逆矩阵作为下一个时间断面所需的初始矩阵;②采用有效措施获得较为精确的Jacobian逆矩阵。本文采用输入、输出校正法[20,21,22]通过一次数值差分得到初始Jacobian逆矩阵,如图2中虚线框所示。测试表明,此方法对通信量和计算量影响均较小,却可大幅提高协调求解收敛速度。

3 算例测试

3.1 测试系统参数

本文的测试环境同文献[20]。测试系统为IEEE 14节点系统、IEEE 39节点系统和某实际电网(下文称HB系统),系统的参数如表1所示。

下文将讨论测试系统测试算法的准确性、收敛速度和仿真效率。用于测试的算法包括原JFNG算法和逆Broyden拟Newton法,2种算法均含有输入、输出环节。3个系统的量测数据是在准确的全网潮流断面基础上,增减一定大小的随机误差而得到的,随机误差小于原数据5%的称为良断面,否则为劣断面。各系统的仿真都是平启动。

3.2 准确性测试

以3个系统为例对算法进行准确性测试。将计算所得的电压幅值、相角、有功、无功分别与集中计算的相应的量进行对比,IEEE 14和IEEE 39系统结果如图3、图4所示,HB系统比较结果见附录A。

从3个系统的误差曲线可知,本文算法很好地保证了计算准确性。3个系统最大误差见表2。

从以上的测试可以看出,在3个规模相差很大的系统中,分布式状态估计都可以取得与集中式算法非常匹配的结果,各项最大误差可以控制在设定的精度以内,从而充分验证了本文算法的正确性和准确性。

3.3 收敛速度测试

对3个系统的集中式状态估计的迭代次数和分布式状态估计的通信次数进行统计,得到表3。

从表3的测试结果可知,JFNG算法可以将通信次数控制在30次以内,而逆Broyden拟Newton法则将其进一步压缩到7次以内,效率很高。

同时注意到,尽管系统规模相差很大,但分布式计算的通信次数基本恒定,表现出良好的稳定性。

取HB系统优、劣2个断面(良断面随机误差设为2%,劣断面随机误差设为10%),进行对比测试,结果表明,虽然集中计算的迭代次数变化很大,但分布式状态估计的通信次数依然保持不变,体现出很好的抗扰动能力。

3.4 耗时分析

集中式状态估计的耗时主要用于求取修正向量,数据准备过程(包括生成H矩阵、不平衡量等)与之相比可以忽略。由文献[23]可知,总耗时T0为:

$Τ_{0} = α n^{1 + η} (13)$

式中:α为常数,与计算环境相关;n为系统节点个数;η为系数,采用Newton法计算时一般取0.4。

在分布式计算中,总计算时间等于协调时间、通信时间和分区计算时间之和[20]。不考虑通信耗时和协调耗时(因为协调时间很短),分布式计算总耗时T1为:

$Τ_{1} = m t_{a r e a} = m α (\frac{n}{k})^{1 + η} (14)$

式中:m为通信次数;k为分区数目。

由式(13)、式(14)可知本文算法的加速比Ψ为:

$Ψ = \frac{Τ_{0}}{Τ_{1}} = \frac{k^{1 + η}}{m} (15)$

由式(15)可以得到不同分区及通信次数时分布式状态估计的加速比,如图5所示。图5中临界点表示加速比为1的状态,对于2分区～6分区的系统而言,当实际的通信次数分别不大于2,4,6,9,12时,系统的加速比就能大于1。

进一步可以得到3个系统的加速比,见表4。

通过以上的分析和测试可知,系统的分区数目越多,则本文算法获得的加速比越大。

4 分析与讨论

本文对分布式状态估计进行建模,对具有普适性的无约束优化问题进行分解协调,从而获得边界量的严格迭代求解公式。需要说明的是,在本文的分解协调算法中,分区状态估计算法只需满足简单的数据接口,其具体实现算法并无限制,可以采用最小二乘估计、加权最小绝对值估计、合格率为最优等状态估计方法。更进一步,本文的分解协调方法也可以用于设计其他可用无约束优化问题建模的电力系统仿真分析的分布式算法。

通过前文的论述可知,分布式与集中式状态估计的结果一致的充分条件是两者KKT条件方程完全等价。因此,分析不同分布式状态估计算法是否具有同样的正确性,只需比较它们的KKT条件方程是否等价。例如,虽然计算过程不同,但本文所提出的算法与文献[13]方法的KKT条件方程是等价的,其所得结果必然相同。相关推导和分析可参见附录A。

5 结语

本文提出了一种新的分布式状态估计算法。该算法具有严格的数学基础,只需要交换边界信息即可获取全网一致收敛的状态估计结果,能够兼容不同的软硬件平台和分区状态估计程序。本算法具有较高收敛速度和鲁棒性,不仅协调计算所需通信次数较少,而且计算性能与系统规模无明显关系。因此,该算法具有较强的实用性,适用于广域网络环境中的互联大电网在线数据整合应用。

附录见本刊网络版(http://aeps.sgepri.sgcc.com.cn/aeps/ch/index.aspx)。

快速分解状态估计篇3

作为电力系统能量管理系统(EMS)核心部分的电力系统状态估计是系统运行、控制和安全评估等方面的基础。当电力系统由于漏测量或其他原因造成不可观察时,系统会出现病态或接近病态,其状态估计结果的数值稳定性将受到很大影响。

电力系统状态估计的基本加权最小二乘法、快速分解法和基于量测变换的状态估计算法[1]在病态条件下可能无法给出估计值。文献[2]以Tank和Hopfield神经网络为基础建立了一种由主从网络构成的电力系统状态估计神经网络模型,用以摆脱病态问题限制。文献[3]提出一种基于分块QR分解的状态估计方法,把虚拟测量处理为等式约束从而避免了由于权因子分散导致的数值病态问题。在每次迭代中,通过对两个分块矩阵的QR分解和一个稀疏三角线性方程组的求解实现系数矩阵的三角分解以保证分解的数值稳定性。文献[4]运用奇异值分解方法进行状态估计,采用节点注入电流相量量测和节点电压相量量测使量测方程线性化,只需进行一次奇异值分解便可得到状态估计结果。文献[5]运用奇异值分解方法进行谐波状态估计。文献[6]将改进的粒子群进化算法应用到状态估计中,使加权最小二乘法的收敛性得到了改善。

无论是基于矩阵分解的方法[3,5],还是基于人工智能算法[2,6],在解决病态问题过程中都以牺牲计算时间为代价。文献[4]的方法尽管保留了奇异值分解不需要进行可观察性分析和解决病态问题的优点,且能缩短程序运行时间,但无法计及节点注入功率量测和支路功率量测。

本文针对电力系统状态估计的病态问题和算法的计算效率进行研究。首先,在直角坐标系下运用等效电流量测变换技术[7]处理节点注入功率量测和支路功率量测,把信息矩阵转换成常数矩阵;利用每一次迭代得到的节点电压修正等效电流量测值。然后,运用奇异值分解方法进行状态估计。由于在迭代过程中只需对信息矩阵进行一次奇异值分解,节省了程序运行时间。仿真算例验证了本文方法在计算时间上的优势和对病态问题良好的处理能力。

1 等效电流量测变换

在用基本加权最小二乘法进行状态估计时,每次迭代过程中状态估计迭代方程组雅可比矩阵的元素都要重新形成,算法计算效率较低。为提高计算效率,本文采用了一种直角坐标系下的等效电流量测变换方法[7]。

1.1 直角坐标形式等效电流量测变换

网络中量测配置通常采用节点注入功率量测iPmea、iQmea,支路功率量测Pijmea、Qijmea、Pjimea、Qjimea,电压量测Uimea。取节点电压实部和虚部为状态量,将节点注入功率量测和支路功率量测变换为直角坐标形式等效电流量测为[7]

式中,ei、fi表示每次状态估计迭代后的节点电压实部和虚部的估计值。

近年来相量量测装置(PMU)在电力系统中逐步得到应用[8,9,10,11],使得电压相量量测得以实现。通过PMU电压相量量测可提高系统可观察性。因此,本文在上述SCADA量测基础上考虑了PMU电压相量量测,其等效变换为

1.2 直角坐标形式线性化量测方程

直角坐标下每次迭代后节点注入电流估计值Iiest、支路电流估计值Ieijst和Iejist、电压估计值Uiest在忽略对地导纳支路后表达式为

式中:Gij和Bij表示节点导纳矩阵中电导和电纳;gij和bij表示支路电导和电纳;N表示节点总数。

由式(1)~式(8)可以看出,经过等效电流量测变换后的量测方程为一个线性方程组,其系数矩阵为一个常数矩阵,可表示为

式中:z表示量测矢量;x表示状态变量矢量;H表示量测系数矩阵;v表示量测误差矢量。

2 基于奇异值分解的状态估计

电力系统状态能够被表征的必要条件是它的可观察性。一般情况下系统只有是可观察的才能进行状态估计。但是,用奇异值分解方法进行估计时可以不需要对系统进行可观测性分析,同时也可以解决病态问题[4]。

2.1 奇异值分解

奇异值分解是一种重要的正交矩阵分解方法,具有强大的数值稳定性。

奇异值分解[12]是指对于任意矩阵A(m×n),存在列正交矩阵U(m×n)和正交矩阵V(n×n)使得

式中,S=diag(α1,α2,…,αn),iα为矩阵ATA第i个特征值的非负平方根值。若A为n阶方阵,A-1为

2.2 基于奇异值分解的加权最小二乘状态估计

经等效电流量测变换后的加权最小二乘状态估计目标函数可表示为

式中:m为量测数;n为状态变量数;R-1表示系统量测的权重矩阵,R=diag(σ12,σ22,…,σm2)。

要使目标函数为最小的条件为

由式(14)可得状态变量的估计值。状态变量估计值的矩阵形式可表示为

令B=HTR-1H为量测方程的信息矩阵。对信息矩阵进行奇异值分解可得

若S矩阵对角元素都大于0,则可以认为系统状态是可观的[4]。

由于量测矢量z中等效电流量测和电压量测均使用了节点电压估计值,因此在得到新的节点电压估计值后,应该用它们修正量测矢量。

新的节点电压估计值由式(18)计算。

由式(18)迭代求解可得各节点电压实部和虚部的状态估计值。由于信息矩阵B是一个常数矩阵,在迭代过程中只需对其进行一次奇异值分解。

3 算法的计算步骤

本文基于奇异值分解的状态估计算法步骤如下:

1)设定状态变量初值、允许误差ε和最大迭代次数kmax,令迭代次数k=0。

2)生成测量矩阵H和信息矩阵B。

3)对信息矩阵B进行奇异值分解,求得B-1。

4)将状态变量估计值代入式(1)~式(4),计算修正的量测矢量z(k)。

5)用式(18)计算。

6)如果或k=kmax,停止计算并输出结果;否则k=k+1,转4)。

4 算例分析

本文用三种方法对IEEE33节点测试系统[13]进行计算。方法1为本文方法;方法2为运用奇异值分解进行状态估计[5],求解时每次迭代都进行奇异值分解;方法3为基本加权最小二乘法[14]。

量测数据由测试系统潮流计算结果叠加量测的随机误差得到。量测的随机误差按文献[14]方法生成。量测配置方案如表1所示。方案1为状态完全可观;方案2存在状态不可观区域。

用上述三种方法对表1的量测配置方案进行状态估计。三种方法状态估计结果如图1和图2所示。

计算时,通过对方法1和方法2的奇异值分解中得到的矩阵S的对角元素观察可进行可观察性分析。计算结果中,量测配置方案1的矩阵S所有对角元素均大于0,说明该配置方案下的系统状态是可观的;量测配置方案2的矩阵S对角元素存在0元素,说明该配置方案下系统状态是不完全可观的。

由图1可见,三种方法对于状态完全可观系统(量测配置方案1)的状态估计结果均能收敛且结果十分接近,表明了本文方法的正确性。由图2可见,对于不完全可观系统(量测配置方案2),即系统存在病态问题时,只有采用奇异值分解技术的方法1和方法2能进行求解,加权最小二乘法(方法3)不能给出状态估计结果。而且,方法1和方法2在对状态可观区域给出较准确的数值解基础上,还能给出状态不可观区域。由图2可知29号和30号节点为状态不可观节点,这说明奇异值分解不仅能对病态问题进行很好的求解且能辨识出状态不可观区域。

考虑到PMU电压相量量测精度相对较低[15],可能对状态估计结果起到负作用,本文考察了电压相量幅值和相角量测误差的标准差提高为文献[14]方法两倍后的状态估计结果,如图3和图4所示。

由图3和图4可见,电压相量幅值和相角量测误差的标准差提高为文献[14]方法两倍后,本文方法依然能得到可行的估计结果;仍能进行可观察性分析且辨识出状态不可观区域。

表2列出了三种方法在未使用稀疏技术情况下对于状态完全可观系统(量测配置方案1)进行状态估计所需的时间和迭代次数。计算时,允许误差取为ε=10-4;最大的迭代次数kmax=10;程序运行计算机的CPU为intel(R)Core(TM)i3 2.93 GHz,内存为2 GB。

由表2可见,在相同迭代次数(4次)情况下,方法2比方法3多用了0.018 s。这是由于方法2每次迭代时进行奇异值分解求解状态修正量,比基本加权最小二乘状态估计在计算速度上存在劣势。方法1(本文方法)虽然多迭代了一次,但由于计算时只需进行一次奇异值分解提高了每次迭代的计算效率,因而方法1比方法2和方法3分别少用了0.052 s和0.034 s,计算速度提高了49.5%和37.4%,显示出计算速度的优势。

5 结论

快速分解状态估计篇4

在无线通信系统中, 由于复杂电磁环境和各种噪声干扰的影响, 信号在无线信道中传输常会引起畸变或延迟, 所以在接收端对信号进行检测、解调时, 往往要先进行信道估计, 以期对接收信号进行幅度和相位校正, 使接收信号成为原始信号的最佳估计。作为4G的核心技术之一, MIMO-OFDM技术在对抗多径效应、频率选择性衰落以及提高频谱利用率方面, 显示了极大的优越性[1]。同样地, 对MIMO-OFDM无线通信系统进行合理的信道估计[2], 有助于降低系统的误码率及优化系统的功率资源。传统的信道估计方法利用插入导频符号或者训练序列的方法达到信道估计的目的, 但是周期性的发送导频符号或训练序列, 造成了系统整体的频谱利用率不高[3]。盲信道估计方法仅利用接收数据的统计相关特性, 使频谱利用率得到了很大程度提高, 但由于算法复杂度较高, 应用受到一定限制[4,5]。

目前, 子空间方法在盲信道估计的研究中应用较多[6], 文献[7]提出了一种基于子空间方法的OFDM系统盲信道估计方法, 但它没有将其引入到MIMO-OFDM系统中。文献[8]讨论了带有循环前缀的MIMO-OFDM系统子空间盲信道估计方法, 在如何利用虚拟载波方面没有提及。本文拟将文献[7]中的方法引入到MIMO-OFDM系统中来, 并针对盲信道估计算法复杂度高的特点, 提出一种基于系统分解的快速算法, 力求在算法时效性和复杂性方面取得一个好的折中。

1MIMO-OFDM系统数据传输模型

对于MIMO-OFDM系统中基于VC的子空间方法, 模型如图1所示。

系统有M个发射天线和N个接收天线, 总的信道带宽被分成Q个正交的子载波, 其中标号k0～k0+P-1用于传输数据, 剩余的Q-P就是虚载波VC。在传输时间段n, 经过线性调制器调制的数据块表示为:

$d (n) = (d^{Τ} (n, 0), d^{Τ} (n, 1), \dots, d^{Τ} (n, Ρ - 1))^{Τ}$

式中:d (n, k) = (d1 (n, k) , d2 (n, k) , …, dM (n, k) ) T, 每一块均构成一个OFDM块, 考虑J (J>1) 个OFDM块, 得到d= (dT (n) , dT (n-1) , …, dT (n-J+1) ) T。通过OFDM调制, 得到的数据表示为s (n, k) = (s1 (n, k) , s2 (n, k) , …, sM (n, k) ) T, s (n) = (sT (n, 0) , sT (n, 1) , …, sT (n, Q-1) ) T, s= (sT (n) , sT (n-1) …, sT (n-J+1) ) T, 定义 $W (k) = (Ι_{Μ} ⨂ W_{Q}^{- k k_{0}}, \dots, Ι_{Μ} ⨂ W_{Q}^{- k (k_{0} + Ρ - 1)}) ‚ W_{Q} = e^{- j 2 π / Q} ‚ W = (W (0), W (1), \dots, W (Q - 1))^{Τ} ‚ \tilde{W} = Ι_{J} ⨂ W (⨂$ 表示Kronecker内积) , 则可以得到s和d的关系如下:

$s = \tilde{W} d (1)$

M个发射天线和N个接收天线之间的频率选择性信道建模为M×N个相互独立的线性时不变FIR滤波器, 这M×N个信道的长度上界假定为L, 第l (0≤l≤L) 阶信道系数如下:

$h_{l} = [\begin{matrix} h_{11} (l) & h_{21} (l) & \dots & h_{Μ 1} (l) \\ h_{12} (l) & h_{21} (l) & \dots & h_{Μ 1} (l) \\ ⋮ & ⋮ & ⋱ & ⋮ \\ h_{1 Ν} (l) & h_{21} (l) & \dots & h_{Μ Ν} (l) \end{matrix}] (2)$

式中:hij (l) 表示第i个发射天线和第j个接收天线之间的第l阶信道系数, i=1, 2, …, M;j=1, 2, …, N。

定义信道矩阵为:

$Η = [\begin{matrix} h_{0} & h_{1} & \dots & h_{L} & 0 & \dots & 0 \\ 0 & h_{0} & h_{1} & \dots & h_{L} & \dots & 0 \\ ⋮ & ⋮ & ⋮ & ⋱ & ⋱ & ⋱ & ⋮ \\ 0 & 0 & \dots & h_{0} & h_{1} & \dots & h_{L} \end{matrix}] (3)$

定义如下结构的接收数据:r (n, k) = (r1 (n, k) , r2 (n, k) , …, rN (n, k) ) T, r (n) = (rT (n, 0) , rT (n, 1) , …, rT (n, Q-1) ) T, r= (rT (n-J+1) (LN+1∶QN) , …, rT (n-1) , rT (n) ) T, 接收数据可以表示为:

$r = Η s + n = Η \tilde{W} d + n = A d + n (4)$

式中:加性噪声 (AWGN) 向量n和r具有相同的结构, $A = Η \tilde{W}$ 。由文献[9]知, 当A满足列满秩条件时, H可被辨识。

2VC子空间盲估计

接收数据向量r的自相关矩阵Rr可以表示为:

$R_{r} = E {r r^{Η}} (5)$

因为噪声和信号源数据之间是相互独立的, 所以Rr又可以表示为:

$R_{r} = A R_{d} A^{Η} + R_{b} (6)$

式中:Rd=E{ddH}, Rb=E{bbH}, 都假定为满秩矩阵, 因为噪声为AWGN, 所以Rb=σ2bIJNQ-NL。为了比较准确的给出Rr, 估计时需要多个数据块进行统计平均。

对Rr做特征值分解:

$R_{r} = U d i a g (λ_{1}, \dots, λ_{J Μ Ρ}, λ_{J Μ Ρ + 1}, \dots, λ_{J Ν Q - Ν L}) U^{Η} (7)$

式中:U的列向量为特征向量;λi为特征值, 并且λ1≥λ2≥…≥λJMP≥λJMP+1≥…≥λJNQ-NL, 对应于λJMP+1, …, λJNQ-NL的特征向量表示为G1, …, Gg, 令G=[G1, …, Gg], 其中g=JNQ-NL-JMP, G张成ARdAH的零空间, 所以可以得到理想条件下的正交关系:

$G_{i}^{Η} A = 0, i = 1, 2, \dots, g (8)$

对于实际的信道估计, 就可以通过优化下面的二次代价函数来实现:

$\hat{h} = \underset{∥ h ∥ = 1}{\arg \min} {\sum_{i = 1}^{g} G_{i}^{Η} A A^{Η} G_{i}} (9)$

3基于系统分解的快速算法

在宽带MIMO-OFDM系统中, 子载波的数目比较大。以上述基于子空间分解的信道盲估计方法为例, 需要对相关矩阵做奇异值分解, 对一个N列矩阵进行奇异值分解需要O (N3) 的矩阵基本操作, 这样当子载波的数目较大时奇异值分解的计算量是非常大的, 这会给应用该信道算法带来很大的困难。考虑到实际的宽带MIMO-OFDM系统中, 多径信道的功率延迟分布持续时间即信道长度, 大大短于一个OFDM符号的持续时间, 当子信道之间的取样间隔大于信道的长度时, 子信道之间是不相关的。这样就可以将整个宽带系统分解成若干个子系统, 每个子系统分别进行信道盲估计, 然后再进行平均处理。这样可以有效地降低的维数, 大大降低奇异值分解带来的计算复杂度, 提高运算速度。这对信道盲估计的初始化和信道跟踪都是极为有意义的。

考虑将JNQ-NL分解成M组, 那么式 (4) 可以写成:

$[\begin{matrix} r^{(1)} \\ r^{(2)} \\ ⋮ \\ r^{(Μ)} \end{matrix}] = [\begin{matrix} Η^{(1)} \\ Η^{(2)} \\ ⋮ \\ Η^{(Μ)} \end{matrix}] s + [\begin{matrix} b^{(1)} \\ b^{(2)} \\ ⋮ \\ b^{(Μ)} \end{matrix}] (10)$

这样每个分块的子系统中数据的长度变为 (JNQ-NL) /M, 此时奇异值分解的计算量变为原来的1/M2。但是此时忽略了分块子系统之间的相关性, 因此需要分别将信道估计出来之后再作平均降噪 (如式 (11) 所示) , 这样既可以降低计算的复杂度又可以避免估计性能的恶化。

$Η = m e a n (Η^{(1)}, Η^{(2)}, \dots, Η^{(Μ)}) (11)$

4仿真实验分析

在MIMO-OFDM系统中, 仿真参数设置如下:M=2, N=3, Q=16, P=12, k0=2。对不同的i, n和k, di (n, k) 是相互独立的, 均匀分布于{1 -1 j -j}, 所以Rd满足满秩条件。信道长度L=2, 假定在每次信道盲估计的间隔期间, 信道是时不变的, 仿真的信道为:

$\begin{array}{l} h_{11} = [1.523 1 + 2.320 7 i - 0.419 2 + 0.171 2 i 0.277 8 + 0.069 2 i] \\ h_{21} = [1.362 6 + 2.758 3 i - 0.284 2 + 0.144 9 i 0.248 5 + 0.391 9 i] \\ h_{12} = [1.201 3 + 2.502 7 i - 0.185 2 + 0.170 6 i 0.317 9 + 0.108 6 i] \\ h_{22} = [1.436 1 + 2.888 3 i - 0.351 4 + 0.267 1 i 0.382 7 + 0.100 9 i] \\ h_{13} = [1.519 1 + 2.882 1 i - 0.273 3 + 0.363 6 i 0.209 4 + 0.350 3 i] \\ h_{23} = [1.109 8 + 2.674 9 i - 0.222 4 + 0.154 6 i 0.352 1 + 0.294 9 i] \end{array}$

为了测量信道估计方法的性能, 给出一个常用的测量指标:

均方根误差 (RMSE) , 其定义为:

$R Μ S E = \frac{1}{∥ h ∥} \sqrt{\frac{1}{D (L + 1)} \sum_{i = 1}^{D} ∥ \hat{h}^{(i)} - h ∥^{2}} (12)$

式中:D表示仿真的次数; $\hat{h}$ (i) 表示第i次估计的结果。

图2给出了当信噪比为10 dB时, RMSE伴随OFDM数据块变化的关系曲线。从图中可以看出, 采用系统分解快速算法后, RMSE指标与没有进行分解的信道估计性能略微有所降低, 但是随着数据块的增多, 二者的变化关系趋于一致。

图3所示为当信噪比变化时RMSE的对应关系曲线, 这里选取数据块个数为300。在图中, 经过系统分解快速算法进行信道估计后, RMSE比原有算法增大了约0.005 dB, 但是从降低系统算法的复杂度方面来看, 这样的误差是可以接受的。

在图4中, 比较了MIMO-OFDM无线通信系统经过系统分解后的信道估计和原算法信道估计下的系统误码率 (误比特率) 。当信噪比较小时, 系统分解信道估计后的误码率大于原算法信道估计的误码率, 当信噪比增大时, 二者的误码率差异显著减小。

通过以上仿真结果说明, 虽然基于系统分解的信道估计算法较原始VC子空间盲信道估计算法在RMSE和误码率方面性能略微有所下降, 但是它对降低系统的算法复杂度具有很大的贡献作用。这比较适合于一些需要对信道进行快速估计的情形中。

5结语

本文在文献[7]的基础上, 将VC子空间盲信道估计的方法引入到MIMO-OFDM系统中, 给出了系统中的数据传输和接收原理, 并对该情形下的信道估计算法进行了推导。针对VC子空间盲信道估计方法需要对接收数据的相关矩阵进行奇异值分解致使算法复杂度较高的特点, 提出了一种基于系统分解的VC子空间盲信道估计快速算法。仿真实验表明, 算法可明显降低原算法的复杂度, RMSE和误码率性能虽然略微有所降低, 但这些都是可以接受的, 并且可以通过增加信噪比或提高数据块数量的方法得到克服。该算法在现实无线通信环境中的时效性方面来考虑, 是可以加以利用的。

摘要：研究MIMO-OFDM系统中的盲信道估计问题。将OFDM系统中的子空间盲信道估计方法引入到MIMO-OFDM系统中, 建立了带有虚拟载波的系统数据传输模型。针对MIMO-OFDM系统VC子空间盲信道估计算法复杂度高的特点, 提出了一种该算法下基于系统分解的快速方法。仿真实验结果表明, 算法在RMSE性能上稍差于原算法, 但却显著降低了信道估计算法的复杂度, 对快速信道估计应用方面具有一定的指导意义。

关键词：MIMO-OFDM,子空间,盲信道估计,虚拟载波,系统分解

参考文献

[1]许光斌, 周围.4G中的MI MO-OFDM技术[J].信息通信, 2007 (1) :22-25.

[2]SHIN Changyong, ROBERT W, EDWARD J.Blindchannel esti mation for MI MO-OFDM systems[J].IEEETrans.on Vehicular Technology, 2007, 56 (2) :670-685.

[3]周鹏, 赵春明, 盛彬.MI MO-OFDM系统中基于导频辅助的信道估计[J].电子与信息学报, 2007, 29 (1) :133-137.

[4]张玲, 张贤达.MI MO-OFDM系统的盲信道估计算法综述[J].电子学报, 2007, 35 (6) :1-6.

[5]黄辉, 纪纯妹, 倪豆.MI MO-MC-CDMA系统中信道估计算法的比较[J].中国新通信, 2009 (2) :39-42.

[6]杨维, 杨薇, 王新生.STBC MC-DS-CDMA系统基于子空间信道盲估计的多用户检测[J].电子与信息学报, 2008 (2) :379-382.

[7]季伟, 吕明.基于子空间方法的OFDM系统盲信道估计[J].大众科技, 2009 (6) :67-69.

[8]佟卫华.采用循环前缀的OFDM信道盲估计算法[J].电子科技大学学报, 2007, 38 (3) :370-373.

[9]BAI W, HE C.Blind channel esti mation in MI MO-OFDMsystems[C]//Proceeding of Global TelecommunicationsConference.[S.l.]:IEEE Press, 2002:317-321.

[10]何飞, 蒋冬初, 李梦醒.基于MI MO-OFDM系统的信道估计算法分析[J].现代电子技术, 2008, 31 (21) :19-22.