Dynamic Simulation of High-Speed Induction Motor

In the drives of high-speed devices, such as a blood centrifuge, dynamic states also play an important role in terms of the time and quality of the tests performed. The article presents the application of modified equations resulting from the mathematical model of an induction motor to model dynamic phenomena during motor start-up, both with mains supply and with frequency start-up. The applied solution considers the phenomenon of current displacement in the rotor bar and the phenomenon of saturation. The comparison of the obtained results with the experiment shows that the method is sufficiently accurate. The obtained results can also be extended to higher power machines and to modeling other dynamic states.


Introduction
Correct modeling of phenomena in dynamic states is significant in driving with frequent rotational speed changes, especially in higher frequency drives. For this purpose, many works use the classic or modified form of the machine equivalent circuit [1][2][3]. In most cases, equations describing an induction motor's dynamics based on the Park transform are solved in the Matlab/Simulink environment [4][5][6][7]. Such models are effective in several situations. However, from work [8], the sources of inaccuracies were indicated that lie in the failure to consider non-linear phenomena, such as saturation or skin effect, as well as the occurrence of higher harmonics generated in the motor as a result of core slotting and winding distribution on the stator and rotor circuits, and motor core losses. Therefore, many studies propose different approaches to refine the model of these phenomena in an induction motor in dynamic states. In the works [9][10][11][12][13], various methods of considering the phenomenon of saturation and losses in the motor in the dynamic state were proposed. However, it should be emphasized that not all the above-mentioned phenomena have a significant impact on the dynamic waveforms of the motor. Therefore, both core losses and higher harmonics can be practically neglected in these calculations. Simultaneously, the omission of the influence of the magnetic circuit's saturation, especially on the leakage inductance of the machine windings, may lead to serious errors, especially when calculating the electromagnetic torque stroke. A separate problem is considering the phenomenon of current displacement, which has a powerful impact, especially on the motor starting runs. A multi-circuit representation of rotor bars is often used to account for this phenomenon, but generally, only bars with the simplest shapes, i.e., rectangular or trapezoidal, are concerned. There are also few studies on applying the field-circuit method using the finite element method for modeling dynamic phenomena [14][15][16].
A centrifuge is a device that rotates objects around a fixed axis, thus generating a significant force perpendicular to the axis of rotation (centrifugal force). The blood centrifuge operation is based on the principle of sedimentation, in which the acceleration along a radius causes the denser fractions to move outwards. At the same time, the less dense elements move toward the axis of rotation. In the blood placed in test tubes, the radial force causes the deposition of denser particles at the bottom of the test tube [17][18][19]. Blood in an adult (4.5 to 6 L) usually makes up 7-8% of the human body weight. Blood is a highly specialized tissue consisting of over 4000 different types of components. Four of the most important are red blood cells, white blood cells, platelets, and plasma. Usually, the purpose of centrifugation is to obtain plasma, which contains biomarkers and substances indicative of the physiological state of the body and can be used to diagnose various diseases or conditions of the body. To improve the accuracy of the diagnostic test, high-purity blood plasma is required. The purity of the plasma depends on the rotational speed. The World Health Organization recommends a protocol applying 2000-3000 g for at least 15 min to prepare cell-free plasma. It is claimed that the centrifugation force and spinning time could be reciprocally adjusted. The rapid centrifugation technique (at 13,000× g for 3 min) for routine coagulation testing can be used safely. This can benefit patients in critical care settings and those on outpatient oral anticoagulant therapy. Due to the desire to shorten the centrifugation time and the influence of the start-up and braking processes on the quality of the results, it is advisable to study the drive motor dynamics.
The article aimed to simulate the dynamic waveforms during the frequency start-up of the induction motor to drive the plasma centrifuge using the circuit model and analyze the influence of the frequency converter control parameters on the start-up dynamics. For this purpose, a mathematical model expressed in a rectangular coordinate system located on a complex plane was used, considering non-linear phenomena occurring in the motor, such as the phenomenon of saturation of the magnetic flux circuit and current displacement in the rotor bar. The presented method makes it possible to consider the skin effect by mapping the bar with three equivalent circuits. Very high compliance of the calculation results with the measurement results for both the direct motor start-up and the frequency start has been demonstrated.
The presented considerations can also be used to study the dynamic states of induction machines in other applications, for example, in electric vehicles.

