CN102801158B - Method for calculating time-lag electric power system eigenvalue and discriminating stability based on Pade approximation - Google Patents

Abstract

The invention discloses a method for calculating time-lag electric power system eigenvalues and discriminating stability based on Pade approximation. The method uses Pade approximation to make wide area feedback time-lag approximate to a rational polynomial. Though connecting with an electric power system without time-lag and a wide area damping controller, a linearized model of the time-lag electric power system is established. Finally, part of characteristic roots of the time-lag system is directly obtained according to a system state matrix, and then time-lag stability of the system is determined. According to an eigenvalue calculating result of a four-machine two-region example system, the method can relatively accurately solve part of the eigenvalues and eigenvectors which are corresponding to dynamic elements in the time-lag system. The method can conveniently and accurately determine time-lag stability of the system, correctly solve the number of the eigenvalues which are corresponding to a time-lag link and calculate precision which is related to order of the rational polynomial. In addition, the method is advantaged by small calculated amount, and short calculating time.

Description

Time-lag power system characteristic value based on Pade is approximate is calculated and stability distinguishing method

Technical field

The present invention relates to a kind of time-lag power system characteristic value and calculate and stability distinguishing method, relate in particular to a kind of based on the Pade approximate calculating of time-lag power system characteristic value and stability distinguishing method.

Background technology

Along with the expansion of electric power system scale, interregional low-frequency oscillation is just becoming the bottleneck of restriction electrical network transmittability, and existing damping controller (being mainly according to the power system stabilizer, PSS of one machine infinity bus system design) can not address this problem well.Basic reason is: (1) can not directly utilize relative merit angle and relative angle speed to form closed-loop control.Although it is the most direct and effective adopting relative merit angle and relative angle speed to realize damping control, lacks for a long time necessary method of measurement, can only adopt other indirect variablees to replace, and causes controlling poor effect.(2) be limited to local local message.Adopt local measurement to form FEEDBACK CONTROL, can not reflect well inter-area oscillation mode, although cause control system energy damping local oscillations pattern, be difficult to effectively suppress inter-area oscillation mode.

Along with the fast development of information processing and the communication technology, synchronized phasor technology and WAMS (Wide Area Measurement System, WAMS), new means are provided to monitoring, the analysis and control of electric power system, for interconnected network damping control has brought new opportunity.Synchronous phasor measurement unit can synchronous acquisition characterizes nearly all variable of operation of power networks state, and the most key is built-in potential, rotor angle, the angular speed that it can measure generator, busbar voltage phase place etc. and the closely-related amount of low-frequency oscillation.Under the support of high-speed communicating network (as electric power data broadband networks), each phasor measurement unit collection with time target data can be delivered to data center station with less time delay, complete synchronous processing and analysis, form WAMS.

Along with the continuous maturation of WAMS and perfect, the closed-loop control that is applied to large electrical network must be one of direction of Future Power System control development.WAMS can obtain relative merit angle and the angular speed between unit in certain time delay, and provide global information to the damping controller of dispersed placement, make to overcome existing damping controller and can only utilize local measurement to form this inherent shortcoming of FEEDBACK CONTROL, and effectively the low-frequency oscillation of the local and interval two kinds of patterns of inhibition becomes possibility.

But Wide-area Measurement Information is transmitted and processed in communication network, there is larger time lag.Time lag is a kind of major incentive [the Wu H X that causes the inefficacy of system control law, operation conditions deterioration and system unstability, Tsakalis K S, Heydt G T.Evaluation of time delay effects to wide-area power system stabilizer design[J] .IEEE Transactions on Power Systems, 2004,19 (4): 1935-1941.], while utilizing Wide-area Measurement Information to carry out the closed-loop control of electric power system, must take into account the impact of time lag.

