Singular Perturbation Theory-Based Qualitative Dynamics Investigation of Flywheel Energy Storage System in Discharge Mode

An investigation on qualitative dynamics in a voltage-current dual-loop controlled flywheel energy storage system (FESS) operating in discharge mode is presented in this paper, providing novel insights into the effect of two-timescale characteristics on the safety and stability of energy transmission of FESS. Based on singular perturbation theory, a two-timescale approach is proposed to separate the FESS into the fast and slow subsystems. Stability analysis of the transient fixed points confirms the effects of systemic parameters on FESS’s dynamics and indicates that the FESS shifts from the spiking state to the quiescent state when the slow variable crosses the bifurcation point of the fast subsystem. Mechanism analysis reveals that the root cause of the qualitative dynamics is the voltage instability of the FESS. Moreover, the feasibility boundaries of key parameters are derived, and application requirements of the proposed approach are also discussed, guiding the extension of the approach to engineering applications and solving the dynamics analysis problem to some extent at a theoretical analysis level. Constant voltage discharge experiment is performed based on the FESS test bench built in Key Laboratory of Smart Grid of Ministry of Education, Tianjin University, which validates the theoretical results.


Introduction
Permanent magnet (PM) brushless dc motor (BLDCM) controlled FESS, with the advantages of high density, low maintenance, long lifetime, and good compactness, has become a new trend for energy storage, applied more and more in the uninterruptible power supply, rail transportation, and smart grids [1,2]. As shown in Figure 1, the FESS mainly consists of two self-contained parts, that is, the mechanical part (flywheel motor system) and the electrical part (power drive system). In general, the flywheel motor is controlled by power electronic circuits. Due to the existence of the intrinsic nonlinearity, various nonlinear dynamics occur during the operation of FESS when the system state changes, which has an influence on the safety and stability of energy transmission.
So far, rotor dynamic problems of FESS's mechanical part have been seen as a cause of decreased rotor dynamic performance and reduced stability [3][4][5][6][7][8][9]. However, few researchers have revealed that, as a strongly coupled system, FESS has more complex nonlinear dynamics due to the interaction and difference between the mechanical part and the electrical part. This complex nonlinear dynamics has much more direct influence on the safety and stability of energy transmission and thus affects the safety and stability of FESS. Zhang et al. have investigated the nonlinear dynamics of FESS from the viewpoint of the interaction between the mechanical part and the electrical part [10]. Turning to the difference of the two parts we find that the electrical variables have significantly faster dynamics than the mechanical variables. As such, the dynamics of the voltage and current are faster than that of the rotate speed of flywheel, making the FESS be a typical two-timescale system [11]. It is well known that bursting phenomenon is observed when a slow variable controls the fast dynamics in some two-timescale systems such as neuronal systems and biological systems [12,13]. Bursting is a state of switching between the spiking state (SP) and the quiescent state (QS). Generally, QS indicates all the variables are at rest or exhibit small amplitude oscillations, In this paper, we aim to present evidences that as FESS operating in discharge mode, a small change in parameter values around the bifurcation points of FESS's fast system will lead to qualitative dynamics of the full-system, and investigate the effect of two-timescale characteristics on such dynamics, which is similar to nonrecurrent bursting. The analysis of dynamical systems with two timescales is a subject whose history interweaves three different viewpoints: nonstandard analysis [14], classical asymptotics [15], and geometric singular perturbation theory [16]. The first two methods lead to relatively large errors, and the geometric singular perturbation method is used to get the analytical solution of simple multiple timescale systems. The two-timescale approach proposed by Rinzel [17] is the classic approach to deeply investigate the two-timescale bifurcation dynamics, which gives a full description of the steady state, and periodic solution set of the fast subsystem, reflecting the global bifurcation structure of the fast subsystem with the slow variables treated as parameters [18]. First, stability analysis of the transient fixed points is proposed to study the bifurcation set of the fast subsystem, showing that as the slow variable varies, the fast subsystem loses stability from the originally stable state to Hopf bifurcation, and the dynamical evolution of the fullsystem is close in accordance with that of the fast subsystem. Not only reveling that the operation state of the FESS shifts when the slow variable crosses the bifurcation point of the fast subsystem, but also giving a way to predict the occurrence and evolution of qualitative dynamics of FESS in discharge mode. Then, the bifurcation mechanism analysis of the fast subsystem is proposed, offering an intuitive explanation of the origin of the nonrecurrent dynamics of the full-system. Furthermore, the feasibility regions are shown and provide instructions to parameters setting of FESS. Finally, the application requirements of the proposed approach are also discussed, guiding the extension of the approach to dynamics analysis of other electromechanical coupling systems. This paper is organized as follows. In Section 2, the normalized dynamic model of FESS is established and numerical simulations have been taken. In Section 3, the two-timescale approach based on singular perturbation theory is proposed and applied. A brief analysis of the application of the proposed approach is shown in Section 4. Also, the observed instability phenomena are observed experimentally, as presented in Section 5. Finally, Section 6 concludes this paper.

