SYNCHRONOUS M A C H I N E M O D E L S F O R SIMULATION OF INDUCTION M O T O R TRANSIENTS by R I C K Y P. K. H U N G B.A.Sc. (EE), The University of British Columbia, 1993 A THESIS SUBMITTED IN P A R T I A L F U L F I L L M E N T OF T H E REQUIREMENTS F O R T H E D E G R E E OF M A S T E R OF APPLIED SCIENCE in T H E F A C U L T Y OF G R A D U A T E STUDIES Department of Electrical Engineering We accept this thesis as conforming to the required standard T H E UNIVERSITY OF BRITISH C O L U M B I A April 1995 © Ricky P. K. Hung, 1995 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department of E L E C T R I C A L EN€.lNEER\N&i The University of British Columbia Vancouver, Canada Date ftflW 3 1 . M& DE-6 (2/88) A B S T R A C T The induction motor is the most common motor type used in industry. Analyzing its performance is either concerned with steady-state behavior, or the transient behavior during start-up and during system disturbances. Due to the advancement in power electronics, induction motors are now commonly used in adjustable-speed drives. The harmonic currents generated from these drives may result in voltage distortions and potential resonance problems. Computer simulations of induction motor performance are therefore more important than ever before. Although many simulation programs are available for induction motor studies, the Electromagnetic Transients Program (EMTP) is probably better than other simulation packages, because it contains detailed models of many other power system components, such as transformers, transmission lines and thyristors, which must be considered in certain types of motor simulations. This thesis presents a method of using the existing EMTP synchronous machine model to simulate induction motor transients. This approach allows simulations for both machine types using essentially the same program code. The data conversion algorithm, which converts standard motor specifications into equivalent circuit parameters, is explained. The steady-state initialization of an induction motor is discussed as well. The start-up of an induction motor is simulated using the proposed model. These simulation results are then compared with an independent program which uses a 4th-order Runge-Kutta solution method. It is observed that the results from the proposed model agree very well with the ones from the Runge-Kutta solution method. T A B L E O F C O N T E N T S Abstract ii Table of Contents . i i i List of Tables v List of Figures vi Acknowledgment vii 1. INTRODUCTION 1 2. SIMILARITY BETWEEN THE SYNCHRONOUS MACHINE AND INDUCTION MACHINE MODELS 4 3. DATA CONVERSION 7 3.1 Introduction. 7 3.2 Conversion Algorithm 9 4. STEADY-STATE REPRESENTATION AND INITIALIZATION 19 4.1 Introduction 19 4.2 The Thevenin Equivalent Based Method 20 4.3 The Compensation Based Method 25 4.4 Initialization of Stator and Rotor Currents • 30 5. TRANSIENT SOLUTION 33 5.1 Introduction : 33 5.2 The Electrical Equations 34 5.3 The Mechanical Equations 39 5.4 Saturation in the Leakage Path ...41 iii 6. OPTIMIZATION ALGORITHM FOR FINDING THE BEST PIECEWISE LINEAR X - i CHARACTERISTIC 44 6.1 Introduction , 44 6.2 Objective Function 44 6.3 Method for Finding the Optimal Solution 47 6.4 Examples , '. 51 7. CASE STUDY 54 8 . CONCLUSIONS 62 REFERENCES 64 APPENDIX A : 65 iv L I S T O F T A B L E S Table 4.1 Machine Data 22 Table 4.2 Machine Parameters of the Motor Network 28 Table 4.3 Transmission Line Data 28 Table 4.4 Transformer Data 28 Table 4.5 Operating Slips Obtained by Using the Compensation Based Method 29 Table 7.1 Specifications of the 11 000 HP Induction Motor 54 V L I S T O F F I G U R E S Figure 2.1 Direct Axis Equivalent Circuit of a Synchronous Machine 4 Figure 2.2 Equivalent Circuit of a Double-Cage Induction Motor 5 Figure 3.1 Simplified Equivalent Circuit of a Double-Cage Induction Motor 8 Figure 3.2 Variation of Leakage Reactance with Current 12 Figure 3.3 Equivalent Circuit of a Double-Cage Induction Motor with Saturable and Unsaturable Leakage Reactances 14 Figure 3.4 Variation of T m a x with I s a t 17 Figure 3.5 Variation of T m a x with m 17 Figure 3.6 Flowchart for the Data Conversion Algorithm 18 Figure 4.1 Definitions of the Parameters in Table 4.1 22 Figure 4.2 Motor Network used to Validate the Compensation Based Method 27 Figure 4.3 Definitions of the Current Phasors 31 Figure 5.1 Saturation Characteristic and its Piecewise Linear-Representation 42 Figure 6.1 X - i Characteristic and its Piecewise Linear Representation 45 Figure 6.2 Comparison of the Saturation Characteristic and its Piecewise Linear Representation ( I s a t =1.5 p.u.) 52 Figure 6.3 Comparison of the Saturation Characteristic and its Piecewise Linear Representation ( I s a t = 3.0 p.u. ) 53 Figure 7.1 Comparison of the Continuous Saturation Curve and its Piecewise Linear Representation 56 Figure 7.2 Simulation Results for an Induction Motor During Start-up 57 Figure 7.3 Steady-State Behavior of an Induction Motor 60 vi A C K N O W L E D G M E N T I would like to express my gratitude to my supervisor Df. H. W. Dommel for his guidance and advice throughout the course of my research. Also, I would like to thank for his patience in reading the original drafts of this thesis. I would like to thank my parents for their care, nurture, and support during the many years I have been in school. I would also like to thank my friends in both the Power and the Communication groups for providing a pleasant working environment during my stay at U.B.C. Last, but not the least, I like to thank B.C. Hydro and Power Authority and the Natural Sciences and Engineering Research Council of. Canada for their financial support in this research project. vii j C H A P T E R 1 I N T R O D U C T I O N The induction motor is the most common motor type used in industry. In some situations, its transient behavior may require detailed analysis..For example, the starting of an induction motor can cause a voltage depression which may significantly reduce the accelerating torque of the motor. The reduction in the accelerating torque can lengthen the starting time of the motor, which may cause an overheating problem, or in severe cases, may fail to accelerate the motor to its rated speed, and cause the motor to draw large currents continuously. This large current could damage the machine windings pennanently. Voltage depression caused by.motor starting may also slow down other motors on the same bus. If the voltage drop is severe, these motors will decelerate significantly and may even stall. The decelerating motors, on the other hand, draw more (reactive) current, causing the bus voltage to be reduced further. In extreme situations, a voltage collapse will occur. Furthermore, due to the advancement in power electronics, induction motors are now commonly used in adjustable-speed drives. These drives can generate harmonic currents and cause voltage distortions. Moreover, these harmonic currents may excite the resonant circuit formed by the supply network. Thus, the analysis of induction motors together with the power electronic circuits becomes important as well. The above examples show that it may be necessary to analyze the machine performance and its impact on the supply network in detail in certain situations. Computer 1 simulation programs are best suited for these types of studies. There are special-purpose programs for induction motor studies, which can be used to analyze certain types of motor transients. The Electromagnetic Transients Program (EMTP) can also be used to study induction motor behavior. It has the advantage that it contains detailed models of many other power system components, such as transmission lines, transformers, and thyristors, etc., which may all affect the transient performance of an induction motor. In this thesis, a method of simulating the induction motor transients using the EMTP synchronous machine model is discussed. One advantage of this approach is that less code is required for the induction motor. This makes program maintenance easier because only one set of code rather than two must be maintained. If the models for a.c. machines are improved, only one set of code has to be modified. This is somewhat similar to the approach taken for the "universal machine" in the DCG/EPRI and ATP versions of the EMTP [1], which can be used for the simulation of induction machines as well as synchronous machines. The remainder of this thesis is divided into seven chapters. In Chapter 2, the similarity between the synchronous machine and the induction machine models is discussed. Moreover, the method of converting a synchronous machine model into an induction machine model is described. A data conversion algorithm which converts standard motor specifications into equivalent circuit parameters is explained in Chapter 3. In Chapter 4, methods for finding the steady-state solution of networks with induction motors are presented. The details of implementing a transient analysis program for an induction motor using the 4th-order Runge-Kutta solution method are discussed in Chapter 5. In Chapter 6, 2 an optimization algorithm for finding the best piecewise linear representation of the continuous saturation.characteristic is described. Finally, simulation results of a motor start-up and the general conclusions of this thesis are presented in Chapters 7 and 8, respectively. C H A P T E R 2 S I M I L A R I T Y B E T W E E N T H E S Y N C H R O N O U S M A C H I N E A N D I N D U C T I O N M A C H I N E M O D E L S In the EMTP, the machine equations are solved in the d-q-o reference frame, which is a reference frame attached to the rotor (d = direct axis of rotor, q = quadrature axis of rotor, o = zero sequence). The d-axis equivalent circuit of the synchronous machine model is shown in Figure 2.1. Figure 2.1 Direct Axis Equivalent Circuit of a Synchronous Machine In Figure 2.1, by setting Vf = 0 and replacing subscripts d, f, and D, respectively, with q, g, and Q, the q-axis equivalent circuit is obtained. The D and Q windings are used to represent the damper bar effects, whereas the g-winding is used to model the eddy current effects in a non-salient pole rotor. The equivalent circuit of a double-cage induction motor is shown in Figure 2.2. Subscripts 1 and 2 are used for the two equivalent circuits on the rotor. 4 R, Ls t -JYYYY Lr r e m Figure 2.2 Equivalent Circuit of a Double-Cage Induction Motor Since there is practically no saliency in an induction motor, the equivalent circuit shown in Figure 2.2 is valid for both d- and q- axes. Comparing Figures 2.1 and 2.2, it can be seen that the synchronous machine model is almost identical with that of the induction motor model. In fact, if the field winding of the synchronous machine model is short-. circuited, the synchronous machine model will become the induction motor model. The flux and current relationships for both machine types are of the same form. For the synchronous machine, the relationship is X, — MDF Lff MFD [f _MDD LDD _ JD. (2.1) and for the induction machine, it is = kn .. k + k + k, L+L kn • k+LM L2+k+k, h (2.2) 5 From Equations (2.1) and (2.2), the two sets of machine inductances can be related as follows: (2.3) (2.4) (2.5) (2.6) (2.7) Since the induction machine parameters are identical for both d- and q- axes, the q-axis parameters can be obtained by replacing subscripts d, f, and D, respectively, with q, g, and Q in Equations (2.3) - (2.7). The conversions with Equations (2.3) - (2.7), together with setting V f = 0, convert the synchronous machine model into an induction machine model. L f f = L l + L r + L m L D D = L 2 + L r + L m M t D = L r + L m M d f = M d D = L m 6 CHAPTER 3 DATA CONVERSION 3.1 INTRODUCTION In order to represent an. induction motor accurately in the EMTP, a detailed equivalent circuit of the machine is required. This equivalent circuit should be general enough to model a wide variety of motors, such as motors with deep-bar rotors, double-cage rotors and single-cage rotors. It has been shown earlier [2,3] that the equivalent circuit of a double-cage induction motor, as shown in Figure 2.2, is adequate for these purposes. Even though a model with one circuit on the rotor can be used to represent a single-cage motor, the double-cage model is preferred because the second circuit can be used to represent the frequency-dependent eddy current effects. In many situations, magnetic saturation has negligible effects on the motor's performance [3]. However, there are cases in which saturation effects do play an important role. In Figure 2.2, three inductances can be influenced by magnetic saturation, namely, the mutual inductance (L m ) , the stator leakage inductance (L s), and the mutual rotor inductance (L r). These inductances have different degrees of saturation depending on the magnitudes of currents flowing through them. For example, during start-up, large inrush currents will flow into the motor, causing L s and L r to saturate, but because of the low rotor resistances, very little inrush current will flow through L m , and hence, L m will remain unsaturated. In fact, the saturation effects of L m are insignificant in most types of studies; nonlinear effects are therefore only considered for L s and L r . An exception is the analysis of reclosing transients 7 of capacitor-compensated motors [4,5], where saturation effects in L m are also important; for such cases, the model discussed here would have to be modified. If the very small inductance of the outer cage is neglected, the simplified equivalent circuit of Figure 3.1 is obtained. This equivalent circuit is considered to be adequate for most system studies, and is therefore used as the induction motor model in the data conversion. The EMTP itself accepts L 1 ; which happens to be zero with the data conversion scheme presented in this thesis. One problem in using this equivalent circuit is that the required parameters may not be known to the user. Even if the manufacturer supplies the user with equivalent circuit parameters, these parameters may not accurately reflect the starting and the full load characteristics of the motor. However, through iterative algorithms, a fairly accurate equivalent circuit can be obtained using the standard motor specifications supplied by the manufacturer. A data conversion algorithm which converts motor specifications into equivalent circuit parameters is therefore developed. Since the manufacturer's data are specified for rated steady-state operating conditions, reactances rather than inductances are used in Figure 3.1 and in the remainder of this chapter. ] r e m xr 7 Y Y Y \ Figure 3.1 Simplified Equivalent Circuit of a Double-Cage Induction Motor 8 3.2 CONVERSION ALGORITHM The data conversion algorithm presented here is based on the work done by Rogers and Shirmohammadi [6], with some modifications. The following data are required as input: • rated apparent power of the motor or rated (active) power on the shaft • rated three-phase voltage (line-to-line) • rated efficiency r\ • rated power factor cos 0 • rated slip sr • starting current at rated voltage I s t (in p.u.) • reduced voltage at which another starting current is available, V r e d (in p.u.) • starting current at reduced voltage, I r e d (in p.u.) • starting torque (in p.u. with respect to full load torque) • breakdown torque (in p.u. with respect to full load torque) • saturation threshold current I s a t (in p.u.) In [6], it is suggested that I s a t = 2.0 p.u. is a good approximation for many induction motors. Hence, I s a t has a default value of 2.0 p.u. In the DCG/EPRI version of the EMTP, V r e ( j has a default value of 0.8 p.u., whereas I r e c i has a default value of 0.78 times the starting current at rated voltage [1]. ' 9 Since core loss, friction loss and windage loss are not part of the equivalent circuit, the "effective" efficiency T|' must be increased. Assuming that these losses contribute 25% of the total losses in a motor [6], the effective efficiency is i f = 0.25 +0.75-ri (3.1) The equivalent circuit parameters may be approximated as follows: R s. = c o s e - [ l - i y 7 ( l - s r ) ] (3.2) Let R r be the equivalent rotor resistance. Then R r = s r - r i ' / [ ( l -s r)-cos0.] • (3.3) X m = r i ' / [ ( l - s r ) - s i n e ] (3.4) During start-up, the equivalent rotor resistance R s t can be approximated by • Rst •= T s t - r i ' - c o s e / [ ( I s t ) 2 - ( l -s r)] (3.5) where T s t is the starting torque and I s t is the starting current at rated voltage. 10 After R r and R s t are determined, the rotor resistances can be calculated as R l = R s t ' (1 + m 2 ) - R r ' m 2 (3.6) R 2 = R 1 - R r / ( R 1 - R r ) (3.7) In Equation (3.6), m is the design ratio. Its normal range is between 0.5 and 1.5 [7]. It is shown in [6] that m is around 1.0 for normal double-cage rotors, and between 0.5 and 0.6 for deep bar rotors. In this conversion procedure, m is initially set to unity and is later adjusted until the breakdown torque requirement is met. After R i and R 2 are determined, X 2 can be calculated as To take the variation of the leakage reactance into account according to Figure 3.2, the starting currents at two different voltage levels are needed. The total leakage reactance at rated voltage starting (X t l s) can be calculated as X 2 = ( R 1 + R 2 . ) / m (3.8) (3.9) 11 whereas the one at reduced voltage starting (X t i 2 ) can be calculated as X ill f V V , y red v I R E D J (3.10) The variation of the leakage reactance with current can be represented by a describing function DF [6]: DF •= 1.0 for I < I sat DF = (2/ 7t) • [a + 0.5 • sin (2- a) ] for I > I sat where a = sin"1 ( I s a t /1) (3.11) Figure 3.2 Variation of Leakage Reactance with Current With this function, the total leakage reactance X t l can be expressed as the sum of an unsaturable part X t 0 and a saturable part X t s , that is, 12 X t l = X t 0 + DF • X t s (3.12) Using Equations (3.9), (3.10), and (3.12), X t s and X t 0 can be calculated as follows: t s ~ (DF2 - DFS) • (3-13)-_ X t I s D F 2 - X t l 2 D F s DF2 - DFS (3-14) where DF S and D F 2 are obtained from Equation (3.11) using I s t and I r e d , respectively. If the total saturable and unsaturable leakage reactances are divided equally between the stator and rotor, the reactances can be formulated as follows: X S 0 = X t 0 / 2 (3.15) X r o = X s o - R r - ( R 1 / R 2 ) - m / ( m 2 + 1) (3.16) X s s = X t s / 2 (3.17) ".; X r s = X l s / 2 (3.18) At this stage, initial estimates for all circuit parameters are found. Furthermore, both X s and X r in Figure 3.1 now consist of two parts, namely, the saturable and the unsaturable parts as shown in Figure 3.3. 13 - c J_JYYY\J^TTL Vi, Zi, J i l l J Figure 3.3 Equivalent Circuit of a Double-Cage Induction Motor .with Saturable and Unsaturable Leakage Reactances The accuracy of X m can be improved to match the manufacturer's data more closely, by using an iterative method involving the power factor [6]. In this method, the induction motor is represented by a transient Thevenin equivalent circuit, which consists of a transient voltage source V = V ' d + j V ' q behind a transient impedance Z ' s = R s + j X ' s . It can be shown that . where X s u is the unsaturated leakage reactance of the stator ( X s u = X s 0 + X s s ) . Using Equation (3.20), the ratio of the imaginary to real parts of V is X ' s = X s u . [ l + X m / ( X m + X s u ) J (3.19) V = ( 1 - R s • cos 9 - X ' s • sin 9 ) + j (- X ' s • cos 0 + R s • sin 0 ) (3.20) Vq -X;-cos0 + /gJ-sin0 T ^ " ~ l - ^ - c o s 0 - X ; - s i n 0 (3.21) 14 Moreover, a new variable X s is defined as x =(Rr/sr)-(cx)s0-R,) smO-X' <3'22> where X s = X m + X s u . After X ' s , IL, and X s are determined, new estimates for R r and X m are: d R = S'r- X s - ( R s - y f + X ' s ) 7? Y ^ ( 3 - 2 3 ) X m - y/Xs • ( X s - X's) (3.24) The new values of R r and X m are then compared to the old ones. If the differences are larger than an acceptable tolerance, parameters R|, R2, X 2 , X r 0 , X ' s , . ^L, X s , R r, and X m are recalculated. Typically, 5 steps are needed for these iterative improvements. Referring to Figure 3.3, the per-unit torque at slip s is T = ( R e ( V i n - I 1 * ) - I I 1 l 2 - R s ) / T r a t e (3.25) where T r a t e is the rated machine torque, defined as follow: T r a t e = T i ' . c o s e / ( l - s r ) (3.26) 15 Care must be exercised in determining I| used in Equation (3.25). Initially 11| I and 1I 2 I are predicted so that the DF(s) for both X s s and X r s can be approximated. After the DF(s) have been estimated, the input impedance Z i n can be determined, which, in turn, leads to new values of I 1^ I and I I 2 I. Ij and I 2 are. iterated in this manner until the changes become negligibly small. After 1^ has been determined in this way, the torque at slip s can then be calculated using Equation (3.25). The breakdown torque for this set of machine parameters is then determined. This actual breakdown torque is compared to the rated breakdown torque. If there is a disagreement, the design ratio m is adjusted [7] and parameters R 1 ; R 2 , X 2 , and X r 0 are recalculated. This process is repeated until the actual and rated breakdown torques are reasonably close. In [6], the breakdown torque requirement is met by adjusting the saturation threshold current I s a t . However, the breakdown torque value is fairly insensitive to changes in I s a t , as shown in Figure 3.4. Thus, after the circuit parameters are calculated to satisfy the starting and the full load characteristics, the breakdown torque value can only be adjusted slightly. As a result, the breakdown torque requirement cannot always be met using this approach. The breakdown torque value is very sensitive to changes in m, however, as shown in Figure 3.5. Modifying m is therefore an effective way of meeting the breakdown torque requirement. Moreover, since m is usually unavailable from the manufacturer, it is better to determine it using the rated breakdown torque, rather than assigning a typical value to it. 16 3.506 Saturation Threshold Current (Isat) in p. Figure 3.4 Variation of T with I, Design Ratio m Figure 3.5 Variation of T with m 17 To conclude this chapter, the data conversion algorithm is summarized in a flowchart shown in Figure 3.6. Read Data . Calculate R s , Rr, Xm, Xsu , X s s , Xrs, Xsu . Figure 3.6. Flowchart for the Data Conversion Algorithm 18 C H A P T E R 4 S T E A D Y - S T A T E REPRESENTATION AND INITIALIZATION 4.1 INTRODUCTION For balanced steady-state operation, the three-phase induction motor can be represented as three identical impedances Z p o s in the three phases. There are several methods for obtaining Z p o s , but they all require that the operating slip s0 be known a priori. The operating slip s0 will be different from the rated slip, however, if the actual and rated terminal voltages are noticeably different. One way to obtain Zp 0 S is to find the input impedance of the equivalent circuit in Figure 3.3. Saturation of the leakage reactances can be neglected, since the current will be close to its rated value under steady-state conditions. Another method of finding Z p o s is to use the following equation [8]: Zpos = R s + J Laa " J ws M i [Rr] + J so ws [Lrr] l " 1 J so [Lra] C4-1) where cos is the synchronous speed of the supply network and L, 'm L. 'aa = L s + L, L, m 0 [L ]=[L L ] IK> 0 R2 [LJ = L\ + Lr+L, L2 + L+L (4.2) 19 In practice, s0 is not known unless the motor is started from rest (s0 = 1, speed = 0), but it can be found iteratively with a given torque-speed characteristic of the load. In steady state, the motor is running at a constant speed. The net torque acting on the rotor is therefore zero, Tmotor-" Tioad ~ 0 (4.3) Since both T m o t o r and T l o a d are functions of slip, Equation (4.3) can be used to find the operating slip s0. T | o a d will be known as a function of speed or slip, but the machine torque-slip characteristic Tm o t o r(s) is strongly influenced by the machine terminal voltage, which in turn depends on Z p o s and s0. Hence, iterative algorithms must be used to find s0. In this section, two methods of finding s0 will be discussed, namely, the Thevenin equivalent based method and the compensation based method. 4.2 THE THEVENIN EQUIVALENT BASED METHOD In the Thevenin equivalent based method, the network to which the machines are connected, is represented by a Thevenin equivalent circuit. The operating slips of the machines are then determined iteratively using this equivalent circuit. The advantage of this method is that it can be implemented easily as a separated routine, thus, very minor modifications are needed in the existing EMTP code. This method, however, works well only in the cases where the motors are connected to the same bus. In spite of this limitation, this method is quite useful in obtaining the operating slip(s) in single motor analysis and in 20 motor-starting studies (i.e., the study of the influence of motor starting on other motors on the same bus). This solution method is summarized as follows: Obtain the positive sequence Thevenin equivalent circuit of the network to which the motors are connected. Set S ( / 0 ^ = s r a t e for each motor. Calculate the input impedance Z p o s for each motor. Calculate the bus voltage, using the voltage divider equation. Find the operating slip sQ(1) for each motor with the voltage obtained from step 4, by solving Equation (4.3) with Newton'smethod. Check whether I sJ 1 ) - s0^1"1) I is less than a given tolerance for each motor. If yes, stop. , Otherwise, assign s0(1) = s0^1 for each motor and go to step 3. A program was written to validate this solution method. This program contains four sets of machine data and allows the user to obtain the operating slip of each machine. The user is required to enter the number of machines on the bus and the Thevenin voltage and impedance of the supply network. The machine data and their definitions are provided, respectively, in Table 4.1 and Figure 4.1. 21 Table 4.1 Machine Data f i l l Motor #1 Motor #2 Motor #3 Motor #4 R l o.07 a 0.25 a o.i9i a 0.076 a R2 o.05 a o.2 a 0.0707 ^ 0.062 a X i 0.2 Q. 1.2 a 0.75398 Q 0.195 a x 2 0.2 £2 1.1.Q 0.75398 Q 0.195 a X m 6.5 a 35.0 a 16.8892 Q 6.386 a # of Poles 8 4 . 6 8 Rated Slip 0.04 0.02222 0.016667 0.03 Load Torque T= 15.467 co T = 3.08-10"3co2 T = 2.4415 co T = 0.11073 co2 Note: rated line-to-line voltage = 460 V R] x i | /TYY\ R 2 / s Figure 4.1 Definitions of the Parameters in Table 4.1 X 2 22 Four cases have been studied using this program; the results are summarized as follows Case #1 Rto = 0.0 £L X t h = 0.0 Q. tli Number of machines studied: 4 Number of iterations required: 1 Operating slip of motor #1 Operating slip of motor #2 Operating slip of motor #3 Operating, slip of motor #4 0.040000 0.022220 0.016667 0.030000 Net torque on motor #1 Net torque on motor #2 Net torque on motor #3 Net torque on motor #4 0.00450378 N-m -0.00006898 N-m -0.00208136 N-m 0.00085635 N-m Case #2 Rth = 0.0 Q, X d l = 0.02 Q. Number of machines studied: 2 Number of iterations required: 3 Operating slip of motor #1: 0.040814 Operating slip of motor #2: 0.022620 Net torque on motor #1: 0.00837990 N • m Net torque on motor #2: 0.00062572 N • m 23 Case #3 R t h = 0.0 Q, Xfi, = 0.02 Q, Number of machines studied: 4 Number of iterations required: 4 Operating slip of motor #1 Operating slip of motor #2 Operating slip of motor #3 Operating slip of motor #4 Net torque on motor #1 Net torque on motor #2 Net torque on motor #3 Net torque on motor #4 0.041580 0.022993 0.017363 0.030986 0.00149221 N-m 0.00010812 N-m 0.00029827 N-m 0.00098647 N-m Case #4 R t h = 0.0 Q,, Xfu = 0.06 a Number of machines studied: 4 Number of iterations required: 6 Operating slip of motor #1 Operating slip of motor #2 Operating slip of motor #3 Operating slip of motor #4 Net torque on motor #1 Net torque on motor #2 Net torque on motor #3 Net torque on motor #4 0.045712 0.024970 0.019215 0.033478 0.00346631 N-m 0.00025674 N-m 0.00073395 N-m 0.00227984 N-m In Case #1, since R^ = X U l = 0.0, the rated voltage is directly applied to the machines; hence, the machines are operating at their rated slips, and the program terminated in one iteration. The situations in Cases #2 and #3 are identical except for the number of 24 machines involved. It is observed that the number of iterations increases with the number of machines being studied; the increase is minor, however. Case #4 shows that this solution method converges to a correct solution even if a large Thevenin reactance is used. All these results suggest that the Thevenin equivalent based method is an efficient way to obtain the operating slips of motors connected to the same bus. 4.3 THE COMPENSA TION BASED METHOD The operating slips of induction motors can also be obtained using the compensation based method. The advantage of this method is that it can be used to determine the operating slips of motors at different locations; however, its implementation requires more changes in the existing EMTP code, compared to the Thevenin equivalent based method. The steady-state solution in the EMTP is based on solving the linear nodal equations [Y] • [V] = [I] . , (4.4) This solution can be modified for steady-state initialization of networks with induction motors as follows: 1) Build the admittance matrix [Y] with Z p o s , assuming that all induction motors are operating at their rated slips, or s0(-°-) = s r a t e for each motor. 2) Triangularize the admittance matrix [Y]. 25 3 ) Find the voltages at each node by solving Equation ( 4 . 4 ) with downward operations i and backward substitutions, using the elements of the triangularized admittance matrix. ' I • 1 4 ) Find the operating slip sJ-1) for each motor with the voltages from step 3 , by solving Equation ( 4 . 3 ) with Newton's method. i 5) Check whether I s J v k V k ^ 1 ^ 1= voltage at node k (at i t h iteration) | Y p 0 S k ( s r a t e ) = input admittance of the motor at node k at slip s r a t e Y p o s (s0^1 -*) = input admittance of the motor at node k at slip sG^1 ^ 7 ) Assign s0\1"1) = s0^1 ^ for each motor. Go to step 3 . In step 6,jthe correction term to the k t h element of [I] compensates for the fact that the motor at node k is not operating at its rated slip. One can achieve the same effect by actually updating' the [Y] matrix; this method, however, requires re-triangularization of [Y] in each iterative step, while the algorithm described here requires only one triangularization of [Y] . Only forward operations and backward substitutions have to be performed in re-26 ! • calculating the voltage vector in each iterative step. For one induction motor, the solution 1 , usually converged within 2-3. iterations. For more induction motors, the number of iterations may increase. j In order tjo validate this solution method, a practical motor network is studied [9]. i This network consists of five interconnected motors, as shown in Figure 4.2. The machine ! parameters and the relevant data of the network are provided in Tables 4.2, 4.3, and 4.4. i Figure 4.2 Motor Network used to Validate the Compensation Based Method i I 27 Table 4.2 Machine Parameters of the Motor Network Based on 345 kV, 100 MVA • M l , M3, M4 M2, M5 Srate 0.0179 0.0083 R2 0.408 p.u. 0.633 p.u. Ri 0.204 p.u. '2.25 p.u. X i , X2 2.01 p.u. 7.0 p.u. 50.0 p.u. 157.5 p.u. # of Poles 2 6 Load Torque for M l Load Torque for M2 Load Torque for M3 Load Torque for M4 Load Torque for M5 •T= 100 + 27.48737 co T = 750 + 78.44743 co T = 7.49709 -10"2 co2 T= 1250 + 6.58522 -10"2 co T = 0.677782 co2 Refer to Figure 4.1 for definitions of the parameters. Table 4.3 Transmission Line Data Based on 345 kV, 100 MVA Series Resistance Series Reactance Sh unt Admittance Line 1 - 4 0.001 p.u. 0.0055 p.u. 0.0008 p.u. Line 5 0.002 p.u. 0.0110 p.u. 0.0016 p.u. Table 4.4 Transformer Data Based on 345 kV, 100 MVA Transformer Reactance Transformers 1, 3, 4 0.3773 p.u. Transformers 2, 5 1.2799 p.u. Note : The transformer reactances are assumed to be 45% of the corresponding machine rating. 28 A program was written to determine the operating slips of these machines using the compensation based method. The results of this program are provided in Table 4.5. In order to test the convergence of this iterative algorithm, the transformer reactances are increased to 2.5 times their nominal values. It is observed that an extra iteration is needed to obtain the operating slips after the transformer reactances have been increased. These results indicate that the compensation based method is efficient and robust. Table 4.5 Operating Slips Obtained by Using the Compensation Based Method Rated Slip Operating Slip nominal transformer • reactance Operating Slip 2.5 times nominal transformer reactance Motor #1 0.017900 0.01827689 0.01888605 Motor #2 .0.008300 0.01045600 0.01091092 Motor #3 0.017900 0.01827924 0.01887615 Motor #4 0.017900 0.01828497 0.01888536 Motor #5 0.008300 0.01043090 0.01087720 Number of iterations for the case of nominal transformer reactance: 3 Number of iterations for the case of 2.5 times nominal transformer reactance: 4 29 4.4 INITIALIZATION OF STATOR AND ROTOR CURRENTS After the operating slip s0 is determined, the stator currents in phases a, b, and c can be calculated. Suppose that the stator currents are as follows: ia(t) = 111 cos(cost + a) ib(t) = 111 cos(cost + a -120°) ic(t) = 111 cos(cost + a - 240°) (4.5) Then, the corresponding currents in the d,q,o- reference frame are [3 id(t) = J— 111 sin(s0 cos t + a - 5) iq(t) = ^ 111 cos(s0 cos t + a - 8) io(t)=0 (4.6) where 8 is the angle between the position of the quadrature axis and the real axis. For induction motors, 8 can be set to any arbitrary value. After 8 is chosen, initial value for p can be determined as follows: (3(0) = 8 + | (4.7) 30 The d- and q- axis stator currents can be represented by a complex phasor Idq = ^ l / l * " - ' ) (4.8) with the understanding that / ?(0 = R e { 4 ^ ' } i,(0 = I m { 4 « ' w } (4.9) As can be seen from Figure 4.3, the current phasors IrUq and Ir2dq can be calculated using current division techniques after Idq has been determined. The real and imaginary parts of these phasors, respectively, represent the corresponding q- and d- axis rotor currents. 'dq - * -c J T Y Y V / Y Y Y V rldq R] / S Xo R 2 / s Figure 4.3 Definitions of the Current Phasors 31 Another way to obtain the current phasors IrUq and Ir2dq is to use the following equation [8] : lr\dq = ~ {IKl+J s0 G>, [LJ T 1 j s0 cos [ L J / dq '(4.10) where [Rr], [L r r] and [L r a] are defined in Equation (4.2). In summary, the steady-state solution of an induction machine consists of two parts: determination of the operating slip and initialization of the machine currents. The operating slip is determined using either the Thevenin equivalent based method or the compensation based method. After the operating slip is found, the stator currents can be determined. These currents are then used to calculate the rotor currents using either the current division techniques or Equation (4.10). 32 C H A P T E R 5 T R A N S I E N T S O L U T I O N 5.7 INTRODUCTION The minor modifications in the synchronous machine model to make it behave as an induction motor model have been implemented in UBC's EMTP version MicroTran®. The actual program changes were minor [10] and were done by Dr. H. W. Dommel because he knows the details of the EMTP code better. In order to validate the EMTP results, an independent simulation program, based on the 4th-order Runge-Kutta solution method, was written. The purpose of this program is to study the induction motor behavior during start-up. The supply network is represented by a Thevenin equivalent circuit, i.e., a voltage source, V, behind a positive sequence inductance, L e x t . The behavior of an induction machine is governed by two sets of equations, namely, the electrical equations and the mechanical equations. These two sets of equations are not independent from each other; in fact, they are closely related to each other through the following three variables: speed, angle position, and the net torque acting on the rotor. In this chapter, the details of solving these two sets of equations using the 4th-order Runge-Kutta method will be discussed. 33 5.2. THE ELECTRICAL EQUATIONS . The voltage and current relationships of an induction motor are governed by the following voltage equations: K b c ] = " [ R ] [iabJ " -T t^abc] (5.1) where [ v ] = [va, v b, v c, v f, v g, v D , vQ] [ i ] - [ ia> ib> 1c> xf> 1 g ' 1 D' XQ T [ R ] = diag [ R A , R A , R A , R F , R G , R D , RQ] (subscript ' a ' for armature). The generator convention is used in Equation (5.1), i.e., positive currents are going out from the machine temiinals. In order to simplify the calculations, the voltage equations are solved in the d-q-o reference frame, which is a reference frame attached to the rotor. The phase quantities are transformed to d-q-o quantities through the transformation matrix [T]"1. That is, [v d q o] = IT]"1 K b J [idq0] = [T]"1 BabJ (5.2) 34 where [T] -1 cos (i cos( (5 -120°) cos(£ +120°) 0 0 0 0 sin /5 sin(j8 - 120°)' sin(j8 +120°) 0 0 .0 0 1 41 I V2 l 41 0 0 0 0 0 0 0 1 0 0 0 0 0 0 .0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 (5:3) and where (3 is the angular position of the rotor relative to the stator, measured in electrical radians. In this transformation, the stator variables in phase quantities are projected onto the rotating reference frame, while the rotor variables are unaffected. Moreover, this transf ormation matrix is orthogonal, that is [T] = ( [TT1 ) T (5.4) To include the inductance of the supply network, an extra term is added to the voltage equations, d • d [vabc] = " [ L exJ — frabJ " [RBabJ " ~ [^abJ dt dt (5-5) where [L e x t ] = diag[L e x t, L e x t , L e x t , 0, 0, 0, 0]. 35 Expressing the a, b, c - quantities as d, q, o- quantities, Equation (5.5) becomes [T] [v d q o] = - [L e x t ] - { [T] [i d q o] } - [R] [T] [i d q o] - - { [T] [X d q o ] } (5.6) dt dt [v d q o] = - [T]"1 [L e x t ] [T] [T]"1 { - [ T ] [i d q o] + [T] - [ i d q o ] } - [T]"1 [R] [T] [i d q o] dt dt - [T]"1 { - [T] [Xd ] + [T] - [Xd ] } (5.7) dt 4 dt H Since [T]"1 — [T] = [co], where dt 0 CO 0 0 0 0 0 - CO 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 (5.8) and [Xdqo] = [L d q o ] [i d q o] (5.9) 36 where [Ldqo ] 0 0 M df 0 M dD 0 0 h 0 0 0 M qQ 0 0 0 0 0 0 M df 0 0 L ff 0 M fD 0 0 M 0 r 0 M Af qg 0 0 0 M fD 0 0 0 M qQ 0 0 M gQ 0 (5.10) with - - L s + L m L Q - L S Lff - L g g - L i + L r + L m L D D - L Q Q - L 2 + L R + L M M f D = M g Q = L r + L M Mdf = M d D = M q g = M q Q = L m (5.11a) (5.11b) (5.11c) (5.lid) (5.11e) (5.111) the voltage equations can be expressed as [v d q o] = - [T] [L e x t ] [T] { [co] [i d q o] + - [ i d q o ] } - [R] [i d q o] { M [Ldq 0] [idqo] + [Ldqo] — [idqol } (5.12) 37 Let [ M j.^m - ^ t L^-m. Then [v d q J = - { [M] [co] + [R] + [co] [L d q o ] }. [i d q o] - { [M] + [L d q o ] } - [i d ] (5.13) dt Define [A] = [M][co] + [R] + [co][Ldqo]. and [B] = [M] + [L d q o ] then [v d q o] = - [ A ] [ i d q J - [ B ] - [ i d q o ] (5.14) dt - [i d q o] = - [B]"1 [A] [i d q o] - [B] 4 [v d q o] (5.15) dt Equation (5.15) is a set of first order ordinary differential equations, which can be solved using 4th-order Runge-Kutta method. In each time step, CO and (3 must be known in order to evaluate the right-hand side of Equation (5.15). These two variables, however, are not known until the mechanical equations, which also depend on [i d q o], are solved. As a result, at the beginning of each time step, (3 and CO are first predicted using Equations (5.16) and (5.17), respectively. 38 p(t) = P(t - At) + 0.5 At [ co(t - At) + co(t) ] (5.16) and co(t) = 2 co(t - At) - co(t - 2At) .(5.17) After P(t) and co(t) are predicted, [idqo(t)] is then determined, which, in turn, enables the mechanical equations to be solved. Upon solving the mechanical equations, new values of P(t) and co(t) are obtained. These values are compared with the ones predicted earlier. If they are reasonably close, the program will proceed to the next time step. Otherwise, [idq0(t)] are re-calculated using the new values of P and co. Alternatively, the mechanical equations [Equations (5.18) and (5.19)], could have been added to the electrical equations [Equation (5.15)] to form a system of nine ordinary differential equations. These equations can then be solved simultaneously without using the prediction formulas. However, the method described here resembles the EMTP solution method more closely, and is therefore used to solve the two sets of differential equations. 5.3 THE MECHANICAL EQUATIONS The mechanical variables of an induction motor are governed by the following differential equations: dco J = T net (5.18) dt and 39 dP • , • 0 3 (5.19) dt where J is the moment of inertia of the rotor and mechanical load, and T n e t is the net torque acting on it. Using the generator convention, T n e t is defined as Tnet = T m e c h " ^elec (5-20) where T m e c h is the mechanical torque of the mechanical load and T e l e c is the electrical torque defined as: T e l e c = V q - V d ( 5 ' 2 1 ) Equations (5.18) and (5.19) can be solved using the trapezoidal integration rule: co(t) = — [ T n e t(t) + T n e t ( t - At) ] + co(t - At) (5.22) 2 7 and A t p(t) = — [ ** °> °] (5.25) and 41 L*** / [ Id " !q ] . (6.5) segment 5: A = { [ ^ i m a x ) - W ] - i +[^Id)'Imax - ^ I m a x ) - Id ]} / [ I m a x - Id ] (6.6) 45 Define M(IVI2) = I2 - I, and C ( / 1 5 / 2 ) h - h Then Equations (6.2) - (6.6) become segment 1: A = L • i segment 2: A = M(Im,I n) • i + C(Im,I n) segment 3: A = M(I 0,I p) • i + C(I0,Ip) segment 4: A = M(Iq,Id) • i + C(Iq,Id) segment 5: A = M(I d ,I m a x) • i + C(I d,Im a x) 46 Using Equations (6.9) -formulated as follows: (6.13), the total area between the exact and the fitted curves can be sat a n m b The goal of the optimization is to find optimal values for independent variables I a, I b, I c, I d, I m , I n, I 0, I p, and lq, such that the objective function A becomes a minimum. 6.3 METHOD FOR FINDING THE OPTIMAL SOLUTION At the minimum point of A, the partial derivatives of A with respect to all independent variables must be zero. That is, d A = L • I a - M(Im,I„) • I a - C(Im,I n) = 0 (6.15) d ia d A - = [ M(Im,I n) - M(I0,Ip) ] • I b + C(Im,I n) - C(I0,Ip) = 0 (6.16) — = [ M(I0,Ip) - M(I q,I d) ] • I c + C(I0,Ip) - C(Iq,Id) = 0 (6.17) d A = 0.5 • ( 2 I q 2 - I c 2 - I d 2 ) • — M(Iq,Id) + ( 2 I q - I e - Id) • — C(Iq,Id) - I d • M(I q,I d) - C(I qJ d) - 0.5 • ( I m a x 2 - I d 2) • — M(I d , I m a x ) dId + I d • M(I d , I m a x ) - ( I m a x - Id) • — C(I d , I m a x ) + C(I d , I m a x ) = 0 • (6.18) 0.5 • ( 2 I m 2 - I a 2 - 2 I n 2 + I b 2) • - i - M(Im,I n) + 2 I m • M(Im ,I n) d I di + ( 2 I m - I a - 2In + Ib) • — C(Im,I n) + 2 C(I m J n ) - 2 A,(Im) di = 0 (6.19) — = 0.5 • ( 2 I m 2 - I a 2 - 2 I n 2 + I b 2) • — M(Im,I n) - 2 I n • M(Im,I n) d I. dl + ( 2 I m - I a - 2In + Ib) • — C(Im,In) - 2 C(Im,I n) + 2 Xan) di . = 0 (6.20) 48 — = 0.5 • ( 2 I 0 2 - I b 2 - 2 I 2 + I c 2) • — M(I0,I ) + 2 I 0 • M(I0,I ) d I * . d I + ( 2 I 0 - I b - 21 + Ic) • — C(I0,Ip) + 2 C(I0,Ip) - 2 X(l0) * di 9 A = 0.5 • ( 2 I 0 2 - I b 2 - 2 I p 2 + I c 2) • — M(I0,Ip) - 2 I p • M(I0,Ip) 3IP . « *i, + ( 2 I 0 - I b - 21 + Ic) • — C(I0,Ip) - 2 C(I0,Ip) + 2 MI p) v • dlp d A ^ o o 5 2 _ T 2 _ T 2\ " - = 0.5 • ( 2 I 2 - I c 2 - L/ ) • — M(I Id) + 2 1 - M(I Id) / 4 5 / 4 i i + ( 2 I q - I c - Id) • j - C(Iq,Id) + 2 C(Iq,Id) - 2 X(Iq) (6.21) (6.22) (6.23) This system of nonlinear equations can be solved with the Newton-Raphson method. The procedure for using this method can be summarized as follows: 1) Setn = 0. 2) Estimate the initial current vector T>n> = [Ia, I b, I c, I d, I m , I n, I 0, I p, I q ] T 3) Evaluate the derivative vector F ( n ) using I (n), where 49 F<»> = [ Fj, F 2 , F 3 , F 4 , F 5 , F 6 , F 7 , F 8 , F 9 ] T 3 A 3 A d A 3 A . 3 A 3 A 3 h 3 !„• 3 lB 3 ld 3 Im 3 /„ 3 A 3 A 3 A d I 4) Evaluate the Jacobian Matrix J ( n ) using I ( n ) . where 0 0 0 UL dF6 dF6 0 0 0 0 0 0 dF2 dh 0 0 d F7 dlb dlc dFs dFs dh 0 d L 0 0 d F-x d F, d h d Id d F4 d F4 dh 0 0 d F7 0 0 0 0 d FQ d-F0 d Ft d Fx dlm dln d Fn d F0 d F0 d Fn dh 0 0 d L d L dh 0 0 0 0 dh . 0 0 d Fc d F< dlm dln dF< d Ft dh 0 0 0 dh 0 0 0 d F7 0 0 0 0 d F7 dh dlp dFj d F& dh dlp 0 0 0 0 dh dlp d F, d F^ d F-x dh dh d FA d L 0 0 0 0 d F9 d I 9 J and where the elements of J are defined in Appendix A. 50 . 5) Estimate the correction current vector AI ( n ) by solving the linear system j(n) . AIM = F^ n ) 6) Check whether the absolute value of each element in AI ( n ) is less than a given tolerance. If yes, stop. 7) Otherwise, adjust the current vector l(n) as follows: |(n+l) _ j(n) . ^j(n) 8) Increment n by 1 and go to step 3. 6.4. Examples Two examples are given to illustrate this optimization algorithm. Example #1 L = 9.08840-IO'5 H, I b a s e = 1137.565 A, I s a t = 1.5.p.u., I m a x = 15.0 p.u. Using the above data, the optimization program yields the following results: L i = 9.08840 IO"5 H la = 1853.75 A L 2 = 1.80852 IO"5 H Ib = 2898.83 A L 3 = 3.18004 •10'6H Ic = 5214.53 A L 4 = 3.63107 •10"7H Id = 10077.05 A L 5 = 8.84719 •10-8-H 51 where L 1 ; L 5 are the slopes (inductances) of the linear segments, and I a, I d are the current points at which the piecewise linear curve changes from one segment to another. A comparison between the actual nonlinear characteristic and its piecewise linear representation is shown in Figure 6.2. 0.2 Current in A Figure 6.2 Comparison of the Saturation Characteristic and its Piecewise Linear Representation ( I s a t =1.5 p.u.) Example #2 L = 9.08840 - IO"5 H, I b a s e = 1137.565 A, I s a t = 3.0 p.u., I m a x = 15.0 p.u. Using the above data, the optimization program yields the following results: L i = 9.08840 IO"5 H Ia = 3639.39 A L 2 = 2.43524 IO'5 H Ib = 5087.09 A L 3 = 6.50988 10"6H Ic = 7773.70 A L 4 = 1.41785 10"6H Id = 12036.04 A L 5 = 5.38722 10"7H 52 Again, Lj , L 5 are the slopes (inductances) of the linear segments, and I a, I d are the current points at which the piecewise linear curve changes from one segment to another. A comparison of the actual nonlinear characteristic and its piecewise linear representation is shown in Figure 6.3. Current in A Figure 6.3 Comparison of the Saturation Characteristic and its Piecewise Linear Representation ( I s a t = 3.0 p.u.) Figures 6.2 and 6.3 show that the piecewise linear curves approximate the nonlinear characteristics very well. 53 CHAPTER 7 CASE STUDY To show the usefulness of the proposed induction motor model, the start-up of a large induction motor is simulated. The supply network of the motor is represented as a Thevenin equivalent circuit, with a voltage of 6797.33 V (RMS, Ime-to-line), behind an inductance of L p o s = 0.5305 mH. The machine specifications are listed in Table 7.1. Table 7.1 Specifications of the 11 000 HP Induction Motor T61 Line-to-Line Voltage: 6600V Full Load Specifications: Efficiency = 0.985, Power Factor = 0.906, Rated Slip = 0.00622. Starting Specifications: at V = 1.0 p.u., I = 8.0 p.u., at V = 0.758 p.u., I = 6.03 p.u., Starting Torque = 1.457 p.u. Other Information: Maximum Torque = 3.5 p.u., Number of Poles = 4, Isat = 2 - ° P - u -Moment of Inertia = 50 590 lb-ft2, Mechanical Load : T m e c h = 1.21*co2, where co is the machine speed in rad/sec. 54 Using the data conversion program described in Chapter 3, the following equivalent circuit parameters (in p.u.) are obtained: R s = 4.586 •IO"3, Xso - 6.009 -IO"2, Xss = 3.616-IO"3 3.094 •10°, Xro = 5.229 -IO"2, Xrs = 3.616 -IO"3 Ri = 2.485 •IO"2, R 2 = 8.756 -IO"3, X 2 = 6.054-IO"2 With the above data and I m a x = 15.0 p.u., the optimization program (described in Chapter 6) yields a piecewise linear curve with the following parameters: L 1 = 9.08840 •10"5H • la = 2454.10 A L 2 = 2.03003 •IO5 H Ib = 3671.75 A L 3 = 4.22453 •10"6H Ic = 6183.43 A L 4 = 6.33421 •10-?H Id = 10869.55 A L 5 = 1.86161 •107 H L l 5 L 5 are the slopes (inductances) of the linear segments, and I a, I d are the current points at which the piecewise linear curve changes from one segment to another. The actual characteristic and its piecewise linear representation are shown in Figure 7.1. It is seen that the piecewise linear representation approximates the continuous characteristic very well. 55 Current in A Figure 7.1 Comparison of the Continuous Saturation Curve and its Piecewise Linear Representation The simulation results for the 11 000 HP induction motor during start-up are shown in Figure 7.2. x 10" 1 i 1 1 1 1 1 r J J I I 1 1 1 2 4 6 8 10 12 14 Time in sec Figure 7.2(a) Current Envelope 56 x 10 1.5 E 1 OJ 0.5 -0.5 •1.5 1 { A. „cu 1 v 1 1 • 0 2 4 6 8 10 12 14 Time in sec Figure 7.2(b) Torque Curve 2000 1500 Q_ en •- 1000 "O QJ dJ Q. CO 500 0 -i r ~i~ r yy yy yy yy yy yy y y ^ . . yy ' yy yy yy •* O y _1 L_ 4 6 8 10 12 Time in sec Figure 7.2(c) Speed Curve 14 Figure 7.2 Simulation Results for an Induction Motor During Start-Up Nonlinear Model - Linear Model 57 The current envelope in Figure 7.2(a) shows the large inrush currents which exist during a motor start-up. The amplitudes of these currents remain practically unchanged until the motor has reached its rated speed. In addition, d.c. offset currents are present in the beginning of the start-up process. The torque curve, on the other hand, shows that the motor torque has large oscillations immediately after the motor is energized. Comparing Figures 7.2(a) and 7.2(b), it is seen that the oscillatory torque decays with the d.c. offset currents. Finally, the speed curve shows that the rotor speed climbs up steadily to its rated value. A small overshoot is observed before the rotor speed settles down to its rated value. The results of Figure 7.2 have been obtained with MicroTran, and have been verified with an independent program, which uses a 4th-order Runge-Kutta solution method, as discussed in Chapter 5. Figure 7.2 shows that the linear model overestimates the starting time of the motor, thus giving a pessimistic result. The saturable inductances are much smaller than their unsaturable counterparts ( L s o / L s s ~ L r o / L r s ~ 15.0); hence, the saturation effects are not very pronounced in this motor, as illustrated by the differences between the linear and nonlinear models. The saturation effects may be more noticeable in other motors. To show the correctness of the proposed steady-state initialization methods, the behavior of this motor is simulated using the steady-state solution as initial conditions. Using the Thevenin equivalent based method, the operating slip of this motor is found to be 0.005906. Using 8 = 15° ( 0.2618 in rad), and the procedures outlined in Section 4.4, the following steady-state solution is obtained: 58 I D = 934.9506 A I 0 = 0.0 A I „ = 297.8980 A I Q - 829.3784 A I Q =-975.2451 A I F = -126.7777 A I D = -393.6692 A co = 374.7645 rad/s (electrical quantity) (3(0) = 1.83260.rad (electrical quantity) This steady-state solution is then used as initial conditions for the transient simulation. The simulation results are shown in Figure 7.3. < cz /rf dld , dld - 2 ^ -C( I q , I d ) + 2 Lj • ^ - M ( I d , I m a x ) - 2 I d • -^-M(I q , I d ) - M(I q,I d) °'yrf ^ ^ r f d2 d + M(I d , I m a x ) - 0.5 ( I m a x 2 - I d 2 ) • — - M(I d , I m a x ) + 2 - - C ( I d ; I m a x ) d1 - ( W - y * 777 COW max/ (A15) d h dFA d „ „ „ d2 2 I q • — M ( I q , I d ) + 0.5 ( 2 I q 2 - I c 2 - I d 2 ) • — — M(I q J d ) d Iq dld dlqdld . + 2 ^ " C(Iq,Id) + ( 2 Iq - I c - Id) • -j—-C(Iq,Id) - I d -^-M(I q ,I d ) C(Iq,Id) (A16) d ' d h dPX ^ " ' dlm dh dF2 dh ' = dlm (A17) (A18) 67 7 7 = °-5 ( 2 lm " Ia2 " 2 I n 2 + I b 2) 77MOWIJ + 4 I m " f M(Im,I n) + ( 2 I m - I a - 2 I n + Ib) T^CfJWIo) + 4 - 7 - C(Im,I n) - 2 r(Im) + 2 M(Im,I n) (A19) 5 7. 2 - T 2 - n T 2 - T" * - ^ — M ^ ^ n ) + 2 Im 7 M ( I m , I n ) di di d I = 0 . 5 ( 2 I m 2 - I a 2 - 2 I n 2 + I b 2) d d + ( 2 I m - I a - 2 I n + Ib) T 7 T - C ( I m , I n ) + 2 — C(Im,I n) 2 - 7 - C(Im,I n) - 2 I n - 7 - M(Im,I n) C/ J d 1 (A20) < ^ 6 dFi = dln dF6 Ms di di (A21) (A22) (A23) m n — = 0.5 ( 2 I m 2 - I a 2 - 2 I n 2 + I b 2 ) -^~2M(lm,ln) - 4 I n f - M ( I m , I n ) 5 5 + ( 2 I m - I a - 2 I n + Ib) — C ( I m , I n ) - 4 77 C(Im,I n) + 2 V(I n) ol. • dl - 2 M(Im,I n) 5 F 7 5 F 2 (A24) (A25) 68 9F7 dF3 (A26) ^ - 0.5 ( 2 I 0 2 - I b 2 - 2 I p 2 + I c 2 ) f7M(I 0 , I p ) + 4 I 0 M(I 0,I p) + ( 2 I 0 - I b - 2 I p + I c) — C(I0,Ip) + 4 — C(I0,Ip) - 2 V(I 0) + 2 M(I0,Ip) (A27) d F l =0.5(2 I 0 2 - I b 2 - 2 I p 2 + I c 2) M(I0,I ) + 2 IG — M ( I 0 , I p ) 5 / dl d l „ d l o + ( 2 I 0 - I b - 2 I p + Ic) 7—--C(I0,Ip) + 2 J- C(I0,Ip) 2 C(I 0,I p)" 2 I p MOoJp) (A28) 5F 8 d F2 ^ " ^ 8 d F3 c p dF-, = d IP (A29) (A30) (A31) d FQ = 0.5 ( 2 I 0 2 - I b 2 - 2 I p 2 + I c 2 ) — M ( I 0 , I ) - 4 I p — M(I0,Ip) + ( 2 I 0 - I b - 2 I p + I c) — - C(I0,Ip) - 4 ±- C(I0,Ip) + 2 X'(Ip) 5 / „ 2M(I 0,I p) (A32) 69 d F9 dF3 3'c d I « dF9 dF4 d F0 d2 d 7 = 0.5 ( 2 I q 2 - V - l d ) — 2 M(Iq,Id) + 4 I q — M(I q,I d) d2 where (A33) (A34) + ( 2 I q - I c - I d ) — C(Iq,Id) + 4 — C(Iq,Id) + 2 M(I q,I d) - 2 X'(Iq) (A35) — Mdi.12) = [ Ml 2 ) " Mil) - ^ ' d l ) ' ( h - h ) ] / ( h - II) 2 (A36) — Mai,I2) = [ Mil) MI 2) + ^'(I2) • ( I 2 - Il )] / ( I 2 - II) 2 (A37) ^ 7 2 — C(lhI2) = I 2 • { X'di) • [ I 2 - Ii ] - MI 2) + Mil) } / ( I 2 - Ii) 2 (A38) I s a t A = - I s a t / ( R - i 2 ) B = - I s a t . R / i 2 + I s a t 3 / ( R . i 4 ) R = [ l - ( I s a t / i ) 2 ] 0 5 (A49) X" = L { 2 • DF'(i) + DF"(i) • i } DF"(i) = 0 ' • if i< I s a t DF"(i) = (2/7t)-[C + D] if i > I s a t C = 2 . I s a t / ( R . i 3 ) + I s a t 3 / ( R 3 . i 5 ) D = 2 - I s a t . R / i 3 - I s a t 5 / ( R 3 - i 7 ) - 5 - I s a t 3 / ( R . i 5 ) R = [ l - ( I s a t / i ) 2 ] 0 - 5 (A50) 72 *