CN105514978A - 一种混合整数二次规划形式的抗差状态估计方法 - Google Patents

Abstract

本发明属于电力系统自动化调度技术领域，尤其涉及一种混合整数二次规划形式的抗差状态估计方法，包括：将节点i的电压向量采用直角坐标形式表示，得到节点i的注入有功功率、注入无功功率；用等效的三个两绕组变压器来表示三绕组变压器，则将电力系统网络中所有支路用统一的π型支路表示；根据π型支路的等值电路图，得到节点i到节点j的支路有功功率和无功功率以及节点i的电压幅值量测方程，并得到混合整数二次规划形式的抗差状态估计模型；利用WLS状态估计法得到MIQCP抗差状态估计模型的初值；对MIQCP抗差状态估计模型进行求解，找出并消除不良量测值；对修正后的量测值再次利用WLS状态估计法进行计算。

Description

一种混合整数二次规划形式的抗差状态估计方法

技术领域

本发明属于电力系统自动化调度技术领域，尤其涉及一种混合整数二次规划形式的抗差状态估计方法(RobustStateEstimationUsingMixedIntegerProgrammingWithQuadraticConstraints,MIQCP)。

背景技术

电力系统状态估计可实现对整个电力系统的全面、实时和精确感知，在此基础上，调度人员可实现对整个电力系统的精确预测、精当决策和精准控制。现在国内外的每一个大型调度中心基本都安装了状态估计器，状态估计已成为电网安全运行的基石。自1970国外学者首次提出状态估计以来，人们对状态估计的研究和应用已经有40多年的历史了，这期间涌现出了各种各样的状态估计方法。

目前，在国内外应用最为广泛的状态估计是加权最小二乘法(Weightedleastsquares,WLS)。WLS模型简洁，求解容易，但是其抗差性很差。为了增强抗差性，一般有两种方法。第一种是在WLS估计之后加入不良数据辨识环节，例如最大正则化残差检验法(Largestnormalresidual，LNR)或估计辨识方法等；另一种是采用抗差状态估计方法。目前，国内外学者已经提出的抗差状态估计方法(Robuststateestimation,RSE)包括加权最小绝对值估计(Weightedleastabsolutevalue,WLAV)、非二次准则法(QL、QC等)、以合格率最大为目标的状态估计(Maximumnormalmeasurementrate,MNMR)以及指数型目标函数状态估计(Maximumexponentialsquare,MES)等。但是这些抗差状态估计方法普遍存在计算效率不够高的特点，从而在一定程度上影响了它们在实际系统中的应用。

发明内容

为了有效提高抗差状态估计的计算效率和抗差性，保证全局最优解，本发明提出了一种混合整数二次规划形式的抗差状态估计方法，包括：

步骤1、将节点i的电压向量采用直角坐标形式表示，得到节点i的注入有功功率、注入无功功率；

步骤2、用等效的三个两绕组变压器来表示三绕组变压器，则将电力系统网络中所有支路用统一的π型支路表示，π型支路的上方左右两端分别为节点i到节点j，π型支路的下方左右两端分别接地；

步骤3、根据π型支路的等值电路图，得到节点i到节点j的支路有功功率和无功功率以及节点i的电压幅值量测方程，并将电压幅值量测方程转化为二次形式；

步骤4、基于节点i的注入有功功率、注入无功功率、电压幅值量测方程得到混合整数二次规划形式的抗差状态估计模型；

步骤5、利用WLS状态估计法得到MIQCP抗差状态估计模型的初值；

步骤6、对MIQCP抗差状态估计模型进行求解，找出并消除不良量测值；

步骤7、对修正后的量测值再次利用WLS状态估计法进行计算。

所述步骤1中节点i的电压向量表示为：

${\stackrel{\·}{v}}_i=e_i+{\mathrm{jf}}_i---\left(1\right)$

e_i和f_i分别表示的实部和虚部；

节点i的注入有功功率和注入无功功率分别表示为：

$P_i=e_i{\Σ}_{j=1}^n(G_{ij}{e}_j-B_{ij}{f}_j)+f_i{\Σ}_{j=1}^n(G_{ij}{f}_j+B_{ij}{e}_j)---\left(2\right)$