Modelling and Two-Timescale Characteristics of FESS
As shown in Figure 2, the modelling of FESS includes the flywheel motor system (the flywheel rotor driven by BLDCM) and the power drive system (the electrical subsystem and feedback control subsystem). The flywheel motor system is designed with a full bridge (IGBTs) at its output electrical terminals and DC-DC converter at the dc link. While diodes perform uncontrolled rectification, the DC-DC converter adjusts the voltage of the BLDCM in order to make it suitable for the load. Electronic commutation is achieved using a microprocessor-based controller with a Hall-effect position and a current sensor as input to generate gating signals for IGBTs.

Modelling of FESS.
Before modeling the FESS, five assumptions of the flywheel motor system are described as follows: (a) the saturation of the core is neglected; (b) the losses of eddy and hysteresis are ignored; (c) the distribution of air gap is uniform; (d) the self-inductance and mutual inductance among the windings are independent of the position of the rotor; (e) ignore the commutation process. The physical structure of the flywheel motor system is shown in Figure 3(a), while the schematic diagram is shown in Figure 3(b). It has been assumed that the phase resistance , the self-inductance , and mutual inductance of all the windings are equal. Assuming further that there is no change in the rotor reluctance with angle, hence, the circuit equations of the three windings in phase variables are   where , , and are phase voltages, , , and are phase currents, and , , and are the induced back EMFs. For + + = 0 and denoting − as , then the coupled circuit of the stator windings in terms of the machine electrical constants can be derived in Figure 3 . (2) According to the switching pattern (commutation function), only two phases are active at the same time, while the third one is silent (Figure 3(d)) [19]. For example, the voltage equation during the operation of phase -can be written as where , V , , = 2 , and = 2 denote line back-EMF constant, friction coefficient, moment of inertia, line resistance, and line inductance, respectively. Considering the equation of motion for BLDCM and treating the susceptibility and flux as constant, the differential equations for flywheel motor system are derived as where denotes torque constant and the angular velocity , the phase current , and the motor voltage 2 are state variables.  As shown in Figure 4, the electrical subsystem include a bidirectional DC-DC converter, which consists of a transfer inductor , capacitors 1 and 2 , IGBT 1 with diode 1 , and IGBT 2 with diode 2 and a three-phase full bridge converter which consists of IGBT-diode 1−6 . Mechanical energy is converted into electrical energy through the diodes of the three-phase full bridge converter when electrical energy is transmitted to the load via bidirectional DC-DC converter and regulated as load demands. During the discharge mode, the switching of and determines the two different modes, which are Mode 1 ( < < + ) when 2 is on and 1 is off and Mode 2 ( + < < ( + 1) ) when 2 is off and 1 is on. A dualloop proportional-integral scheme is applied to regulate the duty cycle of switches, of which the PI coefficients are ( V ,    normalized on the basis of = √ 1 , 1 ( ) = / ref , Then, the dimensionless switched dynamical equations of the present FESS are where , 2.2. Dynamic Characteristics. The operation principle of FESS is concerned with two processes: charge and discharge. In the charge mode, the flywheel motor system is driven by power grid and the electric energy is stored in the form of mechanical energy. Once the unit receives a signal of discharge, the flywheel rotor starts to decelerate and drives the BLDCM to generate electricity. This paper reports that the parametric regions have a significant effect on the twotimescale characteristics of FESS and the discharge performance degrades due to nonlinear dynamics when FESS's systemic parameters fall into a certain area.
Here, based on the exact state (5), a series of numerical simulations are carried out to make an initial evaluation of the possible dynamics. As shown in Table 3, the physical system parameters are set according to the actual FESS, and the PI control parameters of the feedback control subsystem are carefully designed to maintain closed-loop performance of the power drive system in spite of varying conditions. For energy storage unit as FESS, the stability in conjunction with the discharge targets is the primary consideration, so we restrict our attention to the case when ref   along with the decrease of the mechanical variables, which will be proved related to the qualitative dynamics in the following section.

Singular Perturbation Theory-Based Qualitative Dynamics Analysis
In this section, based on the singular perturbation theory, the two-timescale approach is proposed to separate the fullsystem into the fast and slow subsystems, providing a way for analyzing the interaction of the two-timescale dynamics.
Treating the slow variable as constant, stability analysis of the transient fixed points of the full-system is proposed to describe the evolution of dynamics of the full process.

Two-Timescale Model.
Averaged model approach is a common method to analyze the physical mechanism of converters which neglects the switching details but focuses on the envelope of the dynamical motion. The following analysis will adopt this approach. Set the duty ratio as . Set all the derivatives to zero, and we get where 1 ( ), 2 ( 1 ), 3 ( 2 ), 4 ( ), 5 (Ω), 6 ( 1 ), and 7 ( 2 ) are average values of 1 ( ), 2 ( 1 ), 3 ( 2 ), 4 ( ), 5 ( ), 6 ( 1 ), and 7 ( 2 ) during a duty cycle, ramp is the independent sawtooth peak voltage of the control subsystem, 13  proves the existence of state variables with two timescales. The full-system can be divided into fast subsystem and slow subsystem, that is, the fast subsystem which contains ( 1 ( ), 2 ( 1 ), 3 ( 2 ), 4 ( ), 6 ( 1 ), 7 ( 2 )) and the slow subsystem which contains 5 ( ). From (7) we can see 4 reflects the coupling effect between 5 ( ) and 4 ( ) and also the coupling effect between the fast and slow subsystems. Define 4 as coupling coefficient; the influence on dynamics of FESS from 4 will be studied later.

Hopf Bifurcation Set for Fast Subsystem.
To study the mechanism of the full-system's dynamics, the concept of transient where is the eigenvalue, I is the identity matrix, and J( ) is the Jacobian matrix at [ 1 , 2 , 3 , 4 , 6 , 7 ]. Parameters are set as Table 3, and loci of the eigenvalues are shown in Figure 8. When all the eigenvalues are in the left halfside of complex plane, the system is stable. When a couple of complex conjugate eigenvalues simultaneously cross the imaginary axis, Hopf bifurcation occurs [20,21]. From Figures 8(a)-8(c), as 5 changes from 0.9 to 0.1, the eigenvalues ( 3 , 4 , 5 , 6 ) stay on the left half-side of complex plane, while a pair of conjugate complex eigenvalues ( 1 , 2 ) are firstly on the right half-side and cross the imaginary axis as 5 reaches 0.4903, 0.4269, and 0.3774, respectively.

Studies on Qualitative
Dynamics. It would be imperative to know how the influence on the stability of the transient fixed points from the slow variable is reflected in the fullsystem. The study of the internal relations between the properties of the transient fixed points and the full-system's dynamics can predict the occurrence and evolution of nonlinear dynamics of FESS in the discharge mode. Among all the fast variables, the inductor current serves as a link between the access system and the flywheel motor system; thus, the stability and dynamic characteristics of influence the energy transmission process a lot. In this part we focus mainly on the dynamics of 1 ( ) against 5 ( ). Parameters are set as Table 3, and obviously in Figure 9, the transient fixed points curves (ES 1 , ES 2 , and ES 3 ) are L-form curves with various 5 ( ), which is divided in a stable part (dashed line) and an unstable part (solid line). Hopf bifurcation points (H 1 , H 2 , and H 3 ) are the joints of the two parts.
Considering the different kinds of equilibria, the stable node represents the quiescent state (QS), which indicates Hopf bifurcation point Stable fixed points Unstable fixed points all the variables are at rest or exhibit small amplitude oscillations. The stable limit cycle surrounding the unstable focus represents the spiking state (SP), which indicates variables may behave in large amplitude oscillations [22].  (7)). Between H 1 and A 1 , the difference between the fast and slow subsystems on timescales causes large amplitude oscillations around ES 1 from the beginning of the discharge process, leading the system to SP (from A 1 to H 1 ). The repetitive oscillation stays until the trajectory meets H 1 , at which Hopf bifurcation of the transient fixed point takes place. Using the center manifold theory, the curvature coefficient at the Hopf bifurcation point is less than 0, which proved the existence of supercritical bifurcation. Therefore, the amplitudes of the oscillations decrease gradually after the trajectory passing by H 1 and SP settles down to QS (from H 1 to B 1 ). The coupling strength of the fast and slow subsystems causes another type of small amplitude oscillations around ES 1 in QS. Then the phase trajectory reaches B 1 , at which the full-system becomes stable. When the FESS is in QS, the slow subsystem only influences the position of the transient fixed points but does not affect the dynamics of the full-system. For the reason that H 1 is the unique Hopf bifurcation point to join SP and QS, there is only one state conversion in a discharge cycle. The above process completes one period of the nonrecurrent qualitative dynamics.

