# 1. Introduction

Currently, the greenhouse effect and the energy crisis are seriously threatening the sustainable development of the human society. To tackle the intractable issue effectively, it is imperative to find new energy resources to replace the traditional fossil fuels. The good news is that more and more countries in the world are increasingly concerned about the development and utilization of renewable energy resources. As one of the renewable energy sources with the most mature technology and the largest application scale, wind power has been gaining more and more attention in virtue of its great contributions in reducing environmental pollution and adjusting the energy structure. Comparing with the traditional power generation technologies, the wind power is characterized with a nature of uncertainty and volatility. With a high penetration of wind power, how to systematically and effectively assess the impacts of intermittent wind power on power system dynamics, especially on the small signal stability (SSS) becomes one of the most challenging cutting-edge issues to be solved in exploring and exploiting large-scale wind power integration.

So far, much existing work corresponding to the impacts of wind power on the dynamic behavior of the power system mainly focused on the doubly-fed induction generator (DFIG) [1] which is currently widely used in wind farms. F. Mei et al [2] described the modeling and small signal analysis of a grid connected doubly-fed induction generator. The DFIG electromechanical models were analyzed by the eigenvalues and participation factor. F. Wu et al [3] derived the small signal stability model based on the DFIG and its control system. A conclusion with the improvement of damping characteristics when integrating wind farm of DFIG type was drawn. V. Vittal et al [4] proposed a comprehensive framework based on eigenvalue sensitivity to examine the effect of penetration of DFIGs on a large system. The case studies were carried out under the simulation environment DSATools developed by PowerTech Labs, Inc. Shi L. et al [5] discussed the impacts of wind power intermittency on power system small signal stability. The Monte Carlo technique was applied to reveal how the grid-connected wind farm of DFIG type affects the small signal stability of an existing power system. Comparing with DFIG, the direct-driven permanent magnet synchronous generator (PMSG) is considered as an important alternative in the future development of wind power due to its higher generation efficiency, gearless structure, lower maintenance cost and higher operation reliability. Espen H. et al [6] investigated the impacts of wind power integration in Norway on the damping of inter-area mode oscillations. Several types of wind generators involving direct drive synchronous generator were tested during simulations. Wu F. et al [7] presented a model of direct-drive permanent magnet generator (DDPMG) and its controllers. The corresponding small signal stability was carried out to show that the controller can improve the small signal stability of the wind turbine system. Kang L. et al [8] studied the impacts of grid-connected wind farm of direct-driven PMSG type on power system small signal stability. A simplified practical model of PMSG was derived. The IEEE 3-generator-9-bus system was tested to conduct the eigenvalue analysis. T. Knuppel et al [9] assessed the impact of full-load converter interfaced wind turbines on small signal stability. The effects of selected wind turbine and wind farm control parameters were studied based on sensitivity analysis. S. Q. Bu et al [10] proposed a method which can directly calculate the probabilistic density function of critical eigenvalues to investigate the impacts of wind generation uncertainty on small signal stability. H. Huang et al [11] built a small signal simulation model under MATLAB/SIMULINK environment to study the stability of grid-connected wind farm equipped with direct drive permanent magnet synchronous generator. H. Z. Huang et al [12] investigated the effects of plug-in electric vehicles and wind farm of DFIG type on power system small signal stability via applying quasi-Monte Carlo method. Compared to the conventional Monte Carlo simulation, more reliable results with shorter computational efforts can be obtained. M. Jafarian et al [13] developed a method to study the interaction of the dynamics of wind farm of DFIG type with the synchronous generators. The proposed method was based on the sensitivity of synchronous generators electromechanical eigenvalue with respect to variations in Jacobian matrix. It should be noted that the existing achievements in exploring and exploiting the impacts of wind power integration on power system SSS mainly focus on the deterministic evaluation, i.e. the nature of wind power uncertainty and volatility is rarely involved in evaluating the effects of wind generation integration on the dynamic behavior of power systems.

In this paper, a simplified practical model of PMSG and its associated controllers is derived. The Monte Carlo method based on the Weibull distribution of wind speed is applied in the solutions of IEEE 16-generator-68-bus test system as benchmark to reveal the impacts of gridconnected wind farm of PMSG type on power system SSS incorporating wind generation uncertainty and volatility.