$Q_i=f_i{\Σ}_{j=1}^n\left(G_{ij}{e}_j-B_{ij}{f}_j\right)-e_i{\Σ}_{j=1}^n\left(G_{ij}{f}_j+B_{ij}{e}_j\right)---\left(3\right)$

节点导纳矩阵中的对应元素G_ij为节点i到节点j的电导，B_ij为节点i到节点j的电纳；n为电力系统网络中所有支路的数量；e_j和f_j分别表示为节点j的电压向量的实部和虚部。

所述π型支路的串联导纳r_ij+jx_ij为串联阻抗值；b_c为支路的接地电纳，对于变压器支路，b_c＝0；k为理想变压器的变比，对于普通支路，k＝1；

所述π型支路的等值电路图中节点i到节点j的串联导纳y_ij＝g_ij+jb_ij，节点i的接地导纳为g_si+jb_si，节点j的接地导纳为g_sj+jb_sj；节点i到节点j的电导g_ij＝g_s/k，节点i到节点j的电纳b_ij＝b_s/k，节点i到地的电导g_si＝(1-k)g_s/k²，节点i到地的电纳b_si＝(1-k)b_s/k²+b_c/2，节点j到地的电导g_sj＝(k-1)g_s/k，节点j到地的电导b_sj＝(k-1)b_s/k+b_c/2。

所述节点i到节点j的支路有功功率和无功功率分别表示为：

P_ij＝(e_i ²+f_i ²)(g_ij+g_si)-(e_ie_j+f_if_j)g_ij+(e_if_j-e_jf_i)b_ij(4)

Q_ij＝-(e_i ²+f_i ²)(b_ij+b_si)+(e_ie_j+f_if_j)b_ij+(e_if_j-e_jf_i)g_ij(5)

节点i的电压幅值量测方程表示为

$v_i=\sqrt{{e_i}^2+{f_i}^2}---\left(6\right)$

将式(6)转化为二次形式得

v_i ²＝e_i ²+f_i ²(7)

式中：v_i为节点i的电压幅值。

所述抗差状态估计模型为

$Min.J(x,b)={\Σ}_{i=1}^m{b}_i---\left(8\right)$

$s.t.h_i\left(x\right)\≤z_i+t_i^{+}+{\mathrm{Mb}}_i---\left(9\right)$

$h_i\left(x\right)\≥z_i-t_i^{-}-{\mathrm{Mb}}_i,i=1,2,...,m---\left(10\right)$

设z∈R^m为量测矢量，包括注入有功功率量测P_i和注入无功功率量测Q_i、支路有功功率量测P_ij和支路无功功率量测Q_ij，以及电压幅值量测的平方v_i ²，z的第i维为z_i，m为量测量的总个数；x∈Rⁿ为直角坐标形式的状态矢量，其第i维为h∈R^m为二次量测方程，h的第i维为h_i(x)，h_i(x)即为上述P_i、Q_i、P_ij、Q_ij或v_i ²其中之一的表达式；和分别为第i个量测值的上不确定度和下不确定度，可取其对应量测量z_i的1/100；M为足够大的正数，可取M＝1000；b_i为第i个量测值对应的0-1变量，对于不良量测值，b_i＝1，否则b_i＝0，b＝[b₁b₂…b_m]^T。

所述步骤6采用LocalSolver软件进行处理。

本发明的有益效果在于：通过采用直角坐标形式的MIQCP抗差状态估计方法在估计过程中能够保证全局最优解，并具有很高的计算效率和良好的抗差性，非常适宜于实际工程应用。

附图说明

图1为本发明的抗差状态估计方法流程图。

图2为网络中所有支路(包括普通线路和变压器支路)的统一π型支路；

图3为网络支路的二端口π形等值电路图。

具体实施方式

下面结合附图，对实施例作详细说明。

本发明提出了一种混合整数二次规划形式的抗差状态估计方法，如图1所示，包括：

步骤1、将节点i的电压向量采用直角坐标形式表示，得到节点i的注入有功功率、注入无功功率；

步骤2、如图2所示，用等效的三个两绕组变压器来表示三绕组变压器，则将电力系统网络中所有支路用统一的π型支路表示，π型支路的上方左右两端分别为节点i到节点j，π型支路的下方左右两端分别接地；