The Method of Calculating the Dynamic Characteristics of Induction Motors
The basis for calculating the electromechanical behavior of an induction motor in transient operating states is the mathematical model of an induction motor in the form: where [L]-square matrix of windings inductance L ik , dimension [n × n], for i = k selfinductance of the windings and for i = k-mutual inductance of the windings.
To eliminate the dependence of coefficients of the mutual inductance of the stator and rotor on the position of the rotor about the stator, the three-phase coordinate system was transformed into a rectangular coordinate system (x, y, 0) rotating in space with the speed of the rotating field ω s , in which the stator and rotor are stationary concerning each other. In contrast, the movement of the rotor concerning the stator is accounted for by introducing an additional electromotive force of the rotation. In this case, the equation for the electrical coordinates takes the form: This model was placed on a complex plane, which means that the x-axis coincides with the real axis, and the y-axis with the imaginary axis, and then takes the form [20]: where the vectors of voltages and currents have the form: The individual components of the current and voltage vectors are complex numbers; ω k -speed of the coordinate system, which in the (x, y, 0) system is equal to the synchronous speed ω s , [1] = diag(1, 1, . . . , 1)-unit matrix, [K] = diag (0, 1, . . . , 1), N-number of equivalent circuits of the rotor bar.
Disregarding non-linear phenomena in the motor, such as the phenomenon of current displacement in the rotor bars and saturation of the magnetic circuit, we obtain a system of equations for electrical coordinates with dimension n = 4. The phenomenon of current displacement in the bars of the rotor cage can be considered by introducing the resistance and reactance of the bar depending on the value of the current frequency in the rotor bar, i.e., the slip, which is calculated for a bar of any shape using the elementary conductors' method and solving the linear system of equations, the parameters of which in the subsequent calculation steps are different values.
Such mapping of the bar impedance in dynamic calculations allows the accurate transient waveforms of the currents in the windings of the machine, but it results in obtaining significantly overestimated values of the transient electromagnetic torque [10]. This is because when calculating the moment in the transient state, the actual coupling of the magnetic flux of the stator with the fluxes of the individual layers of the rotor bar is not considered, which due to the phase shift of the currents in the individual layers of the bar are weaker than those determined at omitting this shift.
To take this phenomenon into account, it would be necessary to represent the rotor bar with a very large number of elementary conductors (N > 100), which would lead to the need to solve a system of differential equations with a very large dimension n and would significantly extend the time needed to solve the system thus obtained.
To reduce this model to several equivalent rotor circuits, a method based on the solution of the Riccati equation [21] was used, which determines the operator impedance of the rotor bar as a function of the spatial variable x along with the bar height. In this method, in a rotor bar of any shape with a unit length, we separate the layer dx at the x height as close as possible to the magnetic field lines ( Figure 1). Since the boundary lines of this layer are close to each other, we can assume that the current density in the dx layer is constant. The separated layer dx has cross-sectional areas dS x , electrical conductivity dγ x , and magnetic conductivity coefficient dλ x , where where c x -length of the midline of the field in the dx layer.
where x c -length of the midline of the field in the dx layer. Based on this equivalent diagram, we can determine the unit impedance of the bar part with the height x + dx, and as a result, derive the generalized Riccati equation as: where ( ) Z x is the equivalent impedance of the rotor bar of height x, the functions Ф1(x) and Ф2(x) are functions dependent on the resistance and inductance of the bar layer with the height dx, k 2 = j ω μ γ = ν, ω-current angular frequency in the rotor bar, μ-magnetic permeability of the bar material.
Assuming that for a homogeneous material γx = γ, the functions Ф1 (x) and Ф2 (x) are determined by the dependencies: The general solution of the Riccati equation can be obtained by expanding into a power series convergent around the point ν = 0, as a result of which we obtain the equivalent bar impedance in the form: This expression can be decomposed into a continued fraction, and as a result, a system of differential equations can be derived in the form: Based on this equivalent diagram, we can determine the unit impedance of the bar part with the height x + dx, and as a result, derive the generalized Riccati equation as: where Z(x) is the equivalent impedance of the rotor bar of height x, the functions Φ 1 (x) and Φ 2 (x) are functions dependent on the resistance and inductance of the bar layer with the height dx, k 2 = j ω µ γ = ν, ω-current angular frequency in the rotor bar, µ-magnetic permeability of the bar material.
Assuming that for a homogeneous material γ x = γ, the functions Φ 1 (x) and Φ 2 (x) are determined by the dependencies: The general solution of the Riccati equation can be obtained by expanding into a power series convergent around the point ν = 0, as a result of which we obtain the equivalent bar impedance in the form: This expression can be decomposed into a continued fraction, and as a result, a system of differential equations can be derived in the form: where in: with R and X are the rotor bar resistance and reactance assuming a uniform current density distribution.
The presented model is general and can be used for bars of any shape, for which the coefficients k 1 , k 2 for the calculated functions η 0 , η 1 , η 2 , η 3 have been determined. Simultaneously, a sufficient accuracy of calculations is obtained for three equivalent rotor bar circuits, the resistance and reactance of which are determined by the coefficients k 1 and k 2 , respectively. Ultimately, the equivalent diagram of an induction motor has the form shown in Figure 2. The presented model is general and can be used for bars of any shape, for which the coefficients 2 1 k , k for the calculated functions have been determined. Simultaneously, a sufficient accuracy of calculations is obtained for three equivalent rotor bar circuits, the resistance and reactance of which are determined by the coefficients k1 and k2, respectively. Ultimately, the equivalent diagram of an induction motor has the form shown in Figure 2.  [20]; Rs, Xs-resistance and leakage reactance of the stator winding, Xm-magnetizing reactance, Ra, Xa-resistance and leakage reactance of the bar part in which there is no skin effect, Rpr, Xpr-resistance and leakage reactance of the bar part surrounded by a ferromagnetic core, in which the skin effect occurs (these parameters are determined during design calculations). Reproduced from [20], Monograph of the Lodz University of Technology Publishing House: Lodz, Poland, 2011.
This approach makes it possible to consider the phenomenon of skin effect in rotor bars of any shape, which can be built by combining basic shapes for which the functions Ф1(x) and Ф2(x) can be easily determined, such as a trapezoidal rectangle or a lower and upper semicircle, with maintaining the continuity of functions on the border of component shapes.
The presented model has been implemented in the own DYNF program, enabling the determination of dynamic waveforms of an induction motor both when powered from the mains and an inverter. It also considers the influence of the saturation phenomenon on the motor winding leakage reactance by introducing a saturation factor depending on the resultant current flow of the motor windings for each step of the calculation.
It is assumed that only the parts of the core near the air gap are saturated. Thus, the reactance parts Xs and Xa corresponding to the stator and rotor slots opening reactance, differential reactance from higher harmonics of the field, and reactance from the skew of the rotor slots as a result of this phenomenon are reduced [12].
The presented model can be used to calculate various dynamic states, such as starting, overrun, reconnection, or motor reversal.  [20]; R s , X s -resistance and leakage reactance of the stator winding, X mmagnetizing reactance, R a , X a -resistance and leakage reactance of the bar part in which there is no skin effect, R pr , X pr -resistance and leakage reactance of the bar part surrounded by a ferromagnetic core, in which the skin effect occurs (these parameters are determined during design calculations). Reproduced from [20], Monograph of the Lodz University of Technology Publishing House: Lodz, Poland, 2011.
This approach makes it possible to consider the phenomenon of skin effect in rotor bars of any shape, which can be built by combining basic shapes for which the functions Φ 1 (x) and Φ 2 (x) can be easily determined, such as a trapezoidal rectangle or a lower and upper semicircle, with maintaining the continuity of functions on the border of component shapes.
The presented model has been implemented in the own DYNF program, enabling the determination of dynamic waveforms of an induction motor both when powered from the mains and an inverter. It also considers the influence of the saturation phenomenon on the motor winding leakage reactance by introducing a saturation factor depending on the resultant current flow of the motor windings for each step of the calculation.
It is assumed that only the parts of the core near the air gap are saturated. Thus, the reactance parts X s and X a corresponding to the stator and rotor slots opening reactance, differential reactance from higher harmonics of the field, and reactance from the skew of the rotor slots as a result of this phenomenon are reduced [12].
The presented model can be used to calculate various dynamic states, such as starting, overrun, reconnection, or motor reversal.

