Cross-Identification of Synchronous Generator ... - IEEE Xplore

776

IEEE TRANSACTIONS ON ENERGY CONVERSION, VOL. 26, NO. 3, SEPTEMBER 2011

Cross-Identification of Synchronous Generator Parameters From RTDR Test Time-Domain Analytical Responses René Wamkeue, Senior Member, IEEE, Christian Jolette, Augustin B. Mpanda Mabwe, Senior Member, IEEE, and Innocent Kamwa, Fellow, IEEE

Abstract—In a recent work, the authors computed generator parameters from a load-rejection test based on what they called hybrid state model of the synchronous machine. The state matrices were presented in compact numerical form and may prove tedious to implement. Such a model is a complicated function of the unknown machine parameter vector to be computed which is not always easy to comprehend in an identification process for engineering applications. In this paper, therefore, the decoupled property of the so-called hybrid model along with the application of the complete well-known solution of the linear control systems theory is used to derive time-variant analytical waveforms of the phase voltages and the field current following a generator tripping (load rejection) and an open stator field short-circuit tests in terms of the generator parameters. In the cross-identification approach, the field short-circuit test is used to compute the generator d-axis parameters while the q-axis parameters are obtained from the load-rejection test q-axis data. The proposed identification technique is successfully applied for the parameter identification of a 4-pole 1.5-kVA, 208-V, and 60-Hz saturated laboratory synchronous generator. Index Terms—Cross-identification, field short-circuit, load rejection, synchronous generator (SG), time-domain analytical waveforms, validation.

I. INTRODUCTION ODAY’S environmental priorities of clean energy around the world are encouraging small private gas, wind, hydroelectric, and other types of generating plant to be connected to the power system rather than the building of new high-voltage lines. With the increase of future energy demand, this will result

T

Manuscript received August 11, 2010; revised November 4, 2010, December 4, 2010 and January 22, 2011; accepted February 21, 2011. Date of publication May 12, 2011; date of current version August 19, 2011. This work was supported in part by the Natural Sciences and Engineering Research Council of Canada. Paper no. TEC-00323-2010. R. Wamkeue is with Université du Québec en Abitibi-Témiscamingue, Rouyn-Noranda, QC J9X 5E4, Canada, and also with the Laval University, Laval, QC G1K 7P4, Canada (e-mail: [email protected]). C. Jolette is with Université du Québec en Abitibi-Témiscamingue, RouynNoranda, QC J9X 5E4, Canada (e-mail: [email protected]). A. B. Mpanda Mabwe is with Ecole Supérieure d’Ingénieurs en Electronique et Electrotechnique, Amiens, France (e-mail: [email protected]). I. Kamwa is with Hydro-Québec-IREQ, Varennes, QC J3X 1S1, Canada, and also with the Laval University, Laval, QC G1K 7P4, Canada (e-mail: [email protected]). Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org. Digital Object Identifier 10.1109/TEC.2011.2140320

in the rapid growth of present distribution networks and their operation near stability limits. With such a significant integration of small generating plants in today’s power system, it becomes more and more imperative to adapt conventional programs with exact generator parameter values for appropriate stability analysis, control, and system fault prediction in order to account for the future impacts of modern power networks. In fact, the adequate operational behavior of the power system mainly depends on the accuracy of its generator dynamic parameters. The extraction of the generator parameters from data generally involves three basic steps [1]–[3]. 1) A generator dataset obtained from a given test; 2) a generator model structure that can reliably predict the test designed in 1); and 3) the identification method that can determine the optimal set of parameters guided by the data. As reported in [1], selecting the appropriate system model structure 2) is without any doubt the most important and, the same time, the most difficult stage of the identification procedure. In fact, since the beginning of the 1990s, powerful time-domain identification tools have been proposed to compute synchronous machine (SM) parameters [2]–[6]. However, in many of them, less interest is given to the model derivation than identification techniques. This might justify in part why the voltage-controlled state model structure based on short-circuit tests (for large disturbance tests) or small disturbances in the field voltage has been intensively used in last centuries. As a solution to this generator modeling shortcoming, authors have recently successfully applied a so-called state hybrid model to compute the generator parameters [7]. Though state-space model structure is widely used in identification practices with an undoubted success, it is a complicated function of the unknown machine parameter vector to be computed, not always easy to comprehend in engineering applications [8]. In fact, identification procedures based on the system state model are tedious, computationally intensive, and in general requires a minimal system identification skill [8]. For a long time, the classical and standardized shortcircuit phase current analytical waveform has been widely used to graphically compute the generator transient reactances and short-circuit time constants [9]–[11]. With availability of today various system identification commercial tools, a

0885-8969/$26.00 © 2011 IEEE

WAMKEUE et al.: CROSS-IDENTIFICATION OF SYNCHRONOUS GENERATOR PARAMETERS