In electric power system linearisation differential-algebraic equation characteristic of correspondence equation after consideration of time-delay, scalar time lag is converted into exponential term, thereby, characteristic equation is a transcendental equation, for this reason, in the time asking for the time lag stability margin of electric power system, conventionally adopt the method for functional transformation, as Rekasius conversion [Rifat S, Nejat O.A novel stability study on multiple time-delay systems (MTDS) using the root clustering paradigm[C] .Proceedings of the American Control Conference, Boston, MA, 2004, 5422-5427.], Lambert-W function [Yi S, Nelson P W, Ulsoy A G.Time-delay systems:analysis and control suing the Lambert W function[M] .Singapore:World Scientific Publishing Company, 2010], SCF method [Chen J, Gu G, Nett C N.A new method for computing delay margins for stability of linear delay systems[J] .Systems and Control Letters, 1995, 26 (2): 107-117.] convert surmounting item, avoid the difficulty of direct solution characteristic equation, there is conversion complexity and the large deficiency of amount of calculation in these class methods.In addition, based on LMI (Linear Matrix Inequlity, the stable adequacy criterion of Delay-Dependent LMI) [Yu Li. robust control---the LMI method [M] of exerting oneself. Beijing: publishing house of Tsing-Hua University, 2002.], also be widely used in time lag stability analysis and the time lag upper limit solves [Jiang Quanyuan, Zhang Pengxiang, Cao Yijia. take into account the wide area FACTS damping control [J] of feedback signal time-delay. Proceedings of the CSEE, 2006, 26 (7): 82-88.Jiang Quanyuan, Zhang Pengxiang, Cao Yijia.Wide-area FACTS damping control in consideration of feedback signal's time delays[J] .Proceedings of the CSEE, 2006, 26 (7): 82-88.], these class methods combine to reduce amount of calculation with system model depression of order conventionally, but still there is intrinsic conservative shortcoming.

In wide area damp control research, Pade is approximate is a kind of conventional Time Delay processing method.Approach Time Delay by rational polynominal, and then can utilize easily classical and modern control theory design wide area damp control law [Wu H X, Ni H, Heyde G T.The impact of time delay on control design in power systems[C] .Proceedings of IEEE Power Engineering Society Winter Meeting, New York, USA, 2002, 1511-1516. king becomes mountain, Shi Jie. consider the Design of Power System Stabilizer [J] of Time Delay Impact. Proceedings of the CSEE, 2007, 27 (10): 1-6.Wang Chengshan, Shi Jie.PSS designing with consideration of time delay impact[J] .Proceedings of the CSEE, 2007, 27 (10): 1-6. stone a word used in person's names, Wang Chengshan. consider the H ∞ damping controller [J] of Wide-area Measurement Information time delay influence. Proceedings of the CSEE, 2008, 28 (1): 30-34.Shi Jie, Wang Chengshan.Design of H ∞ controller for damping interarea oscillations with consideration of the delay of remote signal[J] .Proceedings of the CSEE, 2008, 28 (1): 30-34. Hu Zhijian, Zhao Yishu. take into account the interconnected electric power system robust stability control [J] of WAMS time lag. Proceedings of the CSEE, 2010, 30 (19): 37-43.Hu Zhijian, Zhao Yishu.Robust stability control of power systems based on WAMS with signal transmission delays[J] .Proceedings of the CSEE, 2010, 30 (19): 37-43.], and by the validity of time-domain-simulation access control device.

Eigenvalues analysis method is the effective method substantially and of microvariations stability analysis.Document [Jia Hongjie, Chen Jianhua, Yu Xiaodan. the impact [J] of Time Delay on Power System microvariations stability. Automation of Electric Systems, 2006,30 (10): 5-8,17.

Jia Hongjie, Chen Jianhua, Yu Xiaodan.Impact of time delay on power system small signal stability[J] .Automation of Electric Power Systems, 2006, 20 (5): 5-8, 17.] and [Yuan Ye, Cheng Lin, Sun Yuanzhang, Deng. the time-delay analysis of wide area damp control and time lag compensation design [J] thereof. Automation of Electric Systems, 2006, 30 (14): 6-9.Yuan Ye, Cheng Lin, Sun Yuanzhang, et al.Effect of delay input on wide-area damping control and design of compensation[J] .Automation of Electric Power Systems, 2006, 30 (14): 6-9.] analyze the impact of time lag on Study of Power System Small Disturbance Stability by the characteristic value of calculating one machine infinity bus system and equivalent two machine systems respectively, but do not provide detailed characteristic value computational methods.Document [Jia Hongjie, the remaining Power System Delay stability margin [J] of knowing under red .2 kind physical constraint. Automation of Electric Systems, 2008, 32 (9): 7-10, 19.Jia Hongjie, Yu Xiaodan.Method of determining power system delay margins with considering two practical constraints[J] .Automation of Electric Power Systems, 2008, 32 (9): 7-10, 19.] proposed a kind of directly, have on efficient search complex plane on specific border (as the imaginary axis, characteristic value real part or damping ratio equal given constant) method of the critical eigenvalue of time-lag system.The accurate characteristic value that the method can obtain, but the amount of calculation of search procedure is larger.

Eigenvalues analysis method has formed comparative maturity and perfect theory, and obtains extensive use in the practice of power industry.If can propose the characteristic value computational methods of time-lag system, and then continue to use the thinking of classical Eigenvalues analysis and theoretical frame and analyze the microvariations stability of time-lag power system, no matter for improving and the abundant microvariations stability analysis theory based on characteristic value, still promote the engineering application of wide area damp control, all will have great importance and be worth.Based on this thought, the present invention proposes a kind of based on the Pade approximate calculating of time-lag power system characteristic value and stability distinguishing method, state space by modeling and embedding Time Delay is expressed, and then can directly utilize routine or sparse eigenvalue method to try to achieve the partial feature value of system, and then the time lag stability of judgement system.For Si Jiliang district example system, by with discretization characteristic value method for solving [Engelborghs K, Roose D.On stability of LMS methods and characteristic roots of delay differential equations[J] .SIAM Journal on Numerical Analysis, 2003,40 (2): 629-650.] contrast of result of calculation, has verified correctness and the validity of the inventive method.

The time lag that wide area measurement system signal produces in transmission and processing procedure, makes electric power system become a time-lag system.Annexation between time-lag power system each several part as shown in Figure 2.

Without time-lag power system model

If describe the differential-Algebraic Equation set of electric power system during without wide area damping control be:

$\left\{\begin{array}{c}\stackrel{\·}{x}=f(x,y)\\ 0=g(x,y)\end{array}\right.---\left(1\right)$

In formula, f is for describing the dynamic differential equation of element, and g is network equation, and x is system state variables, and y is system algebraically variable (node voltage),

for the differential of system state variables.

Be output as u without time-lag power system, it feeds back the input as damping controller by wide area.If f _ffor describing the function of contact between u and (x, y), u=f _f(x, y).Y _cfor the output of wide area damping control, and as the control inputs without time-lag power system.At steady operation point (x ₀, y ₀) formula (1) and u being carried out to linearisation, can obtain:

$\left\{\begin{array}{c}\mathrm{\&Delta;}\stackrel{\&CenterDot;}{x}=\mathrm{A\&Delta;x}+\mathrm{B\&Delta;y}+\mathrm{E\&Delta;}{\mathrm{<figure-callout\; id="0"\; label="y\; c"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">y</figure-callout>}}_c\\ 0=\mathrm{C\&Delta;x}+\mathrm{D\&Delta;y}\\ \mathrm{\&Delta;u}=K_1\mathrm{\&Delta;x}+{\mathrm{<figure-callout\; id="2"\; label="K"\; filenames="HDA00001956071700012.png,HDA00001956071700022.png"\; state="\{\{state\}\}">K</figure-callout>}}_2\mathrm{\&Delta;y}\end{array}\right.---\left(2\right)$

In formula, A, B, C, D are respectively as differential equation f and network equation g are with respect to the partial derivative of system state variables x and algebraically variable y, $A=\frac{\&PartialD;f}{\&PartialD;x},$ $B=\frac{\&PartialD;f}{\&PartialD;y},$ $C=\frac{\&PartialD;g}{\&PartialD;x},$ $D=\frac{\&PartialD;g}{\&PartialD;y},$ ${\mathrm{<figure-callout\; id="1"\; label="K"\; filenames="HDA00001956071700012.png,HDA00001956071700031.png"\; state="\{\{state\}\}">K</figure-callout>}}_1=\frac{\&PartialD;f_f}{\&PartialD;x},$ ${\mathrm{<figure-callout\; id="2"\; label="K"\; filenames="HDA00001956071700012.png,HDA00001956071700022.png"\; state="\{\{state\}\}">K</figure-callout>}}_2=\frac{\&PartialD;f_f}{\&PartialD;y},.$ E is the input matrix of system, and the nonzero element in E has characterized additional wide area damping control and added the annexation between control appliance (as automatic excitation adjustor of generator, FACTS equipment etc.).The row at nonzero element place in E, in corresponding control appliance the output variable of amplifying element in the position without in time-lag power system state variable x, the multiplication factor that nonzero element value is amplifying element and the ratio of time constant.

Wide area damping control state space is expressed

Dynamic and the output of wide area damping control can be represented by following differential-Algebraic Equation set:

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_c=f_c(x_c,y_d)\\ y_c=g_c(x_c,y_d)\end{array}\right.---\left(3\right)$

In formula, f _cfor the differential equation, g _cfor algebraic equation, x _cfor the state variable in wide area damping control, y _dfor the input of wide area damping control.Y _d=ue ^{-τ t}, the input that u is wide area damping control, τ=[τ ₁..., τ _i..., τ _m] ^tfor the vector that the time lag of all wide area damping controls forms, τ _i>0 is the time lag constant of i Time Delay, i=1, and 2 ..., m, m is positive integer, total number of Time Delay in expression system.

The lienarized equation that equation (3) is corresponding is:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}_c=A_c\mathrm{\&Delta;}{x}_c+B_c\mathrm{\&Delta;}{y}_d\\ \mathrm{\&Delta;}{y}_c=C_c\mathrm{\&Delta;}{x}_c+D_c\mathrm{\&Delta;}{y}_d\end{array}\right.---\left(4\right)$

In formula, A _c, B _c, C _c, D _crepresent respectively differential equation f _cwith algebraic equation g _cwith respect to state variable x _cwith algebraically variable y _dpartial derivative, $A_c=\frac{\&PartialD;f_c}{\&PartialD;x_c},$ $B_c=\frac{\&PartialD;f_c}{\&PartialD;y_d},$ $C_c=\frac{{\&PartialD;g}_c}{{\&PartialD;x}_c},$ $D_c=\frac{{\&PartialD;g}_c}{{\&PartialD;y}_c}.$

Time-lag power system model

By the differential equation f of descriptive system dynamic element and the dynamic differential equation f of description wide area damping control _cwrite together, formation matrix equation vector f ',

correspondingly, by the state variable x in the two and x _cwrite together, form state variable vector x ',

so, consider that after wide area Feedback Delays τ, electric power system can be described by following Delay Differential-algebraic equation:

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}^{\&prime;}=f^{\&prime;}(x^{\&prime;}y,x_{\mathrm{\&tau;}1}^{\&prime;},y_{\mathrm{\&tau;}1},x_{\mathrm{\&tau;}2}^{\&prime;},y_{\mathrm{\&tau;}2},...,x_{\mathrm{\&tau;m}}^{\&prime;},y_{\mathrm{\&tau;m}})\\ 0=g(x^{\&prime;},y)\\ 0=g\left({x_{\mathrm{\&tau;}1}^{\&prime;},y}_{\mathrm{\&tau;}1}\right)\\ 0=g(x_{\mathrm{\&tau;}2}^{\&prime;},y_{\mathrm{\&tau;}2})\\ .\\ .\\ .\\ 0=g(x_{\mathrm{\&tau;m}}^{\&prime;},y_{\mathrm{\&tau;m}})\end{array}\right.---\left(5\right)$

In formula, [x ' _{τ i}, y _{τ i}]=[x ' (t-τ _i), y (t-τ _i)] be time lag state variable and algebraically variable, i is positive integer.

Steady operation point (x ' ₀, y ₀) locate (5) to carry out linearisation, can obtain system-wide inearized model:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;}={\mathrm{<figure-callout\; id="0"\; label="A"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">A</figure-callout>}}_0^{\&prime;}\mathrm{\&Delta;}{x}^{\&prime;}{\mathrm{<figure-callout\; id="0"\; label="B"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">B</figure-callout>}}_0^{\&prime;}\mathrm{\&Delta;y}+\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}(A_{\mathrm{\&tau;i}}^{\&prime;}\mathrm{\&Delta;}{x}_{\mathrm{\&tau;i}}^{\&prime;}+B_{\mathrm{\&tau;i}}^{\&prime;}\mathrm{\&Delta;}{y}_{\mathrm{\&tau;i}})\\ 0={\mathrm{<figure-callout\; id="0"\; label="C"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">C</figure-callout>}}_0^{\&prime;}\mathrm{\&Delta;}{x}^{\&prime;}+{\mathrm{<figure-callout\; id="0"\; label="D"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">D</figure-callout>}}_0^{\&prime;}\mathrm{\&Delta;y}\\ 0=C_{\mathrm{\&tau;}1}^{\&prime;}\mathrm{\&Delta;}{x}_{\mathrm{\&tau;}1}^{\&prime;}+D_{\mathrm{\&tau;}1}^{\&prime;}\mathrm{\&Delta;}{y}_{\mathrm{\&tau;}1}\\ .\\ .\\ .\\ 0=C_{\mathrm{\&tau;m}}^{\&prime;}\mathrm{\&Delta;}{x}_{\mathrm{\&tau;m}}^{\&prime;}+D_{\mathrm{\&tau;m}}^{\&prime;}\mathrm{\&Delta;}{y}_{\mathrm{\&tau;m}}\end{array}\right.----\left(6\right)$

In formula: Α ₀', B ' ₀, C ' ₀, D ' ₀represent respectively differential equation f ' and the algebraic equation g partial derivative to state variable x ' and algebraically variable y,

Α ' _{τ i}, B ' _{τ i}, C ' _{τ i}, D ' _{τ i}represent that respectively differential equation f ' and algebraic equation g are to time lag state state variable x ' _{τ i}with algebraically variable y _{τ i}partial derivative

Time-lag power system characteristic value Solve problems

Work as D ₀' and D ' _{τ i}(i=1,2 ..., when m) nonsingular, cancellation algebraically variable, equation (6) can be reduced to:

$\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;}={\tilde{A}}_0^{\&prime;}\mathrm{\&Delta;}{x}^{\&prime;}+\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\tilde{A}}_{\mathrm{\&tau;i}}^{\&prime;}\mathrm{\&Delta;}{x}_{\mathrm{\&tau;i}}^{\&prime;}---\left(7\right)$

In formula,

$\begin{array}{c}{\tilde{A}}_0^{\&prime;}=A_0^{\&prime;}-{\mathrm{<figure-callout\; id="0"\; label="B"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">B</figure-callout>}}_0^{\&prime;}{\left({\mathrm{<figure-callout\; id="0"\; label="D"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">D</figure-callout>}}_0^{\&prime;}\right)}^{-1}{\mathrm{<figure-callout\; id="0"\; label="C"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">C</figure-callout>}}_0^{\&prime;}\\ {\tilde{A}}_{\mathrm{\&tau;i}}^{\&prime;}=A_{\mathrm{\&tau;i}}^{\&prime;}-B_{\mathrm{\&tau;i}}^{\&prime;}{\left(D_{\mathrm{\&tau;i}}^{\&prime;}\right)}^{-1}{C}_{\mathrm{\&tau;i}}^{\&prime;}\end{array}---\left(8\right)$

The characteristic equation of the linearized system that formula (7) represents is:

$\det (\mathrm{\&lambda;I}-{\tilde{A}}_0^{\&prime;}-\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}{\tilde{A}}_{\mathrm{\&tau;i}}^{\&prime;}{e}^{-\mathrm{\&lambda;}{\mathrm{\&tau;}}_i})=0---\left(9\right)$

In formula, the characteristic value that λ is system.

The Theory of Stability of Dynamic System with Delays point out [Liao Xiaoxin. the Theory of Stability of dynamical system and application [M]. Beijing: National Defense Industry Press, 2000.], if the All Eigenvalues of linearized system (7) all has negative real part, time-lag system (5) steady operation point (x ' ₀, y ₀) to locate be that microvariations are stable; Otherwise if at least there is a characteristic value with positive real part, system is that microvariations are unsettled at this some place.

But, due to time lag and exponential term

existence, the characteristic equation (9) of linearized system is transcendental equation, it has infinite multiple solution.Therefore, obtain the critical eigenvalue of system by this equation of direct solution and the time lag stability of judgement system becomes very difficult.

Chinese patent 201010123345.8, proposed a kind of directly, have on efficient search complex plane a method of the critical eigenvalue of (as the imaginary axis, characteristic value real part or damping ratio equal given constant) time-lag system on specific border.The accurate characteristic value that the method can obtain, but the amount of calculation of search procedure is larger.

Chinese patent 200810151217.7,200910070255.4,200910070254.X differentiate the small signal stability of time-lag system based on LMI (Linear Matrix Inequlity, LMI) method.In the time of the stability of judgement system, there is intrinsic conservative in these class methods.

Summary of the invention

Object of the present invention is exactly in order to address the above problem, providing a kind of calculates and stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, it does not need iteration just can solve exactly the partial feature value relevant to dynamic element in time-lag system and characteristic vector with search, the characteristic value that utilization calculates can directly judge the small signal stability of system, there is amount of calculation little, computing time is few, does not have the advantage of any conservative.

To achieve these goals, the present invention adopts following technical scheme:

Time-lag power system characteristic value based on Pade is approximate is calculated and a stability distinguishing method, and concrete steps are as follows:

Step 1: specify wide area Feedback Delays τ _i;

Step 2: carry out Pade and be similar to, obtain the approximate rational polynominal of Time Delay;

Step 3: rational polynominal is converted into state-space expression;

Step 4: carry out equilibrating processing;

Step 5: by Time Delay be connected without time-lag power system, wide area damping control, obtain the inearized model of time-lag power system;

Step 6: the characteristic value of utilizing QR method or sparse eigenvalue Algorithm for Solving system;

Step 7: the directly time lag stability of judgement system.

In described step 2, in Laplace domain, i Time Delay can be expressed as

t is the time, and s is frequency.Pade is approximate is that a kind of [l, k] rank rational polynominal of utilizing is approached

method:

$e^{-{\mathrm{\&tau;}}_{i}s}\&ap;R\left(s\right)=\frac{N_l\left(s\right)}{N_k\left(s\right)}=\frac{a_0+a_1{\mathrm{\&tau;}}_{i}s+...+a_l{\left({\mathrm{\&tau;}}_{i}s\right)}^{\mathrm{<figure-callout\; id="0"\; label="l\; b"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">l</figure-callout>}}}{b_{}}---\left(10\right)$

In formula, l and k are positive integer.Coefficient a _iand b _ifor real number, can be obtained by following formula:

$a_j={(-1)}^j\frac{(l+k-j)!l!}{j!(l-j)!}---\left(11\right)$

$b_j=\frac{(l+k-j)!k!}{j!(k-j)!}---\left(12\right)$

Described step 2 is in Pade is approximate, and exponent number l and k are larger, rational polynominal R (s) more close to

under normal circumstances, get l=k; Time lag is less, R in frequency domain (s) with

the interval that phase place is consistent is larger, and frequency band is wider.

In described step 3, by the molecule of formula (10) and denominator simultaneously divided by

and then utilize transfer function problem of implementation general approach [Liu leopard. modern control theory [M]. Beijing: China Machine Press, 2000.], rational polynominal based on the approximate Time Delay obtaining of Pade is approached to R (s) and be converted to state space expression of equal value, for i Time Delay, represent take variable incremental form can control standard type as:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}_{\mathrm{di}}={\tilde{A}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{B}}_d\mathrm{\&Delta;}{u}_i\\ \mathrm{\&Delta;}{y}_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{D}}_{\mathrm{di}}\mathrm{\&Delta;}{u}_i\end{array}\right.---\left(13\right)$

In formula (13), coefficient matrix

expression be:

K is the approximate exponent number of Pade.

Compared with Small Time Lag τ _iin the situation of higher-order number, between the coefficient in the approximate transfer function R (s) obtaining of Pade, and the coefficient matrices A of the state-space expression being converted to by transfer function _diand C _dinonzero element between, on the order of magnitude, differ greatly, adopt linear transformation to make the coefficient matrix balance more of the state-space expression of Time Delay:

$A_{\mathrm{di}}=T^{-1}{\tilde{A}}_{\mathrm{di}}T,B_{\mathrm{di}}=T^{-1}{\tilde{B}}_{\mathrm{di}},C_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}T,D_{\mathrm{di}}={\tilde{D}}_{\mathrm{di}}---\left(15\right)$

In formula, T is diagonal transformation matrix; Hence one can see that, and equilibrating coefficient matrix after treatment still meets can control standard type;

Described step 5: by the state-space expression of the Time Delay obtaining in step 3, be connected with the inearized model without time-lag power system, the inearized model of wide area damping control, thereby set up the inearized model of the Study of Power System Small Disturbance Stability analysis that comprises Time Delay:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;\&prime;}=A^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}+B^{\&prime;\&prime;}\mathrm{\&Delta;<figure-callout\; id="0"\; label="y"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">y</figure-callout>}\\ 0=C^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}+D^{\&prime;\&prime;}\mathrm{\&Delta;y}\end{array}\right.---\left(16\right)$

" be system state variables, y is system algebraically variable (node voltage), and A ", B ", C ", D " are coefficient matrix in formula, x;

Suppose the state variable x of i Time Delay _dstate variable x with the wide area damping control being attached thereto _cbe arranged in without after time-lag power system state variable x, $\mathrm{\&Delta;}{x}^{\&prime;\&prime;}={\left[\begin{array}{ccccccccc}{\mathrm{\&Delta;x}}^T& .& .& .& {\mathrm{\&Delta;x}}_{\mathrm{di}}^T& {\mathrm{\&Delta;x}}_{\mathrm{ci}}^T& .& .& .\end{array}\right]}^T;$ If be A without the coefficient matrix of time-lag power system, B, C, D; The coefficient matrix of wide area damping control is A _c, B _c, C _c, D _c, the coefficient matrix of wide area feedback signal is respectively K _1iand K _2i; So coefficient matrices A in formula (16) ", B ", C ", D " specifically can be expressed as:

In formula (17), coefficient matrices A _di, B _di, C _di, D _disatisfiedly can control standard type, A is a point Block diagonal matrix, B, C is piecemeal sparse matrix, D is 2 × 2 piecemeal sparse matrix [Du Zhengchun, Liu Wei, Fang Wanliang, Deng. the sparse realization [J] that in analysis on Small Disturbance Stability, a kind of critical eigenvalue calculates. Proceedings of the CSEE, 2005, 25 (2): 17-21.Du Zhengchun, Liu Wei, Fang Wanliang, et al.A sparse method for the calculation of critical eigenvalue in small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (2): 17-21.], E is sparse matrix, and Aci, Bci, Cci, Dci are relevant with the concrete structure of wide area damping control,

In described formula (17), each coefficient matrix has following features:

(1) " be still a point Block diagonal matrix if using the coefficient matrix of the dynamic element relating in wide area feedback signal, Time Delay, wide area damping control state variable as sub-block, an A; " only on Time Delay and wide area damping control state variable are expert at, increased null matrix or a small amount of several row nonzero element, B " is still Sparse Array to B;

(2) C " has only increased null matrix, D " with just the same without the corresponding coefficient matrix of time-lag power system in Time Delay and wide area damping control state variable column.

In summary, the lienarized equation coefficient matrices A of the electric power system that comprises Time Delay ", B ", C ", D ", with without time-lag power system at linearisation equation coefficient matrix A, B, C, D, there is identical sparsity structure.

Described step 6: as D " when nonsingular, cancellation algebraically variable y, formula (16) can be reduced to:

$\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;\&prime;}={\tilde{A}}^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}---\left(18\right)$