This paper is organized as follows. Section 2 details the dynamic model of PMSG. In Section 3 a framework for the SSS evaluation incorporating wind power uncertainty based on Monte Carlo simulation is discussed. The case studies are given in Section 4. Finally, the conclusions are summarized in Section 5.

# 2. PMSG Modeling

In this section, a simplified practical electromechanical (EM) model of wind farm equipped with PMSG is derived elaborately. A general configuration of a PMSG including the power electronic converters and controllers is depicted in Fig. 1. The PMSG is connected to the power grid through a back-to-back power electronic converter.

**Fig. 1.**Typical PMSG configuration

## 2.1 PMSG model

In order to build the dynamic mathematical model of PMSG in d-q axis, the following assumptions are given during analysis: (i) Magnetic saturation effects are negligible; (ii) Eddy current loss and magnetic hysteresis loss are negligible; (iii) One-mass or lumped model is applied for the wind turbine drive train systems; (iv) Stator transients and stator resistance are negligible. With these assumptions mentioned above, and by applying the rotorflux-oriented control technique [7] in modeling PMSG, we have

where uds, uqs, ids and iqs are stator voltage and stator current in d-q axis, respectively; Ld and Lq are the sta tor inductance in d-q axis, respectively; ψf is the mag net flux linkage; ωe is generator electrical speed.

The active and reactive power outputs of the wind turbine generator (WTG) are given by

The corresponding swing equation is as

where H is the inertia constant; D is the damping coefficient; Tm and Te are the mechanical and generator electromagnetic torque, respectively; δ is defined as the angle in which q-axis leads x-axis. ωB is the speed base.

## 2.2 Converter model

The converter system consists of generator side converter, grid side converter, DC link and converter control. In normal operation, there is no reactive power exchange between generator side and grid side converter. Only active power is delivered by the generator side converter to the grid side converter. With the power losses of generator side and grid side converter neglected, the model of the DC link can be derived as

where C is the capacitor of the DC link; UDC is the voltage of C; udg, uqg, idg and iqg are grid side converter voltage and current feeding the grid in d-q axis.

Fig. 2 illustrates the control block diagram of generator side converter, which is utilized to control the active and reactive power of PMSG.

**Fig. 2.**Control block diagram of generator side converter

Based on the generator-voltage-oriented control strategy [8], the corresponding state equations can be derived below

where KP1, KP2, KP3, KP4, KI1, KI2, KI3 and KI4 are the PI controller gain and time constants respectively. x1, x2 ,x3, x4 are the introduced intermediate variables; Pref, Qref are the active power and reactive power control reference values, respectively; idsref, iqsref are the specified d-axis and q-axis stator current reference values, respectively.

Fig. 3 illustrates the control block diagram of grid side converter, which is utilized to maintain the DC link voltage and control the reactive power exchange between the generator side and grid side converter.

**Fig. 3.**Control block diagram of grid side converter

Based on the grid-voltage-oriented control strategy [8], the corresponding state equations can be derived below

where KP5, KP6, KP7, KP8, KI5, KI6, KI7 and KI8 are the PI controller gain and time constants respectively. x5, x6, x7, x8 are the introduced intermediate variables; Qgref is the reference value of the reactive power of the grid; idgref, iqgref are the specified d-axis and q-axis grid current reference values, respectively.

Eqs.(1-23) constitute the 11th order simplified practical electromechanical model of a PMSG system.

## 2.3 Wind farm model

It is known that a large wind farm consists of tens or hundreds of WTGs. It is different for each WTG to capture the wind energy from the incoming wind due to the wake effect [14]. In general, the details how to take wake effect into account for the whole wind farm are negligible during the system impact studies. The corresponding simplifications and assumptions are imperative.

In this paper, one equivalent single WTG model [15] representing the entire wind farm is employed to implement the modeling of wind farm.

## 2.4 Wind turbine model

The following functions approximation [1] with maximum power point tracking (MPPT) control strategy is widely used to model wind turbine.

where vw is the wind speed; Pm is the power output of the wind turbine; ρ is the air density (kg/m3); A is the area covered by the wind turbine rotor (m2); Cpmax is the optimal power coefficient; β is the blade pitch angle; λopt is the optimal tip-speed ratio; vcut-in and vcut-off are the cut-in and cut-off wind speed, respectively; vrated is the rated wind speed at which the power output of wind turbine will be the rated power.

