- Email: [email protected]
Limit cycle analysis of a hydroelectric system: A new approach
Electric Power Systems Research, 11 (1986) 49 - 58
49
Limit Cycle Analysis of a Hydroelectric System: a New Approach
H. C. CHANG*, C. T. PAN t, C. C. WEI* and C. L. HUANG*
*Department of Electrical Engineering, National Cheng Kung University, Tainan, 700 (Taiwan) tDepartment of Electrical Engineering, National Tsing Hua University, Hsinchu, 300 (Taiwan) (Received July 14, 1986)
SUMMARY
In this paper, a systematic method that predicts limit cycle oscillations in a hydroelectric system with governor deadband on the parameter plane is presented. A new technique for constructing limit cycle loci on this plane is also proposed. This approach is based on the characteristic equation as well as the describing function approximation and avoids the use of derivation and separation o f stability equations. For each family o f solution loci in the parameter plane, only two integra. tions are needed. The stability o f the limit cycle can also be determined easily from the plane using the gradient approach. Finally, the proposed method is employed to investigate the influence that several variations in the modeling o f a hydro system may have on the existence o f limit cycles. The results from this study will shed some light on the selection o f automatic generation control (AGC) parameters as a remedy for this phenomenon.
I. INTRODUCTION
It has been shown [ 1 - 4 ] that speedgovernor deadband tends to produce continuing oscillations, that is, limit cycles, and hence has a destabilizing effect on the dynamic performance of automatic generation control (AGC). Since the existence of the limit cycle may incur unnecessary control action and degrade inadvertent interchange and time error control [5], it is most urgent that this phenomenon is analyzed thoroughly so that some strategy may be used to avoid or reduce the influence of this effect. 0378-7796/86/$3.50
Over recent decades, papers have become available [1- 9] which deal with limit cycle analysis of power systems. In early works [1, 2], the exact analog simulation technique was used to arrive at the conditions for the existence of limit cycles. This procedure works well for specific cases but is awkward for parametric studies. Besides, time-domain simulations done by analog, digital or hybrid methods are expensive. Afterwards, many studies [5- 8] were carried out to analyze limit cycles in power systems using the describing function technique together with a Nyquist plot. For example, Wu and Dea [8] successfully used this method to study limit cycles on a classic single-area power system. Although this method may serve as an analytical tool, it cannot provide a sufficiently complete picture about the influence of system parameters on the limit cycles. Very recent research by Liaw et al. [9] presents a new approach based on the stability equation method [10] and the parameter plane method to predict the existence of limit cycles. However, derivation and separation of the stability equations are, though possible, laborious and time consuming, and hence this restricts its applications. In this paper, a new technique based on the describing function approximation and characteristic equations is applied to study the existence of limit cycles on the parameter plane. From the family of characteristic curves in this plane, the loci of the limit cycles and the influence which changes in the values of the parameters has on them can be clearly seen. The stability of the limit cycle is checked by a gradient approach with almost no additional computations. Significantly, the same approach can also be extended to multi© Elsevier Sequoia/Printed in The Netherlands
50
area systems as well as control systems with multiple nonlinearities. Finally, some remedies to control limit cycling are also suggested.
Pref
RATE
PILOT VQLUE LIMIT
[P,r-H_~_
GQTE
POSITION
SE R VO M O TO R
LIMIT
PEP~IQMEHT
~
SPEEDDROOP
2. SYSTEM MODEL INVESTIGATED DQSHPOT
The transfer function block diagram of a hydroelectric power system with governor backlash (deadband) is shown in Fig. 1. This system model under investigation is the same as that used in ref. 5 with the exception of the hydro governor model. The parameter values within typical ranges used in this system model correspond to those of refs. 5, 11 and 12 and are included in Table 1.
Fig. 2. Detailed representation hydraulic speed-governing system.
detailed transfer function of the hydro governor can be closely approximated by the simpler one: (1 +
TRS)/R
TR TuS2/o + [Tg + TR(O
QGC
BQCKLQSN GOVERNOR HYDRO TURB[tie
PO~ER SYSTEM
J
Fig. 1. Transfer function block diagram of a hydroelectric system.
TABLE 1 Typical values and ranges of various constants Parameters
Typical
Range
TR Tg Tw
5.0 0.2 1.0 0.04 0.31 0.05 0.25 0.75 10.6 T r = 5.0 T w
2.5 - 25.0 0.2-0.4 0.5 - 5.0 0.03 - 0.06 0.2 - 1.0 0.05 - 0.167 0.05 - 0.5 0.75 - 1.50 5.3 - 10.6
o 6 R KI D M
Since there are significant differences in the literature concerning the modeling of the hydro governor, the use of a specific approximate model must be explained. Figure 2 shows the detailed transfer function model [11, 12] of a hydro governor. By eliminating nonlinearities such as rate limits, position limits and the small time constant of the pivot value Tp (typically less than 0.05 s), the
of a mechanical
+
6)]S/a
+
(1)
1
Note that if the permanent speed regulation R is expressed in p.u. Hz/p.u. M W instead of Hz/p.u. M W , the symbol R [5, 9] corresponds to o [11, 12]. Hence, R and a are used interchangeably throughout this paper. A comparison of gains and phase margins of the detailed approximate governor model and t h a t u s e d b y P a n t a t o n e a n d P i e g z a [ 5 ] is s h o w n in F i g . 3. I t is o b v i o u s t h a t t h e d e t a i l e d and approximate models show rather good agreement over a range of frequency from 0 . 0 1 r a d s -1 ( o r 6 2 8 . 3 s) t o 0 . 5 t a d s -1 ( o r 1 2 . 6 s), w h i c h e n c o m p a s s e s t h e r a n g e o f frequency of the limit cycles (between 30 s a n d 9 0 s) o b s e r v e d b y E w a r t [ 1 3 ] . H o w e v e r , the governor model used by Pantalone shows less g o o d a g r e e m e n t w i t h t h e d e t a i l e d m o d e l . 40
:6O 1"0
30
120 lOO 80 ~>-
2o z
60
B 10
~o
}
20 ~ 0 -20 -lO OOl
002
0.05
Ol
02
0.5
lO
2o
5.0
10.
£REQt~NCY IN RADIANS/SECONDS
Fig. 3. B o d e diagram o f various hydro governor transfer functions: - - , detailed m o d e l ; - - - - - , approximate m o d e l ; . . . . , Pantalone's m o d e l .
51
Therefore, to simplify the sophisticated hydro governor model, it is more reasonable to use the approximate model.
or, equivalently,
FR(A, co, ~, ~) = 0
(4)
FI (A, co, ~, H) = 0 3. LIMIT CYCLE PREDICTION BY USING THE PARAMETER PLANE METHOD
In the past, the parameter plane method has been generalized and applied to quantitative and qualitative analysis of a wide variety of nonlinear phenomena. The most attractive merit of this m e t h o d is that it provides information about the effects on system overall behavior on changing the operating conditions and parameters. As such, the parameter plane method is applied in this paper for limit cycle analysis of a power system. Now, consider the hydroelectric system with governor deadband (DB) as shown in Fig. 1. As limit cycling is of primary concern in this study, the signal prior to the governor backlash is assumed sinusoidal, that is, X = A sin(cot). The nonlinear element of governor backlash is then characterized by the sinusoidal input describing function (SIDF) [14] as
N(A) =
2j
f Y[A sin(cot)] e-J°~t d(cot)
1 [(2~b- ¢2)1/2(1 --~b) + lr _
¢) + j(¢2 _
2¢)1
H = H(A, ¢o)
(5)
Since eqn. (4) is nonlinear in the unknowns and H, the explicit solution in the form of eqn. (5) may n o t be possible. To overcome this difficulty, a more general and efficient alternative is to construct a family of constant A curves and constant co curves on the parameter plane. From the intersecting points of both curves, the limit cycle can easily be predicted. A new technique for constructing limit cycle loci on this plane is presented in the following.
3.1. Construction of constant A curves
FR(w, a(w), H(w)) = o
ff cos-l(1
Ol(A, co)
o~ =
Each constant A curve is constructed by fixing the A value and considering co as an independent variable. Therefore, eqn. (4) can be written explicitly as
0
_
where FR and F1 are the real and imaginary parts of the characteristic equation respectively. If A and co are considered as independent variables and ~ and /3 as two unknowns, then eqn. (4) may be solved for ~ and H as
(2)
where
F,(w,
(6)
=0
Differentiating eqn. (6) with respect to co yields
¢ = DB/A A typical value of 0.06% for governor backlash is considered in this paper. In this case, the condition for limit cycling with AP d = 0 is the system characteristic equation equal to zero, that is,
F(S)=I +N(A)GGGTGv(I+
~)=O
(3)
To analyze the limit cycle on a specific ~-~ parameter plane, substitution of S = jco into eqn. (3) yields
F(A, co, ~, [3) = FR(A, co, ~, [3) + jR, (A, co, o~, [3) = 0
dFR
~F R d~ -
dw dFi dco
~a
()F R d~ +
dco
OFt d a -
0o~ dco
0F R
0H d w
+
()FI dH +
0H dco
0co 8F I
+"
0co
-0
(7) =0
Equation (7) can be rearranged as
8co [J]
=_
Ld-gwJ
(8)
J
where the Jacobian matrix
52
[J] = ~ F I
3F[
J of
3.2. C o n s t r u c t i o n
c o n s t a n t co c u r v e s
Considering CO as a constant and A an independent variable, eqn. (4) can be expressed as F R ( A , a(A), r(A)) = 0
(9)
F I ( A , c~(A ), r ( A )) = 0
Differentiating eqn. (9) with respect to A yields dR R -
dA
~F R d a 3o~
3FR dr
+
dA
3r
+
dA
3FR
-0
Now, we reduce the problem of limit cycle prediction to solving ordinary differential equations (ODEs) such as eqn. {12). If, to a given pair (A0, COo), there corresponds a pair (a0, r0) such that the characteristic equation F ( S ) = 0 is satisfied, then the constant A0 and constant COo curves can be constructed using an ODE solver. However, as stated before, the explicit form given by eqn. (5) may n o t be possible; but the sufficient conditions for the existence of a solution can be verified by the parameter mapping theorem [15]. In this paper, if the explicit form is not available, the linear h o m o t o p y m e t h o d is applied to solve for (a0, G0) with given (A0, co0) during the initialization process. For completeness, the h o m o t o p y m e t h o d is briefly reviewed here.
3A
(lO) dFi
3FI d a
~FI dr
-
dA
+
3FI +
~a dA
3r
dA
4. EVALUATION OF THE INITIAL CONDITIONS USING THE HOMOTOPY METHOD
=0
~A
Similarly, eqn. (10) can be rewritten as
=--
~A
(11)
[J] [-~-~]
L~AA
It is obvious that a sufficient condition for the inverse transformation of eqn. (5) is that the Jacobian [J] of eqns. (8) and (11) does not vanish. In general, in studying the inverse transformation, it is intentional to avoid points at which the Jacobian of the transformation is equal to zero. Thus, eqns. (8) and (11) can be written in the form
It has been shown that eqn. (3) is the condition for the existence of a limit cycle. It follows that the problem of finding the limit cycles is reduced to solving the complex nonlinear algebraic equation. Since eqn. (3) is rather nonlinear and a good initial guess is not available, the c o m m o n l y used N e w t o n Raphson m e t h o d generally fails to find the solution. Therefore, a recently developed h o m o t o p y m e t h o d which basically possesses the global convergence characteristics is applied to find the initial pair (a0, rio) with given (A0, COo). During 1960, Davidenko [16] produced a new idea to solve the equation F: R n ~ R n,
/
/
i Ai
J
L -J (12)
It should be pointed out that during the actual solution process it is n o t necessary to separate eqn. (3) into the real form as eqn. (4). Each term on the right-hand side of eqn. (12) can be evaluated directly from the characteristic equation.
F(X) = 0
(13)
First, we can start with an easy equation whose solution is trivial, for example, V(X) = X--a
(14)
= 0
Then, we can construct a h o m o t o p y H, for example, H:R n× [0,1]-+R n H ( X , t) = (1 -- t ) ( X - - a) + t F ( X )
(15)
to ' h o m o t o p y ' G ( X ) = 0 to F ( X ) = 0. It is obvious from eqn. (15) that t = O,
H ( X , O) = X - - a
t = 1,
H(X,
(16) 1) =
F(X)
53
Suppose that X ( t ) exists such that H ( X ( t ) , t) = 0 for all 0 ~< t ~< 1. Then the solution of eqn. (13) can be found as follows. Differentiate the following identity with respect to t:
FR(~, fl, A, co) = 0
H ( X , t) = (1 -- t ) ( X - - a ) + t F ( X ) = 0
F,(cz, fl, A, co) = 0
we obtain dH aH d X dt
-
aX dt
aH +
--
at
-- 0
(17)
(18)
which implies \ ax ]
aoL
at
(19) X(0)
= a
Hence, we can obtain the solution by using any O D E solver to integrate eqn. (19) up to t = 1. Unfortunately, there is a serious problem in Davidenko's method, since the assumption X =X(t) m a y not be true. To overcome this difficulty, Li [16] used the homotopy H ( X ( s ) , t(s)) = (1 -- t ) ( X - - a ) + t F ( X ) = 0
(20) where s is the parameter representing the arc length. It has also been proved theoretically that for a randomly chosen initial guess value the proposed algorithm will eventually lead to the solution with probability one. As a result, the following equation can be used to find the solution:
aH
a~
The detailed numerical technique of their method can be found from ref. 16 and will be n o t repeated here.
air
aco
aw
o~
as
a~
~F~
aA
aA
a/~
as
a/3
< 0
(23)
where the arguments of FR and F1 are omitted for simplicity. Note that the first determinant of eqn. (23) has been evaluated in the construction of the constant A and constant co curves as shown in eqn. (12). The second determinant can be either calculated by the algebraic method or evaluated by inspecting the gradient directions of the characteristic curves at the points where they cross each other. Now, we define a third coordinate % which is perpendicular to both c~ and ~, by Uv = Us × U~
(24)
where Us, U# and Uv are the unit vectors along the ~, fl and ~ axis, respectively. Then, the second determinant can be c o m p u t e d from
acoacoI
all" ~ --
(22)
If there is a limit cycle at (s0, rio) in the q-fl plane, with frequency coo and amplitude A0, then the limit cycle is stable if and only if
a~
dX dt
T h e o r e m 1. If F ( S ) is the linearized characteristic equation with parameters ~ and fl, the stability equation can be written as
Z~_
as aA
aA
= U~.(Vco XVA)
(25)
as where aco •
aco
aA VA= ~Us+
aA ~U~ op
V co = 5. S T A B I L I T Y A N A L Y S I S O F T H E LIMIT C Y C L E S
The stability properties of limit cycles can be understood by inspecting the gradient directions of the constant A and constant co curves in the parameter plane. The theoretical background is based on the following theorem [17].
(26)
Since only the sign of the determinant is of interest, it suffices to determine the stability of the limit cycle by inspecting the directions of VA and Vco. Therefore, almost no addi-
54
tional calculations are needed in determining the stability properties. 6. E F F E C T S O F P A R A M E T E R THE LIMIT CYCLES
VARIATIONS
ON
To investigate the influence of parameter variations on the existence of the limit cycle and gain a deeper insight into its characteristics, different parameter planes showing limit cycling have been constructed by the present method. Although most of the results are consistent with those concluded by Wu and Pantalone [5, 8], such studies can provide qualitative and quantitative analyses of the limit cycles. This study will also help to clarify an apparent ambiguity in these works [5, 8]. It should be noted that only stable limit cycle loci are plotted on the plane. Also, unless stated otherwise, the periods up to 90 s which were observed in the real power system [13] are considered. Each plot o f the parameter plane is conducted b y varying parameters examined within a typical working range and keeping all other parameters fixed at the typical values shown in Table 1. In all plots the locus of broken lines represents constant period (or constant frequency) curves while solid lines represent constant amplitude curves. Since the permanent speed droop regulation R and the AGC integral controller gain Ki are t w o important feedback gains for power system operation and control, special emphasis is placed on the analysis of these parameters. Figure 4 shows the effect of I(I 0.60
~
-
%
~ '.o-~o~j'
t I 001I
o;::[ i/,,,?/
varying K~ and R on the limit cycle characteristics. It is significant that there exists an asymptotically stable region, which was not reported in the literature and may be meaningful to power engineers. In this region, a proper set of (KI,R) can be chosen to construct a stable system without limit cycles. For example, keeping KI at a typical value of 0.25, the limit cycling condition begins to vanish when R is increased to about 0.043. Thus, when the range of K1 is small (say K~ <~ 0.36), decreasing R will make a limit cycle more likely to happen. This appears to be in agreement with the conclusions reached by Wu and Pantalone [ 5 , 8 ] . However, for larger values of KI (say KI/> 0.36) the effect of varying R is negligible. Similarly, the effect of varying K~ can be seen by keeping R fixed. It is obvious that for larger values of R (say R ~> 0.024), increasing KI will also make a limit cycle more likely to happen. On the other hand, for small values of R (say R ~< 0.024), the effect of R dominates that of KI. In this region, a larger value of K1 will result in a larger amplitude and frequency of the limit cycle. Since the describing function is just an approximate one, it is necessary to confirm the prediction of a limit cycle by digital simulations. Table 2 shows the results of a comparison between the predicted values and steady state values resulting from digital simulations on some of the sample points (P1 - P8) shown in Fig. 4. From Table 2, it is shown that the predicted values agree rather well with simulation results from an engineering standpoint. Some simulated trajectories corresponding to P3 in Fig. 4 are depicted in Fig. 5. Next, the gradient approach based on Theorem 1 is applied to determine the stability of the limit cycle. For example, in the limit cycle region (refer to Q in Fig. 4), one has AL = U~.(V00 × VA) > 0 and
o, til?tl/
.........................
0"000.00
0.08
0.02 0.0#+ 0.06
010 0.12 0.1#+ 0.16 0.10 R Fig. 4. L i m i t c y c l e loci o n t h e KI-R p a r a m e t e r p l a n e for DB = 0.06%.
[J[ = - - 1 . 4 2 0 < 0 Therefore, the limit cycle is stable. In fact, all the limit cycles plotted in this plane were tested to be stable. The effect of decreasing governor deadband can be shown b y comparing Figs. 4 and 6. The frequencies of the limit cycles are almost
55 TABLE 2 Comparison between the calculated and simulated results Sample point
Value o f parameters
Amplitude and period
Calculated results
Simulated results
Error (%)
P1
K I = 0.3626 R = 0.06955 K I = 0.3805 R = 0.05477 K I = 0.3981 R = 0.03544 K I = 0.4795 R = 0.02430 K I = 0.5348 R = 0.1003 K I = 0.5196 R = 0.1489 K! = 0.5491 R = 0.1507 K I = 0.5774 R = 0.1642
A T A T A T A T A T A T A T A T
0.0003800 60 0.0004000 52.0 0.0004200 46.0 0.0004400 42.0 0.0004600 44.0 0.0004800 44.0 0.0005000 42.0 0.0005300 40.0
0.0003783 58.0 0.0003950 51.7 0.0004148 45.5 0.0004343 41.5 0.0004391 45.6 0.0004534 46.2 0.0004715 44.0 0.0004991 41.8
0.45 3.33 1.25 0.58 1.24 1.09 1.30 1.19 4.54 3.64 5.54 5.00 5.70 4.76 5.83 4.50
P2 P3 P4 P5 P6 P7 P8
(X10-3)
(XlO-3)
S-O
5.0 4-0
4.0 A=0.0004148
Tffi45.5 ApdffiO.O005
3.0
3-0
X 2.0
Y 2.0
1.0
i-0
0.0
0.0
-I.0
(a)
-I.0 0
30
60
90
120
150 180 T(SEc)
210
240
270
300
(b)
30
60
90
120
1S0 180 T(S~C)
210
240
270
300
(XlO-3) 1.0
(XIO-L 0.4
0.8
0.0
0.6
/
-0.4 Y 0.4
AF -0.8
.//
/
/
/
0.2 -1-2 ,
0.0
-1.6
(c)
30
60
90
120
150 180 T(SEc)
210
240
270
300
-0.2 - - / -0.6 -0.4
(d)
~.0.0006 .,~:-- . . . . . . . . . . . -0.2
0.0
0.2
0.4
0,6
0.8
1.0
1.2
1.4 (XIO"3)
x
Fig. 5. Results of computer simulations at limit cycle P3 for a 0.0005 p.u. load increase: trajectories o f (a) the signal at input of governor backlash; (b) the signal at ouput of governor backlash; (c) the frequency deviation in p.u.; (d) governor backlash (deadband) nonlinearity.
unchanged. However, the amplitudes of l i m i t c y c l e s are d i r e c t l y p r o p o r t i o n a l t o magnitude of the deadband. This suggests the straightforward way of relieving
the the that the
d e t r i m e n t a l e f f e c t s c a u s e d b y l i m i t c y c l i n g is t o m a k e t h e d e a d b a n d as s m a l l as p o s s i b l e . Unfortunately, only modest success has been a c h i e v e d o v e r t h e last d e c a d e s [ 1 8 ] .
56
time constant TR, and the governor response time Tg. From Fig. 8, it is shown that decreasing TR tends to enhance limit cycling, which is similar to the effect of decreasing 6. On the other hand, decreasing Tg will decrease the likelihood of limit cycling. However, the effect of varying Tg is rather small.
060 0.5~ 0.48
0.~2 0.36 Kx 0.30 0.2~
6.0
0.1~
5L,
0.12
¢.8
0.06 0.00
¢2
0.00
002
00~,
0.06
0.08
0.10 R
0.12
0.1L,
0.16
0.18
3.6
020
Fig. 6. Limit cycle loci on the KI-R parameter plane for DB = 0.04%•
TR
3.0 2z~ 18 12
Another aspect examined was the influence of the hydro governor model. As shown in Fig. 2, the temporary speed droop 5 and the permanent speed droop R are two possible adjustable parameters in the hydro governor. Variation of the parameters R and 5 yield the results of Fig. 7. It is seen that there also exists an asymptotically stable region. As 5 is increased, the limit cycling condition eventually vanishes. The effect of increasing will result in smaller magnitude and lower frequency of the limit cycle, which is totally contrary to the effect of increasing K~ shown in Fig. 4. Since the effect of varying R has been discussed previously, it is omitted here. Two other time constants to be examined in the detailed governor model are the dashpot
0.6 0.0
,
,
000
[
] CONSTANT
(
} CONSTANT
,
0.05
AMPLITUDE PERIOD ,
0.10
O1S
020
0.25 030 Tg
,
CURVE CU~E
,
035
,
O.t,O
A
,
0.t, 5
050
Fig. 8. Limit cycle loci on the TR-Tgparameter plane for DB = 0.06%.
Next, the influence of the hydro turbine will be considered. Figure 9 shows the effect of varying the water starting time Tw on the KI-Tw parameter plane. As Tw is increased, it tends to enhance the limit cycling. Since the value of TR is closely related with that of Tw in practical applications, typically TR/ Tw = 5.0, it is more reasonable to consider TR and Tw at the same time. From previous
0.70 [
063
(
]
CONSTANT AMPLITUDE CURVE ) CONSTANT PERIOD C U R V E
0-1
.~.~ O.L~9 ~ , Off6
..........................
t 0.1L~
~
0.12 ,"
0.07 0.00
/ /
Fig. for
v /
//
: / / / / / /~ /
/
l
I
I
'
0.000 7.
2,,\',:i,1
0.00 ~
"
Limit
cycle
DB = 0.06%.
loci
on
the
~-R
parameter
plane
'
L ,
O0 0.2 Ot~ 05 08
0.016 0.032 O.OA,8 O.06Z+ 0.080 0.096 0.112 0.128 0.11,1, 0.160 R Fig. for
9. Limit
cycle
DB = 0.06%.
loci
on
~
,
1.0 1.2 Tw the
,
,
•
,
1.¢ 15
K I - T w parameter
.
i ....... 18 2.0 plane
57
discussions, it can be concluded that, as the ratio of TR/T. is decreased, limit cycling is more likely to happen. Finally, the influence o f varying the parameters of the power system model, that is, the damping torque coefficient D and the effective rotary inertia of area M, will be explored. It has been reported that the limit cycling condition can be suppressed by either a decrease in M or an increase in D [5]. In contrast, Wu and Dea [8] concluded that the effect of decreasing M will enhance the limit cycling condition. Such contradictory reports can be clarified by referring to Fig. 10. For larger values of D (say D I> 0.4), a decrease in M will indeed enhance the limit cycling as observed by Wu and Dea. However, for smaller values o f D (say D ~< 0.4), the effect is opposite. Therefore, the effect of varying M not only is determined by itself but also depends on the value of D. Similarly, the effect of varying D cannot solely be determined by itself. For larger values of M (say M/> 3.6), an increase in D can indeed suppress the limit cycling, which is consistent with the conclusions of Pantalone and Piegza. However, for small values of M (say M ~< 3.6) the conclusion is completely contrary.
~
1.2
[ ] CONSTANT AMPLITUDE CURVE ( } CONSTANT PERIOD CURVE
\ II
1.0
~it%--
ASYMPTOTICALLY STABLE REGION
D0.8 0.6 O.h 0.2 0.0
0.O
. . . . . . 1.2 2.6
3.6
4.8
6.0 H
7.2
8.6.
9.6
10.8
REFERENCES 1 C. Concordia, L. K. Kirchmayer and E. A. Szymanski, Effect of speed governor deadband on tie-line power and frequency control performance, A I E E Trans., 76 (1957) 429 - 435. 2 L. K. Kirchmayer, Economic Control of Interconnected Systems, Wiley, N e w York, 1959. 3 S. C. Tripathy, G. S. Hope and O. P. Malik, Optimisation of load-frequency control parameters for power systems with reheat steam turbines and governor deadband nonlinearities, Proc. Inst. Electr. Eng., Part C, 129 (1982) 10 15. 4 S. C. Tripathy, T. S. Bhatti, C. S. Jha, O. P. Malik and G. S. Hope, Sampled data automatic generation control analysis with reheat steam turbines and governor dead-band effects, I E E E Trans., PAS-103 ( 1 9 8 4 ) 1 0 4 5 - 1 0 5 1 .
1.6 1.~*
presented. A new technique for constructing limit cycle loci on this plane is also proposed. This new technique is based completely on the characteristic equation and avoids derivation and separation of stability equations. Therefore, the same approach can also be extended to multiarea systems as well as control systems with multiple nonlinearities. Furthermore, the gradient approach was employed to investigate the stability of the limit cycles. Although most of the results are consistent with previous works, such studies can give deeper insight into limit cycle phenomena. In particular, from the quantitative analysis of limit cycles some apparently conflicting conclusions of previous works can be clarified and reconciled. Results from this study also shed some light on the selection of system parameters. It is possible to choose proper sets of parameters such that the limit cycling condition can be eliminated.
12.0
Fig. 10. L i m i t cycle loci o n t h e D - M p a r a m e t e r p l a n e for DB = 0.06%.
5 D. K. P a n t a l o n e a n d D. M. Piegza, L i m i t cycle analysis o f h y d r o e l e c t r i c systems, IEEE Trans., PAS-IO0 ( 1 9 8 1 ) 6 2 9 - 638. 6 J. M. Bailey a n d G. F. Pierce, E f f e c t s o f n o n linearities o n g o v e r n o r p e r f o r m a n c e at B U L L R U N s t e a m p l a n t , Proc. 7th Annu. South-eastern
Syrup. 7. C O N C L U S I O N S
In this paper, a systematic method that predicts limit cycle oscillations in a hydroelectric system on the parameter plane is
on
System
Theory,
Alabama,
1975,
pp. 137 - 138. 7 J. M. Bailey a n d G. F. Pierce, Backlash a n d r a t e s a t u r a t i o n effects o n g o v e r n o r p e r f o r m a n c e at t h e B U L L R U N s t e a m p l a n t , Proc. 1EEE South-
eastern Region 3 Conf. on Engineering in a Changing Economy, South Carolina, 1976, pp. 10
-
12.
58 8 F. F. Wu and V. S. Dea, Describing function analysis of automatic generation control system with governor deadband, Electr. Power Syst. Res., 1 (1977/78) 113 - 116. 9 C. M. Liaw, C. T. Pan and K. W. Han, Limit cycle analysis of power systems with governor deadband nonlinearities, in preparation. 10 C. H. Ai and K. W. Han, Stability analysis of a nuclear reactor control system with multiple transport lags and asymmetrical nonlinearities, IEEE Trans., NS-22 (1975) 2103 - 2112. 1 IEEE Committee Report, Dynamic models for steam and hydro turbines in power system studies, IEEE Trans., PAS-92 (1973) 1904 - 1916. 12 D. G. Ramey and J. W. Skooglund, Detailed hydrogovernor representation for system stability studies, IEEE Trans., PAS-89 (1970) 106 - 111. 13 D. Ewart, Automatic generation control: per-
formance under normal conditions, Proc. Engi-
neering Foundation Conf. on Systems Engineering for Power: Status and Prospects, Henniker, New Hampshire, August 17- 22, 1975, pp. 1 - 14. 14 A. Gelb and W. E. Van der Velde, MultipleInput Describing Function and Nonlinear S y s t e m Design, McGraw-Hill, New York, 1968. 15 D. D. Siljak, Nonlinear Systems, Wiley, New York, 1969. 16 T. Y. Li, Lectures on the numerical method of finding solutions of a system of nonlinear equations, Seminar o f Numerical Analysis, Chia~Tung Univ., Taiwan, 1983, Lecture 1 - 4. 17 Y. T. Tsay and K. W. Han, A gradient approach to the determination of stability of limit cycles, J. Franklin Inst., 300 (1975) 391 - 418. 18 C. Concordia, Discussion of ref. 8.
49
Limit Cycle Analysis of a Hydroelectric System: a New Approach
H. C. CHANG*, C. T. PAN t, C. C. WEI* and C. L. HUANG*
*Department of Electrical Engineering, National Cheng Kung University, Tainan, 700 (Taiwan) tDepartment of Electrical Engineering, National Tsing Hua University, Hsinchu, 300 (Taiwan) (Received July 14, 1986)
SUMMARY
In this paper, a systematic method that predicts limit cycle oscillations in a hydroelectric system with governor deadband on the parameter plane is presented. A new technique for constructing limit cycle loci on this plane is also proposed. This approach is based on the characteristic equation as well as the describing function approximation and avoids the use of derivation and separation o f stability equations. For each family o f solution loci in the parameter plane, only two integra. tions are needed. The stability o f the limit cycle can also be determined easily from the plane using the gradient approach. Finally, the proposed method is employed to investigate the influence that several variations in the modeling o f a hydro system may have on the existence o f limit cycles. The results from this study will shed some light on the selection o f automatic generation control (AGC) parameters as a remedy for this phenomenon.
I. INTRODUCTION
It has been shown [ 1 - 4 ] that speedgovernor deadband tends to produce continuing oscillations, that is, limit cycles, and hence has a destabilizing effect on the dynamic performance of automatic generation control (AGC). Since the existence of the limit cycle may incur unnecessary control action and degrade inadvertent interchange and time error control [5], it is most urgent that this phenomenon is analyzed thoroughly so that some strategy may be used to avoid or reduce the influence of this effect. 0378-7796/86/$3.50
Over recent decades, papers have become available [1- 9] which deal with limit cycle analysis of power systems. In early works [1, 2], the exact analog simulation technique was used to arrive at the conditions for the existence of limit cycles. This procedure works well for specific cases but is awkward for parametric studies. Besides, time-domain simulations done by analog, digital or hybrid methods are expensive. Afterwards, many studies [5- 8] were carried out to analyze limit cycles in power systems using the describing function technique together with a Nyquist plot. For example, Wu and Dea [8] successfully used this method to study limit cycles on a classic single-area power system. Although this method may serve as an analytical tool, it cannot provide a sufficiently complete picture about the influence of system parameters on the limit cycles. Very recent research by Liaw et al. [9] presents a new approach based on the stability equation method [10] and the parameter plane method to predict the existence of limit cycles. However, derivation and separation of the stability equations are, though possible, laborious and time consuming, and hence this restricts its applications. In this paper, a new technique based on the describing function approximation and characteristic equations is applied to study the existence of limit cycles on the parameter plane. From the family of characteristic curves in this plane, the loci of the limit cycles and the influence which changes in the values of the parameters has on them can be clearly seen. The stability of the limit cycle is checked by a gradient approach with almost no additional computations. Significantly, the same approach can also be extended to multi© Elsevier Sequoia/Printed in The Netherlands
50
area systems as well as control systems with multiple nonlinearities. Finally, some remedies to control limit cycling are also suggested.
Pref
RATE
PILOT VQLUE LIMIT
[P,r-H_~_
GQTE
POSITION
SE R VO M O TO R
LIMIT
PEP~IQMEHT
~
SPEEDDROOP
2. SYSTEM MODEL INVESTIGATED DQSHPOT
The transfer function block diagram of a hydroelectric power system with governor backlash (deadband) is shown in Fig. 1. This system model under investigation is the same as that used in ref. 5 with the exception of the hydro governor model. The parameter values within typical ranges used in this system model correspond to those of refs. 5, 11 and 12 and are included in Table 1.
Fig. 2. Detailed representation hydraulic speed-governing system.
detailed transfer function of the hydro governor can be closely approximated by the simpler one: (1 +
TRS)/R
TR TuS2/o + [Tg + TR(O
QGC
BQCKLQSN GOVERNOR HYDRO TURB[tie
PO~ER SYSTEM
J
Fig. 1. Transfer function block diagram of a hydroelectric system.
TABLE 1 Typical values and ranges of various constants Parameters
Typical
Range
TR Tg Tw
5.0 0.2 1.0 0.04 0.31 0.05 0.25 0.75 10.6 T r = 5.0 T w
2.5 - 25.0 0.2-0.4 0.5 - 5.0 0.03 - 0.06 0.2 - 1.0 0.05 - 0.167 0.05 - 0.5 0.75 - 1.50 5.3 - 10.6
o 6 R KI D M
Since there are significant differences in the literature concerning the modeling of the hydro governor, the use of a specific approximate model must be explained. Figure 2 shows the detailed transfer function model [11, 12] of a hydro governor. By eliminating nonlinearities such as rate limits, position limits and the small time constant of the pivot value Tp (typically less than 0.05 s), the
of a mechanical
+
6)]S/a
+
(1)
1
Note that if the permanent speed regulation R is expressed in p.u. Hz/p.u. M W instead of Hz/p.u. M W , the symbol R [5, 9] corresponds to o [11, 12]. Hence, R and a are used interchangeably throughout this paper. A comparison of gains and phase margins of the detailed approximate governor model and t h a t u s e d b y P a n t a t o n e a n d P i e g z a [ 5 ] is s h o w n in F i g . 3. I t is o b v i o u s t h a t t h e d e t a i l e d and approximate models show rather good agreement over a range of frequency from 0 . 0 1 r a d s -1 ( o r 6 2 8 . 3 s) t o 0 . 5 t a d s -1 ( o r 1 2 . 6 s), w h i c h e n c o m p a s s e s t h e r a n g e o f frequency of the limit cycles (between 30 s a n d 9 0 s) o b s e r v e d b y E w a r t [ 1 3 ] . H o w e v e r , the governor model used by Pantalone shows less g o o d a g r e e m e n t w i t h t h e d e t a i l e d m o d e l . 40
:6O 1"0
30
120 lOO 80 ~>-
2o z
60
B 10
~o
}
20 ~ 0 -20 -lO OOl
002
0.05
Ol
02
0.5
lO
2o
5.0
10.
£REQt~NCY IN RADIANS/SECONDS
Fig. 3. B o d e diagram o f various hydro governor transfer functions: - - , detailed m o d e l ; - - - - - , approximate m o d e l ; . . . . , Pantalone's m o d e l .
51
Therefore, to simplify the sophisticated hydro governor model, it is more reasonable to use the approximate model.
or, equivalently,
FR(A, co, ~, ~) = 0
(4)
FI (A, co, ~, H) = 0 3. LIMIT CYCLE PREDICTION BY USING THE PARAMETER PLANE METHOD
In the past, the parameter plane method has been generalized and applied to quantitative and qualitative analysis of a wide variety of nonlinear phenomena. The most attractive merit of this m e t h o d is that it provides information about the effects on system overall behavior on changing the operating conditions and parameters. As such, the parameter plane method is applied in this paper for limit cycle analysis of a power system. Now, consider the hydroelectric system with governor deadband (DB) as shown in Fig. 1. As limit cycling is of primary concern in this study, the signal prior to the governor backlash is assumed sinusoidal, that is, X = A sin(cot). The nonlinear element of governor backlash is then characterized by the sinusoidal input describing function (SIDF) [14] as
N(A) =
2j
f Y[A sin(cot)] e-J°~t d(cot)
1 [(2~b- ¢2)1/2(1 --~b) + lr _
¢) + j(¢2 _
2¢)1
H = H(A, ¢o)
(5)
Since eqn. (4) is nonlinear in the unknowns and H, the explicit solution in the form of eqn. (5) may n o t be possible. To overcome this difficulty, a more general and efficient alternative is to construct a family of constant A curves and constant co curves on the parameter plane. From the intersecting points of both curves, the limit cycle can easily be predicted. A new technique for constructing limit cycle loci on this plane is presented in the following.
3.1. Construction of constant A curves
FR(w, a(w), H(w)) = o
ff cos-l(1
Ol(A, co)
o~ =
Each constant A curve is constructed by fixing the A value and considering co as an independent variable. Therefore, eqn. (4) can be written explicitly as
0
_
where FR and F1 are the real and imaginary parts of the characteristic equation respectively. If A and co are considered as independent variables and ~ and /3 as two unknowns, then eqn. (4) may be solved for ~ and H as
(2)
where
F,(w,
(6)
=0
Differentiating eqn. (6) with respect to co yields
¢ = DB/A A typical value of 0.06% for governor backlash is considered in this paper. In this case, the condition for limit cycling with AP d = 0 is the system characteristic equation equal to zero, that is,
F(S)=I +N(A)GGGTGv(I+
~)=O
(3)
To analyze the limit cycle on a specific ~-~ parameter plane, substitution of S = jco into eqn. (3) yields
F(A, co, ~, [3) = FR(A, co, ~, [3) + jR, (A, co, o~, [3) = 0
dFR
~F R d~ -
dw dFi dco
~a
()F R d~ +
dco
OFt d a -
0o~ dco
0F R
0H d w
+
()FI dH +
0H dco
0co 8F I
+"
0co
-0
(7) =0
Equation (7) can be rearranged as
8co [J]
=_
Ld-gwJ
(8)
J
where the Jacobian matrix
52
[J] = ~ F I
3F[
J of
3.2. C o n s t r u c t i o n
c o n s t a n t co c u r v e s
Considering CO as a constant and A an independent variable, eqn. (4) can be expressed as F R ( A , a(A), r(A)) = 0
(9)
F I ( A , c~(A ), r ( A )) = 0
Differentiating eqn. (9) with respect to A yields dR R -
dA
~F R d a 3o~
3FR dr
+
dA
3r
+
dA
3FR
-0
Now, we reduce the problem of limit cycle prediction to solving ordinary differential equations (ODEs) such as eqn. {12). If, to a given pair (A0, COo), there corresponds a pair (a0, r0) such that the characteristic equation F ( S ) = 0 is satisfied, then the constant A0 and constant COo curves can be constructed using an ODE solver. However, as stated before, the explicit form given by eqn. (5) may n o t be possible; but the sufficient conditions for the existence of a solution can be verified by the parameter mapping theorem [15]. In this paper, if the explicit form is not available, the linear h o m o t o p y m e t h o d is applied to solve for (a0, G0) with given (A0, co0) during the initialization process. For completeness, the h o m o t o p y m e t h o d is briefly reviewed here.
3A
(lO) dFi
3FI d a
~FI dr
-
dA
+
3FI +
~a dA
3r
dA
4. EVALUATION OF THE INITIAL CONDITIONS USING THE HOMOTOPY METHOD
=0
~A
Similarly, eqn. (10) can be rewritten as
=--
~A
(11)
[J] [-~-~]
L~AA
It is obvious that a sufficient condition for the inverse transformation of eqn. (5) is that the Jacobian [J] of eqns. (8) and (11) does not vanish. In general, in studying the inverse transformation, it is intentional to avoid points at which the Jacobian of the transformation is equal to zero. Thus, eqns. (8) and (11) can be written in the form
It has been shown that eqn. (3) is the condition for the existence of a limit cycle. It follows that the problem of finding the limit cycles is reduced to solving the complex nonlinear algebraic equation. Since eqn. (3) is rather nonlinear and a good initial guess is not available, the c o m m o n l y used N e w t o n Raphson m e t h o d generally fails to find the solution. Therefore, a recently developed h o m o t o p y m e t h o d which basically possesses the global convergence characteristics is applied to find the initial pair (a0, rio) with given (A0, COo). During 1960, Davidenko [16] produced a new idea to solve the equation F: R n ~ R n,
/
/
i Ai
J
L -J (12)
It should be pointed out that during the actual solution process it is n o t necessary to separate eqn. (3) into the real form as eqn. (4). Each term on the right-hand side of eqn. (12) can be evaluated directly from the characteristic equation.
F(X) = 0
(13)
First, we can start with an easy equation whose solution is trivial, for example, V(X) = X--a
(14)
= 0
Then, we can construct a h o m o t o p y H, for example, H:R n× [0,1]-+R n H ( X , t) = (1 -- t ) ( X - - a) + t F ( X )
(15)
to ' h o m o t o p y ' G ( X ) = 0 to F ( X ) = 0. It is obvious from eqn. (15) that t = O,
H ( X , O) = X - - a
t = 1,
H(X,
(16) 1) =
F(X)
53
Suppose that X ( t ) exists such that H ( X ( t ) , t) = 0 for all 0 ~< t ~< 1. Then the solution of eqn. (13) can be found as follows. Differentiate the following identity with respect to t:
FR(~, fl, A, co) = 0
H ( X , t) = (1 -- t ) ( X - - a ) + t F ( X ) = 0
F,(cz, fl, A, co) = 0
we obtain dH aH d X dt
-
aX dt
aH +
--
at
-- 0
(17)
(18)
which implies \ ax ]
aoL
at
(19) X(0)
= a
Hence, we can obtain the solution by using any O D E solver to integrate eqn. (19) up to t = 1. Unfortunately, there is a serious problem in Davidenko's method, since the assumption X =X(t) m a y not be true. To overcome this difficulty, Li [16] used the homotopy H ( X ( s ) , t(s)) = (1 -- t ) ( X - - a ) + t F ( X ) = 0
(20) where s is the parameter representing the arc length. It has also been proved theoretically that for a randomly chosen initial guess value the proposed algorithm will eventually lead to the solution with probability one. As a result, the following equation can be used to find the solution:
aH
a~
The detailed numerical technique of their method can be found from ref. 16 and will be n o t repeated here.
air
aco
aw
o~
as
a~
~F~
aA
aA
a/~
as
a/3
< 0
(23)
where the arguments of FR and F1 are omitted for simplicity. Note that the first determinant of eqn. (23) has been evaluated in the construction of the constant A and constant co curves as shown in eqn. (12). The second determinant can be either calculated by the algebraic method or evaluated by inspecting the gradient directions of the characteristic curves at the points where they cross each other. Now, we define a third coordinate % which is perpendicular to both c~ and ~, by Uv = Us × U~
(24)
where Us, U# and Uv are the unit vectors along the ~, fl and ~ axis, respectively. Then, the second determinant can be c o m p u t e d from
acoacoI
all" ~ --
(22)
If there is a limit cycle at (s0, rio) in the q-fl plane, with frequency coo and amplitude A0, then the limit cycle is stable if and only if
a~
dX dt
T h e o r e m 1. If F ( S ) is the linearized characteristic equation with parameters ~ and fl, the stability equation can be written as
Z~_
as aA
aA
= U~.(Vco XVA)
(25)
as where aco •
aco
aA VA= ~Us+
aA ~U~ op
V co = 5. S T A B I L I T Y A N A L Y S I S O F T H E LIMIT C Y C L E S
The stability properties of limit cycles can be understood by inspecting the gradient directions of the constant A and constant co curves in the parameter plane. The theoretical background is based on the following theorem [17].
(26)
Since only the sign of the determinant is of interest, it suffices to determine the stability of the limit cycle by inspecting the directions of VA and Vco. Therefore, almost no addi-
54
tional calculations are needed in determining the stability properties. 6. E F F E C T S O F P A R A M E T E R THE LIMIT CYCLES
VARIATIONS
ON
To investigate the influence of parameter variations on the existence of the limit cycle and gain a deeper insight into its characteristics, different parameter planes showing limit cycling have been constructed by the present method. Although most of the results are consistent with those concluded by Wu and Pantalone [5, 8], such studies can provide qualitative and quantitative analyses of the limit cycles. This study will also help to clarify an apparent ambiguity in these works [5, 8]. It should be noted that only stable limit cycle loci are plotted on the plane. Also, unless stated otherwise, the periods up to 90 s which were observed in the real power system [13] are considered. Each plot o f the parameter plane is conducted b y varying parameters examined within a typical working range and keeping all other parameters fixed at the typical values shown in Table 1. In all plots the locus of broken lines represents constant period (or constant frequency) curves while solid lines represent constant amplitude curves. Since the permanent speed droop regulation R and the AGC integral controller gain Ki are t w o important feedback gains for power system operation and control, special emphasis is placed on the analysis of these parameters. Figure 4 shows the effect of I(I 0.60
~
-
%
~ '.o-~o~j'
t I 001I
o;::[ i/,,,?/
varying K~ and R on the limit cycle characteristics. It is significant that there exists an asymptotically stable region, which was not reported in the literature and may be meaningful to power engineers. In this region, a proper set of (KI,R) can be chosen to construct a stable system without limit cycles. For example, keeping KI at a typical value of 0.25, the limit cycling condition begins to vanish when R is increased to about 0.043. Thus, when the range of K1 is small (say K~ <~ 0.36), decreasing R will make a limit cycle more likely to happen. This appears to be in agreement with the conclusions reached by Wu and Pantalone [ 5 , 8 ] . However, for larger values of KI (say KI/> 0.36) the effect of varying R is negligible. Similarly, the effect of varying K~ can be seen by keeping R fixed. It is obvious that for larger values of R (say R ~> 0.024), increasing KI will also make a limit cycle more likely to happen. On the other hand, for small values of R (say R ~< 0.024), the effect of R dominates that of KI. In this region, a larger value of K1 will result in a larger amplitude and frequency of the limit cycle. Since the describing function is just an approximate one, it is necessary to confirm the prediction of a limit cycle by digital simulations. Table 2 shows the results of a comparison between the predicted values and steady state values resulting from digital simulations on some of the sample points (P1 - P8) shown in Fig. 4. From Table 2, it is shown that the predicted values agree rather well with simulation results from an engineering standpoint. Some simulated trajectories corresponding to P3 in Fig. 4 are depicted in Fig. 5. Next, the gradient approach based on Theorem 1 is applied to determine the stability of the limit cycle. For example, in the limit cycle region (refer to Q in Fig. 4), one has AL = U~.(V00 × VA) > 0 and
o, til?tl/
.........................
0"000.00
0.08
0.02 0.0#+ 0.06
010 0.12 0.1#+ 0.16 0.10 R Fig. 4. L i m i t c y c l e loci o n t h e KI-R p a r a m e t e r p l a n e for DB = 0.06%.
[J[ = - - 1 . 4 2 0 < 0 Therefore, the limit cycle is stable. In fact, all the limit cycles plotted in this plane were tested to be stable. The effect of decreasing governor deadband can be shown b y comparing Figs. 4 and 6. The frequencies of the limit cycles are almost
55 TABLE 2 Comparison between the calculated and simulated results Sample point
Value o f parameters
Amplitude and period
Calculated results
Simulated results
Error (%)
P1
K I = 0.3626 R = 0.06955 K I = 0.3805 R = 0.05477 K I = 0.3981 R = 0.03544 K I = 0.4795 R = 0.02430 K I = 0.5348 R = 0.1003 K I = 0.5196 R = 0.1489 K! = 0.5491 R = 0.1507 K I = 0.5774 R = 0.1642
A T A T A T A T A T A T A T A T
0.0003800 60 0.0004000 52.0 0.0004200 46.0 0.0004400 42.0 0.0004600 44.0 0.0004800 44.0 0.0005000 42.0 0.0005300 40.0
0.0003783 58.0 0.0003950 51.7 0.0004148 45.5 0.0004343 41.5 0.0004391 45.6 0.0004534 46.2 0.0004715 44.0 0.0004991 41.8
0.45 3.33 1.25 0.58 1.24 1.09 1.30 1.19 4.54 3.64 5.54 5.00 5.70 4.76 5.83 4.50
P2 P3 P4 P5 P6 P7 P8
(X10-3)
(XlO-3)
S-O
5.0 4-0
4.0 A=0.0004148
Tffi45.5 ApdffiO.O005
3.0
3-0
X 2.0
Y 2.0
1.0
i-0
0.0
0.0
-I.0
(a)
-I.0 0
30
60
90
120
150 180 T(SEc)
210
240
270
300
(b)
30
60
90
120
1S0 180 T(S~C)
210
240
270
300
(XlO-3) 1.0
(XIO-L 0.4
0.8
0.0
0.6
/
-0.4 Y 0.4
AF -0.8
.//
/
/
/
0.2 -1-2 ,
0.0
-1.6
(c)
30
60
90
120
150 180 T(SEc)
210
240
270
300
-0.2 - - / -0.6 -0.4
(d)
~.0.0006 .,~:-- . . . . . . . . . . . -0.2
0.0
0.2
0.4
0,6
0.8
1.0
1.2
1.4 (XIO"3)
x
Fig. 5. Results of computer simulations at limit cycle P3 for a 0.0005 p.u. load increase: trajectories o f (a) the signal at input of governor backlash; (b) the signal at ouput of governor backlash; (c) the frequency deviation in p.u.; (d) governor backlash (deadband) nonlinearity.
unchanged. However, the amplitudes of l i m i t c y c l e s are d i r e c t l y p r o p o r t i o n a l t o magnitude of the deadband. This suggests the straightforward way of relieving
the the that the
d e t r i m e n t a l e f f e c t s c a u s e d b y l i m i t c y c l i n g is t o m a k e t h e d e a d b a n d as s m a l l as p o s s i b l e . Unfortunately, only modest success has been a c h i e v e d o v e r t h e last d e c a d e s [ 1 8 ] .
56
time constant TR, and the governor response time Tg. From Fig. 8, it is shown that decreasing TR tends to enhance limit cycling, which is similar to the effect of decreasing 6. On the other hand, decreasing Tg will decrease the likelihood of limit cycling. However, the effect of varying Tg is rather small.
060 0.5~ 0.48
0.~2 0.36 Kx 0.30 0.2~
6.0
0.1~
5L,
0.12
¢.8
0.06 0.00
¢2
0.00
002
00~,
0.06
0.08
0.10 R
0.12
0.1L,
0.16
0.18
3.6
020
Fig. 6. Limit cycle loci on the KI-R parameter plane for DB = 0.04%•
TR
3.0 2z~ 18 12
Another aspect examined was the influence of the hydro governor model. As shown in Fig. 2, the temporary speed droop 5 and the permanent speed droop R are two possible adjustable parameters in the hydro governor. Variation of the parameters R and 5 yield the results of Fig. 7. It is seen that there also exists an asymptotically stable region. As 5 is increased, the limit cycling condition eventually vanishes. The effect of increasing will result in smaller magnitude and lower frequency of the limit cycle, which is totally contrary to the effect of increasing K~ shown in Fig. 4. Since the effect of varying R has been discussed previously, it is omitted here. Two other time constants to be examined in the detailed governor model are the dashpot
0.6 0.0
,
,
000
[
] CONSTANT
(
} CONSTANT
,
0.05
AMPLITUDE PERIOD ,
0.10
O1S
020
0.25 030 Tg
,
CURVE CU~E
,
035
,
O.t,O
A
,
0.t, 5
050
Fig. 8. Limit cycle loci on the TR-Tgparameter plane for DB = 0.06%.
Next, the influence of the hydro turbine will be considered. Figure 9 shows the effect of varying the water starting time Tw on the KI-Tw parameter plane. As Tw is increased, it tends to enhance the limit cycling. Since the value of TR is closely related with that of Tw in practical applications, typically TR/ Tw = 5.0, it is more reasonable to consider TR and Tw at the same time. From previous
0.70 [
063
(
]
CONSTANT AMPLITUDE CURVE ) CONSTANT PERIOD C U R V E
0-1
.~.~ O.L~9 ~ , Off6
..........................
t 0.1L~
~
0.12 ,"
0.07 0.00
/ /
Fig. for
v /
//
: / / / / / /~ /
/
l
I
I
'
0.000 7.
2,,\',:i,1
0.00 ~
"
Limit
cycle
DB = 0.06%.
loci
on
the
~-R
parameter
plane
'
L ,
O0 0.2 Ot~ 05 08
0.016 0.032 O.OA,8 O.06Z+ 0.080 0.096 0.112 0.128 0.11,1, 0.160 R Fig. for
9. Limit
cycle
DB = 0.06%.
loci
on
~
,
1.0 1.2 Tw the
,
,
•
,
1.¢ 15
K I - T w parameter
.
i ....... 18 2.0 plane
57
discussions, it can be concluded that, as the ratio of TR/T. is decreased, limit cycling is more likely to happen. Finally, the influence o f varying the parameters of the power system model, that is, the damping torque coefficient D and the effective rotary inertia of area M, will be explored. It has been reported that the limit cycling condition can be suppressed by either a decrease in M or an increase in D [5]. In contrast, Wu and Dea [8] concluded that the effect of decreasing M will enhance the limit cycling condition. Such contradictory reports can be clarified by referring to Fig. 10. For larger values of D (say D I> 0.4), a decrease in M will indeed enhance the limit cycling as observed by Wu and Dea. However, for smaller values o f D (say D ~< 0.4), the effect is opposite. Therefore, the effect of varying M not only is determined by itself but also depends on the value of D. Similarly, the effect of varying D cannot solely be determined by itself. For larger values of M (say M/> 3.6), an increase in D can indeed suppress the limit cycling, which is consistent with the conclusions of Pantalone and Piegza. However, for small values of M (say M ~< 3.6) the conclusion is completely contrary.
~
1.2
[ ] CONSTANT AMPLITUDE CURVE ( } CONSTANT PERIOD CURVE
\ II
1.0
~it%--
ASYMPTOTICALLY STABLE REGION
D0.8 0.6 O.h 0.2 0.0
0.O
. . . . . . 1.2 2.6
3.6
4.8
6.0 H
7.2
8.6.
9.6
10.8
REFERENCES 1 C. Concordia, L. K. Kirchmayer and E. A. Szymanski, Effect of speed governor deadband on tie-line power and frequency control performance, A I E E Trans., 76 (1957) 429 - 435. 2 L. K. Kirchmayer, Economic Control of Interconnected Systems, Wiley, N e w York, 1959. 3 S. C. Tripathy, G. S. Hope and O. P. Malik, Optimisation of load-frequency control parameters for power systems with reheat steam turbines and governor deadband nonlinearities, Proc. Inst. Electr. Eng., Part C, 129 (1982) 10 15. 4 S. C. Tripathy, T. S. Bhatti, C. S. Jha, O. P. Malik and G. S. Hope, Sampled data automatic generation control analysis with reheat steam turbines and governor dead-band effects, I E E E Trans., PAS-103 ( 1 9 8 4 ) 1 0 4 5 - 1 0 5 1 .
1.6 1.~*
presented. A new technique for constructing limit cycle loci on this plane is also proposed. This new technique is based completely on the characteristic equation and avoids derivation and separation of stability equations. Therefore, the same approach can also be extended to multiarea systems as well as control systems with multiple nonlinearities. Furthermore, the gradient approach was employed to investigate the stability of the limit cycles. Although most of the results are consistent with previous works, such studies can give deeper insight into limit cycle phenomena. In particular, from the quantitative analysis of limit cycles some apparently conflicting conclusions of previous works can be clarified and reconciled. Results from this study also shed some light on the selection of system parameters. It is possible to choose proper sets of parameters such that the limit cycling condition can be eliminated.
12.0
Fig. 10. L i m i t cycle loci o n t h e D - M p a r a m e t e r p l a n e for DB = 0.06%.
5 D. K. P a n t a l o n e a n d D. M. Piegza, L i m i t cycle analysis o f h y d r o e l e c t r i c systems, IEEE Trans., PAS-IO0 ( 1 9 8 1 ) 6 2 9 - 638. 6 J. M. Bailey a n d G. F. Pierce, E f f e c t s o f n o n linearities o n g o v e r n o r p e r f o r m a n c e at B U L L R U N s t e a m p l a n t , Proc. 7th Annu. South-eastern
Syrup. 7. C O N C L U S I O N S
In this paper, a systematic method that predicts limit cycle oscillations in a hydroelectric system on the parameter plane is
on
System
Theory,
Alabama,
1975,
pp. 137 - 138. 7 J. M. Bailey a n d G. F. Pierce, Backlash a n d r a t e s a t u r a t i o n effects o n g o v e r n o r p e r f o r m a n c e at t h e B U L L R U N s t e a m p l a n t , Proc. 1EEE South-
eastern Region 3 Conf. on Engineering in a Changing Economy, South Carolina, 1976, pp. 10
-
12.
58 8 F. F. Wu and V. S. Dea, Describing function analysis of automatic generation control system with governor deadband, Electr. Power Syst. Res., 1 (1977/78) 113 - 116. 9 C. M. Liaw, C. T. Pan and K. W. Han, Limit cycle analysis of power systems with governor deadband nonlinearities, in preparation. 10 C. H. Ai and K. W. Han, Stability analysis of a nuclear reactor control system with multiple transport lags and asymmetrical nonlinearities, IEEE Trans., NS-22 (1975) 2103 - 2112. 1 IEEE Committee Report, Dynamic models for steam and hydro turbines in power system studies, IEEE Trans., PAS-92 (1973) 1904 - 1916. 12 D. G. Ramey and J. W. Skooglund, Detailed hydrogovernor representation for system stability studies, IEEE Trans., PAS-89 (1970) 106 - 111. 13 D. Ewart, Automatic generation control: per-
formance under normal conditions, Proc. Engi-
neering Foundation Conf. on Systems Engineering for Power: Status and Prospects, Henniker, New Hampshire, August 17- 22, 1975, pp. 1 - 14. 14 A. Gelb and W. E. Van der Velde, MultipleInput Describing Function and Nonlinear S y s t e m Design, McGraw-Hill, New York, 1968. 15 D. D. Siljak, Nonlinear Systems, Wiley, New York, 1969. 16 T. Y. Li, Lectures on the numerical method of finding solutions of a system of nonlinear equations, Seminar o f Numerical Analysis, Chia~Tung Univ., Taiwan, 1983, Lecture 1 - 4. 17 Y. T. Tsay and K. W. Han, A gradient approach to the determination of stability of limit cycles, J. Franklin Inst., 300 (1975) 391 - 418. 18 C. Concordia, Discussion of ref. 8.