In formula,

for the state matrix based on the approximate time-lag system obtaining of Pade.

By calculating

characteristic value, just can obtain the partial feature value of time-lag power system.In addition, because the electric power system that comprises Time Delay has identical sparse characteristic with conventional without the coefficient matrix of time-lag power system, sparse treatment technology [the Du Zhengchun using in the time calculating on a large scale without time-lag power system part critical eigenvalue, Liu Wei, Fang Wanliang, Deng. the sparse realization [J] that in analysis on Small Disturbance Stability, a kind of critical eigenvalue calculates. Proceedings of the CSEE, 2005, 25 (2): 17-21.Du Zhengchun, Liu Wei, Fang Wanliang, et al.A sparse method for the calculation of critical eigenvalue in small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (2): 17-21.] and computational methods [Du Zhengchun, Liu Wei, Fang Wanliang, Deng. the critical eigenvalue in the analysis on Small Disturbance Stability based on Jacobi-Davidson method calculates [J]. Proceedings of the CSEE, 2005, 25 (14): 19-24.Du Zhengchun, Liu Wei, Fang Wanliang, et al.The application of the Jacobi-Davidson method to the calculation of critical eigenvalues in the small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (14): 19-24.], stand good in the partial feature value of calculating time-lag power system.

Beneficial effect of the present invention: utilizing approximate wide area Feedback Delays is approached of Pade is a rational polynominal, by with without being connected between time-lag power system and wide area damping control, set up time-lag power system inearized model, finally directly obtain the Partial Feature root of time-lag system according to system mode matrix; The method can solve the partial feature value relevant to dynamic element in time-lag system and characteristic vector more exactly, correctly solves characteristic value number and the computational accuracy relevant to Time Delay, relevant with the exponent number of rational polynominal.