Eq.(24) denotes the relationship between wind speed and the mechanical power extracted from the wind.

# 3. Monte Carlo Simulation Based SSS Evaluation Incorporating Wind Power Volatility

## 3.1 SSS evaluation incorporating wind farm of PMSG type

SSS is defined as the ability to maintain synchronism when subjected to small disturbances [16]. In accordance with the SSS theory, for a dynamic autonomous power system described by a set of non-linear ordinary differential and algebraic equations, the corresponding state matrix with appropriate linearization treatment can by derived as

where Λ is called as the state matrix; A, B, C and D are the Jacobian matrices of the test system.

The corresponding eigenvalues and eigenvectors of Λ can denote the stability of the test system at a specific operating condition.

For a power system with integration of wind farm equipped with PMSG, the deduction of state matrix can be briefly described as follows. According to the aforementioned derived 11th order simplified model of PMSG system, the state variables are ωe, δ, x1, x2, x3, x4, UDC, x5, x6, x7 and x8, respectively; the algebraic variables are uds, uqs, udg, uqg, ids, iqs, idg and iqg, respectively. Then the corresponding Jacobian matrices Apmsg, Bpmsg, Cpmsg and Dpmsg can be derived as

Where

By transforming the algebraic variables in d-q coordinate system to the synchronous rotating coordinate system (also called as x-y coordinate system), we have

where AWF=Apmsg; BWF=BpmsgTYWFWF; CWF=-T-1Cpmsg; DWF=I4×4-T-1DpmsgTYWFWF

Assume that a power system consists of n generators (wind farm involved as well as traditional synchronous generators) and N buses with m state variables. By eliminating ΔUxy descried in (27) and (28), finally we can obtain the system state matrix as given below

where A = [AG]m×m ; B=[BG 0]m×2N ;

## 3.2 Monte carlo based SSS evaluation incorporating wind power volatility

In order to simulate and study the impacts of wind power uncertainty and volatility on SSS more accurately and systematically, the Monte Carlo simulation technique is applied during analysis. In our work, the following Weibull distribution [17] is employed to describe the uncertainties of wind speed vw, the source of causing wind power uncertainty and volatility.

where k and c denote the shape and scale parameter, respectively. In this paper, the Weibull distribution with k=2 and c=12m/s as shown in Fig. 4 are selected as the probability distribution of regional wind speed. The sample size n is set to be 6000.

**Fig. 4.**Weibull distribution of wind speed (k=2, c=12m/s) and frequency distribution based on Monte Carlo method

The main procedures of implementing SSS evaluation with wind farm involved based on Monte Carlo simulation are as follows