步骤3、如图3所示，根据π型支路的等值电路图，得到节点i到节点j的支路有功功率和无功功率以及节点i的电压幅值量测方程，并将电压幅值量测方程转化为二次形式；

步骤4、基于节点i的注入有功功率、注入无功功率、电压幅值量测方程得到混合整数二次规划形式的抗差状态估计模型；

步骤5、利用WLS状态估计法得到MIQCP抗差状态估计模型的初值；

步骤6、利用LocalSolver对MIQCP抗差状态估计模型进行求解，找出并消除不良量测值；

步骤7、对修正后的量测值再次利用WLS状态估计法进行计算。

所述步骤1中节点i的电压向量表示为：

${\stackrel{\·}{v}}_i=e_i+{\mathrm{jf}}_i---\left(1\right)$

e_i和f_i分别表示的实部和虚部；

节点i的注入有功功率和注入无功功率分别表示为：

$P_i=e_i{\Σ}_{j=1}^n(G_{ij}{e}_j-B_{ij}{f}_j)+f_i{\Σ}_{j=1}^n(G_{ij}{f}_j+B_{ij}{e}_j)---\left(2\right)$

$Q_i=f_i{\Σ}_{j=1}^n(G_{ij}{e}_j-B_{ij}{f}_j)-e_i{\Σ}_{j=1}^n(G_{ij}{f}_j+B_{ij}{e}_j)---\left(3\right)$

节点导纳矩阵中的对应元素G_ij为节点i到节点j的电导，B_ij为节点i到节点j的电纳；n为电力系统网络中所有支路的数量；e_j和f_j分别表示为节点j的电压向量的实部和虚部。

所述π型支路的串联导纳r_ij+jx_ij为串联阻抗值；b_c为支路的接地电纳，对于变压器支路，b_c＝0；k为理想变压器的变比，对于普通支路，k＝1。

所述π型支路的等值电路图中节点i到节点j的串联导纳y_ij＝g_ij+jb_ij，节点i的接地导纳为g_si+jb_si，节点j的接地导纳为g_sj+jb_sj；节点i到节点j的电导g_ij＝g_s/k，节点i到节点j的电纳b_ij＝b_s/k，节点i到地的电导g_si＝(1-k)g_s/k²，节点i到地的电纳b_si＝(1-k)b_s/k²+b_c/2，节点j到地的电导g_sj＝(k-1)g_s/k，节点j到地的电导b_sj＝(k-1)b_s/k+b_c/2。

所述节点i到节点j的支路有功功率和无功功率分别表示为：

P_ij＝(e_i ²+f_i ²)(g_ij+g_si)-(e_ie_j+f_if_j)g_ij+(e_if_j-e_jf_i)b_ij(4)

Q_ij＝-(e_i ²+f_i ²)(b_ij+b_si)+(e_ie_j+f_if_j)b_ij+(e_if_j-e_jf_i)g_ij(5)

节点i的电压幅值量测方程表示为

$v_i=\sqrt{{e_i}^2+{f_i}^2}---\left(6\right)$

将式(6)转化为二次形式得

v_i ²＝e_i ²+f_i ²(7)

式中：v_i为节点i的电压幅值。

所述抗差状态估计模型为

$Min.J(x,b)={\Σ}_{i=1}^m{b}_i---\left(8\right)$

$s.t.h_i\left(x\right)\≤z_i+t_i^{+}+{\mathrm{Mb}}_i---\left(9\right)$

$h_i\left(x\right)\≥z_i-t_i^{-}-{\mathrm{Mb}}_i,i=1,2,...,m---\left(10\right)$

设z∈R^m为量测矢量，包括注入有功功率量测P_i和注入无功功率量测Q_i、支路有功功率量测P_ij和支路无功功率量测Q_ij，以及电压幅值量测的平方v_i ²，z的第i维为z_i，m为量测量的总个数；x∈Rⁿ为直角坐标形式的状态矢量，其第i维为h∈R^m为二次量测方程，h的第i维为h_i(x)，h_i(x)即为上述P_i、Q_i、P_ij、Q_ij或v_i ²其中之一的表达式；为第i个量测值的上/下不确定度，可取其对应量测量z_i的1/100；M为足够大的正数，可取M＝1000；b_i为第i个量测值对应的0-1变量，对于不良量测值，b_i＝1，否则b_i＝0，b＝[b₁b₂…b_m]^T。