Utilize the approximate rational polynominal of Pade to approach Time Delay, and then calculate the partial feature value that contains wide-area communication time-lag power system, thereby directly judge the small signal stability of electric power system.

Accompanying drawing explanation

Fig. 1 is overall flow figure of the present invention;

Fig. 2 is time-lag power system schematic diagram;

Tu3Wei Liang district four machine system diagrams;

Fig. 4 is the partial feature value that real part is greater than-50;

Fig. 5 is the partial feature value that real part is greater than-10.

Wherein, 1. without time-lag power system, 2.e ^{-s τ}, 3. wide area damping control.

Embodiment

Below in conjunction with accompanying drawing and embodiment, the invention will be further described.

As shown in Figure 1, a kind of based on the Pade approximate calculating of time-lag power system characteristic value and stability distinguishing method, concrete steps are as follows:

Step 1: specify wide area Feedback Delays τ _iwith the approximate exponent number k of Pade;

Step 2: carry out Pade and be similar to, obtain the approximate rational polynominal of Time Delay;

Step 3: rational polynominal is converted into state-space expression;

Step 4: carry out equilibrating processing;

Step 5: by Time Delay be connected without time-lag power system, wide area damping control, obtain the inearized model of time-lag power system;

Step 6: the characteristic value of utilizing QR method or sparse eigenvalue Algorithm for Solving system;