(i) Set the total sample size for Monte Carlo simulation;(ii) Generate a random wind speed vw in terms of frequency distribution of wind speed, and verify if wind farm can be connected to power grid in accordance with the grid-connected condition (vcut-in

The flow chart of SSS analysis incorporating wind power based on Monte Carlo method is given in Fig. 5.

**Fig. 5.**Flow chart of SSS analysis incorporating wind power based on Monte Carlo method

# 4. Application Example

The 16-generator-68-bus test system [18] as shown in Fig. 6 is employed to assess the effects of wind generation uncertainty and volatility on SSS. In this system, generator G13 (slack bus) is modeled as the classic model; the remaining generators are modeled as the 4th order model with a simplified 3rd order excitation system model. Two simulation scenarios are designed below for the simulation studies: (i) In scenario 1, the wind farm equipped with PMSGs substitutes for the traditional synchronous generator; (ii) In scenario 2, the wind farm is connected to a load bus directly. Furthermore, the wind farm equipped with DFIGs modeled as a 7th order simplified model [5] is also introduced for comparisons during simulations. The details of DFIG modeling and corresponding parameters can be found in [5]. The parameters of PMSG are given in Table 1. The optimized parameters of PI controllers shown in Fig. 2 and Fig. 3 are given in Table 2 and Table 3 with respect to different scenarios. All simulations are executed under the MATLABTM environment.

**Fig. 6.**IEEE 16-generator-68-bus test system

**Table 1.**Parameters of PMSG

**Table 2.**Parameters of PI controllers applied in scenario 1

**Table 3.**Parameters of PI controllers applied in scenario 2

## 4.1 Scenario 1: Replacing synchronous generator

Firstly, the wind farms equipped with PMSGs and DFIGs substitute for all the synchronous generators except the reference machine G13, respectively. We found after the initial simulations that when the wind farm of PMSG type substitutes for the traditional synchronous generators G14 and G16, the test system keeps stable. When substituting for the remaining synchronous generators, the test system is unstable. When the wind farm of DFIG type substitutes for the synchronous generators G2-G4, G7 and G12, the test system is small signal stable. Furthermore, the simulation results show that whatever the wind farm of PMSG type or DFIG type is integrated into the power grid, they all do not participate in the EM oscillations. Table 4 gives the obtained EM oscillatory modes and their properties when G16 is replaced by the wind farms of PMSG type and DFIG type, respectively.

**Table 4.**The EM oscillatory properties when substituting for G16 in scenario 1

It can be seen from Table 4 that when G16 is replaced by the wind farm of PMSG type, there exist 15 EM oscillatory modes consisting of 13 local modes (0.7-2.5Hz) and 2 interarea modes (0.2-0.7Hz). When G16 is replaced by the wind farm of DFIG type, there exist 14 EM oscillatory modes consisting of 13 local modes and 1 interarea mode.

Table 5 shows the most relevant generators and the corresponding mode shapes when G16 is replaced by PMSG and DFIG, respectively.

**Table 5.**The most relevant generators and mode shapes when substituting for G16 in scenario 1

From Table 5, when G16 is replaced by the wind farm of PMSG type, the mode shape of EM15 represents a kind of special characteristics. In this EM mode, the magnitudes and angles of right eigenvectors with respect to the generator electrical speed bear relatively small difference, and no typical oscillations among generators can be found.

In the following studies, G14 is replaced by the wind farm consisting of 357 5MW PMSGs to specifically investigate the impacts of wind generation uncertainty and volatility on SSS via applying Monte Carlo method. Fig. 4 illustrates the Weibull frequency distribution histogram based on Monte Carlo with 6000 sample size.

With combination of (24), we can obtain the frequency distribution of power output of wind farm as shown in Fig. 7. From Fig. 7, two concentrations of probability masses can be found in the distribution. The left concentration with zero power output corresponds to the situation in which the wind farm is cut off with a probability of 6.22% (373/6000). The right concentration with the maximum power output corresponds to the situation in which the wind farm generates the rated power with a probability of 39.63% (2378/6000).

**Fig. 7.**Frequency distribution of wind farm power output

The distribution of calculated eigenvalues within the specific dimensions: real axis [-100,100] and imaginary axis [-200,200] based on the pre-set samples is given in Fig. 8. In Fig. 8, in accordance with the statistical results, the stable probability of the test is 49.28% (2957/6000). Fig. 9 shows the distributions of picked out total 15 EM oscillatory modes on the complex plane. It should be noted that if the wind farm is cut off (373 samples), the EM oscillatory mode EM14 does not exist. It can be seen from Fig. 9 that in the 5627 simulations with wind farm connected, except EM15, the remaining EM oscillatory modes EM1-EM14 fall into the left half complex plane (stable area). The stable probability for EM15 is 67.63% (4058/6000). In these two situations, the changes of EM1-EM4, EM6-EM9 and EM11-EM12 are relatively small.

**Fig. 8.**Distribution of eigenvalues within specific dimensions in scenario 1

**Fig. 9.**Distributions of EM oscillatory modes in scenario 1

Table 6 lists the statistical features of all EM oscillatory modes. From Table 6, it can be seen that the standard deviations in frequency with respect to EM5, EM10, EM13, EM14 and EM15 are relatively large, which denote that these oscillatory models are liable to be influenced by the wind speed.

**Table 6.*** Only for the situation of wind farm connected

## 4.2 Scenario 2: Grid connected directly

In this scenario, a wind farm consisting of 80 with single 5MW capacity PMSGs and a wind farm consisting of 200 with single 2MW capacity DFIGs are successively connected to the non-generator bus (bus 1-52), respectively. We found from the simulation results that when the wind farm of DFIG type is connected to bus 42, there exists an EM oscillatory mode with positive real part 0.0013±j3.5. When connecting to the remaining non-generator buses, the test system is small signal stable. For the wind farm of PMSG type, whichever non-generator bus is connected, the test system is small signal unstable. Similarly, the wind farm whatever is equipped with PMSG type or DFIG type does not participate in any EM oscillations. Table 7 gives the solved EM oscillatory modes and their properties when the wind farms of PMSG type and DFIG type are connected to bus 40 as shown in Fig. 6, respectively.

**Table 7.**The EM oscillatory properties when connecting to bus 40 in scenario 2

From Table 7, it can be found that when the wind farm of DFIG type is connected to bus 40, there exist 14 EM oscillatory modes, which is same as the number given in scenario 1. When the wind farm of PMSG type is connected to bus 40, there exist 16 EM oscillatory modes. Especially, there are 3 EM modes with negative damping ration, i.e. the integration of wind farm of PMSG type leads to a worsening of system damping. The most relevant generators pertinent to the EM modes and the mode shapes are summarized in Table 8.

**Table 8.**The most relevant generators and mode shapes when connecting to bus 40 in scenario 2

From Table 8, compared with the results given in scenario 1, the relatively big changes in the most relevant generators and the mode shapes happened when the wind farm equipped with PMSG or DFIG is connected to a specific non-generator bus. In summary, in accordance with the aforementioned simulation results, it is more suitable for the wind farm of DFIG type than the wind farm of PMSG type to be connected to the power grid directly.

In the succedent simulation studies, we mainly focus on the impacts of wind generation uncertainty and volatility on SSS via applying Monte Carlo method when the wind farm of PMSG type is connected to bus 40. The sample size for Monte Carlo simulation is set to be 5000.

Similarly, for the frequency distribution of wind farm power output, the left concentration with respect to the wind farm cut off is happened in a probability of 6.4% (320/5000); the right concentration with respect to the rated power output of wind farm is happened in a probability of 40.12% (2006/5000).

The distribution of calculated eigenvalues within the specific dimensions: real axis [-10, 10] and imaginary axis [-10,10] based on the pre-set samples is given in Fig. 10. It can be seen from Fig. 10 that some eigenvalues are located in the right half complex plane, i.e. the small signal stability is unstable.

**Fig. 10.**Distribution of eigenvalues within specific dimensions in scenario 2

According to the simulation analysis, total 16 EM oscillatory modes are identified. Fig. 11 illustrates the distributions of these EM oscillatory modes. Similarly, if the wind farm is cut off (320 samples), the EM oscillatory mode EM16 does not exist. In this situation, there exist only 15 EM oscillatory modes, which all fall into the left half complex plane. When the wind farm is connected (total 4680 samples) to the power grid, EM14 and EM16 fall into the right half complex plane with the probability of 72.74% (3637/5000) and 82.58% (4129/5000), respectively. The remaining EM oscillatory modes are all located in the left half complex plane. In these two situations, the changes of EM1-EM4 and EM7-EM11 are relatively small. The change of EM16 is relatively remarkable.

**Fig. 11.**Distributions of EM oscillatory modes in scenario 2

The characteristics of EM oscillatory modes involving damping ratio, oscillation frequency etc., are given in Table 9. From Table 9, it can be seen that these EM oscillatory modes from EM12 to EM16 are liable to be influenced by the wind speed. It should be noted that for EM6 and EM15, there exists a certain uncertainty for the most relevant generator. For example, EM6 is most relevant to G10 with a probability of 36.62%, and to G8 with a probability of 64.38%. The similar phenomenon is happened to EM15.

**Table 9.*** Only for the situation of wind farm connected

# 5. Conclusion

With large scale integration of wind power, how to evaluate the effects of wind power uncertainty and volatility on power system dynamics, especially on small signal stability is posing great challenges for the electrical engineers and scientists. In this paper, based on the modeling of wind farm equipped with PMSG, two scenarios involving the replacement of traditional synchronous generator with wind farm and direct connection of wind farm are designed and studied for the IEEE 16-generator-68-bus test system. From the simulation results, it can be found that there exist two concentrations for the distribution of wind farm power output due to the uncertainties of wind speed. And the probability of keeping test system stable is determined by some eigenvalues or wind speed conditions. Furthermore, the statistical analysis for the EM oscillatory modes further shows the effects of wind power uncertainty on power system small signal stability.