Mechanism Analysis Based on Homotopy Method.
According to the analysis above, the transient fixed points curves obtained possesses Hopf bifurcation points at which the full-system can be divided into a stable part and an unstable part, and qualitative dynamics is closely bound up with the properties of the transient fixed points. Therefore, the mechanism analysis of the bifurcation of the transient fixed points can give an intuitive explanation of the origin of complex oscillations of the full-system. A state-to-eigenvalue correspondence can be set up to reveal the physical mechanism of the qualitative dynamics by tracing the changing trend of the eigenvalues. The homotopy method [24][25][26] can be applied to link eigenvalues of the Jacobian matrix J( ) of a dynamic model to the corresponding state variables through the following homotopy relation: where F( ) = diag[ 11 , 22 , 33 , 44 , 55 , 66 ] is the diagonal Jacobian and is the homotopy parameter. When varies in the interval [0, 1] and the difference between each two adjacent values is sufficiently small, homotopy method takes the trajectory of eigenvalues of H( ) as a continuous path. Following the paths from = 0 to = 1, the correspondence between eigenvalues and state variables of the fast subsystem is established. Parameters are set as Table 3, and Figure 12 shows the trace of the sorted eigenvalues by making as abscissa and the real part of eigenvalues of H( ) as vertical coordinate at ref = 10 V, 5 = 0.9. More traces are calculated and show the same correspondence, which is shown in Table 4. From the analysis above, Hopf bifurcation occurs when the root inducement of the qualitative dynamics. Therefore, dynamic voltage stability control scheme plays a key role in stabilizing the FESS, which will be the focus of our further work. A detailed look into the influencing factors on the qualitative dynamics is taken further by presenting the boundaries of SP and QS in terms of 2 and 1 . It shows that the region of QS gets smaller with the increase of 1 , which represents the load level ( Figure 13).