Step 7: the directly time lag stability of judgement system.

In described step 2, in Laplace domain, i Time Delay

can be expressed as

pade is approximate is that a kind of [l, k] rank rational polynominal of utilizing is approached

method:

$e^{-{\mathrm{\&tau;}}_{i}s}\&ap;R\left(s\right)=\frac{N_l\left(s\right)}{N_k\left(s\right)}=\frac{a_0+a_1{\mathrm{\&tau;}}_{i}s+...+a_l{\left({\mathrm{\&tau;}}_{i}s\right)}^{\mathrm{<figure-callout\; id="0"\; label="l\; b"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">l</figure-callout>}}}{b_{}}---\left(19\right)$

In formula, coefficient ai and bi can be obtained by following formula:

$a_j={(-1)}^j\frac{(l+k-j)!l!}{j!(l-j)!}---\left(20\right)$

$b_j=\frac{(l+k-j)!k!}{j!(k-j)!}---\left(21\right)$

Described step 2 is in Pade is approximate, and exponent number l and k are larger, rational polynominal R (s) more close to under normal circumstances, get l=k; Time lag is less, R in frequency domain (s) with

the interval that phase place is consistent is larger, and frequency band is wider.

In described step 3, by the molecule of formula (10) and denominator simultaneously divided by

and then utilize transfer function problem of implementation general approach [Liu leopard. modern control theory [M]. Beijing: China Machine Press, 2000.], rational polynominal based on the approximate Time Delay obtaining of Pade is approached to R (s) and be converted to state space expression of equal value, for i Time Delay, represent take variable incremental form can control standard type as:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}_{\mathrm{di}}={\tilde{A}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{B}}_d\mathrm{\&Delta;}{u}_i\\ \mathrm{\&Delta;}{y}_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{D}}_{\mathrm{di}}\mathrm{\&Delta;}{u}_i\end{array}\right.---\left(22\right)$

In formula (13),

Compared with Small Time Lag τ _iin the situation of higher-order number, between the coefficient in the approximate transfer function R (s) obtaining of Pade, and the coefficient matrices A of the state-space expression being converted to by transfer function _diand C _dinonzero element between, on the order of magnitude, differ greatly, adopt linear transformation to make the coefficient matrix balance more of the state-space expression of Time Delay:

$A_{\mathrm{di}}=T^{-1}{\tilde{A}}_{\mathrm{di}}T,B_{\mathrm{di}}=T^{-1}{\tilde{B}}_{\mathrm{di}},C_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}T,D_{\mathrm{di}}={\tilde{D}}_{\mathrm{di}}---\left(24\right)$

