CN105958470A - 一种基于双线性凸优化理论电力系统双线性抗差估计方法 - Google Patents

Abstract

本发明公布了一种基于双线性凸优化理论的电力系统双线性抗差估计方法，本发明首先引入双线性理论，将非线性量测方程转化为两个阶段的线性量测方程；然后计及量测粗差的稀疏特性，将抗差估计转化为两个阶段的严格凸优化问题；每个阶段首先基于ADMM辨识稀疏的量测粗差，剔除量测中的粗差后采用WLS求解，保留了WLS的优点。IEEE标准系统以及国内实际电网的测试结果表明，由于双线性理论的引入，本发明提出的方法计算效率高于传统的WLS估计器，而ADMM技术很好地辨识了稀疏的量测粗差，使得本发明提出的方法的估计精度也优于传统的抗差估计器。

Description

一种基于双线性凸优化理论电力系统双线性抗差估计方法

技术领域

发明涉及一种基于双线性凸优化理论电力系统双线性抗差估计方法，属于电力系统监测、分析和控制技术领域。

背景技术

状态估计根据遥测的生数据估计出电力系统实时的运行状态，基于状态估计的结果，能量管理系统(energy management system，EMS)进行一系列后续分析计算，因而状态估计是EMS至关重要的一部分。传统的加权最小二乘(weighted least estimation，WLS)估计在量测噪声服从严格的高斯分布时，能够高效地估计出系统最佳的状态。然而由于量测仪表的老化、数据的远距离传输、甚至是人为恶意注入的坏数据，不可避免地使WLS的估计结果受不良数据(或量测粗差)的影响，从而偏离实际的真值。

由于能够抑制量测粗差对估计精度的影响，抗差估计引起了国内外学者的广泛研究，其中以加权最小绝对值(weighted least absolute values，WLAV)、非二次准则、最小中位数(least median of squares，LMS)为主。此外，基于指数型目标函数、最大合格率的状态估计也是抗差估计领域的新方法。抗差估计器以增加计算复杂度为代价，提高了状态估计的精度，然而相比于WLS估计器，较低的计算效率一定程度上也限制了其在工程实践中的应用。

在交流电力网络模型中，目前状态估计主要依赖于数据采集与监控系统(supervisory control and data acquisition,SCADA)提供的非线性量测，这使得状态估计本质上是非线性(非凸)优化问题，常用的方法是以高斯－牛顿法近似线性化迭代求解，这种处理方法可能会有以下不足：1)对初值敏感；2)易陷入局部最优解；3)收敛性难以保证。

综合以上叙述可知，目前电力系统状态估计主要面临两方面的困难：1)量测量与状态量之间的非线性关系，使得状态估计等价于求解非凸优化问题；2)量测粗差的存在对WLS估计精度有较大影响，传统的抗差估计方法虽能抑制量测粗差的影响，但计算效率偏低。为此，本发明首先引入提出的双线性理论，通过变量代换，将非线性量测方程转化为两个阶段的线性量测方程；然后计及量测粗差的稀疏特性，以凸l₁范数正则化稀疏的量测粗差向量，将抗差估计转化为两个阶段的严格凸优化问题；每个阶段首先基于交替方向乘子算法(alternating direction method of multipliers，ADMM)辨识稀疏的量测粗差，剔除量测中的粗差后可认为量测噪声服从高斯分布，因而采用WLS求解，保留了WLS的优点。IEEE30、118节点标准系统以及我国2个实际省网系统的测试结果表明，由于双线性理论的引入，本发明提出的方法计算效率高于传统的WLS估计器，而ADMM技术很好地辨识了稀疏的量测粗差，使得本发明提出的方法的估计精度也优于传统的抗差估计器。

发明内容：

本发明所要解决的技术问题是针对现有技术存在的不足而提供一种基于双线性凸优化理论电力系统双线性抗差估计方法。

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

本发明为一种一种基于双线性凸优化理论电力系统双线性抗差估计方法，其特征在于所述方法是在计算机中依次按以下步骤实现的：

1)获取电力网络的参数信息，包括：输电线路的首端节点和末端节点编号、支路π型等效电路的电阻、电抗、对地并联电导、电纳以及变压器变比和阻抗。

2)获取电力系统的量测参数，包括：母线电压幅值、节点注入有功、节点注入无功、支路有功、支路无功。