Feasibility Regions
Analysis. This part will apply the twotimescale approach to derive the feasibility regions of FESS in the discharge mode and then provide instructions to parameters setting of FESS. From the analysis above we know that the feasibility regions of the fast subsystem dominate that of the full-system. Therefore, the FESS is stable when all  Figure 14. Figure 14 and , all of which clearly illustrate the effect of those sensitive parameters on the feasibility regions. The space in front of the critical surface corresponds to stable operation and the space behind it corresponds to unstable operation. The results can be used as instructions to the parameters setting of the access unit of FESS itself and constraints to improve the safety and stability of FESS and the power system.

Application Requirements.
Obviously, the applicability and rationality of the proposed two-timescale approach with transient fixed points analysis mainly depend on the existence of state variables with two timescales, which is not to be considered as an exact criterion but as a guideline [27]. From model (7) we can see 7 = − ref √ 1 /( ref ) and 8 = − V √ 1 / represent the change rates of 5 ( ); thus, dominates the difference between the fast and slow variables on timescales. For typical electromechanical coupling systems, the two-timescale characteristics are ubiquitous but in degree. As is shown in the foregoing analysis, the proposed approach does well in predicting the qualitative dynamics when the magnitude difference of fast and slow variables is about 100 times. Set to 0.5, 0.25, and 0.1 times of its original value; other parameters are set as Table 3; then the corresponding magnitude difference will be 50, 25, and 10 times. Figure 15 shows the stroboscopic phase trajectory of the full-system in the phase plane 5 ( ) − 1 ( ) as well as the transient fixed points bifurcation diagram at ref = 10 V. At = 0.5 ⋅ in Figure 15(a), H 4 divides the phase trajectory of the full-system into two qualitatively different parts, which agrees very well with ES 4 . At = 0.25 ⋅ in Figure 15(b), ES 5 basically corresponds to the dynamics of the full-system. At = 0.1 ⋅ in Figure 15(c), ES 6 has deviated from the correct equilibrium position.
Considering that the two-timescale approach based on singular perturbation theory is a kind of model reduction method, Ghorbel and Spong [28] have given out the condition for the reduction of multitimescales system model; the equilibrium of the fast system must be close to that of the full-system. We can see that though the SP and QS caused by Hopf bifurcation are especially apparent in qualitative dynamics, they can still become weakening or even disappear as the difference of state variables on timescales diminishes. And meanwhile, the deviation between the transient fixed points curves and the stroboscopic phase trajectory of the full-system becomes larger, which reflects the inapplicability of the proposed approach. Therefore, the applicability and rationality of the proposed approach we concern here mainly refers to the bifurcation characteristics of the fast subsystem can actually reflect what extent of the full-system's dynamical evolution. From the numerical simulations in Figure 15 we can see the proposed approach is applicable when the magnitude difference of state variables is bigger than 25 times; more simulations have been done and gave out the same conclusion that when the magnitude difference is not big enough, the proposed approach is not applicable.