In formula, T is diagonal transformation matrix; Hence one can see that, and equilibrating coefficient matrix after treatment still meets can control standard type;

Described step 5: by the state-space expression of the Time Delay obtaining in step 3, be connected with the inearized model without time-lag power system, the inearized model of wide area damping control, thereby set up the inearized model of the Study of Power System Small Disturbance Stability analysis that comprises Time Delay:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;\&prime;}=A^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}+B^{\&prime;\&prime;}\mathrm{\&Delta;<figure-callout\; id="0"\; label="y"\; filenames="HDA00001956071700031.png"\; state="\{\{state\}\}">y</figure-callout>}\\ 0=C^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}+D^{\&prime;\&prime;}\mathrm{\&Delta;y}\end{array}\right.---\left(25\right)$

" be system state variables, y is system algebraically variable (node voltage), and A ", B ", C ", D " are coefficient matrix in formula, x;

Suppose the state variable x of i Time Delay _dstate variable x with the wide area damping control being attached thereto _cbe arranged in without after time-lag power system state variable x, $\mathrm{\&Delta;}{x}^{\&prime;\&prime;}={\left[\begin{array}{ccccccccc}{\mathrm{\&Delta;x}}^T& .& .& .& {\mathrm{\&Delta;x}}_{\mathrm{di}}^T& {\mathrm{\&Delta;x}}_{\mathrm{ci}}^T& .& .& .\end{array}\right]}^T;$ If be A without the coefficient matrix of time-lag power system, B, C, D; The coefficient matrix of wide area damping control is A _c, B _c, C _c, D _c, the coefficient matrix of wide area feedback signal is respectively K _1iand K _2i; So coefficient matrices A in formula (16) ", B ", C ", D " specifically can be expressed as:

In formula (17), coefficient matrices A _di, B _di, C _di, D _disatisfiedly can control standard type, A is a point Block diagonal matrix, B, C is piecemeal sparse matrix, D is 2 × 2 piecemeal sparse matrix [Du Zhengchun, Liu Wei, Fang Wanliang, Deng. the sparse realization [J] that in analysis on Small Disturbance Stability, a kind of critical eigenvalue calculates. Proceedings of the CSEE, 2005, 25 (2): 17-21.Du Zhengchun, Liu Wei, Fang Wanliang, et al.A sparse method for the calculation of critical eigenvalue in small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (2): 17-21.], E is sparse matrix, and Aci, Bci, Cci, Dci are relevant with the concrete structure of wide area damping control,

In described formula (17), each coefficient matrix has following features:

(1) " be still a point Block diagonal matrix if using the coefficient matrix of the dynamic element relating in wide area feedback signal, Time Delay, wide area damping control state variable as sub-block, an A; " only on Time Delay and wide area damping control state variable are expert at, increased null matrix or a small amount of several row nonzero element, B " is still Sparse Array to B;

(2) C " has only increased null matrix, D " with just the same without the corresponding coefficient matrix of time-lag power system in Time Delay and wide area damping control state variable column.

In summary, the lienarized equation coefficient matrices A of the electric power system that comprises Time Delay ", B ", C ", D ", with without time-lag power system at linearisation equation coefficient matrix A, B, C, D, there is identical sparsity structure.

Described step 6: as D " when nonsingular, cancellation algebraically variable y, formula (16) can be reduced to:

$\mathrm{\&Delta;}{\stackrel{\&CenterDot;}{x}}^{\&prime;\&prime;}={\tilde{A}}^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\&prime;\&prime;}---\left(27\right)$

In formula,

for the state matrix based on the approximate time-lag system obtaining of Pade.

By calculating

characteristic value, just can obtain the partial feature value of time-lag power system.In addition, because the electric power system that comprises Time Delay has identical sparse characteristic with conventional without the coefficient matrix of time-lag power system, sparse treatment technology [the Du Zhengchun using in the time calculating on a large scale without time-lag power system part critical eigenvalue, Liu Wei, Fang Wanliang, Deng. the sparse realization [J] that in analysis on Small Disturbance Stability, a kind of critical eigenvalue calculates. Proceedings of the CSEE, 2005, 25 (2): 17-21.Du Zhengchun, Liu Wei, Fang Wanliang, et al.A sparse method for the calculation of critical eigenvalue in small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (2): 17-21.] and computational methods [Du Zhengchun, Liu Wei, Fang Wanliang, Deng. the critical eigenvalue in the analysis on Small Disturbance Stability based on Jacobi-Davidson method calculates [J]. Proceedings of the CSEE, 2005, 25 (14): 19-24.Du Zhengchun, Liu Wei, Fang Wanliang, et al.The application of the Jacobi-Davidson method to the calculation of critical eigenvalues in the small signal stability analysis[J] .Proceedings of the CSEE, 2005, 25 (14): 19-24.], stand good in the partial feature value of calculating time-lag power system.

The time lag that wide area measurement system signal produces in transmission and processing procedure, makes electric power system become a time-lag system.Annexation between time-lag power system each several part as shown in Figure 2, without time-lag power system and e ^{-s τ}connect e ^{-s τ}be connected with wide area damping control.

Liang district four electro-mechanical force systems are example as shown in Figure 3, correctness and the validity of the time-lag power system characteristic value computational methods that checking the present invention proposes, the parameter of system refers to [Kunder P.Power System Stabilityand Control[M] .New York:McGraw-Hill, 1994.].

At generator G1 and G3, local power system stabilizer, PSS (Power System Stabilizer is installed, PSS) on basis, consider the wide area PSS of installing take G1 and G3 relative rotation speed deviation delta ω 13 as feedback signal on G1, further improve the damping level of system.

The partial feature value that the present invention calculates using the time-lag system characteristic value method for solving based on discretization is as exact solution, and propose as checking the present invention based on the approximate characteristic value computational methods validity of Pade and the benchmark of accuracy.