time-domain curve fitting approach along with the standardized phase short-circuit current analytical waveform was developed in [12]. In [13], a modern identification procedure based on the generator phase voltage analytical waveform used as a model predictor was also applied to compute the machine transient reactances and open-circuit time constants. The great advantage of identification procedures using state models is that they provide equivalent circuits parameters of the generator. In fact, modern stability analysis based on the eigenvalues determination, vector control techniques, and performances predictions of the generator is more suitable with equivalent circuits’ parameters than time constants. Accordingly, it is then easy from test book formulas to compute generator various transient reactances and time constants from equivalent circuit parameters [14], [15]. From the aforementioned analysis, it will be, therefore, interesting to derive an identification procedure that could combine advantages of the estimation based on the state-model approach (provide equivalent circuit parameters) and those of the analytical system responses method (easy to comprehend, suitable, and time saving). This summarizes the motivation of this paper: to propose simple and suitable two-step approach to compute the generator-equivalent circuit parameters based on field shortcircuit and load rejection time-variant analytical responses. As defined in IEEE Std-115 [9], the field short-circuit test with stator on open circuit is performed to compute d-axis opencircuit time constants. As established in this paper, running time-domain response (RTDR) tests such as the load rejection and field short-circuit tests have the same transfer function. Accordingly, in the cross-identification method proposed, the field short test is used for the first time in this paper to derive generator d-axis-equivalent circuit parameters while the q-axis parameters are obtained from the load-rejection test. The main objectives of the paper are: 1) to organize the hybrid state model recently developed by authors into decoupled d- and q-axes state subsystems [7]; 2) to derive the time-variant analytical responses during the field short-circuit and load-rejection tests in terms of machine-equivalent circuit parameters using the wellknown complete solution of the state equation of the linear control systems’ theory and the MATLAB symbolic toolbox; 3) to present a numerical procedure to perform the field short circuit and the q-axis load-rejection tests of the synchronous generator (SG) from analytical formula previously derived; 4) to perform experimental data of the aforementioned two RTDR tests for the identification process [16], [17]; 5) to provide the suitable prediction-error cross-identification procedure based on the asymptotic weighted least-squares estimator (AWLSE) and a Newton-type finite difference constrained optimization algorithm; 6) to prove the effectiveness and suitability of the propose approach in the estimation of the SG-equivalent circuits parameters; 7) to assess more the reliability of the estimated parameters in wide range of applications with cross validations of

Fig. 1.

777

(a) d-axis and (b) q-axis-equivalent circuits of the SG.

the optimal obtained model with an arbitrary axis loadrejection test. II. HYBRID MODEL OF THE SG The electrical and mechanical equations of the SM (see Fig. 1) are defined by (1)–(4), where ωn is the nominal angular frequency (in rad/s) and ωm is the rotor speed (in p.u.). Flux equations ⎡

ψd

⎢ ψq ⎢ ⎢ ⎢ ψf ⎢ ⎣ ψD ψQ

⎡

⎤

−xdd

⎥ ⎢ 0 ⎥ ⎢ ⎥ ⎢ ⎥ = ⎢ −xf d ⎥ ⎢ ⎦ ⎣ −xD d 0

0

xdf

xdD

−xq q 0

0 xf f

0 xf D

0 −xQ q

xD f 0

xD D 0

⎤⎡

0

id

⎥ ⎢ iq ⎥⎢ ⎥⎢ ⎥ ⎢ if ⎥⎢ 0 ⎦ ⎣ iD xQ Q iQ xq Q 0

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (1)

or equivalently

Ψs Ψr

−Xs = t −Xsr

Xsr Xr

Is . Ir

(2)

Voltage equations ⎤ ⎡ −ra vd ⎢ vq ⎥ ⎢ 0 ⎢ ⎥ ⎢ ⎢ ⎥ ⎢ ⎢ vf ⎥ = ⎢ 0 ⎢ ⎥ ⎢ ⎣ 0 ⎦ ⎣ 0 ⎡

0

0 −ra

0 0

0 0

0 0

rf 0

0 rD

0 ⎤

0

0

0 ⎡ ψd ⎢ ψq 1 d ⎢ ⎢ + ⎢ ψf ωn dt ⎢ ⎣ ψD ψQ

⎥ ⎥ ⎥ ⎥ + ωm ⎥ ⎦

⎤⎡ id 0 ⎢ ⎥ 0 ⎥ ⎢ iq ⎥⎢ 0 ⎥ ⎢ if ⎥⎢ 0 ⎦ ⎣ iD ⎡

rQ

⎤

−ψq ⎢ ψd ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ 0 ⎥ ⎢ ⎥ ⎣ 0 ⎦ 0

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

iQ

(3)

778


or equivalently −Rs O2,3 Is Vs 1 d I5 = + ωn dt Vr O3,2 Rr Ir ωm Ξ O2,3 Ψs + . O3,2 O3,3 Ψr

respectively

p (Ψr ) =

Ad O1,2

y = Cd

(4)

Motion equations

Ur d O2,1 Bd O2,1 Ψr d + Aq O1,2 Bq Ur q Ψr q

Ψr d Ur d Cq D D d q + (13) Ψr q Ur q

y = yd + yq = if

1 d (ωm ) = [(Tm − Te ) − d0 ωm ] dt 2H d (δ) = ωn (ωm − 1) dt

vd

vq

t

.

(14)

A. d-Axis Decoupled Hybrid Model of the SM (5)

The d-axis subsystem derived from (13) and (14) is given by p (Ψr d ) = Ad Ψr d + Bd Ur d ;

where Ψs = [ ψd

t

ψq ] ;

Is = [ id

t

iq ] ;

V s = [ vd

t

vq ]

Ψr d0

yd = Cd Ψr d + Dd Ur d

(6)

(15)

where Ψr = [ ψf ψD ψQ ]t ; 0 −1 Ξ= . 1 0

V r = [ vf

0 0 ]t ;

ψD ]t ;

Ur d = [ vf

id ]t ;

Ψr d0 = [ ψf 0 yd = [ if

ψD 0 ]t vq ]t

0

(16)

(7)