3)定义网络参数集合A＝{网络拓扑，量测类型，系统参数}，状态估计计算集合T＝{B^*，C^*，((B^*)^TB^*)^-1和((C^*)^TC^*)^-1的分解因子表，G_b，S}。若集合A不变，则直接读取上个的集合T信息，若集合A变动，则重新计算并保存集合T。

4)以ADMM辨识一阶段的稀疏量测粗差，剔除粗差后采用WLS求解一阶段中间变量y。

5)中间变量非线性变量，即u＝f(y)，并求解中间非线性变换的雅克比矩阵。

6)计算而阶段等级权重，以ADMM辨识二阶段稀疏量测粗差，剔除量测粗差后采用WLS求解。

7)重新计算二阶段量测变量，以及非线性雅克比矩阵。

8)以WLS重复二阶段线性状态估计。

9)输出状态量以及量测量的估计值，结束。

附图说明：

图1：本发明方法流程图。

图2：以IEEE30为标准测试节点，本发明方法与WLAV、SHGM算法S^V、S^θ的pdf分布图。

图3：以IEEE118为标准测试节点，本发明方法与WLAV、SHGM算法S^V、S^θ的pdf分布图。

具体实施方式：

下面结合附图对发明的技术流程进行详细说明：

1双线性状态估计模型

双线性理论利用变量代换的思想，将电力系统非线性状态估计转化为两个阶段的分步线性状态估计问题，且两个线性状态估计之间包含一步变量的非线性变换。

1.1一阶段线性状态估计

对于连接母线i与母线j的每条支路，定义如下变量：

K_ij＝V_iV_jcosθ_ij

L_ij＝V_iV_jsinθ_ij

式中：V_i、V_j分别为母线i、j的电压幅值，θ_i、θ_j分别为母线i、j的电压相角，θ_ij＝θ_i-θ_j。

对于系统中的每条母线，定义电压幅值平方为新的变量：

假定系统包含N条母线，T条支路，则一阶段线性状态估计引入N+2T维状态量y：

y＝{U_i,K_ij,L_ij}

则SCADA系统提供的量测量与状态量y成如下线性关系：

支路功率量测：

式中：g_ij、b_ij分别为支路π型等效电路的电导、电纳，g_si、b_si分别为母线i侧对地电导、电纳。

节点注入量测：

电压幅值量测：

式中：e为量测误差向量，且假定e服从正态分布；

m维量测向量z与状态量y可表示为如下线性关系：

z＝By+e_z

基于传统的WLS算法求解，y的估计值为：

式中：W_z为量测z的权重矩阵，假设第i个量测标准差为σ_i，则W_z＝diag(1/σ₁ ²,…,1/σ_m ²)。G_b为增益矩阵，且G_b＝B^TWB。

1.2中间变量非线性变换

中间变量的非线性变换为等维数变换，定义如下N维变量α，T维变量α_ij,θ_ij：

α_i＝lnU_i＝2lnV_i

α_ij＝ln(K_ij+L_ij)＝α_i+α_j

令u＝{α_i,α_ij,θ_ij}，则N+2T维变量u与y呈非线性关系

中间变量u的权重矩阵W_u为：

式中：F为f_u在处的雅克比矩阵。

1.3二阶段线性状态估计

定义2N-1维状态量x＝[α θ]^T(参考母线的相角固定为0)，则二阶段状态量x与中间变量u呈如下线性关系：

式中：I为单位阵，A为节点关联矩阵，A_r为不包含参考母线的节点关联矩阵。

基于WLS求解，x的估计值为：

式中：G_c为增益矩阵，且G_c＝C^TW_uC。

2凸优化抗差估计模型

2.1量测噪声的标准化

对于一阶段线性状态估计，令r＝z-By，W_z＝L_z ^TL_z，则WLS算法的优化目标等价于：

式中：v为服从标准正态分布的随机变量，其协方差矩阵为单位阵I；L_z＝diag(1/σ₁,…,1/σ_m)。

将一阶段量测方程标准化，即将量测方程两侧同时乘以L_z，则可得：

z^*＝B^*y+e

式中：z^*＝L_zz，B^*＝L_zB；e_i～N(0,1)，

同样，对于二阶段线性状态估计，其量测权重W_u可按如下的方式分解：

将二阶段量测方程两侧乘以L_u标准化二阶段量测方程，可得：