The present invention is by the exact solution of the partial feature value of the example system that in Matlab tool box DDE-BIFTOOL, p_stabil function calculates, and wherein, the convergence precision of Newton iterative method is taken as 1e-8.

For example system, be similar under the exponent number k of rational polynominal at different Feedback Delays τ and different Pade, that utilizes that the present invention proposes has carried out a large amount of calculating based on the approximate characteristic value computational methods of Pade, and carried out contrast and analysis with the result of the characteristic value computational methods based on discretization, correctness and the validity of the inventive method is only described with some numerical results below.

In the time of Feedback Delays τ=0.30s, utilize characteristic value computational methods based on discretization, based on 5 rank and the approximate characteristic value computational methods of 20 rank Pade, the partial feature value that the real part obtaining is greater than-50, as shown in Figure 4.Detailed distribution in the rectangle frame of Fig. 4 lower right corner, that real part is greater than-10 partial feature value as shown in Figure 5, as shown in Figure 5:

(1) in the time of k=5 and 20, utilize based on the approximate characteristic value computational methods of Pade, can calculate comparatively exactly real part and be greater than the partial feature value in-10 complex plane, known by further analysis, it is corresponding with the state variable of descriptive system dynamic element.

(2) for k=5, be less than at real part in-10 complex plane, based on the approximate characteristic value computational methods of Pade, can obtain several wrong characteristic values, as :-15.06761008 ± 23.47922927i ,-21.53755703 ± 15.83642620i ,-36.68588408 ,-49.61254513; In addition, also have multiple characteristic values not to be calculated, known by further analysis, omit and the wrong characteristic value solving, main corresponding with each rank variable of rational polynominal.

(3) only have in the time that the exponent number of rational polynominal significantly increases to k=20, guarantee the inventive method is greater than the situation that does not occur leaking root, wrong root in-41 complex plane at real part, the partial feature value corresponding with rational polynominal variable also can be tried to achieve comparatively exactly.In the time that exponent number is increased to 30, be greater than at real part in-50 complex plane, it is identical that the present invention and the characteristic value method for solving based on discretization obtain characteristic value.

When table 1 τ=0.3s, the characteristic value obtaining based on discretization and Pade approximate calculation

(4) be similar to rational polynominal exponent number with in accurately asking for the partial feature value corresponding with rational polynominal variable at increase Pade, the precision of the partial feature value corresponding with system dynamic element is also improved significantly.As shown in table 1, for example, in the time of k=5, maximum absolute error max (abs (Δ Re (λ) between the inventive method and the characteristic value that obtains based on discretization method, Δ Im (λ)))=1.042e-4, in the time that k is increased to 7, maximum absolute error is correspondingly reduced to 4.673e-8, and characteristic value computational accuracy is significantly improved.

(5) from the characteristic value of the system calculating, in the time of τ=0.3s, system is small interference stability.

For the exponent number k of the approximate rational polynominal of the deep investigation Pade impact on characteristic value computational accuracy, under Unequal time lag value, the characteristic value obtaining based on the approximate characteristic value computational methods of different rank Pade, and maximum absolute error between the characteristic value obtaining based on discretization method is as shown in table 2.Be not difficult to find out, under given accuracy requires, along with the increase of time lag, the approximate rational polynominal exponent number k of Pade must correspondingly increase.From in table 2, for 4 machine 2 district systems, when τ changes in [0.05,0.5] scope, the approximate characteristic value that just can make the inventive method calculate of 6 rank Pade reaches the computational accuracy of 1e-5.

Under the approximate exponent number of table 2 Unequal time lag and Pade, the computational accuracy of characteristic value

In the time of Feedback Delays τ=0.30, the exact value of right characteristic vector corresponding to components of system as directed characteristic value that the characteristic value computational methods based on discretization obtain, as shown in table 3 2nd ~ 4 row.Can obtain three electromechanic oscillation modes and the corresponding mode of system by further model analysis, as shown in table 3 2nd ~ 4 row.Wherein, λ 1 shows as G1, the G2 vibration with respect to G3, G4; λ 2 shows as G1 vibration with respect to G4 with respect to G2, G3, but the mould value of corresponding G1, G2 modal components is larger; λ 3 shows as G1 vibration with respect to G4 with respect to G2, G3, but the mould value of corresponding G3, G4 modal components is larger.

Utilize in the approximate components of system as directed characteristic value calculating of 5 rank and 7 rank Pade and corresponding right characteristic vector and the corresponding component of each generator speed, as shown in 5th ~ 10 row of table 3.By with the contrast of the exact value of mode, known the present invention can accurately calculate amplitude and the phase place of Oscillatory mode shape.When k=5, the maximum absolute error of mode amplitude is 4e-8, and the maximum absolute error of phase place is 1e-2 °; When k=7, it is 1e-8 that the maximum of amplitude is spent absolutely error; The maximum absolute error of mode phase place is 9e-6 °.

When table 3 τ=0.3s, the mode that Pade approximate calculation obtains and exact value thereof

By reference to the accompanying drawings the specific embodiment of the present invention is described although above-mentioned; but not limiting the scope of the invention; one of ordinary skill in the art should be understood that; on the basis of technical scheme of the present invention, those skilled in the art do not need to pay various modifications that creative work can make or distortion still in protection scope of the present invention.

Claims (7)

1. calculate and a stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, it is characterized in that, concrete steps are as follows:

Step 1: specify wide area Feedback Delays τ _iwith the approximate exponent number k of Pade;

Step 2: carry out Pade and be similar to, obtain the approximate rational polynominal of Time Delay; In described step 2, in Laplace domain, i Time Delay

be expressed as

t is the time, and s is frequency; The value of i is positive integer; Pade is approximate is that a kind of [l, k] rank rational polynominal of utilizing is approached

method:

$e^{-{\mathrm{\&tau;}}_{i}s}\&ap;R\left(s\right)=\frac{N_l\left(s\right)}{N_k\left(s\right)}=\frac{a_0+a_1{\mathrm{\&tau;}}_{i}s+...+a_l{\left({\mathrm{\&tau;}}_{i}s\right)}^l}{b_0+b_1{\mathrm{\&tau;}}_{i}s+...+b_k{\left({\mathrm{\&tau;}}_{i}s\right)}^k}---\left(10\right)$

In formula, l and k are positive integer; Coefficient a _mand b _jobtained by following formula:

$a_m={(-1)}^m\frac{(l+k-m)!l!}{m!(l-m)!}---\left(11\right);$

$b_j=\frac{(l+k-j)!k!}{j!(k-j)!}---\left(12\right);$