Relationships between the reactances of (1) are given in the Appendix. I5 in (4) is the five-order identity matrix. The rotor speed and internal rotor angle of the generator are given by (8) and (9), respectively. P is the active power consumed by the initial load; H is the equivalent of the inertia constant; and d0 is the damping friction. For a poor active power and reactive loads and for a field short-circuit test with stator on open circuit (P ≈ 0), the generator maintains the synchronism ωm = ωm 0 = 1 (p.u.) in (8) and the rotor internal angle in (9) δ(t) = δ0 ≈ 0(rad). Under such operating conditions, the mechanical transients do not have any influence on the model obtained in (10). Accordingly, ωm = ωm 0 = 1 in system output matrices defined in (10)

t 1 1 (P )dζ + ωm 0 = P.t + ωm 0 (8) ωm (t) = 2H 0 2H

t ωn δ(t) = P t2 . ωn (ωm (ς) − 1)dς + δ0 = (9) 4H 0 The so-called hybrid state model of the SM given in (10) is derived from (2) and (4) (see the Appendix for details) d (Ψr ) = A.Ψr + B.U ; dt y = C.Ψr + D.U.

Ψr d = [ ψf

⎡ −r x rf xf D ⎤ f DD a11 a12 ⎢ ⎥ N N Ad = ωn ⎣ = −rD xf f ⎦ rD xf D a21 a22 N N ⎡ rf (−xD D xdf + xdD xf D ) ⎤ 1 ωn ⎢ ⎥ N Bd = ωn ⎣ ⎦= 0 rD (−xdD xf f + xdf xf D ) 0 N θdt = [ ra

rf

rD

xa

xf

xD

xm d

xk f ] .

(19)

B. q-Axis Decoupled Hybrid Model of the SM As previously mentioned with the d-axis decoupled subsystem, the q-axis subsystem is given by Ψr q 0

yq = Cq Ψr q + Dq Ur q

(10)

From system matrices (11)–(12), the hybrid state model (10) can be organized into d- and q- axis decoupled subsystems,

Output matrices of the d-axis subsystem in (15) are given in the Appendix. All state matrices of the d-axis subsystem are computed from the generator d-axis parameters vector (19). θd is the d-axis parameter vector.

p (Ψr q ) = Aq Ψr q + Bq Ur q ;

III. DECOUPLED HYBRID STATE MODEL OF THE GENERATOR

b12 b22

(18)

Ψr 0

The system state matrices defined in (10) can be organized in the following decoupled forms in d- and q-axes: Ad O2,1 Bd O2,1 A= ; B= (11) O1,2 Aq O1,2 Bq C = Cd Cq ; D = [ Dd Dq ] . (12)

(17)

(20)

where Ψr q = [ψQ ] ;

Ψr q 0 = [ψQ 0 ]

0 ]t rQ rQ xq Q Aq = −ωn ; Bq = −ωn xQ Q xQ Q ⎡ ⎤ ⎤ ⎡ 0 0 ⎢ xQ q ⎥ x2Q q ⎥ ⎢ ⎥ ⎢ Cq = ⎣ ωm x ; Dq = ⎢ ωm xq − ⎥ ⎦ QQ xQ Q ⎦ ⎣ 0 0

Ur q = [iq ] ;

