**7. Conclusions**

In this work, we compared a proposed method of estimating the Floquet multiplier with time series data and a calculated Floquet multiplier based on a mathematical model. We described how to construct a Poincaré map in time series data, starting with the definition of the Poincaré map. We stressed the possibility of disagreement between the calculated and estimated Floquet multipliers by projection concept. A linearized form of the Floquet multiplier was provided to estimate the time series data with period.

We found that the estimated Floquet multipliers followed a similar trend as the Floquet multipliers calculated with the monodromy matrix in two power system networks after conducting correlation analysis. A Poincaré map selected by a local peak value searching algorithm provided enough information on an arbitrary system by considering their standard deviation. Thus, the critical multiplier of the estimated value was unity for the stable oscillation. For practical application with noisy measurement, we have conducted Floquet multiplier estimation to actual oscillation incident. Results show that estimated Floquet multiplier is not exactly unity due to load fluctuation or generator excitor response. However, still the estimated Floquet multiplier act as an indicator for feature of sustained oscillation while the value could stand for signification of the oscillation sources. Therefore, the proposed method can successfully function as a period oscillation indicator for time series data acquired from power systems.

**Author Contributions:** N.C. conceived and build up the research methodology, conducted the system simulations, and wrote this paper. B.L. and H.C. supervised the research, improved the system simulation, and made suggestions regarding this research.

**Funding:** This research received no external funding.

**Acknowledgments:** This work was supported by "Human Resources program in Energy Technology" of the Korea Institute of Energy Technology Evaluation and Planning (KETEP)-granted financial resource from the Ministry of Trade, Industry, and Energy, Republic of Korea (no. 20174030201820).

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Appendix A. Mathematical Model of System**

On the base of three-bus voltage collapse dynamics including induction motor at load bus provided by Chiang et al. [17], DFIG model of [18] is replaced to synchronous generator. Specific model of the system is given as follows:

*Appendix A.1. Aerodynamics of DFIG*

$$T\_{\rm m} = \frac{1}{2} \rho \mathbb{C}\_p(\lambda, \theta) \frac{\omega\_s}{\omega\_r} A\_{\rm wt} V\_{\rm wird}^3 \tag{A1}$$

$$\mathcal{C}\_p(\lambda, \theta) = 0.22(\frac{116}{\lambda} - 0.4\theta - 5)e^{-\frac{12.5}{\lambda}} \tag{A2}$$

$$
\lambda = (\frac{1}{\lambda' + 0.08\theta} - \frac{0.035}{\theta^3 + 1})^{-1} \tag{A3}
$$

where:

*ρ* : Air density [kg/m3] *Cp*(*<sup>λ</sup>*, *θ*) : Power coefficient *Awt* : Wind turbine swept area [m2] *Vwind* : Wind speed [m/s] *ωr* : Electrical rotor speed [rad/s] *ωs* : Electrical speed base [rad/s] *λ* : Tip speed ratio of a WTG

*θ* : Pitch angle *Appendix A.2. Differential Equations—DFIG*

$$T\_0' \frac{dE\_{qD}'}{dt} = -(E\_{qD}' + (X\_s - X\_s')I\_{ds})$$

$$+ T\_0' (\omega\_s \frac{X\_m}{\nu} V\_{dr} - (\omega\_s - \omega\_r)E\_{dD}') \tag{A4}$$

*Xr T* 0 *dEdD dt* = −(*EdD* − (*Xs* − *X s*)*Iqs*) *Xm*

$$+T\_0'( -\omega\_s \frac{X\_m}{X\_r} V\_{qr} + (\omega\_s - \omega\_r) E\_{qD}') \tag{A5}$$

$$\frac{2H\_{\rm D}}{\omega\_{\rm s}} \frac{d\omega\_{\rm r}}{dt} = T\_m - E\_{\rm dD}' I\_{\rm ds} - E\_{q\rm D}' I\_{q\rm s} \tag{A6}$$