注意到MIQCP抗差状态估计模型中的二次约束为非半正定二次约束，无法利用CPLEX进行求解，宜用LocalSolver进行求解。而LocalSolver对浮点型变量的初值要求较高，故先利用WLS状态估计法计算得MIQCP抗差状态估计模型的初值。

设定利用IEEE标准系统检验混合整数二次规划形式(MIQCP)的抗差状态估计的性能。试验采用全量测，量测值通过在潮流计算的结果上叠加白噪声(均值为0，标准差为τ)来获得。对于电压量测，取τ＝0.001p.u.；对于功率量测，取τ＝0.02。测试环境为PC机，CPU为Intel(R)Core(TM)i3M370、主频为2.40GHz、内存2.00GB。

1、抗差性能测试

将本发明的MIQCP与其他状态估计器进行比较，来测试MIQCP的抗差性。

在IEEE-300系统上，将支路1-5的电抗缩小为原来的1/10，从而创造了一个杠杆点，同时设置4个一致性不良数据(P_1-2、Q_1-2、P₁、Q₁)。所设置的不良量测值以及量测量的正确值如表1所示。

表1MIQCP对IEEE300系统一致性不良数据的辨识

作为对比，首先用广为应用的WLS进行估计，并用LNR进行不良数据的辨识(简记为WLS+LNR)。首次辨识的结果为：10个量测量的标准化残差大于门槛值(3.0)，这10个量测量被认为是可疑数据；其中标准化残差最大的量测量为P_2-1，删去该量测后重新运行WLS；此时发现P₂的标准化残差最大。以上过程循环4次，4个良好的量测量被LNR误认为是可疑数据而被删去，但真正的不良数据仍然存在。可见，WLS+LNR不能辨识一致性不良数据。

应用MIQCP方法的估计结果如表1所示。可以发现，即使量测量中存在一致性不良数据，MIQCP的估计值与真值也可很好地吻合。在IEEE其他系统的多次试验也表明MIQCP在估计的过程中可以自动抑制不良数据，具有良好的抗差性。

2、计算效率测试

在正常量测条件下分别对四种状态估计器WLS、WLAV、MNMR以及MIQCP进行了测试，其中后三种属于抗差状态估计器。在试验中，WLS采用牛顿法求解，其他三种状态估计采用内点法求解；且MNMR采用两阶段法，即第一阶段进行WLS估计，第二阶段将WLS的估计值作为MNMR估计的初值进行计算。

共进行50次仿真试验，状态估计收敛时的迭代次数以及平均计算耗时如表2所示。由表2可见，在这四种状态估计器中，WLS的计算效率最高；而在后三种抗差状态估计器中，MIQCP的计算效率最高；而且随着系统规模的增大，MIQCP的迭代次数以及计算耗时增长的很缓慢，因而MIQCP适用于实际的大规模系统的估计。

表2四种状态估计器的迭代次数以及计算耗时

此实施例仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。

Claims (7)

1.一种混合整数二次规划形式的抗差状态估计方法，其特征在于，包括：

步骤1、将节点i的电压向量采用直角坐标形式表示，得到节点i的注入有功功率、注入无功功率；

步骤2、用等效的三个两绕组变压器来表示三绕组变压器，则将电力系统网络中所有支路用统一的π型支路表示，π型支路的上方左右两端分别为节点i到节点j，π型支路的下方左右两端分别接地；

步骤3、根据π型支路的等值电路图，得到节点i到节点j的支路有功功率和无功功率以及节点i的电压幅值量测方程，并将电压幅值量测方程转化为二次形式；

步骤4、基于节点i的注入有功功率、注入无功功率、电压幅值量测方程得到混合整数二次规划形式的抗差状态估计模型；

步骤5、利用WLS状态估计法得到MIQCP抗差状态估计模型的初值；

步骤6、对MIQCP抗差状态估计模型进行求解，找出并消除不良量测值；

步骤7、对修正后的量测值再次利用WLS状态估计法进行计算。

