CN103248043A - Power system multi-zone distributed state estimation method based on synchronous phase angle measurement device - Google Patents

Abstract

The invention discloses a power system multi-zone distributed state estimation method based on a synchronous phase angle measurement device, which quickly estimates voltage amplitude values and phase angle values of all nodes of a power grid in real time by utilizing measurement data from a SCADA (supervisory control and data acquisition) system and a PMU (phasor measurement unit), and performs real-time interaction on state information of all zones. According to the method, a geographical characteristic-based non-overlapping decoupling strategy is adopted to decouple a large power system into subsystems, and states of all the subsystems are integrated by a coordination system; and with adoption of the multi-zone distributed state estimation method, estimation of states of the large system turns into estimation of states of a series of partial small zones, so that the computation speed is greatly improved. Real-time accurate measurement information of voltage, phase angle and the like can be offered to the system, and the system is enabled to obtain higher measurement redundancy, so that the accuracy of state estimation can be improved. According to the method, zone interaction management facilitates secure and effective access of new energy, so the development requirement of a future intelligent power grid can be met.

Description

A kind of electric power system multizone distributed state estimation method based on the synchronous phase angle measurement mechanism

Technical field

The invention belongs to power system operation and control technology field, be specifically related to a kind of electric power system multizone distributed state estimation method based on synchronous phase angle measurement mechanism (PMU).

Background technology

Power system state estimation is basis and the core of energy management system (EMS), its effect is that filtering is carried out in the real time information that data collection and supervisory control system (SCADA) provide, with the raising data precision, the interference of debug information, thus obtain a high-quality real-time data base.Along with the continuous expansion of interconnected network scale, the complexity that operation of power networks shows aspect three in space, time, control target increases day by day, and the scheduling controlling of the electrical network of traditional E MS lacks good coordination aspect these three.The world-shaking U.S.A that took place in 2003 adds and has a power failure on a large scale, exposed traditional E MS three-dimensional coordinated and aspect the real time data processing of system, have serious problems, and state estimation is as basis and the core of EMS, bearing the magnanimity metric data of system is being handled in real time, the real-time monitoring system running status, the important task of real time data is provided for more senior EMS software, therefore, Power system state estimation is badly in need of transforming and is developed to adapt to following intelligent grid real-time monitoring system, coordinates each regional requirement of electrical network.

At present, the metric data of Power system state estimation mainly comes from the SCADA system, generally comprise node injecting power, branch power and node voltage amplitude etc., these data transmitted once in per 2 seconds, but because the measurement information of these systems often all is to be sent to the control centre by telemechanical apparatus, and the error of each link makes the precision of quantity of states such as voltage that iterative comes out, phase angle be difficult to be guaranteed in the error of telemechanical apparatus and the transport process.In recent years, progressively be applied in the electric power system based on the PMU of GPS (Global Position System), it have the metric data of collection fast, can measure phase angle information, and metric data ratio of precision SCADA advantages of higher.

On the other hand, present method for estimating state, just at power transmission network, just at power distribution network, method for estimating state (multizone state estimation) research of each Regional Integration of electrical network being got up to carry out unified management is also fewer, and following intelligent grid will be generating, transmission of electricity, the electric power system of distribution and each link of electricity consumption is as object, by the novel power grid control technology of continuous research, and with every technology combination, the intelligent AC of realization from generating to all link information of electricity consumption, thereby systematically optimize electrical production, carry and use, and the multizone Power system state estimation can be realized this goal, it can integrate each regional real-time status information of electrical network, realizes each regional information interchange of electrical network, for the optimization scheduling of intelligent grid and the realization of self-healing function lay the foundation.Therefore, invent a kind of electric power system multizone distributed state estimation method that can effectively merge the PMU technology and have very important and practical meanings, and be conducive to EMS constantly perfect of following intelligent grid.

Summary of the invention

Technical problem to be solved by this invention is to overcome the deficiencies in the prior art, a kind of electric power system multizone distributed state estimation method based on PMU is provided, make this method can be used to the metric data voltage magnitude that estimates each node of electrical network and angle values real-time from SCADA system and PMU, and the disposal ability with bad data, the state information that this method can be regional with each is carried out real-time, interactive in addition, be suitable for following intelligent grid, have very important and practical meanings.

The present invention is for solving its technical problem, and the technology that adopts comprises following steps:

(1) electric power system multizone decoupling zero

Adopt the non-overlapping copies electric power system multizone decoupling zero strategy based on region characteristic among the present invention, particularly, suppose an electric power system, contain the wiring of n bar, and divide into r nonoverlapping zone (or subsystem) S according to ground _i, each S _iN is arranged _iThe bar wiring, and link to each other with connecting line, as shown in Figure 1.After the decoupling zero, all by the local control centre control of oneself, the state of oneself is responsible for estimating in this center, and is connected to the control centre of a coordination by communication line in each zone.

(2) each Local Area Network data reads

In this step, reading of each Local Area Network data comprises network configuration and line impedance, and forms node admittance matrix and branch road-node incidence matrices thus.