y q = [ 0 vd

(21) (22)

(23)


θqt = [ ra

rQ

xa

xQ

xm q ]

(24)

where θq is the q-axis parameter vector. The decoupled hybrid state model clearly illustrates that a given test can be performed separately into d- and q-axis subtests using d- and q-axis subsystems, respectively. In other words, an RTDR test is the superposition of d- and q-axis subtests. When the generator operates at the synchronous speed ωm = ωm 0 = 1 as explained earlier, the phase A and terminal voltages are defined by va (t) = vd cos (ωn t) + vq sin (ωn t) vt (t) =

vd2 (t) + vq2 (t).

(25)

779

Φr d (s) = (sI − Ad )−1 =

adj (sI − Ad ) . det (sI − Ad )

(33)

We can define the d-axis transfer function Hd (s) for running frequency-domain responses introduced for the first time in this paper Yd (s) = Cd Φr d (s)Bd + Dd Hd (s) = Ur d (s) Ψ r d 0 =0 hf f (s) hf d (s) = . (34) hdf (s) hdd (s) The transition matrix is computed as follows:

IV. TIME-VARIANT q-AXIS LOAD-REJECTION TEST RESPONSES To better understand the methodology used, let us begin the derivation with q-axis load-rejection solution whose subsystem is only one order. The q-axis subsystem (20) in Laplace domain leads to the q-axis transfer function (26) for RTDR tests with the transition matrix defined in (27) Yq (s) = Cq Φr q (s)Bq + Dq (26) Hq (s) = Ur q (s) Ψ r q 0 =0 −1 ωn rQ adj (sI − Aq ) Φr q (s) = (sI − Aq )−1 = = s+ det (sI − Aq ) xQ Q

adj(sI − Ad ) det(sI − Ad ) 1 a12 s − a22 = (35) a21 s − a11 (s + p1 ) (s + p2 )

Φr d (s) = (sI − Ad )−1 =

det (sI − Ad ) = s2 − (a11 + a22 ) s + a11 a22 − a12 a21 = (s + p1 ) (s + p2 ) with (a11 + a22 ) − p1 = 2

(27) where adj(M) and det(M) denote, respectively, the adjoin and the determinant of the matrix M. Φr q (s) is the q-axis hybrid subsystem transition matrix. The q-axis control input variable Ur q = [0] for the load-rejection test [18]. Using the complete solution of a state equation of linear systems’ theory applied to the q-axis subsystem of (20) for only homogeneous equation (Ur q = [0], Ψr q 0 = 0) yields −1 ωn rQ ψQ 0 (28) Ψr q (s) = ΨQ (s) = s + xQ Q

ψQ (t) = L−1 {ΨQ (s)} = ψQ 0 e−t/T q 0 xQ q vd (t) = Cq ψQ (t) = − ψQ 0 e−t/T q 0 xQ Q x QQ Tq0 = . ωn rQ

(29) (30) (31)

Here, Tq0 is the generator q-axis subtransient open-circuit time constant. V. TIME-VARIANT FIELD SHORT-CIRCUIT TEST RESPONSES The hybrid d-axis subsystem (15) can be written in Laplace domain as follows: sΨr d (s) − Ψr d0 = Ad Ψr d (s) + Bd Ur d (s) Yd (s) = Cd Ψr d (s) + Dd Ur d (s). The solution of (32) is given by Ψr d (s) = Φr d (s)Ψr d0 + Φr d (s)Bd Ur d (s) Yd (s) = Φr d (s)Ψr d0 + [Cd Φr d (s)Bd + Dd ] Ur d (s)

(32)

p2 =

(a11 + a22 ) + 2

√ Δ √ Δ

(36)

=−

1 Td0

=−

1 Td0

Δ = a211 + a222 − 2a11 a22 + 4a12 a21

(37)

Td0 =

2N 2 ωn rf xD D + rD xf f − (rf xD D + rD xf f ) − 4rf rD N Td0 =

2N . 2 ωn rf xD D + rD xf f + (rf xD D + rD xf f ) − 4rf rD N (38) Td0

Td0

and are the generator transient and subtransient Here, d-axis open-circuit time constants. Prior to the field short-circuit test, the unloaded generator is running at synchronous speed as defined in IEEE Std-115 in [9, subsection 11.10]. The test is performed by suddenly short-circuiting the field voltage vf = 0. The d-axis input vector is then Ur d = [ 0 0 ]t . As previously mentioned for the q-axis load-rejection test, the complete solution of the state model (15) is then reduced to a homogeneous case (Ur d = [ 0 0 ]t ; Ψr d0 = 0). The first equation of (33) becomes Ψr d (s) = Φr d (s)Ψr d0

(39)

780


t where Ψr d (t) = L−1 {Ψr d (s)} = ψf (t) ψD (t) in the time domain. L is the Laplace function. From (39), the d-axis subsystem state variables ψf (t) and ψD (t) yield kf 10 e(−t/T d 0 ) + kf 20 e(−t/T d 0 ) ψf (t) −1 ψf (s) =L = ψD (t) ψD (s) kD 10 e(−t/T d 0 ) + kD 20 e(−t/T d 0 ). (40) The output vector yd (t) following the field short-circuit test is computed from the second equation of (15) using (40) ψf (t) t yd (t) = [ if (t) 0 vq (t) ] = Cd (41) ψD (t) xD D kf 10 − xf D kD 10 if (t) = e−t/T d 0 N xD D kf 20 − xf D kD 20 + (42) e−t/T d 0 N xD D xdf − xdD xf D kf 10 e−t/T d 0 + kf 20 e−t/T d 0 vq (t) = N xdD xf f −xdf xf D + kD 10 e−t/T d 0 +kD 20 e−t/T d 0 . N (43) The field short-circuit test is chosen instead of the d-axis loadrejection because of the simplicity of (43) that will provide an easy identification process. For the d-axis load-rejection test deriving the output vector (41) from (33) involves more complicated time-domain analytical formulations. The experimental phase A voltage is entire on the d-axis. In fact, vd (t) = 0 from (30) and (25) yields va (t) = vq sin (ωn t) vt (t) = vq .

(44)

VI. CROSS-IDENTIFICATION METHOD AND MODEL PARAMETERIZATION A. Identifiability of a Deterministic Linear State Model of the SG The technique used in this paper is resulted from the identifiability property of deterministic linear models. In fact, as specified in [8], time-domain identification of the SG parameter vector θ from the state model (10) is equivalent to the frequency-domain identification of the SG parameter vector from the RTDR transfer function (45) obtained from the superposition principle using (34) and (26). As previously mentioned in Section I, the direct use of the state form (10) or the state-model solution (46) for the identifiability problem quickly involves the solution of very complicated nonlinear equations with tedious error prone derivation. As suggested in [8], this difficulty is avoided in this study by a suitable use of the MATLAB symbolic algebra program to derive (46) [19]. It is the main goal aimed in Sections IV and V where the q-axis load-rejection test and field short-circuit test time-variant analytical responses are derived, respectively, in terms of SG parameters to be estimated. H(s) = C(sI3 − A)−1 B + D

Fig. 2. Field voltage during the field short circuit with stator open: step change in the field voltage.

= [Cd (sI2 − Ad )−1 + Dd ] + [Cq (sI1 − Aq )−1 + Dq ] = Hd (s) + Hq (s)

(45)

y(θ, t) = C(θ)eA (θ )t Ψd0

t + C(θ) eA (θ )(t−τ ) B(θ)U (τ )dτ + D(θ)U (τ ). 0

(46) B. Prediction Model for the Identification Procedure According to the decoupled property of the hybrid model previously demonstrated, the d- and q-axis parameters can be estimated separately; thus, (46) yields y(θ, t) = yd (θd , t) + yq (θq , t) θ = {θd } ∪ {θq } .

(47)

The prediction model vector is such that yq = vd (30) for q-axis parameter estimation, while yd = [ if vq ]t (41)–(42) when daxis parameters are to be estimated. The analytical formulation derived in this paper greatly simplifies the computation of the optimization stage during the identification process. A given predictor is selected if it allows computing the corresponding axis SG parameter vector. The two-step cross-identification method can be accordingly summarized as follows. 1) Step 1: Estimation of d-Axis Parameters: The field shortcircuit test is first used to estimate the d-axis parameter vector θd (19). In fact, since the SG is unloaded during the test, the initial rotor angle is zero and the rotor of the generator is aligned with the d-axis. The field voltage waveform following the test is shown in Fig. 2. The SG becomes rapidly unsaturated during this test. From unsaturated initial conditions (low magnitude terminal voltage), the SG operates under linear conditions. 2) Step 2: Estimation of q-Axis Parameters: The q-axis parameter vector θq (24) is computed from the q-axis loadrejection test data. In this step, the procedure developed in [7] is used to compute the steady-state conditions. A resistive/inductive initial load is chosen. If the terminal voltage waveform magnitude remains small during the load-rejection test, the