2.根据权利要求1所述方法，其特征在于，所述步骤1中节点i的电压向量表示为：

${\stackrel{\·}{v}}_i=e_i+{\mathrm{jf}}_i---\left(1\right)$

e_i和f_i分别表示的实部和虚部；

节点i的注入有功功率和注入无功功率分别表示为：

$P_i=e_i{\Σ}_{j=1}^n(G_{ij}{e}_j-B_{ij}{f}_j)+f_i{\Σ}_{j=1}^n(G_{ij}{f}_j+B_{ij}{e}_j)---\left(2\right)$

$Q_i=f_i{\Σ}_{j=1}^n(G_{ij}{e}_j-B_{ij}{f}_i)-e_i{\Σ}_{j=1}^n(G_{ij}{f}_j+B_{ij}{e}_j)---\left(3\right)$

节点导纳矩阵中的对应元素G_ij为节点i到节点j的电导，B_ij为节点i到节点j的电纳；n为电力系统网络中所有支路的数量；e_j和f_j分别表示为节点j的电压向量的实部和虚部。

3.根据权利要求1所述方法，其特征在于，所述π型支路的串联导纳r_ij+jx_ij为串联阻抗值；b_c为支路的接地电纳，对于变压器支路，b_c＝0；k为理想变压器的变比，对于普通支路，k＝1；

所述π型支路的等值电路图中节点i到节点j的串联导纳y_ij＝g_ij+jb_ij，节点i的接地导纳为g_si+jb_si，节点j的接地导纳为g_sj+jb_sj；节点i到节点j的电导g_ij＝g_s/k，节点i到节点j的电纳b_ij＝b_s/k，节点i到地的电导g_si＝(1-k)g_s/k²，节点i到地的电纳b_si＝(1-k)b_s/k²+b_c/2，节点j到地的电导g_sj＝(k-1)g_s/k，节点j到地的电导b_sj＝(k-1)b_s/k+b_c/2。

4.根据权利要求3所述方法，其特征在于，所述节点i到节点j的支路有功功率和无功功率分别表示为：

P_ij＝(e_i ²+f_i ²)(g_ij+g_si)-(e_ie_j+f_if_j)g_ij+(e_if_j-e_jf_i)b_ij(4)

Q_ij＝-(e_i ²+f_i ²)(b_ij+b_si)+(e_ie_j+f_if_j)b_ij+(e_if_j-e_jf_i)g_ij(5)

节点i的电压幅值量测方程表示为

$v_i=\sqrt{{e_i}^2+{f_i}^2}---\left(6\right)$

将式(6)转化为二次形式得

v_i ²＝e_i ²+f_i ²(7)

式中：v_i为节点i的电压幅值。

5.根据权利要求1所述方法，其特征在于，所述抗差状态估计模型为

$Min.J(x,b)={\Σ}_{i=1}^m{b}_i---\left(8\right)$

$\begin{array}{cc}s.t.& h_i\left(x\right)\≤z_i+t_i^{+}+{\mathrm{Mb}}_i\end{array}---\left(9\right)$

$h_i\left(x\right)\≥z_i-t_i^{-}-{\mathrm{Mb}}_i,i=1,2,\mathrm{...},m---\left(10\right)$

设z∈R^m为量测矢量，包括注入有功功率量测P_i和注入无功功率量测Q_i、支路有功功率量测P_ij和支路无功功率量测Q_ij，以及电压幅值量测的平方v_i ²，z的第i维为z_i，m为量测量的总个数；x∈Rⁿ为直角坐标形式的状态矢量，其第i维为e_i和f_i分别表示节点i的电压向量的实部和虚部；h∈R^m为二次量测方程，h的第i维为h_i(x)，h_i(x)即为上述P_i、Q_i、P_ij、Q_ij或v_i ²其中之一的表达式；和分别为第i个量测值的上不确定度和下不确定度；M为正数；b_i为第i个量测值对应的0-1变量，对于不良量测值，b_i＝1，否则b_i＝0，b＝[b₁b₂…b_m]^T。

6.根据权利要求5所述方法，其特征在于，所述和取其对应量测量z_i的1/100；M取值1000。

7.根据权利要求1所述方法，其特征在于，所述步骤6采用LocalSolver软件进行处理。