Wherein, 0≤m≤l, 0≤j≤k;

Step 3: rational polynominal is converted into state-space expression;

Step 4: carry out equilibrating processing;

Step 5: by Time Delay be connected without time-lag power system, wide area damping control, obtain the inearized model of time-lag power system;

Step 6: the characteristic value of utilizing QR method or sparse eigenvalue Algorithm for Solving system;

Step 7: the directly time lag stability of judgement system.

2. as claimed in claim 1 a kind ofly calculate and stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, is characterized in that, described step 2 is in Pade is approximate, and exponent number l and k are larger, rational polynominal R (s) more close to

under normal circumstances, get l=k; Time lag is less, R in frequency domain (s) with

the interval that phase place is consistent is larger, and frequency band is wider.

3. as claimed in claim 1 a kind ofly calculate and stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, is characterized in that, in described step 3, by the molecule of formula (10) and denominator simultaneously divided by

and then utilize the general approach of transfer function problem of implementation, rational polynominal based on the approximate Time Delay obtaining of Pade is approached to R (s) and be converted to state space expression of equal value, for i Time Delay, represent take variable incremental form can control standard type as:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{.}{x}}_{\mathrm{di}}={\tilde{A}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{B}}_{\mathrm{di}}\mathrm{\&Delta;}{u}_i\\ \mathrm{\&Delta;}{y}_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}\mathrm{\&Delta;}{x}_{\mathrm{di}}+{\tilde{D}}_{\mathrm{di}}\mathrm{\&Delta;}{u}_i\end{array}\right.---\left(13\right)$

In formula (13), coefficient matrix

expression be:

4. as claimed in claim 1 a kind of based on the Pade approximate calculating of time-lag power system characteristic value and stability distinguishing method, it is characterized in that, described step 4, in the situation that counting compared with Small Time Lag τ i and higher-order, between the coefficient in the approximate transfer function R (s) obtaining of Pade, and the coefficient matrices A of the state-space expression being converted to by transfer function _diand C _dinonzero element between, on the order of magnitude, differ greatly, adopt linear transformation to make the coefficient matrix balance more of the state-space expression of Time Delay:

$A_{\mathrm{di}}=T^{-1}{\tilde{A}}_{\mathrm{di}}T,B_{\mathrm{di}}=T^{-1}{\tilde{B}}_{\mathrm{di}},C_{\mathrm{di}}={\tilde{C}}_{\mathrm{di}}T,D_{\mathrm{di}}={\tilde{D}}_{\mathrm{di}}---\left(15\right)$

In formula, T is diagonal transformation matrix; Hence one can see that, and equilibrating coefficient matrix after treatment still meets can control standard type.

5. as claimed in claim 1 a kind of based on the Pade approximate calculating of time-lag power system characteristic value and stability distinguishing method, it is characterized in that, described step 5: by the state-space expression of the Time Delay obtaining in step 3, be connected with the inearized model without time-lag power system, the inearized model of wide area damping control, thereby set up the inearized model of the Study of Power System Small Disturbance Stability analysis that comprises Time Delay:

$\left\{\begin{array}{c}\mathrm{\&Delta;}{\stackrel{.}{x}}^{\text{'}}=A^{\&prime;\&prime;}\mathrm{\&Delta;}{x}^{\text{'}}+B^{\text{'}}\mathrm{\&Delta;y}\\ 0=C^{\text{'}}\mathrm{\&Delta;}{x}^{\text{'}}+D^{\text{'}}\mathrm{\&Delta;y}\end{array}\right.---\left(16\right)$

" be system state variables, y is system algebraically variable, and A ", B ", C ", D " are coefficient matrix in formula, x;

The state variable x of the wide area damping control of supposing the state variable xd of i Time Delay and be attached thereto _cbe arranged in without after time-lag power system state variable x, $\mathrm{\&Delta;}{x}^{\text{'}}={\left[\begin{array}{ccccc}\mathrm{\&Delta;}{x}^T& \&CenterDot;\&CenterDot;\&CenterDot;& \mathrm{\&Delta;}{x}_{\mathrm{di}}^T& \mathrm{\&Delta;}{x}_{\mathrm{ci}}^T& \&CenterDot;\&CenterDot;\&CenterDot;\end{array}\right]}^T;$ If be A without the coefficient matrix of time-lag power system, B, C, D; The coefficient matrix of wide area damping control is A _c, B _c, C _c, D _c, the coefficient matrix of wide area feedback signal is respectively K _1iand K _2i; So coefficient matrices A in formula (16) ", B ", C ", D " are specifically expressed as:

In formula (17), coefficient matrices A _di, B _di, C _di, D _disatisfiedly can control standard type; A is a point Block diagonal matrix, and B, C are piecemeal sparse matrix, and D is 2 × 2 piecemeal sparse matrixes; E is sparse matrix, A _ci, B _ci, C _ci, D _cirelevant with the concrete structure of wide area damping control.

6. as claimed in claim 5 a kind ofly calculate and stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, is characterized in that, in described step 5 formula (17), each coefficient matrix has following features:

(1) " be still a point Block diagonal matrix if using the coefficient matrix of the dynamic element relating in wide area feedback signal, Time Delay, wide area damping control state variable as sub-block, an A; " only on Time Delay and wide area damping control state variable are expert at, increased null matrix or a small amount of several row nonzero element, B " is still Sparse Array to B;

(2) C " has only increased null matrix, D " with just the same without the corresponding coefficient matrix of time-lag power system in Time Delay and wide area damping control state variable column;

In summary, the lienarized equation coefficient matrices A of the electric power system that comprises Time Delay ", B ", C ", D ", with without time-lag power system at linearisation equation coefficient matrix A, B, C, D, there is identical sparsity structure.

7. as claimed in claim 1 a kind ofly calculate and stability distinguishing method based on the approximate time-lag power system characteristic value of Pade, is characterized in that described step 6: as D " when nonsingular, cancellation algebraically variable y, formula (16) is reduced to:

$\mathrm{\&Delta;}{\stackrel{.}{x}}^{\text{'}}={\tilde{A}}^{\text{'}}\mathrm{\&Delta;}{x}^{\text{'}}---\left(18\right)$

In formula,

for the state matrix based on the approximate time-lag system obtaining of Pade;

By calculating

characteristic value, must arrive the partial feature value of time-lag power system; In addition, because the electric power system that comprises Time Delay has identical sparse characteristic with conventional without the coefficient matrix of time-lag power system, the sparse treatment technology and the computational methods that in the time calculating on a large scale without time-lag power system part critical eigenvalue, use, stand good in the partial feature value of calculating time-lag power system.

Images