781

Fig. 3. MATLAB/Simulink diagram of a space-state hybrid model of the SG for cross validation.

saturation can be ignored if the initial load voltage is on the unsaturated region of the unloaded saturation characteristic curve. The d-axis transient voltage data for the estimation process of the q-axis parameters are computed from armature voltage data as follows: 2π vd (t) = va cos (ωn t) + vb cos ωn t − 3 2π . (48) + vc cos ωn t + 3

Fig. 4.

Flowchart of the SG estimation procedure.

Fig. 5.

Diagram field for the short-circuit test principle.

The hybrid d-axis state submodel (13) is applied for the cross-validation of d-axis parameters estimated in step 1. The MATLAB/Simulink diagram in Fig. 3 is used to implement the d-axis submodel (13) for simulations purposes. C. Estimator, Parameterization, and Optimization Algorithm The prediction error approach based on the AWLSE is used in this paper [1], [8]. In fact, assuming that N observations are available, a good model is one that minimizes the discrepancies σ(θ, k) (49) between the model prediction yp (θ, k) and the actual data y(k). This is achieved by minimizing a certain norm of covariance of σ which is defined by (50) x0 = x(0) (49) y(tk ) = yp (θ, tk ) + σ(θ, tk ) V (θ) =

N 1 σ(tk )t W σ(tk ) N

(50)

k =1

min V (θ) subject to G(θ)

For q-axis parameters

⎧ ⎪ ⎨ θq − θq m in ≥ 0 G(θ) = Gq (θ) = θq m ax − θq ≥ 0 ⎪ ⎩ Tq 0 (θ) > 0

(53)

where θm ax and θm in are upper and the lower bounds of the parameter vector θ to be estimated. The identification procedure is defined by the flowchart given in Fig. 4.

(51)

where W is a diagonal weighting matrix which is helpful for scaling purposes. The parameter estimation is done by a nonlinear programming procedure defined as the optimization problem (51), which can be solved starting with an initial guess of parameter θ0 using a Newton-type constrained optimization algorithm [1], [8], [20]. The “fmincon” function of the MATLAB optimization toolbox was used for the present paper [20]. G(θ) is a function of physical estimation constraints defined by (52) and (53) when estimating d- and q-axis parameters, respectively. For d-axis parameters ⎧ ⎪ ⎨ θd − θd m in ≥ 0; θd m ax − θd ≥ 0 Td0 (θ) > 0; Td0 (θ) > 0 G(θ) = Gd (θ) = ⎪ ⎩ Td0 (θ) > Td0 (θ). (52)

VII. IDENTIFICATION RESULTS AND DISCUSSIONS A. Identification Experiments The load-rejection test principle is well defined in [7]. The test loads are laboratory Labvolt modules [7]. For the field shortcircuit test with the stator in open circuit (IEEE Std-115, [9]), the unloaded generator is running at the synchronous speed and the field circuit is suddenly short-circuited (k1 = 1 is closed, k2 = 1 is closed, and k3 = 0 is open, see Fig. 5). For both tests, a reduced steady-state terminal voltage was chosen to keep the model linear. Since the generator does not supply any real power to the load, p0 = 0 p.u. and ωm ≈ ωm 0 = 1 p.u. The machine to be identified is a 208-V, 1.5 -kVA, 4-pole, 60-Hz synchronous turbine generator driven by a dc motor at 1800 r/min. A reduced steady-state terminal voltage has been chosen to keep the model linear. The reciprocal per unit base

782


TABLE I PARAMETER (P.U.) OF SG OBTAINED FROM IEEE STD-115 TESTS

TABLE IV D-AXIS ESTIMATED PARAMETER VALUES OF THE GENERATOR (P.U.)

TABLE II INITIAL PARAMETER (P.U.) FOR THE IDENTIFICATION PROCESS

TABLE III STEADY-STATE VARIABLES (P.U.) FOR THE FIELD SHORT-CIRCUIT TEST OF THE SG

Fig. 6.

Experimental setup for the load-rejection test.

Fig. 7.

Phase A voltage following the field short circuit: model validation.

Fig. 8.

Field current following the field short circuit: model validation.

system defined in [9] is used in this paper. If b0 is the field current for 1 p.u. open-circuit terminal voltage on the air-gap line [9]. Initial parameter values of Table I were computed from the classical method given in IEEE Std-115 [9]. Some sensitive parameters of Table I given in Table II were perturbed to prove the robustness of the identification procedure. The lower and upper bounds of the inequality constraints (52) were θm ax = 10θ0 and θm in = 0.05θ0 where θ0 is the initial parameter vector listed in Table I √ √ Vb = 2VN = 2 × 120 V = 169.7 V; Vdb = Vq b = VN = 120 V Sb = SN = 3VN IN = 1517 VA ≈ 1.5 kVA √ Ib = 2IN = 5.96 A; Idb = Iq b = IN = 4.21 A If b0 = 0.759 A;