u^*＝C^*x+e

式中：u^*＝L_uu，C^*＝L_uC。

对量测噪声标准化后，两个线性阶段的量测z^*、u^*的噪声均服从标准正态分布，因而状态估计的目标函数无需计及权重矩阵。

2.2稀疏量测粗差向量模型

实际量测的误差除了服从正态分布的量测噪声外，还包括量测粗差，其主要产生于量测仪表、通信故障甚至是人为恶意注入的坏数据(对状态估计精度的影响更为严重)。以向量o描述量测粗差，则量测量与状态量的关系，可以更合理地表示为：

z^*＝B^*y+e+o

式中：o(i)当且仅当第i个量测为坏数据时不为0，因而向量o是非常稀疏的。

利用向量o的稀疏特性，一阶段线性状态估计问题等价于求解

式中：||o||₀为o非零元素的个数，λ₀>0为l₀范数的正则化因子。

由于l₀范数的正则化，上式属于NP难求解的优化问题，借鉴压缩传感方面的理论，凸l₁范数可以用作求解稀疏向量的启发式算法，即

式中：λ₁>0为l₁范数的正则化因子。

2.3Lasso优化辨识量测粗差向量

在上式中，假设o已知，则y的估计值可表示为：

因而可由o描述y的最佳估计，则稀疏向量o等价于求解：

式中：S＝I-B^*((B^*)^TB^*)^-1(B^*)^T，为残差灵敏度矩阵。

上式属于经典的Lasso优化问题，本发明采用可高效求解Lasso的分布式凸优化算法ADMM。ADMM本质上是求解含等式约束的凸优化问题，为构造等式约束，引入向量p，将上式转化为：

min f(o)+g(p)

s.t.o-p＝0

式中：

引入拉格朗日乘子u以及惩罚因子ρ＞0，上式的增广拉格朗日函数为：

则上式的迭代求解步骤为：

式中：k为迭代次数；II_a(φ)＝(φ-a)₊-(-φ-a)₊，(φ)₊＝max(φ,0)。

以ADMM算法辨识出稀疏量测粗差o后，理论上可以认为z^*－o的误差服从标准正态分布，因而WLS可以高效地解出y的无偏估计，即采用计算本发明提出的方法与WLS的区别在于，利用ADMM辨识稀疏的量测粗差，抑制了量测粗差对状态估计结果的影响，同时保留了WLS的优点。

2.4参数λ₁和ρ的选择

参数λ₁衡量了l₁范数正则化的程度，因而参数λ₁的选择会影响算法的抗差性能。以ADMM求解Lasso优化问题，λ₁可选取：

λ₁＝Cλ_max

式中：λ_max＝||S(Sz^*)||_∞＝||Sz^*||_∞，C为大于0的常数，本发明选取C＝0.1。

参数ρ的设置不会影响ΑDMM的全局最优解(即不影响算法的估计精度)，但一定程度上会影响算法的收敛性能，本发明选取ρ＝1。

3算例分析

本发明测试的算例包括IEEE30、118节点标准系统以及2个国内实际的省网系统。IEEE标准节点量测数据由严格潮流真值添加随机噪声得到，其中功率量测噪声的标准差为0.01，电压幅值量测噪声的标准差为0.004，不良数据在真值的基础上随机加减[5,30]倍量测标准差，系统的量测冗余度介于3～4之间，随机添加3％不良数据。

不同抗差估计算法的抗差性能主要取决于其优化的目标函数，本发明选取抗差性能较好的WLAV估计器、广义M估计(Schweppe-type generalized M-estimator with Huber psi-function，SHGM)，与本发明提出的方法进行抗差性能比较。

3.1IEEE标准系统估计精度比较

为更好地比较三种抗差估计器的对状态量的估计精度，定义如下指标：

式中：V^ex、θ^ex分别为电压幅值与相角的真值，默认节点1为平衡节点。

随机模拟1000次，S^V、S^θ的概率密度函数(probability density functions，pdf)及均值μ、标准差σ如图2、图3所示。

由图2、图3可知，本发明提出的方法S^V、S^θ的pdf曲线更接近于坐标原点(即均值μ更小)，说明了在不良数据含不同比例、不同随机组合的情况下，本发明提出的方法整体而言有着更高的估计精度；此外，本发明方法的S^V、S^θ的pdf曲线相对更为“瘦高”(即标准差σ更小)，说明了本发明提出的方法在不良数据比例、组合变化的情况下估计精度波动相对更小，因而估计性能更为稳定。