(3) system measurements and PMU configuration

The measurement of system comprises that voltage magnitude measures, gains merit and reactive power flow measures, meritorious and idle injection trend measures, and measurement equation is as follows:

● voltage magnitude: ${V_l}^{\mathrm{mea}}=V_l+e_{V_l}$

● meritorious and reactive power flow:

$P_{\mathrm{lm}}^{\mathrm{mea}}=P_{\mathrm{lm}}+e_{P_{\mathrm{lm}}}$

$Q_{\mathrm{lm}}^{\mathrm{mea}}=Q_{\mathrm{lm}}+e_{Q_{\mathrm{lm}}}$

● meritorious and idle injection trend:

${P_l}^{\mathrm{mea}}=({V_l}^2{g}_l+\underset{m\∈a\left(l\right)}{\mathrm{\Σ}}{P}_{\mathrm{lm}})+\underset{m\∈b\left(l\right)}{\mathrm{\Σ}}{P}_{\mathrm{lm}}+e_{P_l}$

${Q_l}^{\mathrm{mea}}=(-{V_l}^2{b}_l+\underset{m\∈a\left(l\right)}{\mathrm{\Σ}}{Q}_{\mathrm{lm}})+\underset{m\∈b\left(l\right)}{\mathrm{\Σ}}{Q}_{\mathrm{lm}}+e_{Q_l}$

Wherein, injection is meritorious in the interior zone is respectively with reactive power flow measurement functional equation:

$P_{\mathrm{lm}}={V_l}^2(g_{\mathrm{lm}}+g_{\mathrm{slm}})-V_l{V}_m(g_{\mathrm{lm}}\cos ({\mathrm{\θ}}_l-{\mathrm{\θ}}_m)+b_{\mathrm{lm}}\sin ({\mathrm{\θ}}_l-{\mathrm{\θ}}_m))$

$Q_{\mathrm{lm}}=-{V_l}^2(b_{\mathrm{lm}}+b_{\mathrm{slm}})-V_l{V}_m(g_{\mathrm{lm}}\sin ({\mathrm{\θ}}_l-{\mathrm{\θ}}_m)-b_{\mathrm{lm}}\cos ({\mathrm{\θ}}_l-{\mathrm{\θ}}_m));$

● be to belong to regional S _iA wiring;

With

It is respectively the measuring value of wiring l and m place voltage magnitude, meritorious and reactive power flow, meritorious and idle injection trend; V _l, V _mIt is respectively the voltage magnitude of wiring l and m; θ _lAnd θ _mIt is respectively the angle values of wiring l and m;

With

It is respectively the error in measurement of voltage magnitude, meritorious and reactive power flow, meritorious and idle injection trend; g _l+ jb _LmIt is branch road l-m series admittance value; g _Slm+ jb _SlmIt is branch road l-m shunt admittance value; g _l+ jb _lIt is the shunt admittance value that is connected to line l; A(l) for linking l upward and belonging to regional S _iWired set; B (l) is for linking l upward and belonging to regional S _iWired set (i ≠ j);

● the PMU equipping rules

Guarantee the observability of whole system, among the present invention to each subsystem all dispose at least a PMU with and PMU all be configured on the connecting line between subsystem and the subsystem.

(4) the distributed state estimation of multizone

On the decoupling zero of electric power system multizone and system measurements and PMU configure base, this step will adopt least square method to carry out state estimation respectively to all subregion, and multizone state estimation model is as follows:

z _i=h _i(x _i)+e _i,i=1,2,...r

z _c=h _c(x)+e _c

Wherein, z _iBe regional S _iM _i* 1 dimension internal node or the local vector that measures; z _cBe m _c* 1 dimension boundary number direction finding amount; $x_i=\left(\begin{array}{c}{\mathrm{\θ}}_i\\ V_i\end{array}\right)$ Be regional S _iMiddle 2n _iThe local state vector of * 1 dimension comprises n _iIndividual magnitude of voltage and n _iIndividual angle values;

It is system-wide state vector

h _i(), h _c() is the Nonlinear Vector measurement function based on kirchhoff voltage or current law; e _i, e _cBe the gaussian random error vector.

The state estimation that can get all subregion according to model is expressed as:

$z=\left[\begin{array}{c}{z}_i\\ z_c\\ z_p\end{array}\right]=\left[\begin{array}{c}{h}_i\left(x_i\right)\\ h_c\left(x\right)\\ h_p\left(x\right)\end{array}\right]+r$

Wherein, z _iAnd h _i(x _i) represent that respectively remote terminal unit measures the nonlinear function that vector sum measures about remote terminal unit; z _pExpression PMU measures, and comprises that voltage measures and angle measurements mutually; h _p() is the PMU measurement functional vector about system mode; $z_p=h_p\left(x\right)=\left[\begin{array}{c}{h}_{\mathrm{Vp}}\·x\\ h_{\mathrm{Ip}}\left(x\right)\end{array}\right]+r,$ This is because the introducing of PMU has made about the measurement function linearisation of voltage.