Dynamic Time Waveforms of the Induction Motor to Drive the Plasma Centrifuge
The motor for the drive of the plasma centrifuge should meet quite stringent requirements, that is, ensure smooth regulation of the rotational speed in an extensive range and the ability to work with variable load, depending on the work cycle, rotational speed, and weight of the drum suspended on the motor shaft. Therefore, this motor must have the required dynamic properties concerning the motor's acceleration to the highest possible rotational speed and braking in the shortest possible time. The paper presents the results of measurements and calculations of dynamic waveforms in the starting conditions of an energy-saving bipolar induction motor for the drive of the blood centrifuge with a shaft axis height h = 63 mm, with rated power 200 W, powered by 230 V, frequency 300 Hz, with the number of series turns of the stator winding N s = 150 and a core made of M270-35A sheet with a thickness of 0.35 mm. Other design parameters of the motor are given in Table 1.

The Method of Measuring Dynamic Parameters of the Motor for the Centrifuge Drive
The measurements of the start-up parameters were performed with various parameters of the frequency converter supplying the motor. To better observe the acceleration and the start-up time, the motor was loaded with an additional moment of inertia in the form of a disc mounted on the shaft, whose task was to simulate the load of the motor through the mounting elements and a set of containers containing the material to be centrifuged. To determine the instantaneous speed as a function of absolute time, the motor shaft's end was coupled to an angular position encoder with a total number of divisions of 8192 per revolution. Angle readings were taken at regular intervals with sampling frequency of 1600 Hz. The measuring system is shown in Figure 3. In that measuring system, currents and voltages were measured using precise current transducers and a multi-channel analog-to-digital converter of instantaneous currents and voltages. The measurement was controlled by a measuring system running under the control of the LabView program.