Experimental Verification
To verify the analysis in this paper, constant voltage discharge experiment is carried out; the parameter values are set as those in Table 3. The driving motor of the flywheel was solved by BLDCM, which indicates high reliability and the rotational speed being up to 8000 r/min. First charge up the FESS to 60% of its rated speed; then catch the time-variant dynamics under discharge mode at different ref . Set ref = 15 V; the trajectory of evolves with ; when 5 ( ) = 0.60 (60% of rated speed) in Figure 16(a), shows a large-scale oscillation with a frequency of 666.7 Hz and an amplitude of 4 A, implying the FESS is in SP; when 5 ( ) = 0.37 in Figure 16

Conclusion
This paper investigates the qualitative dynamics of the voltage-current dual-loop controlled FESS, which is mainly shown as the fast oscillations of the inductor current and the motor voltage weakens along with the slowdown of the flywheel rotor. By the proposed two-timescale approach based on singular perturbation theory, the state variables are separated into fast and slow variables. First, it is shown that the stability of FESS is closely bound up with that of the transient fixed points. The FESS shifts from SP to QS when the slow variable crosses the bifurcation point on the transient fixed points curve. Larger eigenvalues' real parts determine more intense oscillations in SP, and larger coupling coefficient leads to longer duration of QS. Further analysis shows that the evolution of the full-system's dynamics is dominated by the difference between the slow and fast variables on timescales, whereas the qualitative dynamics are mainly caused by the voltage instability. Moreover, the feasibility regions of the main system parameters are derived, in which stability operation and power transmission quality of FESS can be guaranteed. Finally, an applicability investigation shows that when the difference of state variables on timescales is not big enough, the proposed approach is not applicable. This paper provides insights into the effect of twotimescale characteristics on the safety and stability of energy transmission of FESS. The results can be used as instructions to the parameters setting of the FESS itself and constraints to improve the safety and stability of FESS and the smart grids.