3.2省网系统的测试

为验证本发明提出的方法在更大规模系统上的估计性能，选取2个实际的省网系统(记为省网A和省网B)，其中省网A含736条母线，959条支路，省网B含1518条母线，2034条支路，表1给出了不同估计器下省网A、B的合格率。由表1可知，相比于WLAV、SHGM、WLS估计器，基于本发明提出的方法，省网A、B状态估计的合格率均得到了提高。

3.3计算效率的测试

本小节比较不同估计器的计算效率，本发明提出的方法计算流程如图1，当网络参数集合A不变时，状态估计计算集合T只需读取上个量测断面的数值；而当A变化时，需要重新计算T。不同规模的系统下，各种状态估计器的计算效率如表2所示。

由表2可知，传统的抗差估计器WLAV、SHGM的计算效率明显不如WLS，这也是WLS广泛应用于工程实践最重要的原因。本发明提出的方法，随着系统规模的增大，计算效率上的优势也越明显。特别是对于较大规模的省网系统，即便网络参数集合A发生变化，本发明提出的方法计算效率依然高于WLS。

表1不同估计器下实际系统的合格率

表2不同估计器计算效率比较

Claims (1)

1.一种基于双线性凸优化理论电力系统双线性抗差估计方法，其特征主要在于：

1)首先引入双线性理论，对于连接母线i与母线j的每条支路，定义如下变量：

K_ij＝V_iV_jcosθ_ij

L_ij＝V_iV_jsinθ_ij

式中：V_i、V_j分别为母线i、j的电压幅值，θ_i、θ_j分别为母线i、j的电压相角，θ_ij＝θ_i-θ_j。

对于系统中的每条母线，定义电压幅值平方为新的变量：

U_i＝V_i ²

假定系统包含N条母线，T条支路，则一阶段线性状态估计引入N+2T维状态量y：

y＝{U_i,K_ij,L_ij}

则m维量测向量z与状态量y可表示为如下线性关系：

z＝By+e_z

2)中间变量的非线性变换为等维数变换，定义如下N维变量α，T维变量α_ij,θ_ij：

α_i＝lnU_i＝2lnV_i

α_ij＝ln(K_ij+L_ij)＝α_i+α_j

令u＝{α_i,α_ij,θ_ij}，则N+2T维变量u与y呈非线性关系

3)定义2N-1维状态量x＝[αθ]^T(参考母线的相角固定为0)，则二阶段状态量x与中间变量u呈如下线性关系：

u＝Cx+e_u

式中：I为单位阵，A为节点关联矩阵，A_r为不含参考母线的节点关联矩阵。

4)然后计及量测粗差的稀疏特性，以向量o描述量测粗差，则量测量与状态量的关系，可表示为：

z^*＝B^*y+e+o

5)借鉴压缩传感方面的理论，凸范数可以用作求解稀疏向量的启发式算法，即：

式中：λ₁>0为范数的正则化因子。

6)假设o已知，则y的估计值可表示为：

7)可由o描述y的最佳估计，即

8)步骤7属于经典的Lasso优化问题，本发明采用ADMM求解，为构造等式约束，引入变量p,即：

min f(o)+g(p)

s.t.o-p＝0

式中：g(p)＝λ₁||p||₁。

9)步骤(8)的迭代求解步骤为：

o^k+1＝(S^TS+ρI)-¹(S^TSz*+ρ(p^k-u^k))

＝(S+ρI)-¹(Sz*+ρ(p^k-u^k))

p^k+1＝II_{λ 1/ ρ}(u^k+o^k+1)

u^k+1＝u^k+o^k+1-p^k+1

式中：k为迭代次数；II_a(φ)＝(φ-a)₊-(-φ-a)₊，(φ)₊＝max(φ,0)。

10)以ADMM算法辨识出稀疏量测粗差o后，理论上可以认为z^*－o的误差服从标准正态分布，因而WLS可以高效地解出y的无偏估计，即采用计算本发明提出的方法与WLS的区别在于，利用ADMM辨识稀疏的量测粗差，抑制了量测粗差对状态估计结果的影响，同时保留了WLS的优点。