The Method of Measuring Dynamic Parameters of the Motor for the Centrifuge Drive
The measurements of the start-up parameters were performed with various parameters of the frequency converter supplying the motor. To better observe the acceleration and the start-up time, the motor was loaded with an additional moment of inertia in the form of a disc mounted on the shaft, whose task was to simulate the load of the motor through the mounting elements and a set of containers containing the material to be centrifuged. To determine the instantaneous speed as a function of absolute time, the motor shaft's end was coupled to an angular position encoder with a total number of divisions of 8192 per revolution. Angle readings were taken at regular intervals with sampling frequency of 1600 Hz. The measuring system is shown in Figure 3. In that measuring system, currents and voltages were measured using precise current transducers and a multi-channel analog-to-digital converter of instantaneous currents and voltages. The measurement was controlled by a measuring system running under the control of the LabView program. The application in the LabView environment was used to record the entire start-up run with a fixed time step. A total of 524,288 samples were collected on each channel during 32 s, so the time step was approximately 61 microseconds. Three voltages, 3 currents, and the signal from the shaft position sensor giving one pulse per revolution were recorded using a card with simultaneous sampling. After the registration, the data were saved in a disk file.
Estimated measurement uncertainty type B was for current measurement not worse than 1.2% and torque 2%.
The frequency start-up was carried out while maintaining a constant voltage-frequency ratio, while the supply system allowed for a slight increase in the initial voltage, which led to an increase in the flux value for low frequencies and an increase in the torque in the initial start-up phase, of course with an increase in the rotor and stator current.
To investigate the influence of the frequency converter control parameters affecting the start-up dynamics, a series of starts was performed, changing the two most important parameters of the supply frequency converter settings at the base frequency of 300 Hz and the voltage at this frequency equal to 230 V: initial voltage and frequency rise time. Comparing the dependence of the rotor acceleration time on the frequency rise time at the initial voltage equal to zero and on the initial value of the inverter voltage for the selected frequency rise time of 5 s is shown in Figure 4.
The frequency start-up was carried out while maintaining a constant voltage-frequency ratio, while the supply system allowed for a slight increase in the initial voltage, which led to an increase in the flux value for low frequencies and an increase in the torque in the initial start-up phase, of course with an increase in the rotor and stator current.
To investigate the influence of the frequency converter control parameters affecting the start-up dynamics, a series of starts was performed, changing the two most important parameters of the supply frequency converter settings at the base frequency of 300 Hz and the voltage at this frequency equal to 230 V: initial voltage and frequency rise time.
Comparing the dependence of the rotor acceleration time on the frequency rise time at the initial voltage equal to zero and on the initial value of the inverter voltage for the selected frequency rise time of 5 s is shown in Figure 4.  Based on the diagrams in Figure 4, it is possible to set the parameters of the frequency converter supplying the motor, ensuring the shortest rotor acceleration time while limiting the current surges in the stator windings.  Based on the diagrams in Figure 4, it is possible to set the parameters of the frequency converter supplying the motor, ensuring the shortest rotor acceleration time while limiting the current surges in the stator windings.   Parameters of the equivalent circuit of the investigated induction motor determined during design calculations are shown in Table 2. Parameters of the equivalent circuit of the investigated induction motor determined during design calculations are shown in Table 2. Table 2. Parameters of the equivalent circuit of the investigated induction motor.

