CN102664397B - Electric power system transient stability simulation method based on implicit fine numerical integral - Google Patents

Electric power system transient stability simulation method based on implicit fine numerical integral Download PDF

Info

Publication number
CN102664397B
CN102664397B CN201210079453.9A CN201210079453A CN102664397B CN 102664397 B CN102664397 B CN 102664397B CN 201210079453 A CN201210079453 A CN 201210079453A CN 102664397 B CN102664397 B CN 102664397B
Authority
CN
China
Prior art keywords
state
generator
transient
value
constantly
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201210079453.9A
Other languages
Chinese (zh)
Other versions
CN102664397A (en
Inventor
王建全
赵志奇
吕晓祥
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University ZJU
Original Assignee
Zhejiang University ZJU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University ZJU filed Critical Zhejiang University ZJU
Priority to CN201210079453.9A priority Critical patent/CN102664397B/en
Publication of CN102664397A publication Critical patent/CN102664397A/en
Application granted granted Critical
Publication of CN102664397B publication Critical patent/CN102664397B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Abstract

The invention discloses an electric power system transient stability simulation method based on a fine integral algorithm. Compared to an existing electric power system transient stability value integral method, in the method of the invention, non-linear differential equations describing an electric power system transient process are expressed as a linear part and a nonlinear part. Through reasonably selecting a system matrix of the linear part and using an integrable function to approximate a nonlinear term function, an implicit fine single step integral formula is derived. Through finely solving a state transition matrix corresponding to a linear part system matrix, the integral formula with high precision can be obtained so as to reduce a calculated amount of the transient stability simulation.

Description

A kind of electric power system transient stability emulation mode based on implicit expression Precise Numerical Integration
Technical field
The invention belongs to power system automation technology field, particularly a kind of electric power system transient stability emulation mode based on implicit expression Precise Numerical Integration.
Background technology
Transient stability analysis of power system is one of content basic, the most most crucial during power system analysis is calculated.Along with the development of power system operation control technology, the advanced technologies such as online dynamic security analysis, safety and stability emergency control, prevention and control, intelligent scheduling are progressively applied in electric power system.The precondition that realizes these advanced technologies be can calculate large-scale electrical power system transient stability accomplish fast, accurately, reliable.
Numerical integrating and direct method, and be conventional main method in the middle of transient stability analysis of power system by the hybrid analysis that numerical integration and direct method combine.In the middle of these methods, numerical integration is that electric power system transient stability calculates the most accurately, the most reliable method.The disadvantage of numerical integrating is that amount of calculation is large.Although the operational speed of a computer had had raising than the past, for large-scale electrical power system, computational speed wants to meet the requirement that online dynamic security analysis, prevention and control, emergency control etc. are calculated, and remains far from being enough.
The various elements of electric power system can be showed by Mathematical Modeling, and its transient process can be described by the differential-Algebraic Equation set of following form
Figure 2012100794539100002DEST_PATH_IMAGE002
(1)
In formula,
Figure 2012100794539100002DEST_PATH_IMAGE004
the state variable that represents descriptive system dynamic characteristic in differential equation group, conventionally vector
Figure 401760DEST_PATH_IMAGE004
the quantity of state that comprises the dynamic elements such as generator's power and angle and rotating speed;
Figure 2012100794539100002DEST_PATH_IMAGE006
the operation variable of system in representation algebra equation group, comprises the operation variable relevant to network, conventionally as the amplitude of node voltage and phase place etc.
The general flow that solves transient state process of electric power system with numerical integrating as shown in Figure 1.Its core procedure is the utilization of frame shown in 8.
Figure 2012100794539100002DEST_PATH_IMAGE008
,
Figure DEST_PATH_IMAGE010
solve the solution that the represented differential-Algebraic Equation set of (1) formula obtains next integration step
Figure DEST_PATH_IMAGE012
with
Figure DEST_PATH_IMAGE014
.At present, the common method that solves differential equation group in (1) formula in electric power system numerical simulation field has implicit expression trapezoidal integration, improved Euler method, Runge-Kutta method etc.Implicit expression trapezoidal integration better numerical value stability, but need repeatedly iterative, amount of calculation is large, and the electric power system business calculation procedure of this integration method adopting at present has PSD-BPA, PSASP.Runge-Kutta method, the higher order Taylor method of development is typical explicit integral, although without iteration, but in order to reach certain computational accuracy, also need solve repeatedly differential-Algebraic Equation set in each integration step and guarantee its precision, and its numerical stability is poor, easily causes calculating unsuccessfully.Therefore, generally progressive failure will, according to the truncated error of algorithm, by selecting rational integration step, guarantee convergence.
In order to guarantee stability and the simulation accuracy of algorithm, the integration step of getting to be inversely proportional to the truncated error of algorithm, the truncated error of numerical integration algorithm is less, under same precision requires, integration step
Figure DEST_PATH_IMAGE016
can obtain larger, otherwise integration step
Figure 194267DEST_PATH_IMAGE016
obtain smaller.Conventionally the truncated error of each integration step is less, and amount of calculation is also larger.As the local truncation error of Euler method is
Figure DEST_PATH_IMAGE018
, each integration step only needs to calculate a differential algebraic equations; The local truncation error of improved Euler method is
Figure DEST_PATH_IMAGE020
, each integration step needs to calculate twice differential algebraic equations; The local truncation error of the explicit Runge-Kutta method of quadravalence is
Figure DEST_PATH_IMAGE022
, each integration step needs to calculate four subdifferential algebraic equations.And the local truncation error of implicit expression trapezoidal integration is
Figure 659490DEST_PATH_IMAGE020
, need, through iterative differential-algebraic equation repeatedly, just can be met the solution of required precision.If can not increase the amount of calculation of algorithm when improving algorithm truncated error, can reduce the amount of calculation of whole transient emulation, accelerate computational speed.
At present, the numerical integration method adopting in electric power system transient stability numerical integration emulation mode, all directly adopt the general-purpose algorithm in computational methods theory, as implicit expression trapezoidal integration, improved Euler method, Runge-Kutta method etc., according to the feature of describing the differential equation of transient state process of electric power system, algorithm is not improved.When alternately solving differential algebraic equations, normally obtain whole state variables, ask for operation variable again.The feature of numerical integrating in the present invention based on generator model, by describing the nonlinear differential equation group of transient state process of electric power system, is expressed as linear segment and non-linear partial.By reasonably choosing the sytem matrix of linear segment and approaching nonlinear terms function by an integrable function, derived the meticulous single step integral formula of implicit expression.This algorithm solves the corresponding state-transition matrix of linear segment sytem matrix by becoming more meticulous, make integral formula have higher precision, thereby has reduced the amount of calculation of transient stability emulation.
Summary of the invention
The present invention seeks in order to solve in electric power system transient stability simulation calculation, existing numerical integration method amount of calculation is large, computational speed can not meet the shortcoming of the online calculation requirement of electric power system, Doha Mil's integration and precision integration based on nonlinear equation, proposed a kind of new transient stability numerical integration emulation mode.The method takes full advantage of generator model and other equipment can be expressed as linearity and non-linear partial also can be write as the feature of transfer function, has derived the single step integral formula based on implicit expression Precise integration method.This integration method amount of calculation is less than local truncation error
Figure 117016DEST_PATH_IMAGE020
implicit expression trapezoidal method; The method also takes full advantage of generator amature angle equation can, with the feature of other system element decoupling zero, improve the versatility of algorithm.
The present invention seeks to be achieved through the following technical solutions: a kind of electric power system transient stability emulation mode based on implicit expression precise integration, comprises the following steps:
Step 1: the initial parameter of input system and information, carry out trend and calculate the operation variate-value under steady state condition
Figure 2012100794539100002DEST_PATH_IMAGE024
, comprise the voltage of each generator node
Figure 2012100794539100002DEST_PATH_IMAGE026
, the electric current of injection network
Figure DEST_PATH_IMAGE028
, generator electromagnetic power
Figure DEST_PATH_IMAGE030
, wherein
Figure DEST_PATH_IMAGE032
Figure DEST_PATH_IMAGE034
(
Figure 322868DEST_PATH_IMAGE034
for generator number of units).
Step 2: computing mode variable initial value comprises generator's power and angle
Figure DEST_PATH_IMAGE036
, angular frequency
Figure DEST_PATH_IMAGE038
transient state and time transient potential, excitation and each Dynamic mode state variable initial value of governing system with generator.
Step 3: form the differential equation and the network algebra equation of descriptive system transient process, and carry out network algebra equation factor table and decompose.
Step 4: put transient stability and calculate initial value constantly
Figure DEST_PATH_IMAGE040
.
Step 5: judged whether that fault or operation occur.If nothing, turns to step 8; If have, perform step 6.
Step 6: according to fault or operational circumstances, revise the differential equation and network algebra equation and factor table thereof.
Step 7: solve network algebra equation, obtain
Figure DEST_PATH_IMAGE042
operation variable constantly.
Step 8: calculate
Figure DEST_PATH_IMAGE044
the state variable value of system constantly comprises each generator's power and angle
Figure DEST_PATH_IMAGE046
, angular frequency
Figure DEST_PATH_IMAGE048
with transient state and time transient potential, excitation and each Dynamic mode state variable value of governing system of generator, and operation variate-value comprises the voltage of generator node
Figure DEST_PATH_IMAGE050
, the electric current of injection network
Figure DEST_PATH_IMAGE052
and electromagnetic power
Figure DEST_PATH_IMAGE054
, this step detailed process is as follows:
Step 8.1: the step-length that judges this step
Figure 706182DEST_PATH_IMAGE016
whether identical with previous moment step-length, if identical, jump to step 8.3, if different, carry out the operation of step 8.2.
Step 8.2: obtain state matrix according to the differential equation
Figure DEST_PATH_IMAGE056
and nonlinear terms
Figure DEST_PATH_IMAGE058
linearized expression, utilize step-length
Figure 661631DEST_PATH_IMAGE016
, calculate state-transition matrix
Figure DEST_PATH_IMAGE060
, constant term state-transition matrix after the linearisation of nonlinear terms state
Figure DEST_PATH_IMAGE062
with once state-transition matrix after nonlinear terms linearisation
Figure DEST_PATH_IMAGE064
.
Step 8.3: given
Figure DEST_PATH_IMAGE066
state variable value constantly and operation variate-value, put iterations
Figure DEST_PATH_IMAGE068
.
Step 8.4: estimate
Figure DEST_PATH_IMAGE070
state variable value constantly.
Step 8.5: according to network algebra equation and
Figure 676510DEST_PATH_IMAGE070
state variable value calculates operation variate-value.
Step 8.6: try to achieve according to calculating
Figure 407706DEST_PATH_IMAGE070
running status variate-value, by following formula, solve
Figure 964851DEST_PATH_IMAGE070
state variable value constantly
Figure DEST_PATH_IMAGE072
Wherein, in formula
Figure DEST_PATH_IMAGE074
represent
Figure 843815DEST_PATH_IMAGE070
the generator's power and angle of each generator constantly
Figure 643143DEST_PATH_IMAGE046
, angular frequency , the state vector subvector that forms of each Dynamic mode state variable of transient state and inferior transient potential, excitation and governing system,
Figure DEST_PATH_IMAGE076
for the generator's power and angle of each generator of the moment
Figure DEST_PATH_IMAGE078
, angular frequency
Figure DEST_PATH_IMAGE080
, the vector that forms of the constant term after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system,
Figure DEST_PATH_IMAGE082
for
Figure 211900DEST_PATH_IMAGE066
with
Figure DEST_PATH_IMAGE084
the generator's power and angle of each generator of the moment , angular frequency
Figure DEST_PATH_IMAGE088
, the vector that forms of the linear segment after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system.
Step 8.7: check the maximum electromagnetic power deviate of each generator of iteration twice, if deviation is greater than given accuracy
Figure DEST_PATH_IMAGE090
, order
Figure DEST_PATH_IMAGE092
, return to step 8.5 and continue iteration; Otherwise, perform step 9.
Step 9: judge that whether system is stable, whether the maximal phase of any two generators is greater than a certain set-point to waving merit angle, if perform step 12; Otherwise, perform step 10.
Step 10: simulation time is advanced to a step-length, order
Figure DEST_PATH_IMAGE094
, and with now
Figure DEST_PATH_IMAGE096
value is as next zequin value constantly, value.
Step 11: judge whether to arrive given simulation time in advance
Figure DEST_PATH_IMAGE098
.If
Figure DEST_PATH_IMAGE100
perform step 12, otherwise return to step 5.
Step 12: output result of calculation also finishes to calculate.
In transient stability simulation calculation step 8.2, calculate state-transition matrix
Figure 685487DEST_PATH_IMAGE060
state-transition matrix with constant term after nonlinear terms linearisation with the once state-transition matrix of item after nonlinear terms linearisation
Figure DEST_PATH_IMAGE102
the method of using is Precise integration method.
Beneficial effect of the present invention: the method takes full advantage of the feature of precision integration, has derived Precise integration method single step integral formula, and the local truncation error of this integral formula is low, and amount of calculation is less than local truncation error and is
Figure 372744DEST_PATH_IMAGE020
implicit expression trapezoidal integration.
Accompanying drawing explanation
Fig. 1 is the general flow figure of transient stability numerical solution;
Fig. 2 is the calculation process that transient stability calculates each integration step;
Fig. 3 is IEEE39 node system maximum work angle swing curve;
Fig. 4 is implicit expression trapezoidal integration maximum work angular error calculation;
Fig. 5 is implicit expression Precise integration method maximum work angular error calculation.
Embodiment
Below in conjunction with accompanying drawing, the invention will be further described.
The present invention is based on Doha Mil's integration of nonlinear equation, proposed a kind of new high-precision numerical value emulation method.The difference of the inventive method and traditional numerical integration emulation mode is, by describing the nonlinear differential equation group of transient state process of electric power system, is expressed as linearity and non-linear partial.By the state-transition matrix that solves linear segment that becomes more meticulous with reasonably choose the integrable function of approaching non-linear partial, obtain the Numerical Integral Formulas of degree of precision.This method is better for the adaptability of model, the new established model of interpolation that can aspect, and do not need loaded down with trivial details step.
Below in conjunction with accompanying drawing 1, the invention will be further described.
A kind of electric power system transient stability emulation mode based on implicit expression Precise Numerical Integration that the present invention proposes, comprises the following steps:
Step 1: the initial parameter of input system and information, carry out trend and calculate the operation variate-value under steady state condition
Figure 513875DEST_PATH_IMAGE024
, comprise the voltage of each generator node
Figure 74170DEST_PATH_IMAGE026
, the electric current of injection network
Figure 845816DEST_PATH_IMAGE028
, electromagnetic power
Figure 858772DEST_PATH_IMAGE030
, wherein
Figure 305059DEST_PATH_IMAGE032
(
Figure 724725DEST_PATH_IMAGE034
for generator number of units).
Step 2: computing mode variable initial value comprises each generator's power and angle
Figure 162659DEST_PATH_IMAGE036
, angular frequency
Figure 645593DEST_PATH_IMAGE038
transient state and time transient potential, excitation and each Dynamic mode state variable initial value of governing system with generator.
Step 3: form the differential equation and the network algebra equation of descriptive system transient process, and form factor table.
Step 4: put transient stability and calculate initial value constantly
Figure 684874DEST_PATH_IMAGE040
.
Step 5: judged whether that fault or operation occur.If nothing, turns to step 8; If have, perform step 6.
Step 6: according to fault or operational circumstances, revise the differential equation and network algebra equation and factor table thereof.
Step 7: solve network algebra equation, obtain
Figure 860641DEST_PATH_IMAGE042
operation variable constantly.
Step 8: calculate
Figure 848188DEST_PATH_IMAGE044
the state variable value of system constantly comprises each generator's power and angle
Figure 134813DEST_PATH_IMAGE046
, angular frequency
Figure 196310DEST_PATH_IMAGE048
with transient state and time transient potential, excitation and each Dynamic mode state variable value of governing system of generator, and operation variate-value comprises the voltage of generator node
Figure 44443DEST_PATH_IMAGE050
, the electric current of injection network
Figure 519286DEST_PATH_IMAGE052
and electromagnetic power
Figure 609602DEST_PATH_IMAGE054
, this step detailed process is as follows:
Step 8.1: the step-length that judges this step
Figure 853502DEST_PATH_IMAGE016
whether identical with previous moment step-length, if identical, jump to step 8.3, if different, carry out the operation of step 8.2.
Step 8.2: obtain state matrix according to the differential equation
Figure 43175DEST_PATH_IMAGE056
and nonlinear terms linearized expression, utilize step-length
Figure 132277DEST_PATH_IMAGE016
, calculate state-transition matrix
Figure 496262DEST_PATH_IMAGE060
, constant term state-transition matrix after the linearisation of nonlinear terms state
Figure 919153DEST_PATH_IMAGE062
with once state-transition matrix after nonlinear terms linearisation .
Computing mode transfer matrix
Figure 239856DEST_PATH_IMAGE060
, constant term state-transition matrix after the linearisation of nonlinear terms state with once state-transition matrix after nonlinear terms linearisation
Figure 52140DEST_PATH_IMAGE064
method used is precision integration.
Step 8.3: given
Figure 723292DEST_PATH_IMAGE066
state variable value constantly and operation variate-value, put iterations
Figure 475215DEST_PATH_IMAGE068
.
Step 8.4: estimate
Figure 485896DEST_PATH_IMAGE044
state variable value constantly.
Step 8.5: according to network algebra equation and
Figure 516169DEST_PATH_IMAGE044
state variable value calculates operation variate-value.
According to
Figure DEST_PATH_IMAGE104
solve, wherein
Figure DEST_PATH_IMAGE106
for electric power networks admittance matrix, solving virtual Injection Current first
Figure DEST_PATH_IMAGE108
, obtain time etching system operation variable
Figure DEST_PATH_IMAGE110
, and then obtain the electromagnetic power of each generator .
Step 8.6 is tried to achieve according to calculating
Figure 199271DEST_PATH_IMAGE044
running status variate-value, by following formula, solve
Figure 625311DEST_PATH_IMAGE044
each generator state variables value constantly
Wherein, in formula
Figure 409914DEST_PATH_IMAGE074
represent
Figure 721946DEST_PATH_IMAGE070
the generator's power and angle of each generator constantly
Figure 769537DEST_PATH_IMAGE046
, angular frequency , the state vector subvector that forms of each Dynamic mode state variable of transient state and inferior transient potential, excitation and governing system,
Figure 776118DEST_PATH_IMAGE076
for
Figure 829525DEST_PATH_IMAGE066
the generator's power and angle of each generator of the moment
Figure 731622DEST_PATH_IMAGE078
, angular frequency
Figure 274598DEST_PATH_IMAGE080
, the vector that forms of the constant term after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system,
Figure 133751DEST_PATH_IMAGE082
for
Figure 53165DEST_PATH_IMAGE066
with
Figure 747452DEST_PATH_IMAGE084
the generator's power and angle of each generator of the moment
Figure 461330DEST_PATH_IMAGE086
, angular frequency
Figure 568963DEST_PATH_IMAGE088
, the vector that forms of an item parts after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system.
Step 8.7: check the maximum electromagnetic power deviate of each generator of iteration twice, if deviation is greater than given accuracy
Figure 527954DEST_PATH_IMAGE090
, order
Figure 139064DEST_PATH_IMAGE092
, return to step 8.5 and continue iteration; Otherwise, perform step 9.
Step 9: judge that whether system is stable, whether the maximal phase of any two generators is greater than a certain set-point to waving merit angle, if perform step 12; Otherwise, perform step 10.
Step 10: simulation time is advanced to a step-length, order
Figure 227106DEST_PATH_IMAGE094
, and with now
Figure 556456DEST_PATH_IMAGE096
value as next zequin value constantly
Figure 817673DEST_PATH_IMAGE066
.
Step 11: judge whether to arrive given simulation time in advance
Figure 781825DEST_PATH_IMAGE098
.If
Figure 103085DEST_PATH_IMAGE100
perform step 12, otherwise return to step 5.
Step 12: output result of calculation also finishes to calculate.
Below introduce in detail the detailed process of the inventive method.
Differential equation group (1) mainly comprises the differential equation of describing generating set and other dynamic apparatus dynamic characteristic, and wherein the differential equation of each generating set can be expressed as:
Figure DEST_PATH_IMAGE114
(2)
In formula, the generator's power and angle, angular frequency, prime mover mechanical output, electromagnetic power and the inertia time constant that represent respectively each generator, for system synchronization electric angle speed,
Figure DEST_PATH_IMAGE120
represent Generator Damping coefficient.
Figure DEST_PATH_IMAGE122
for the transient state of each generator and the state vector subvector of each Dynamic mode state variable composition of time transient potential, excitation and governing system,
Figure DEST_PATH_IMAGE124
for with state vector subvector
Figure 732780DEST_PATH_IMAGE122
the functional vector on corresponding differential equation right side.Like this, the state vector of each generating set
Figure DEST_PATH_IMAGE126
can be expressed as:
Figure DEST_PATH_IMAGE128
.
Each generator differential equation group (2) just can further be expressed as:
Figure DEST_PATH_IMAGE130
(3)
In formula,
Figure DEST_PATH_IMAGE132
for the sytem matrix of nonlinear differential equation group linear segment,
Figure DEST_PATH_IMAGE134
for the transient state of generator and the sytem matrix of time transient potential, excitation and each Dynamic mode differential equation group linear segment of governing system;
Figure DEST_PATH_IMAGE136
for Nonlinear System of Equations, the functional vector of non-linear partial.
The sytem matrix of nonlinear differential equation group linear segment
Figure 775517DEST_PATH_IMAGE056
central submatrix
Figure 361219DEST_PATH_IMAGE134
selection rule be: according to the state vector of each generating set
Figure 853381DEST_PATH_IMAGE126
and differential equation group, extract the central generating set state vector of equation group the coefficient of the linear term of each corresponding element, is write as matrix form side by side.
To nonlinear function vector
Figure DEST_PATH_IMAGE138
, with an approximate linear function, represent:
(4)
Wherein for the constant term after linearisation,
Figure DEST_PATH_IMAGE144
once item for after linearisation, is specifically expressed as:
Figure DEST_PATH_IMAGE146
,
Figure DEST_PATH_IMAGE148
,
Wherein, for
Figure 449819DEST_PATH_IMAGE066
the generator's power and angle of each generator of the moment
Figure 847302DEST_PATH_IMAGE078
, angular frequency
Figure 638541DEST_PATH_IMAGE080
, the vector that forms of the constant term after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system,
Figure 284066DEST_PATH_IMAGE082
for with
Figure 350428DEST_PATH_IMAGE084
the generator's power and angle of each generator of the moment
Figure 628963DEST_PATH_IMAGE086
, angular frequency
Figure 573785DEST_PATH_IMAGE088
, the vector that forms of an item parts after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system.
Figure DEST_PATH_IMAGE150
for
Figure 552368DEST_PATH_IMAGE066
the vector that constant term after the transient state of generator and time transient potential, excitation and each Dynamic mode nonlinear terms linearisation of governing system forms constantly,
Figure DEST_PATH_IMAGE152
for
Figure 557233DEST_PATH_IMAGE070
the vector that constant term after the transient state of generator and time transient potential, excitation and each Dynamic mode nonlinear terms linearisation of governing system forms constantly,
Figure DEST_PATH_IMAGE154
for generator inertia time constant,
Figure 618336DEST_PATH_IMAGE120
represent each Generator Damping coefficient,
Figure DEST_PATH_IMAGE156
for generator prime machine mechanical output,
Figure 101270DEST_PATH_IMAGE118
for system synchronization electric angle speed.
How the differential equation for above foundation solves, and is core of the present invention place.For the differential equation , precision integration of the present invention is mainly passed through following formula:
Figure DEST_PATH_IMAGE160
try to achieve.After substitution linearisation
Figure DEST_PATH_IMAGE162
further arrangement can obtain , in formula represent the constant term after linearisation,
Figure DEST_PATH_IMAGE168
represent the once item after linearisation.
After being calculated, generator differential equation substitution precise integration can obtain following formula:
Figure 449206DEST_PATH_IMAGE072
Wherein, in formula represent
Figure 863054DEST_PATH_IMAGE070
the generator's power and angle of each generator constantly
Figure 87362DEST_PATH_IMAGE046
, angular frequency , the state vector subvector that forms of each Dynamic mode state variable of transient state and inferior transient potential, excitation and governing system, for
Figure 799732DEST_PATH_IMAGE066
the generator's power and angle of each generator of the moment
Figure 624468DEST_PATH_IMAGE078
, angular frequency , the vector that forms of the constant term after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system, for
Figure 285760DEST_PATH_IMAGE066
with
Figure 412722DEST_PATH_IMAGE084
the generator's power and angle of each generator of the moment
Figure 511128DEST_PATH_IMAGE086
, angular frequency
Figure 199598DEST_PATH_IMAGE088
, the vector that forms of an item parts after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system.Here
Figure DEST_PATH_IMAGE170
in the corresponding differential equation
Figure DEST_PATH_IMAGE172
,
Figure DEST_PATH_IMAGE174
corresponding
Figure DEST_PATH_IMAGE176
,
Figure DEST_PATH_IMAGE178
corresponding
Figure DEST_PATH_IMAGE180
For this reason, according to the present invention, the calculation procedure of each integration step of a kind of electric power system transient stability emulation mode based on implicit expression precise integration is as follows:
1: the step-length that judges this step
Figure 966478DEST_PATH_IMAGE016
whether identical with previous moment step-length, if identical, jump to 3, if different, carry out the operation of step 2.
2: according to the differential equation, obtain state matrix, utilize step-length after step change
Figure 133017DEST_PATH_IMAGE016
, calculate state-transition matrix
Figure 351509DEST_PATH_IMAGE060
, constant term state-transition matrix after nonlinear terms linearisation
Figure 210881DEST_PATH_IMAGE062
with linear term state-transition matrix after nonlinear terms linearisation
Figure 383498DEST_PATH_IMAGE064
.
3: given state variable value constantly and operation variate-value, put iterations
Figure 692306DEST_PATH_IMAGE068
.
4: estimate
Figure 660262DEST_PATH_IMAGE044
state variable value constantly.
5: according to network algebra equation and
Figure 818711DEST_PATH_IMAGE044
state variable value calculates operation variate-value.
6: according to calculating, try to achieve
Figure 825588DEST_PATH_IMAGE044
running status variate-value, by formula (5), (6) and (7), solve
Figure 18671DEST_PATH_IMAGE044
state variable value constantly.
7: check the maximum electromagnetic power deviate of each generator of iteration twice, if deviation is greater than given accuracy
Figure 219846DEST_PATH_IMAGE090
, order
Figure 803274DEST_PATH_IMAGE092
, return to step 5 and continue iteration; Otherwise this integration step iterative process finishes;
Its calculation process as shown in Figure 2.
Be below an embodiment of the inventive method, with IEEE39 node system, carry out emulation experiment and make embodiment, further illustrate as follows:
The generator model that simulation example adopts is that generator third-order model is variation model, generator auxiliary device is field regulator, load adopts constant-impedance model.For there is three phase short circuit fault in 0 moment on the top of the circuit between No. 34 circuits and node 28 and node 29 in fault, fault is removed when 0.1s.Whole trouble duration is 0.1s.The whole time span of emulation is 3.0s, and system emulation result is for unstable phenomenon occurs.
Fig. 3 be IEEE39 node system maximal phase to merit angle swing curve, step-length is 0.01s, adopts implicit expression trapezoidal integration algorithm, its result is as normal data value.Fig. 4 is implicit expression trapezoidal integration maximum work angular error calculation, and step-length is respectively 0.02s and 0.05s, and Fig. 5 is implicit expression Precise integration method maximum work angular error calculation, and step-length is respectively 0.02s and 0.05s.
By Fig. 4 and Fig. 5 as seen when integration step is got 0.02-0.05s maximum work angle error of the present invention be all no more than 10 degree.And error precision of the present invention will be much smaller than implicit expression trapezoidal integration algorithm under large step-length.And in amount of calculation, the number of times that the present invention solves network algebra equation when step-length is got 0.05s is 539 times, and implicit expression trapezoidal integration solves the number of times of network algebra equation, it is 577 times.The inventive method has been saved the amount of calculation of about 6-7%.

Claims (2)

1. the electric power system transient stability emulation mode based on implicit expression Precise Numerical Integration, is characterized in that the method comprises the steps:
Step 1: the initial parameter of input system and information, carry out trend calculating, obtain the operation variate-value y (0) under steady state condition, comprise the voltage V of each generator node i(0), the electric current I of injection network i(0), electromagnetic power , i=1 wherein, 2 ... N g, N gfor generator number of units;
Step 2: computing mode variable initial value, comprises generator's power and angle δ i(0), angular frequency iand the transient state of generator and time transient potential, excitation and each Dynamic mode state variable initial value of governing system (0);
Step 3: form the differential equation and the network algebra equation of descriptive system transient process, and carry out network algebra equation factor table and decompose;
Step 4: put transient stability and calculate initial value t=0 constantly;
Step 5: judged whether that fault or operation occur; If nothing, turns to step 8; If have, perform step 6;
Step 6: according to fault or operational circumstances, revise the differential equation and network algebra equation and factor table thereof;
Step 7: solve network algebra equation, obtain t 0operation variable constantly;
Step 8: calculate t 0the state variable value of+h system constantly comprises generator's power and angle δ i(t 0+ h), angular frequency i(t 0and operation variate-value comprises the voltage V of generator node+h) and the transient state of generator and time transient potential, excitation and each Dynamic mode state variable value of governing system, i(t 0+ h), the electric current I of injection network i(t 0+ h) and electromagnetic power
Figure FDA0000403558120000012
detailed process is as follows:
Step 8.1: whether the step-length h that judges this step is identical with previous moment step-length, if identical, jumps to step 8.3, if different, carries out the operation of step 8.2;
Step 8.2: the linearized expression that obtains state matrix H and nonlinear terms F (t) according to the differential equation, utilize step-length h, calculate the Tb (h) of the constant term state-transition matrix after state-transition matrix Ta (h), the linearisation of nonlinear terms state and the T1b (h) of the linear term state-transition matrix after nonlinear terms linearisation;
Step 8.3: given t 0state variable value constantly and operation variate-value, put iterations m=0;
Step 8.4: estimate t 0+ h state variable value constantly;
Step 8.5: according to network algebra equation and t 0+ h state variable value calculates operation variate-value;
Step 8.6: the t trying to achieve according to calculating 0the running status variate-value of+h, solves t by following formula 0+ h state variable value constantly
X i (m+1)(t 0+h)=Ta i(h)X i(t 0)+Tb i(h)b i(t 0)+T1b i(h)b1 i (m+1)(t 0+h),
Wherein, the X in formula i (m+1)(t 0+ h) represent t 0the generator's power and angle δ of+h each generator constantly i(t 0+ h), angular frequency i(t 0+ h), the state vector subvector that forms of each Dynamic mode state variable of transient state and inferior transient potential, excitation and governing system, b i(t 0) be t 0the generator's power and angle δ of each generator of the moment i(t 0), angular frequency i(t 0), the vector that forms of the constant term after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system, b1 i (m+1)(t 0+ h) be t 0and t 0+ h is the generator's power and angle δ of each generator constantly i, angular frequency i, the vector that forms of an item parts after transient state and each Dynamic mode nonlinear terms linearisation of inferior transient potential, excitation and governing system;
Step 8.7: check the maximum electromagnetic power deviate of each generator of iteration twice, if deviation is greater than given accuracy ε, make m=m+1, return to step 8.5 and continue iteration; Otherwise, perform step 9;
Step 9: judge that whether system is stable, whether the maximal phase of any two generators is greater than a certain set-point to waving merit angle, if perform step 12; Otherwise, perform step 10;
Step 10: simulation time is advanced to a step-length, make t=t 0+ h, and the t value of usining is now as next zequin value constantly, i.e. t 0value;
Step 11: judge whether to arrive given simulation time T in advance; If t >=T performs step 12, otherwise return to step 5;
Step 12: output result of calculation also finishes to calculate.
2. a kind of electric power system transient stability emulation mode based on implicit expression Precise Numerical Integration according to claim 1, is characterized in that: the method for calculating the state-transition matrix Tb (h) of the constant term after state-transition matrix Ta (h) and nonlinear terms linearisation and the state-transition matrix Tb1 (h) of the linear term after nonlinear terms linearisation use in transient stability calculation procedure 8.2 is Precise integration method.
CN201210079453.9A 2012-03-23 2012-03-23 Electric power system transient stability simulation method based on implicit fine numerical integral Expired - Fee Related CN102664397B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201210079453.9A CN102664397B (en) 2012-03-23 2012-03-23 Electric power system transient stability simulation method based on implicit fine numerical integral

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210079453.9A CN102664397B (en) 2012-03-23 2012-03-23 Electric power system transient stability simulation method based on implicit fine numerical integral

Publications (2)

Publication Number Publication Date
CN102664397A CN102664397A (en) 2012-09-12
CN102664397B true CN102664397B (en) 2014-01-29

Family

ID=46773839

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210079453.9A Expired - Fee Related CN102664397B (en) 2012-03-23 2012-03-23 Electric power system transient stability simulation method based on implicit fine numerical integral

Country Status (1)

Country Link
CN (1) CN102664397B (en)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105468864B (en) * 2015-12-14 2018-05-15 三峡大学 Based on the ultra-high-tension power transmission line electro-magnetic transient numerical computation method for increasing dimension precise integration
CN105867360B (en) * 2016-06-14 2018-05-08 江南大学 A kind of initial value of Mechatronic control system estimates iterative learning fault diagnosis algorithm
CN107706948B (en) * 2017-11-03 2020-12-01 华北电力大学 Multi-dimensional order controlled multi-step Taylor series transient stability analysis method
CN110119523A (en) * 2019-01-29 2019-08-13 中国电力科学研究院有限公司 A kind of network algebra equation solution preprocess method and system based on the generator rotor angle method of equal effect
CN110135031A (en) * 2019-04-30 2019-08-16 东南大学 Electric power system transient stability calculation method based on half implicit runge kutta method
CN111404825A (en) * 2020-03-13 2020-07-10 重庆第二师范学院 Data transmission method, system, computer device and storage medium
CN113221298B (en) * 2021-04-21 2023-02-24 南方电网科学研究院有限责任公司 Method and system for simulating electromechanical transient process

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101699448A (en) * 2009-10-26 2010-04-28 清华大学 Transient stability distributed simulation method of electric power system

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4616206B2 (en) * 2006-04-14 2011-01-19 株式会社日立製作所 Power system stability determination method and apparatus

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101699448A (en) * 2009-10-26 2010-04-28 清华大学 Transient stability distributed simulation method of electric power system

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
JP特开2007-288878A 2007.11.01
基于隐式Taylor级数法的电力系统暂态稳定计算;王宇宾等;《华北电力大学学报》;20050330;第32卷(第02期);第1-6页 *
王宇宾等.基于隐式Taylor级数法的电力系统暂态稳定计算.《华北电力大学学报》.2005,第32卷(第02期),第1-6页.

Also Published As

Publication number Publication date
CN102664397A (en) 2012-09-12

Similar Documents

Publication Publication Date Title
CN102664397B (en) Electric power system transient stability simulation method based on implicit fine numerical integral
CN102609575B (en) Power system transient stability simulating method based on implicit numerical integration
Yao et al. Efficient and robust dynamic simulation of power systems with holomorphic embedding
CN102545263B (en) Power system transient stability simulation method based on explicit numerical integration
CN104156542B (en) It is a kind of based on the active distribution system Simulation of stability method implicitly projected
CN106099922B (en) Transient Voltage Stability in Electric Power System judgment method based on Heuristic Energy Function method
CN106451576B (en) A kind of control method of the electric power electric transformer of single-phase multi output
CN102611380B (en) Online identification method for parameters of double-fed motor
CN103700036A (en) Transient stability projection integral method suitable for multi-time scale of electrical power system
CN104375876B (en) Electromagnetical transient emulation method is immunized in a kind of 0+ errors under input quantity catastrophe
CN103810646A (en) Improved projection integral algorithm based active power distribution system dynamic simulation method
WO2018102720A1 (en) System and method for a fast power network simulator
CN110875600A (en) Approximate analysis model for dynamic frequency response of two-machine equivalent power system
CN112380670A (en) Virtual rotor-based modeling method and system for sectional power supply linear induction motor
Ajala et al. A library of second-order models for synchronous machines
CN105512367A (en) Method for determining critical stable value of grid connected thermal power generating unit primary frequency regulation rotating speed diversity factor
CN109428340A (en) A kind of emulation mode and system of flexible DC transmission device
CN102609576B (en) Power system transient stability simulating method based on estimation-correction numerical integration
Döşoğlu Decoupled power-based sliding mode control modeling enhancement for dynamic stability in doubly-fed induction generator-based wind turbines
CN103117692A (en) Control method of mechanical elastic energy storing permanent magnet motor group under various external disturbances
CN102354332B (en) Method for simplifying relative gain matrix (RGA) calculation in flexible alternating-current/direct-current electricity transmission system
CN107203139A (en) The stabilized control method of multiple-input and multiple-output non-linear differential algebraic subsystem
CN105140957B (en) Electromechanic oscillation mode evaluation method based on wind power plant and photovoltaic plant polymerization model
CN105305429A (en) Three-phase reclosing sequence resetting method for improving transient frequency stability of two-unit electric power system
Germanos Power System Stability Response and Control Using Small Signal Analysis

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20140129

Termination date: 20210323