If b = If b0 x0m d = 0.50 A; Vf b =

Sb = 3055.645 V. If b

B. SG d-Axis Parameter Identification Results and Model Validation Using Field Short-Circuit Test: Step 1 Table III shows initial steady-state conditions before shortcircuiting the field voltage. A sample of 2500 points per channel of the Tektronix numerical scope was selected. The weighted matrix in (50) was W = diag(4, 1) for the variables scaling. The estimated d-axis parameters are given in Table IV. Comparisons of the estimated output responses of the model and experimental data shown in Figs. 7 and 8 attest to the reliability and efficacy of the estimation process.


TABLE V STEADY-STATE VARIABLES (P.U.) FOR THE LOAD-REJECTION TEST OF THE INDUCTIVE/RESISTIVE LOAD (q-AXIS PARAMETER ESTIMATION)

Fig. 9. d-axis voltage following load rejection of resistive/inductive load: model validation.

783

TABLE VII STEADY-STATE VARIABLES (P.U.) FOR THE LOAD-REJECTION TEST OF THE INDUCTIVE/RESISTIVE LOAD (CROSS VALIDATION)

Fig. 10. Field current following the load rejection of inductive/resistive load: estimation model cross validation.

TABLE VI q-AXIS ESTIMATED PARAMETER VALUES OF THE GENERATOR (P.U.)

C. SG q-Axis Parameter Identification Results and Model Validation Using q-Axis Load-Rejection Test: Step 2 For the estimation process, an initial resistive/inductive load is connected to the generator. The initial steady-state conditions are given in Table V. Due to the decoupled property of the hybrid model of the generator, only q-axis variables of Table V are selected for the model prediction. The d-axis voltage data for qaxis parameter estimation are computed from (48) using phase voltage data. The d-axis voltage (30) is chosen as the model predictor for the q-axis parameter identification. The d-axis estimated voltage is compared with the actual d-axis voltage in Fig. 9 for the model validation. The estimated q-axis parameters are given in Table VI. D. SG Model Cross Validation Using Load-Rejection Test For the model cross-validation process, estimated parameters of Tables III and VI along with the complete hybrid state model (10) were used for simulations of another load-rejection test of a resistive/inductive. The initial conditions for cross-validation stage are given in Table VII.

e Fig. 11. Evolution of residuals σ = im f − if (θ e ): estimation model cross validation.

The MATLAB/Simulink diagram of the space-state hybrid model of the SG is given in Figs. 3, 10, and 12 and shows comparisons of the estimated model and measured data for an inductive/resistive load-rejection test, while Figs. 11 and 13 illustrate the discrepancy between the estimated and experimental data. xm and xe denote the measured and estimated values of variable x, respectively. It is well observed again in the cross-validation stage that the estimated model and output data waveforms are very close. It should be noted that in saturated initial conditions (at rated voltage, for example), the load-rejection test transients are influenced by iron saturation. In such a context, the simulation model should use the magnetizing-reactance updating procedure developed in [7].

784


SG stator flux (2) and voltage (4) subequations Ψs = −Xs Is + Xsr Ir

(a6)

1 d (Ψs ) + ωm ΞΨs . Vs = −Rs Is + ωn dt

(a7)

SG rotor flux (2) and voltage (4) subequations t Ψr = −Xsr Is + Xr Ir ;

Vr = Rr Ir +

1 d (Ψr ). ωn dt

(a8)

Derivation of SG hybrid State-Space model [7] From the first equation of (a8) derive the rotor current vector Fig. 12. Phase A voltage following the load rejection of inductive/resistive load: estimation model cross validation.

t Ir = Xr−1 (Ψr + Xsr Is ).

(a9)

From the second equation of (a8) and (a9) derive d t (Ψr ) = −ωn Rr Xr−1 Ψr + ωn Vr − ωn Rr Xr−1 Xsr Is dt t = − ωn Rr Xr−1 Ψr +ωn [1 0 0]t .vf −ωn Rr Xr−1 Xsr Is t .U = − ωn Rr Xr−1 Ψr + ωn [1 0 0]t −ωn Rr Xr−1 Xsr

= A.Ψr + B.U.

(a10)

Derivation of the field current [7] ψf = F1 .Ψr = xf f if + F2 .U + F3 .Ir = xf f if Fig. 13. Evolution of residuals σ = v am − v ae (θe ): estimation model cross validation.

t + F2 .U + F3 .Xr−1 (Ψr + Xsr Is ) t = xf f if + F2 .U + F3 .Xr−1 Ψr + F3 Xr−1 Xsr Is