Parameter Symbol Value
Stator winding resistance at 20 • (Ω)  As shown in Figures 5-8, the applied mathematical model gives excellent compliance of the calculation results with the measurement results. For direct motor start-up when powered by the voltage with the mains frequency (50 Hz), the discrepancy between the measured and calculated first current amplitude is (3.37-3.21)/3.37 × 100 = 4.75%, while the discrepancy between the measured and calculated first torque amplitude is (0.555 − 0.547)/0.547 × 100 = 1.45%. The measured and calculated rotational speed rise time is practically the same. Figures 9 and 10 show the measured and calculated motor start-up waveforms with a frequency start up to a frequency of 300 Hz, with an initial voltage of zero and a frequency rise time of 2 s for the motor with additional inertia mass (J = 0.00454 kg m 2 ). As shown in Figures 5-8, the applied mathematical model gives excellent compliance of the calculation results with the measurement results. For direct motor start-up when powered by the voltage with the mains frequency (50 Hz), the discrepancy between the measured and calculated first current amplitude is (3.37 -3.21)/3.37 × 100 = 4.75%, while the discrepancy between the measured and calculated first torque amplitude is (0.555 − 0.547)/0.547 × 100 = 1.45%. The measured and calculated rotational speed rise time is practically the same. Figures 9 and 10 show the measured and calculated motor start-up waveforms with a frequency start up to a frequency of 300 Hz, with an initial voltage of zero and a frequency rise time of 2 s for the motor with additional inertia mass (J = 0.00454 kg m 2 ).  (a) (b) Figure 9. Measured and calculated waveforms of the instantaneous value of the current in the stator phase winding (the entire course (a) and the initial part (b)) at the frequency start of the motor with the initial voltage equal to zero and the frequency rise time equal to 2 s, for the motor with additional inertia mass (J = 0.00454 kg m 2 ).    The comparison of the measured dynamic electromagnetic torque and the rotational speed at the initial voltage equal to zero and different times of frequency rise is shown in Figure 12 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). The comparison of the measured dynamic electromagnetic torque and the rotational speed at the initial voltage equal to zero and different times of frequency rise is shown in Figure 12 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). The comparison of the measured dynamic electromagnetic torque and the rotational speed at the initial voltage equal to zero and different times of frequency rise is shown in Figure 12 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). As shown in Figures 9-12, increasing the speed rise time with the same initial voltage value reduces the current and starting torque surges and extends the motor starting time. Extending the frequency rise time from 2 to 5 s extends the motor starting time by about 30% but at the same time reduces the maximum value of the starting current by about 35% and the maximum value of the electromagnetic torque by about 15% compared to the values obtained for the time equal to 2 s. Even greater extension of this time (up to 7 s) extends the motor start-up time by about 60%. It reduces the maximum value of the starting current by about 50% and the maximum value of the electromagnetic torque by about 35%. Figure 13 shows the measured motor start-up waveforms at the frequency start with the frequency rise time to 300 Hz, which is equal to 5 s and with different values of the initial voltage for the motor with additional inertia mass (J = 0.00454 kg m 2 ). As shown in Figures 9-12, increasing the speed rise time with the same initial voltage value reduces the current and starting torque surges and extends the motor starting time. Extending the frequency rise time from 2 to 5 s extends the motor starting time by about 30% but at the same time reduces the maximum value of the starting current by about 35% and the maximum value of the electromagnetic torque by about 15% compared to the values obtained for the time equal to 2 s. Even greater extension of this time (up to 7 s) extends the motor start-up time by about 60%. It reduces the maximum value of the starting current by about 50% and the maximum value of the electromagnetic torque by about 35%. Figure 13 shows the measured motor start-up waveforms at the frequency start with the frequency rise time to 300 Hz, which is equal to 5 s and with different values of the initial voltage for the motor with additional inertia mass (J = 0.00454 kg m 2 ). A comparison of the measured dynamic electromagnetic torque and the rotational speed with a frequency rise time of 5 s and a different value of the initial voltage is shown in Figure 14 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). A comparison of the measured dynamic electromagnetic torque and the rotational speed with a frequency rise time of 5 s and a different value of the initial voltage is shown in Figure 14 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). Figure 13. The measured waveforms of the instantaneous value of the current in the stator phase winding at the frequency start-up of the motor with a frequency rise time of 5 s and an initial voltage of 5% of the rated voltage (a) and 15% (b), for the motor with additional inertia mass (J = 0.00454 kg m 2 ).
A comparison of the measured dynamic electromagnetic torque and the rotational speed with a frequency rise time of 5 s and a different value of the initial voltage is shown in Figure 14 for the motor with additional inertia mass (J = 0.00454 kg m 2 ). Increasing the initial voltage value with the same speed rise time increases the current and starting torque surges in the initial starting period and practically does not affect the motor start-up time (Figures 13 and 14).
By analyzing the presented waveforms, it can be concluded that the shortest acceleration time of the motor is obtained with the shortest possible frequency rise time (on the order of 2 s) and with the initial voltage within the range of 0% to 5% of the rated voltage. This is a major requirement of the manufacturers of blood centrifuges. The increase of the initial value of the supply voltage has to be applied only in a higher load on the centrifuge drum. Increasing the initial voltage value with the same speed rise time increases the current and starting torque surges in the initial starting period and practically does not affect the motor start-up time (Figures 13 and 14).
By analyzing the presented waveforms, it can be concluded that the shortest acceleration time of the motor is obtained with the shortest possible frequency rise time (on the order of 2 s) and with the initial voltage within the range of 0% to 5% of the rated voltage. This is a major requirement of the manufacturers of blood centrifuges. The increase of the initial value of the supply voltage has to be applied only in a higher load on the centrifuge drum.

Conclusions
Comparing the simulation results with the measurement shows that the proposed model allows for high simulation accuracy for dynamic states both for the direct start-up at the mains frequency and for the frequency start-up. The obtained results allow for selecting the optimal control method from the point of view of the driven device, especially in the analyzed case of selecting the frequency rise time at frequency start-up and the initial voltage to obtain the optimal start-up time.
The developed model can also be successfully used to model other dynamic states, such as free coasting, reclosing, switching the number of pole pairs, voltage fluctuations and short-circuiting at the motor terminals, reversing, braking by counter-connection or dynamic, both in low and higher power motors.
It enables the improvement of the motor design for static operating states and dynamic states. In most cases, motors are designed to achieve the best efficiency at the selected frequency. The motor's work in a wide frequency range, and even more so the work in dynamic states, poses new challenges for the designer, in which the presented tool can be of significant help.
It should be noted that the models of dynamic phenomena used in control systems are often very simplified, as they do not consider the saturation aspect and changes in the rotor cage parameters caused by the skin effect in the rotor cage. As it is easy to see, the motor parameters during the start-up are significantly changed, so models with constant parameters are burdened with a significant error. Of course, if static control is considered for efficiency optimization, the models used must be extended to include core losses.