*Appendix A.3. Differential Equations—Active and Reactive Power Controller*

$$\frac{d\mathbf{x}\_1}{dt} = K\_{I1}(P\_{ref} - P\_{Gen})\tag{A7}$$

$$\frac{d\mathbf{x\_2}}{dt} = \mathcal{K}\_{l2}(\mathcal{K}\_{P1}(P\_{ref} - P\_{\text{Gen}}) + \mathbf{x\_1} - I\_{qr}) \tag{A8}$$

$$\frac{d\mathbf{x}\_3}{dt} = \mathbf{K}\_{l3} (\mathbf{Q}\_{ref} - \mathbf{Q}\_{Gen}) \tag{A9}$$

$$\frac{d\mathbf{x\_4}}{dt} = K\_{I4}(K\_{P3}(Q\_{ref} - Q\_{Gen}) + \mathbf{x\_3} - I\_{dr})\tag{A10}$$

*Appendix A.4. Differential Equations—Load Dynamics(Induction Motor)*

$$\frac{d\delta\_L}{dt} = \frac{1}{k\_{qw}}(-k\_{qv2}V\_L^2 - k\_{qv}V\_L - Q\_0 - Q\_1 + Q) \tag{A11}$$

$$\begin{split} \frac{dV\_{\rm L}}{dt} &= \frac{1}{Tk\_{q\rm pv}k\_{pv}} (k\_{pw}k\_{q\rm p2}V\_{\rm L}^{2} + (k\_{pw}k\_{q\rm v} - k\_{q\rm p}k\_{pv})) \\ &+ k\_{pw}(Q\_{0} + Q\_{1} - Q) - k\_{q\rm pv}(P\_{0} + P\_{1} - P) \end{split} \tag{A12}$$

*Appendix A.5. Algebraic Equations*

$$-V\_{qr} + Kp\_2\{K\_{P1}(P\_{ref} - P) + \mathbf{x}\_1 - I\_{qr}\} + \mathbf{x}\_2 = 0\tag{A13}$$

$$\mathbf{V}\_{dr} + \mathbf{K}\_{P4} \{ \mathbf{K}\_{P3} (\mathbf{Q}\_{r\mathbf{r}f} - \mathbf{Q}) + \mathbf{x}\_3 - I\_{dr} \} + \mathbf{x}\_4 = \mathbf{0} \tag{A14}$$

$$I - P\_{\rm Ger} + E\_{dD}' I\_{ds} + E\_{qD}' I\_{qs} - R\_s \left( I\_{ds}^2 + I\_{qs}^2 \right) = 0 \tag{A15}$$

$$-Q\_{\rm Gen} + E\_{qD}' I\_{\rm ds} - E\_{\rm dD}' I\_{\rm qs} - X\_{\rm s}' (I\_{\rm ds}^2 + I\_{\rm qs}^2) = 0 \tag{A16}$$

$$I\_r - I\_{dr} + \frac{E\_{qD}^{\prime}}{X\_m} + \frac{X\_m}{X\_r} I\_{ds} = 0 \tag{A17}$$

$$-I\_{qr} - \frac{E\_{dD}^{\prime}}{X\_m} + \frac{X\_m}{X\_r} I\_{qs} = 0\tag{A18}$$

where:

*T* 0: Transient open-circuit time constant

−


*Xr* : Rotor reactance

*Xm* : Mutual reactance

*<sup>E</sup>qD*, *<sup>E</sup>dD*: Transient rotor voltage

*x*1, *x*2, *x*3, *x*4 : Active, reactive power controller

*δL*, *VL* : Load angle and voltage magnitude *P*0, *Q*0 : Constant active and reactive power of the motor

*P*1, *Q*1 : Constant active and reactive power of the load

*kp<sup>ω</sup>*, *kpv* : Constant impedance parameter related to active power with frequency and voltage

*kq<sup>ω</sup>*, *kqv*, *kqv*2 : Constant impedance parameter related to reactive power with frequency and voltage