= xf f if + F2 .U + F3 .Xr−1 Ψr + [0

VIII. CONCLUSION Compartmental model structures such as state deterministic linear models frequently used in modern identification procedures generally involve the solution of complicated and burden nonlinear equations. In this paper, the complete solution of the linear system control theory along with MATLAB symbolic algebra toolbox is used to derive the hybrid state model timedomain outputs responses in terms of SG-equivalent circuit parameters. From the decoupled property of the SG hybrid model established in this paper, the aforementioned time-variant outputs are organized into d- and q-axis estimators for the so-called cross-identification experience based on the open stator field short circuit and the load-rejection tests. The proposed estimation technique appears very suitable and easy to comprehend and to apply for the SG parameter determination. It provides interesting and competitive results compared to methods based on SG state models. APPENDIX

if =

xQ Q = xm q + xQ

(a2)

xf f = xm d + xf + xk f 1 ;

xdf = xf d = xm d

(a3)

xQ q = xq Q = xm q

xD f = xf D = xD D − xD .

(a11)

1 1 [F1 − F3 .Xr−1 ].Ψr + [F4 − F2 ].U xf f xf f

= cf .Ψr + df .U F1 = [ 1

0

F3 = −F2 ;

(a12)

0 ]; F2 = [ 0 F4 = [ 0

−xm d

0]

t F2 Xsr Xr−1 ].

Derivation SG system outputs vd −ψq Vs = ≈ ΞΨs = vq ψd π ⎤ π ⎡ sin − cos − 2π ⎦ ψd 2π =⎣ ψq cos − − sin − 2 2

xq q = xq = xm q + xa (a1)

xD D = xm d + xD + xk f 1 ;

xD d = xdD = xm d ;

= xf f if + F3 .Xr−1 Ψr + (F2 − F4 ).U

(a13) (a14)

(a15)

t )Is Ψs = Xsr Xr−1 Ψr + ε(Xs − Xsr Xr−1 Xsr

Relationships between SG reactances xdd = xd = xm d + xa ;

t F3 Xr−1 Xsr ].U

(a4) (a5)

= cϕ Ψr + dϕ Is = cϕ .Ψr + [ 0 dϕ ].U ⎡

(a16)

⎤ if c df y = ⎣ vd ⎦ = f .Ψr + .U = C.Ψr + D.U. cϕ 0 dϕ vq (a17)


State matrices in terms of SG parameters ⎡ −r x ⎤ rf xf D f DD 0 N N ⎢ ⎥ ⎢ rD xf D ⎥ −rD xf f ⎥ = Ad 0 A = ωn ⎢ ⎢ ⎥ N N O1,2 ⎣ −rQ ⎦ 0 0 xQ Q N = xD D xf f − x2f D ⎡

rf (−xD D xdf + xdD xf D ) N rD (−xdD xf f + xdf xf D ) N 0

1

⎢ ⎢ B = ωn ⎢ ⎢0 ⎣ 0

Bd = O1,2

O2,1 Bq

xD D ⎢ N ⎢ xD D xdf − xdD xf D ⎢ C = Υ⎢ N ⎢ ⎣ 0

c31

c32

⎡ 0

⎢ ⎢ 0 D = Υ⎢ ⎢ ⎣ 0 ⎡ 0

⎢ ⎢ D=⎢ ⎢0 ⎣

kf 20 = (a18)

kD 10 =

⎤

kD 20 =

0

x2Q q − xq xQ Q

0

⎥ ⎥ 0 ⎥ ⎥ −xQ q ⎥ ⎦ xQ Q (a20)

(a21)

⎤ 0

⎥ ⎥ ⎥ ⎥ ⎦

(a22)

⎤ ωm

ωm d22

0 x2Q q xq − xQ Q 0

Dq ]

⎥ ⎥ ⎥ ⎥ ⎦ (a23)

d22 = xD D x2df − xd xD D xf + xf f x2df − 2xdD xdf xf D + xd x2f D Υ=

1 O2,1

(a24) O1,2 ; Ξ

Ξ=

0 1

−1 . 0

Intermediate coefficients of (40) and (41) ωn rf xD D ωn rf xf D ; a12 = a11 = − N N

ωn rD xf D ωn rD xf f ; a22 = − N N √ 1 a12 a11 − a12 + Δ ψf 0 + √ ψD 0 2 Δ √ 1 a12 −a11 + a12 + Δ ψf 0 − √ ψD 0 2 Δ √ 1 a21 −a11 + a12 + Δ ψD 0 + √ ψf 0 2 Δ √ 1 a21 a11 − a12 + Δ ψD 0 − √ ψf 0 . 2 Δ

(a27) (a28) (a29) (a30) (a31)

REFERENCES

0

Cq ]

0

xf d xD D − xD f xD d N