Using least square method to carry out state estimation has:

$\left[G\left(x^k\right)\right]\&CenterDot;{\mathrm{\&Delta;x}}^{k+1}=\left[\begin{array}{c}{H}_i\left(x_i^k\right)\\ H_c\left(x^k\right)\\ H_p\left(x^k\right)\end{array}\right]\&CenterDot;{\left[\begin{array}{ccc}{R}_i& 0& 0\\ 0& {\mathrm{<figure-callout\; id="0"\; label="R\; c"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">R</figure-callout>}}_c& 0\\ 0& 0& R_p\end{array}\right]}^{-1}\&CenterDot;\left[\begin{array}{c}{z}_i-h_i\left({x_i}^k\right)\\ z_c-h_c\left(x^k\right)\\ z_p-h_p\left(x^k\right)\end{array}\right]$

$G\left(x^k\right)={\left[\begin{array}{c}{H}_i\left(x_i^k\right)\\ H_c\left(x^k\right)\\ H_p\left(x^k\right)\end{array}\right]}^T\&CenterDot;{\left[\begin{array}{ccc}{R}_i& 0& 0\\ 0& {\mathrm{<figure-callout\; id="0"\; label="R\; c"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">R</figure-callout>}}_c& 0\\ 0& 0& R_p\end{array}\right]}^{-1}\&CenterDot;\left[\begin{array}{c}{H}_i\left(x_i^k\right)\\ H_c\left(x^k\right)\\ H_p\left(x^k\right)\end{array}\right]$

x ^k+1=x ^k+Δx ^k+1

Wherein k is iterations, can find the solution each regional state estimation according to iterating of above formula.

(5) PMU tuning algorithm

Adopt accurately measuring in real time all subregion state information of system being carried out linearity coordination of PMU among the present invention, coordinate scheme as shown in Figure 2, each detailed zone state information interaction and coordination principle figure such as accompanying drawing 3.

At this moment, because the utilization of PMU metric data, state estimation becomes linearity in the coordination level, and the linear condition estimation model is:

$z_p=\left[\begin{array}{c}{z}_{\mathrm{Vp}}\\ z_{\mathrm{Ip}}\end{array}\right]=\left[\begin{array}{c}I\\ Y\end{array}\right]\&CenterDot;x$

Wherein, z _VpAnd z _IpThe voltage and current that is respectively PMU measures vector; I represents the measurement unit matrix about voltage vector; Y represents the linear matrix about the current measurement vector, is made up of node and branch road admittance; B is the transmission matrix of being set up by I and Y.

Following formula is rewritten into: z=Bx, then use method for estimating state to get

In order to improve the redundancy of system, the present invention measures the state estimation result of each sub regions imports in the coherent system as puppet, so state estimation need be modified to $z_p=\left[\begin{array}{c}{z}^{\&prime;}\\ z_{\mathrm{ps}}\end{array}\right]=\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]\&CenterDot;x,$ And error in measurement matrix R _pNeed be modified to $R^{\&prime;}=\left[\begin{array}{cc}{\mathrm{<figure-callout\; id="0"\; label="R\; p"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">R</figure-callout>}}_p& 0\\ 0& R_{\mathrm{ps}}\end{array}\right],$ Subscript ps represents to measure with puppet a series of vectors of vector correlation; Easily know transition matrix B and I _PsAll be linear, so state estimation do not need iteration just can obtain the result, thereby saved the time of state estimation greatly, and have:

$x=\left(\begin{array}{}\end{array}\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]\&CenterDot;{\left(R^{\&prime;}\right)}^{-1}\&CenterDot;{\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]}^{-1}\&CenterDot;{\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]}^T\&CenterDot;{\left(R^{\&prime;}\right)}^{-1}\&CenterDot;\left[\begin{array}{c}{z}_p\\ I_{\mathrm{ps}}\end{array}\right]\right)$

Thus, the state of each sub regions and borderline region is all estimated.

(6) distributed bad data is handled

This step will be on the basis of Distributed Calculation residual error covariance matrix, can be divided into the processing of subregion bad data and coherent system bad data two stages of processing to the processing of bad data, the present invention adopts normalization residual test method to detect bad data in addition.Concrete enforcement is as follows:

After state estimation, exist residual error to be expressed as between measuring value and the estimated value:

$\left(\begin{array}{c}{r_1}^N\\ M\\ {r_r}^N\\ {r_c}^N\end{array}\right)={\left(\mathrm{diagP}\right)}^{-1/2}\left(\begin{array}{c}{r}_1\\ M\\ r_r\\ r_c\end{array}\right)$

Wherein, r _i=z _i-h _i(x _i), i=1 ..., r;

r _c=z _c-h _c(x)； $P=\mathrm{cov}\left(\begin{array}{c}{r}_1\\ M\\ r_r\\ r_c\end{array}\right)=\left(\begin{array}{cc}{R}_a& 0\\ 0& R_c\end{array}\right)-\left(\begin{array}{cc}{H}_a& 0\\ 0& H_c\end{array}\right)A_1{\left(\begin{array}{cc}{H}_a& 0\\ 0& H_c\end{array}\right)}^T$

And $R_a=\left(\begin{array}{ccc}{\mathrm{<figure-callout\; id="1"\; label="R"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">R</figure-callout>}}_1& 0& 0\\ 0& O& 0\\ 0& 0& R_r\end{array}\right),$ $H_a=\left(\begin{array}{ccc}{H}_1& 0& 0\\ 0& O& 0\\ 0& 0& H_r\end{array}\right),$ $G_a=\left(\begin{array}{ccc}{\mathrm{<figure-callout\; id="1"\; label="G"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">G</figure-callout>}}_1& 0& 0\\ 0& \mathrm{<figure-callout\; id="0"\; label="O"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">O</figure-callout>}& 0\\ 0& 0& G_r\end{array}\right)$

$\left(\begin{array}{cc}{A}_1& A_2^T\\ A_2& -A_3\end{array}\right)={\left(\begin{array}{cc}{G}_a& H_c^T\\ H_c& -R_c\end{array}\right)}^{-1},$ $A_1={G_a}^{-1}-{G_a}^{-1}{H}_c^T{G_c}^{-1}{H}_c{G_a}^{-1}$

Can be got by above formula:

I=1 ..., r,

Wherein have $P_i=\mathrm{cov}\left(r_i\right)=R_i-H_i{G}_i^{-1}{H}_i^T+H_i{G}_i^{-1}\left(H_{\mathrm{ci}}^T{G}_c^{-1}{H}_{\mathrm{ci}}\right)G_i^{-1}{H}_i^T,$ $P_c=\mathrm{cov}\left(r_c\right)=R_c{G}_c^{-1}{R}_c.$ Compute matrix P _iAnd P _cThe effective method of diagonal entry is the compute sparse inverse matrix

With And carry out necessity operation that is defined in two formulas, so the coherent system bad data is handled and subregion bad data two stages of processing are as follows:

1) by coherent system compute sparse inverse matrix

Finish the processing of coherent system bad data and with result transmission in the state estimator of subregion;

2) receiving the sparse inverse matrix that the coherent system biography is come

After carry out subregion state estimation and compute sparse inverse matrix

And then finish the subregion bad data and handle.

(7) condition of convergence is judged.

If all state estimation all restrain, the bad data processing finishes, and then exports the voltage magnitude of each node and the identification result of angle values and output bad data, otherwise changes step 4.

The invention provides a kind of electric power system multizone distributed state estimation method based on PMU, have following advantage than the method for estimating state of present existence:

(1) has the ability of handling interconnected large power system.The present invention adopts the non-overlapping copies decoupling zero strategy based on region characteristic, be several subsystems with a large power system decoupling zero, and integrate the state of each subsystem cleverly by a coherent system, and the utilization of distributed method makes us can obtain the real-time status information of system fast and accurately.

(2) state estimation speed is fast, the precision height.The present invention adopts the multizone distributed state estimation method, and big system state estimation is become a series of local zonules state estimation, has improved computational speed greatly; Reasonably adopt PMU in addition in coherent system, making does not need the state estimation linearisation of coherent system side to iterate and finds the solution, and then improved state estimation speed again as traditional least square state estimation; At last, the introducing of PMU can provide in real time measurement information such as voltage, phase angle accurately for system, and can the assurance system obtain higher measurement redundancy, and then improves the precision of state estimation.

(3) has the disposal ability of distributed bad data.By distributed mode bad data is handled among the present invention, improved the efficient that bad data is handled greatly.

(4) application prospect is good.The present invention carries out multizone with electric power system and handles, and by coherent system in time and other system carry out the real-time, interactive of data message; In addition, regional interactive maintenance is more conducive to new forms of energy and inserts safely and effectively, satisfies following intelligent grid demand for development.

Description of drawings

Fig. 1 is the non-overlapping copies multizone decoupling zero figure that the present invention is based on region characteristic;

Fig. 2 is the linear coordination principle figure of PMU of the present invention;

Fig. 3 is the mutual and coordination principle figure of each zone state information of the present invention;

Fig. 4 is flow chart of the present invention;

Fig. 5 is embodiment of the invention IEEE14 resolution chart;

Fig. 6 is each node phase angle test result;

Fig. 7 is each node voltage amplitude test result;

Embodiment

Below in conjunction with accompanying drawing technical scheme of the present invention is described in detail.

In order to find the solution based on the distributed state estimation problem of the electric power system multizone of PMU, shown in Figure 4 according to flow process, need take following steps:

Step 1: based on the non-overlapping copies multizone electric power system decoupling zero of region characteristic, as shown in Figure 5, IEEE14 node test system is S by decoupling zero ₁, S ₂, S ₃And S ₄, S wherein ₁ Comprise 1,2,5 nodes; S ₂Comprise 3,4,7,8 nodes; S ₃Comprise 9,10,14 nodes; S ₄Comprise 6,11,12,13 nodes.

Step 2: read each regional network data, comprise network configuration and line impedance, and form node admittance matrix and branch road-node incidence matrices.

Measurement configuration and the PMU of step 3:IEEE14 system install.Each regional internal node measures: S ₁{ P _1-2, P _1-5, P _2-5, P ₁, S ₂{ P _3-4, P _4-7, P _7-8, S ₃{ P _6-11, P _6-12, P _6-13, P _12-13, P ₁₂, S ₄{ P _9-10, P _9-14; The border measures and is configured to: { P ₅, P _4-5, P _4-9, P ₃, P _13-14, P ₁₃, P _10-11, P ₁₄; Have 4 PMUs and be installed in node 2,6 respectively, on 7,9, wherein the PMU on 2 nodes can measure 2-4, the voltage of the phase angle of 2-3 and node 2; PMU on 6 nodes measures 6-11,6-12, the magnitude of voltage of the phase angle of 6-13 and node 6; PMU on node 7 measures 7-9,7-4, the magnitude of voltage of the phase angle of 7-8 and node 7; PMU on node 9 measures 9-10, the magnitude of voltage of the phase angle of 9-14 and node 9; All angle values all represent with the number of degrees, and all voltage measured values all are that to have added an average in real calculation of tidal current be zero, and variance is 0.01 gaussian random error; The rough error of σ and-20 σ is added to P respectively in addition+20 _1-5And P ₁₃On.

Step 4: subregion state estimation.Each zone there is following measuring value:

${\mathrm{<figure-callout\; id="1"\; label="z"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">z</figure-callout>}}_1=\left(\begin{array}{c}0.0000\\ 0.0869\\ 0.3532\\ 0.0663\\ 0.2402\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_1\\ P_{1-2}\\ P_{1-5}\\ P_{2-5}\\ {\mathrm{<figure-callout\; id="1"\; label="P"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">P</figure-callout>}}_1\end{array},$ ${\mathrm{<figure-callout\; id="2"\; label="z"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">z</figure-callout>}}_2=\left(\begin{array}{c}0.0000\\ -0.0243\\ 0.356\\ -0.0002\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_3\\ P_{3-4}\\ P_{4-7}\\ P_{7-8}\end{array},$ $z_3=\left(\begin{array}{c}0.0000\\ 0.0099\\ 0.0148\\ 0.0164\\ 0.0016\\ -0.0133\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_6\\ P_{6-11}\\ P_{6-12}\\ P_{6-13}\\ P_{12-13}\\ {\mathrm{<figure-callout\; id="12"\; label="P"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">P</figure-callout>}}_{12}\end{array},$ ${\mathrm{<figure-callout\; id="4"\; label="z"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">z</figure-callout>}}_4=\left(\begin{array}{c}0.0000\\ 0.0028\\ 0.0192\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_9\\ P_{9-10}\\ P_{9-14}\end{array}$

Wherein, the measuring value that is connected with PMU all is actual value, and other measuring value to contain average be zero, variance is 0.01 gaussian random error.Thereby Jacobian matrix is arranged:

${\mathrm{<figure-callout\; id="1"\; label="H"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">H</figure-callout>}}_1=\left(\begin{array}{ccc}\frac{{\mathrm{\&delta;}}_1}{1}& \frac{{\mathrm{\&delta;}}_2}{0}& \frac{{\mathrm{\&delta;}}_5}{0}\\ 1& -1& 0\\ 1& 0& -1\\ 0& 1& -1\\ 2& -1& -1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_1\\ P_{1-2}\\ P_{1-5}\\ P_{2-5}\\ P_1\end{array}$ ${\mathrm{<figure-callout\; id="2"\; label="H"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_2=\left(\begin{array}{cccc}\frac{{\mathrm{\&delta;}}_3}{1}& \frac{{\mathrm{\&delta;}}_4}{0}& \frac{{\mathrm{\&delta;}}_7}{0}& \frac{{\mathrm{\&delta;}}_8}{0}\\ 1& -1& 0& 0\\ 0& 1& -1& 0\\ 0& 0& 1& -1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_3\\ P_{3-4}\\ P_{4-7}\\ P_{7-8}\end{array}$ $H_3=\left(\begin{array}{cccc}\frac{{\mathrm{\&delta;}}_6}{1}& \frac{{\mathrm{\&delta;}}_{11}}{0}& \frac{{\mathrm{\&delta;}}_{12}}{0}& \frac{{\mathrm{\&delta;}}_{13}}{0}\\ 1& -1& 0& 0\\ 1& 0& -1& 0\\ 1& 0& 0& -1\\ 0& 0& 1& -1\\ -1& 0& 2& -1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_6\\ P_{6-11}\\ P_{6-12}\\ P_{6-13}\\ P_{12-13}\\ {\mathrm{<figure-callout\; id="12"\; label="P"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">P</figure-callout>}}_{12}\end{array}$

${\mathrm{<figure-callout\; id="4"\; label="H"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">H</figure-callout>}}_4=\left(\begin{array}{ccc}\frac{{\mathrm{\&delta;}}_9}{1}& \frac{{\mathrm{\&delta;}}_{10}}{0}& \frac{{\mathrm{\&delta;}}_{14}}{0}\\ 1& -1& 0\\ 1& 0& -1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_9\\ P_{9-10}\\ P_{9-14}\end{array}$

Wherein, δ ₁, δ ₃, δ ₆And δ ₉Be respectively regional S ₁, S ₂, S ₃And S ₄The middle puppet measurement amount that adds; The diagonal entry of local covariance matrix is (1/ σ ²)=(1/0.01 ²), thereby gain matrix has:

${\mathrm{<figure-callout\; id="1"\; label="G"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">G</figure-callout>}}_1={10}^4\&times;\left(\begin{array}{ccc}7& -3& -3\\ -3& 3& 0\\ -3& 0& 3\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_1\\ {\mathrm{\&delta;}}_2\\ {\mathrm{\&delta;}}_5\end{array}$ ${\mathrm{<figure-callout\; id="2"\; label="G"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">G</figure-callout>}}_2={10}^4\&times;\left(\begin{array}{cccc}2& -1& 0& 0\\ -1& 2& -1& 0\\ 0& -1& 2& -1\\ 0& 0& -1& 1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_3\\ {\mathrm{\&delta;}}_4\\ {\mathrm{\&delta;}}_7\\ {\mathrm{\&delta;}}_8\end{array}$ $G_3={10}^4\&times;\left(\begin{array}{cccc}5& -1& -3& 0\\ -1& 1& 0& 0\\ -3& 0& 6& -3\\ 0& 0& -3& 3\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_6\\ {\mathrm{\&delta;}}_{11}\\ {\mathrm{\&delta;}}_{12}\\ {\mathrm{\&delta;}}_{13}\end{array}$ ${\mathrm{<figure-callout\; id="4"\; label="G"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">G</figure-callout>}}_4={10}^4\&times;\left(\begin{array}{ccc}3& -1& -1\\ 1& 1& 0\\ 1& 0& 1\end{array}\right)\begin{array}{c}{\mathrm{\&delta;}}_9\\ {\mathrm{\&delta;}}_{10}\\ {\mathrm{\&delta;}}_{14}\end{array}$

Therefore, the basis in each zone

The angle values that calculates the one's respective area is as follows:

${\mathrm{\&theta;}}_1=\left(\begin{array}{c}0.00\\ -4.98\\ -12.60\end{array}\right)\begin{array}{c}\frac{\mathrm{<figure-callout\; id="1"\; label="bus"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">bus</figure-callout>}}{1}\\ 2\\ 5\end{array},$ ${\mathrm{\&theta;}}_2=\left(\begin{array}{c}0.00\\ 1.39\\ -0.65\\ -0.64\end{array}\right)\begin{array}{c}\frac{\mathrm{bus}}{3}\\ 4\\ 7\\ 8\end{array},$ ${\mathrm{\&theta;}}_3=\left(\begin{array}{c}0.00\\ -0.57\\ -0.85\\ -0.94\end{array}\right)\begin{array}{c}\frac{\mathrm{<figure-callout\; id="6"\; label="bus"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">bus</figure-callout>}}{6}\\ 11\\ 12\\ 13\end{array},$ ${\mathrm{\&theta;}}_4=\left(\begin{array}{c}0.00\\ -0.16\\ -1.10\end{array}\right)\begin{array}{c}\frac{\mathrm{bus}}{9}\\ 10\\ 14\end{array}$

Step 5: coherent system state estimation and data interaction.At this moment, because PMU measures each zone is carried out, so correlated state estimates to become linearity in the coordination level, and the linear condition estimation model is:

$z_p=\left[\begin{array}{c}{z}_{\mathrm{Vp}}\\ z_{\mathrm{Ip}}\end{array}\right]=\left[\begin{array}{c}I\\ Y\end{array}\right]\&CenterDot;x$

Wherein, z _VpAnd z _IpThe voltage and current that is respectively PMUs measures vector, and wherein the PMU on 2 nodes can measure 2-4, the voltage of the phase angle of 2-3 and node 2; PMU on 6 nodes measures 6-11,6-12, the magnitude of voltage of the phase angle of 6-13 and node 6; PMU on node 7 measures 7-9,7-4, the magnitude of voltage of the phase angle of 7-8 and node 7; PMU on node 9 measures 9-10, the magnitude of voltage of the phase angle of 9-14 and node 9; I represents the measurement unit matrix about voltage vector; Y represents the linear matrix about the current measurement vector, is made up of node and branch road admittance; B is the transmission matrix of being set up by I and Y.

Because the adding of pseudo-measurement amount, so state estimation need be modified to $z_p=\left[\begin{array}{c}{z}^{\&prime;}\\ z_{\mathrm{ps}}\end{array}\right]=\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]\&CenterDot;x,$ And error in measurement matrix R _pNeed be modified to $R^{\&prime;}=\left[\begin{array}{cc}{\mathrm{<figure-callout\; id="0"\; label="R\; p"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">R</figure-callout>}}_p& 0\\ 0& R_{\mathrm{ps}}\end{array}\right],$ Use least square method to get:

$x=\left(\begin{array}{}\end{array}\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]\&CenterDot;{\left(R^{\&prime;}\right)}^{-1}\&CenterDot;{\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]}^{-1}\&CenterDot;{\left[\begin{array}{c}B\\ I_{\mathrm{ps}}\end{array}\right]}^T\&CenterDot;{\left(R^{\&prime;}\right)}^{-1}\&CenterDot;\left[\begin{array}{c}{z}_p\\ I_{\mathrm{ps}}\end{array}\right]\right)$

After coordinating, have:

${\mathrm{\&theta;}}_1=\left(\begin{array}{c}0.00\\ -5.54\\ -11.60\end{array}\right)\begin{array}{c}\frac{\mathrm{<figure-callout\; id="1"\; label="bus"\; filenames="2013101574437100003DEST\_PATH\_IMAGE001.png,HDA00003130546000012.png"\; state="\{\{state\}\}">bus</figure-callout>}}{1}\\ 2\\ 5\end{array},$ ${\mathrm{\&theta;}}_2=\left(\begin{array}{c}-12.20\\ -15.39\\ -14.35\\ -14.64\end{array}\right)\begin{array}{c}\frac{\mathrm{bus}}{3}\\ 4\\ 7\\ 8\end{array},$ ${\mathrm{\&theta;}}_3=\left(\begin{array}{c}-14.00\\ -14.27\\ -14.45\\ -14.84\end{array}\right)\begin{array}{c}\frac{\mathrm{<figure-callout\; id="6"\; label="bus"\; filenames="HDA00003130546000032.png"\; state="\{\{state\}\}">bus</figure-callout>}}{6}\\ 11\\ 12\\ 13\end{array},$ ${\mathrm{\&theta;}}_4=\left(\begin{array}{c}-14.90\\ -14.76\\ -15.98\end{array}\right)\begin{array}{c}\frac{\mathrm{bus}}{9}\\ 10\\ 14\end{array}$

Step 6: distributed bad data is handled.By I=1 ..., r,

Calculate the normalization maximum residul difference on each intra-zone and border respectively, for the zone

To the zone $2(P_{3-4},{r_i}^N=-3.02);$ To the zone $3(P_{6-13},{r_i}^N=-6.35);$ To the zone $4(P_{9-14},{r_i}^N=-4.20);$ To coherent system

Therefore regional Maximum normalization residual error is arranged, so be recognized as bad data, after removing this bad data, maximum normalization residual error is in the zone

And coherent system

Yi Zhi measures P ₁₃Be that bad data and coordinated system remove; After removing these two bad datas, system is proceeded state estimation, do not find extra bad data.

Step 7: the condition of convergence is judged.If all state estimation all restrain, the bad data processing finishes, and then exports the voltage magnitude of each node and the identification result of angle values and output bad data, otherwise changes step 4.

Each final node phase angle and test result such as Fig. 6 of voltage magnitude, shown in 7, as can be seen from the figure, method of the present invention is than traditional least square method state estimation, and estimated accuracy is higher; The whole state estimation procedure of method of the present invention has only been spent 0.2496s in addition, and estimating speed is fast.It is pointed out that the present invention just adopts a few electric power system of node to elaborate implementation process, under the big and more situation of node number, computational speed of the present invention is fast in electric power system, and the advantage that estimated accuracy is high can be more outstanding.

To sum up, the present invention can the online real-time electric power system that estimates use state, and the estimated accuracy height, and computational speed is fast, can reject bad data timely, guarantees the stable operation of system safety; In addition coherent system can coordinate each regional running status preferably, be one of following intelligent grid energy management center important module, to promoting further developing of intelligent grid significant.

Claims (1)

1. electric power system multizone distributed state estimation method based on the synchronous phase angle measurement mechanism, be used to carry out real-time, interactive from the metric data of SCADA system and the PMU voltage magnitude that estimates each node of electrical network and angle values and the state information that each is regional real-time, comprise following treatment step:

Step 1: based on the non-overlapping copies multizone electric power system decoupling zero of region characteristic

In an electric power system that contains the wiring of n bar, divide into r nonoverlapping regional S according to ground _i, each S _iN is arranged _iThe bar wiring, and link to each other with connecting line; After the decoupling zero, all by the local control centre control of oneself, the state of oneself is responsible for estimating in this center, and is connected to the control centre of a coordination by communication line in each zone;

Step 2: read each regional network data, comprise network configuration and line impedance, and form node admittance matrix and branch road-node incidence matrices;

Step 3: system measurements and PMU configuration;

The measurement of system comprises that voltage magnitude measures, gains merit and reactive power flow measures, meritorious and idle injection trend measures amount

The survey equation is as follows:

● voltage magnitude: ${V_l}^{\mathrm{mea}}=V_l+e_{V_l}$

● meritorious and reactive power flow:

$P_{\mathrm{lm}}^{\mathrm{mea}}=P_{\mathrm{lm}}+e_{P_{\mathrm{lm}}}$

$Q_{\mathrm{lm}}^{\mathrm{mea}}=Q_{\mathrm{lm}}+e_{Q_{\mathrm{lm}}}$

● meritorious and idle injection trend:

${P_l}^{\mathrm{mea}}=({V_l}^2{g}_l+\underset{m\&Element;a\left(l\right)}{\mathrm{\&Sigma;}}{P}_{\mathrm{lm}})+\underset{m\&Element;b\left(l\right)}{\mathrm{\&Sigma;}}{P}_{\mathrm{lm}}+e_{P_l}$

${Q_l}^{\mathrm{mea}}=(-{V_l}^2{b}_l+\underset{m\&Element;a\left(l\right)}{\mathrm{\&Sigma;}}{Q}_{\mathrm{lm}})+\underset{m\&Element;b\left(l\right)}{\mathrm{\&Sigma;}}{Q}_{\mathrm{lm}}+e_{Q_l}$

Wherein, injection is meritorious in the interior zone is respectively with reactive power flow measurement functional equation:

$P_{\mathrm{lm}}={V_l}^2(g_{\mathrm{lm}}+g_{\mathrm{slm}})-V_l{V}_m(g_{\mathrm{lm}}\cos ({\mathrm{\&theta;}}_l-{\mathrm{\&theta;}}_m)+b_{\mathrm{lm}}\sin ({\mathrm{\&theta;}}_l-{\mathrm{\&theta;}}_m))$

$Q_{\mathrm{lm}}=-{V_l}^2(b_{\mathrm{lm}}+b_{\mathrm{slm}})-V_l{V}_m(g_{\mathrm{lm}}\sin ({\mathrm{\&theta;}}_l-{\mathrm{\&theta;}}_m)-b_{\mathrm{lm}}\cos ({\mathrm{\&theta;}}_l-{\mathrm{\&theta;}}_m));$

L belongs to regional S _iA wiring;

With It is respectively the measuring value of wiring l and m place voltage magnitude, meritorious and reactive power flow, meritorious and idle injection trend; V _l, V _mIt is respectively the voltage magnitude of wiring l and m; θ _lAnd θ _mIt is respectively the angle values of wiring l and m;

With

It is respectively the error in measurement of voltage magnitude, meritorious and reactive power flow, meritorious and idle injection trend; g _l+ jb _LmIt is branch road l-m series admittance value; g _Slm+ jb _SlmIt is branch road l-m shunt admittance value; g _l+ jb _lIt is the shunt admittance value that is connected to line l; A(l) for linking l upward and belonging to regional S _iWired set; B (l) is for linking l upward and belonging to regional S _iWired set (i ≠ j);

● the PMU equipping rules

In order to guarantee the observability of whole system, among the present invention to each subsystem all dispose at least a PMU with and PMU all be configured on the connecting line between subsystem and the subsystem

Step 4: subregion state estimation

Adopt least square method to carry out state estimation respectively to all subregion, subregion state estimation model is as follows:

z _i=h _i(x _i)+e _i,i=1,2,...r

z _c=h _c(x)+e _c

Wherein, z _iBe regional S _iM _i* 1 dimension internal node or the local vector that measures; z _cBe m _c* 1 dimension boundary number direction finding amount; $x_i=\left(\begin{array}{c}{\mathrm{\&theta;}}_i\\ V_i\end{array}\right)$ Be regional S _iMiddle 2n _iThe local state vector of * 1 dimension comprises n _iIndividual magnitude of voltage and n _iIndividual angle values;

It is system-wide state vector

h _i(), h _c() is the Nonlinear Vector measurement function based on kirchhoff voltage or current law; e _i, e _cBe the gaussian random error vector;

Step 5:PMU coordinates

Adopt the accurate measurement in real time of PMU to come all subregion state information of system is carried out the linearity coordination, state estimation becomes linearity in the coordination level, and the linear condition estimation model is:

$z_p=\left[\begin{array}{c}{z}_{\mathrm{Vp}}\\ z_{\mathrm{Ip}}\end{array}\right]=\left[\begin{array}{c}I\\ Y\end{array}\right]\&CenterDot;x$

Wherein, z _VpAnd z _IpThe voltage and current that is respectively PMU measures vector; I represents the measurement unit matrix about voltage vector; Y represents the linear matrix about the current measurement vector, is made up of node and branch road admittance;

Step 6: distributed bad data is handled

Adopt the coherent system bad data to handle and two phase process of subregion bad data processing, that is:

1) calculates the sparse inverse matrix of interior zone by coherent system, finish the processing of coherent system bad data and with the result

Be transferred in the state estimator of subregion;

2) receiving the sparse inverse matrix of carrying out the subregion state estimation after coherent system passes the sparse inverse matrix of interior zone of coming and calculating each sub regions, and then finishing the subregion bad data and handle;

Step 7: if all state estimation all restrain, the bad data processing finishes, and then exports the voltage magnitude of each node and the identification result of angle values and output bad data, otherwise changes step 4.

Images