= [ Dd

kf 10 =

⎤

−xf D N xdD xf f − xdf xf D N 0

xf d xD D − xD f xD d N d22

0

⎥ ⎥ ⎥ 0 ⎥ −rQ xq Q ⎦ xQ Q

⎤ 0 xQ q ⎥ ⎥ = [ Cd ωm ⎦ xQ Q 0

−xf D N 0

⎢ N C=⎢ ⎣ 0

O2,1 Aq

( a19)

⎡

⎡ xD D

a21 =

0

785

(a25)

(a26)

[1] L. Ljung, System Identification: Theory for Users, 2nd ed. Englewood Cliffs: NJ: Prentice-Hall, 1999. [2] E. Kyriakides, G. T. Heydt, and V. Vittal, “On line parameter estimation of round rotor synchronous generators including magnetic saturation,” IEEE Trans. Energy Convers., vol. 20, no. 3, pp. 529–537, Sep. 2005. [3] R. Wamkeue, I. Kamwa, and X. Dai-Do, “Line-to-line short-circuit test based maximum likelihood estimation of stability model of large generators,” IEEE Trans. Energy Convers., vol. 14, no. 2, pp. 167–174, Jun. 1999. [4] M. Karrari and O. P. Malik, “Identification of physical parameters of a synchronous generator from online measurements,” IEEE Trans. Energy Convers., vol. 19, no. 2, pp. 407–415, Jun. 2004. [5] H. B. Karayaka and A. Keyhani, “Synchronous generator model identification and parameter estimation from operating data,” IEEE Trans. Energy Convers., vol. 18, no. 1, pp. 121–125, Mar. 2003. [6] J. C. Wang, H. D. Chiang, C. T. Huang, Y. T. Chen, C. L. Chang, and C. Y. Chiou, “On-line measurement-based model parameter estimation for synchronous generators: Solution algorithm and numerical studies,” IEEE Trans. Energy Convers., vol. 9, no. 2, pp. 337–343, Jun. 1994. [7] R. Wamkeue, F. Baetscher, and I. Kamwa, “Hybrid state model based timedomain identification of synchronous machine parameters from saturated load rejection test records,” IEEE Trans. Energy Convers., vol. 23, no. 1, pp. 68–77, Mar. 2008. [8] G. A. F. Seber and C. J. Wild, Nonlinear Regression. Hoboken, NJ: Wiley, 2003. [9] Test Procedures for Synchronous Machine. IEEE/ANSI Standard 115, 1996. [10] I. Kamwa, P. Viarouge, B. Mpanda-Mabwe, and R. Mahfoudi, “Experience with computer-aided graphical analysis of sudden-short-circuit oscillograms of large synchronous machines,” IEEE Trans. Energy Convers., vol. 10, no. 3, pp. 407–414, Sep. 1995. [11] M. Mehmedovic, Z. Rabuzin, and T. Bicanic, “Estimation of parameters for synchronous generator from sudden short-circuit test data considering rotor speed variation,” Automatika, vol. 14, no. 3, pp. 454–359, 2007. [12] J. P. Martin, C. E. Tindall, and D. J. Marrow, “Synchronous machine parameter determination using the short-circuit axis currents,” IEEE Trans. Energy Convers., vol. 14, no. 3, pp. 454–459, Sep. 1999. [13] E. D. C. Bortoni and J. A. Jardini, “Synchronous machines parameter identification using load rejection test data,” IEEE Trans. Energy Convers., vol. 17, no. 2, pp. 242–247, Jun. 2002. [14] P. Krause, O. Wasynczuc, and S. D. Sudhoff, Analysis of Electric Machinery, NJ: IEEE Press, 1994. [15] P. Kundur, Power System Stability and Control. New York: McGrawHill, 1994. [16] Confirmation of Test Methods for Synchronous Machine Dynamic Performance Models, EPRI, EL-5736, 1988. [17] P. A. E. Rusche, G. J. Brock, L. N. Hannet, and J. R. Willis, “Test and simulation of network dynamic response using SSFR and RTDR derived synchronous machine models,” IEEE Trans. Energy Convers., vol. 5, no. 1, pp. 145–155, Mar. 1990. [18] R. Wamkeue, C. Jolette, and I. Kamwa, “Alternative Approaches for linear analysis and prediction of a synchronous generator under partial- and fullload rejection tests,” IET J. Electr. Power Appl., vol. 1, no. 4, pp. 581–590, 2007.

786


[19] Symbolic Math Toolbox for Use With Matlab, The MathWorks, Inc. Natick, MA, 2000. [20] A. Grace, Optimisation Toolbox for Use With Matlab. The MathWorks, Inc., Natic, MA, 2008.

René Wamkeue (S’95–M’98–SM’11) received the Ph.D. degree in electrical ´ engineering from Ecole Polytechnique de Montréal, Montreal, QC, Canada, in 1998. Since 1998, he has been with the Université du Québec en AbitibiTémiscamingue, Rouyn-Noranda, QC, Canada, where he is currently a Full Professor of Electrical Engineering. He is also an Associate Professor of Electrical Engineering at Université Laval, Québec City. His research interests include control, power electronics, modeling and identification of electric machines, power system cogeneration with induction generators, and wind energy conversion systems. Prof. Wamkeue is a member of the IEEE-PES Electric Machine Committee and the current secretary of Working Group 7 for revision of IEEE Std-115. He is a recipient of the 2009 IEEE Power Engineering Society Prize Paper Awards.

Christian Jolette received the B.Eng. degree in electrical engineering from ´ Ecole Polytechnique de Montréal, Montreal, QC, Canada, in 2001 and the M.Sc. degree in electrical engineering from Université du Québec en AbitibiTémiscamingue, Québec, in 2006, where he is currently working toward the Ph.D. degree. His main research interests include electric machine design, modeling, and parameter identification of electric machines and drive system.

Augustin B. Mpanda Mabwe (SM’97) received the Ph.D. degree in electrical engineering from the Faculté Polytechnique de Mons, Mons, Belgium, in 1990. He is currently a Professor and Head of Electrical Engineering Department at Ecole Supérieure d’Ingénieurs en Electronique et Electrotechnique, Amiens, France. Since 2007, he has been a Visiting Professor within the French’South African Institute of Technology, Tshwane University of Technology, Pretoria, South Africa. His research interests include distributed generation, power quality, energy conversion, diagnostic, and parameter identification of electric machines.

Innocent Kamwa (S’83–M’88–SM’98–F’05) received the Ph.D. degree in electrical engineering from Laval University, Laval, QC, Canada, in 1988. Since 1988, he has been with the Hydro-Québec Research Institute, Power System Analysis, Operation and Control, Varennes, QC, where he is currently a Principal Researcher in bulk system dynamic performance. Since 1990, he has been an Associate Professor of Electrical Engineering at Laval University. He has been active for the last 15 years on the IEEE Electric Machinery Committee where, he is presently the Standards Coordinator. Dr. Kamwa is a member of Conference Internationale des Grands Reseaux Electriques and a Registered Professional Engineer. He is a recipient of the 1998, 2003, and 2009 IEEE Power Engineering Society Prize Paper Awards and is currently serving on the Adcom of the IEEE System Dynamic Performance Committee.