**Aggregation of Type-4 Large Wind Farms Based on Admittance Model Order Reduction**

#### **Jaime Martínez-Turégano \*, Salvador Añó-Villalba, Soledad Bernal-Perez and Ramon Blasco-Gimenez**

Instituto de Automática e Informática Industrial, Universitat Politècnica de València, 46022 Valencia, Spain; sanyo@die.upv.es (S.A.-V.); sbernal@die.upv.es (S.B.-P.); r.blasco@ieee.org (R.B.-G.) **\*** Correspondence: jaumartu@upv.es

Received: 10 April 2019; Accepted: 30 April 2019; Published: 7 May 2019

**Abstract:** This paper presents an aggregation technique based on the resolution of a multi-objective optimization problem applied to the admittance model of a wind power plant (WPP). The purpose of the presented aggregation technique is to reduce the order of the wind power plant model in order to accelerate WPP simulation while keeping a very similar control performance for both the simplified and the detailed models. The proposed aggregation technique, based on the admittance model order reduction, ensures the same DC gain, the same gain at the operating band frequency, and the same resonant peak frequency as the detailed admittance model. The proposed aggregation method is validated considering three 400-MW grid-forming Type-4 WPPs connected to a diode rectifier HVDC link. The proposed aggregation technique is compared to two existing aggregation techniques, both in terms of frequency and time response. The detailed and aggregated models have been tested using PSCAD-EMTsimulations, with the proposed aggregated model leading to a 350-fold reduction of the simulation time with respect to the detailed model. Moreover, for the considered scenario, the proposed aggregation technique offers simulation errors that are, at least, three-times smaller than previously-published aggregation techniques.

**Keywords:** off-shore wind farms; wind farm aggregation; admittance model order reduction; HVDC diode rectifiers; grid-forming wind turbines

#### **1. Introduction**

Large off-shore wind power plants (WPPs) with more than one hundred wind turbine generators (WTG) are currently operational or in the planning or approval phases. Simulation studies covering individual WTGs during the design phase are very demanding from the computational point of view. When studying the impact of WPPs in the overall HVDC or HVAC transmission network, very detailed models are not generally required. Therefore, simulation complexity and computing time can be reduced by using aggregated WTG models. However, these aggregated models should provide a faithful representation of the actual WPP dynamics.

This is particularly important if the WPP consists of grid-forming converters, such as when using diode rectifier (DR)-based HVDC stations for the connection of the offshore WPP. The use of DR HVDC stations can significantly reduce capital expenditure (CAPEX) and operational expenditure (OPEX) by increasing the efficiency and robustness of the overall system [1–6]. However, a large number of case studies is required in order to verify the correct integration of DR-connected WPPs.

Therefore, this paper is focused on the procedure to obtain an aggregated WTG model for large Type-4 WPPs that allows the verification of grid-forming WTG controllers for DR HVDC station connection. The presented aggregation technique aims at reducing the computational requirements of such a verification, while achieving better simulation accuracy than current state-of-the-art aggregation techniques.

A large amount of literature has been developed on aggregation techniques for different applications [7,8]. For example, the work in [9] presented an aggregation technique based on power losses, while [10] presented an aggregation technique for grid disturbances studies. The work in [11] included a method for Type-3 and Type-1 WTG aggregation suitable for large power systems simulation. In this case, the WPP grid was aggregated using a short-circuit impedance method. In [12], an aggregated method based on voltage drop assumptions was performed as the most appropriate aggregation method for stability assessment studies. On the other hand, in [13], an aggregation method was proposed in order to represent the whole power system seen from the WTG AC terminals. Finally, an improvement on the aggregation method based on power losses for off-shore WPP with a diode-based HVDC link was proposed in [14].

However, none of the previous methods were based on a model order reduction of the total WPP admittance model. The proposed model order reduction leads to better matching between the dynamic behavior of the aggregated and full WPP models, as the aggregated system keeps the original DC admittance, fundamental frequency admittance, and first resonance peak characteristics of the original system.

The paper is organized as follows. The second section shows the procedure to obtain the full admittance model of a typical WPP, consisting of a number of radially-connected WTG strings, as this is the typical topology used in WPPs. Other WPP topologies (i.e., with different degrees of meshing) can also be aggregated with the proposed method, by suitable modification of the analytical admittance calculation. The third section explains how to apply the proposed aggregation technique. The fourth section describes the system that is used in the fifth section to compare the proposed aggregation technique with other existing aggregation techniques found in the literature. The sixth section compares the aggregated and detailed model of a WPP using PSCADsimulations, considering three different cases. The last section includes the discussion of the proposed aggregation technique compared with other existing aggregation algorithms.

#### **2. Full-Order Admittance System Modeling**

The proposed model aggregation technique was based on the resolution of a multi-objective optimization problem applied to a full-WPP admittance model. This section includes the analytical development of the WPP admittance model. The inputs to the model were the individual WTG voltages, and its output was the WPP current at the point of common coupling (PCC).

Each considered WPP consisted of *n* radial-connected strings, with each string consisting of *mi* cascaded connected WTGs. This topology corresponds to the great majority of existing WPPs, as ring-connected strings are always exploited radially. In any case, it is very easy to modify the proposed analytical model to calculate the overall grid-admittance of other WPP topologies.

Figure 1 shows the WTG-*j* (*j* = 1, 2, ..., *mi*) in the string-*i* (*i* = 1, 2, ..., *n*) of the WPP. *Uij* is the WTG grid-side converter, and (*LTij* and *RTij*) represents the WTG transformer. Several alternatives have been considered to model the array cables, e.g., [15] added parallel L-R branches to a PI-model in order to represent the frequency-dependent impedance of the cable. This model has been used in [16,17] to obtain a state-space model and reduce the order of the cable model. Additionally, a frequency-dependent PI model of a three-core submarine cable was studied in [18].

In the presented case study, PI-sections will be used to represent the cable dynamics (*CLij*, *LLij*, and *RLij*), since the cable length between WTGs was less than 3 km with its resonant peak at frequencies higher than 1000 Hz. The use of PI-models offers a good trade-off between complexity and accuracy.

**Figure 1.** WTG-*j* (*j* = 1, 2, ..., *mi*) in the string-*i* (*i* = 1, 2, ..., *n*) of the wind power plant (WPP).

#### *2.1. String Admittance Model*

The dynamics of the WTG-*j* (*j* = 1, 2, ..., *mi*) and cable shown in Figure 1 can be written as:

$$dI\_{\bar{i}\bar{j}} = R\_{T\bar{i}\bar{j}} \cdot I\_{\bar{i}\bar{j}} + L\_{T\bar{i}\bar{j}} \cdot \frac{dI\_{\bar{i}\bar{j}}}{dt} + V\_{F\bar{i}\bar{j}} \tag{1}$$

$$I\_{Fi,j-1} + I\_{ij} - I\_{Fij} = \left(\frac{C\_{L,j-1}}{2} + \frac{C\_{Lij}}{2}\right) \cdot \frac{dV\_{Fij}}{dt} \tag{2}$$

$$V\_{Fij} = R\_{Lij} \cdot I\_{Fij} + L\_{Lij} \cdot \frac{dI\_{Fij}}{dt} + V\_{Fi,j+1} \tag{3}$$

where *Iij* is the current from WTG-*j*, *VFij* is the voltage at the secondary side of the transformer, and *IFij* is the current through the cable. Note that for *j* = 1 (leftmost WTG), the current *IFi*,*j*−<sup>1</sup> did not exist, and for *j* = *mi* (rightmost WTG), the voltage *VFi*,*j*+<sup>1</sup> was the PCC voltage (*VF*), which was considered a disturbance.

From Equations (1)–(3), the state-space admittance model of string-*i* (*i* = 1, 2, ..., *n*) can be calculated as:

$$\dot{\mathbf{x}}\_{st,i} = A\_{st,i} \cdot \mathbf{x}\_{st,i} + B\_{st,i} \cdot \mathbf{u}\_{st,i} + \mathbf{W}\_{st,i} \cdot V\_{\mathbf{F}} \tag{4}$$

$$y\_{st,i} = \mathbb{C}\_{st,i} \cdot x\_{st,i} \tag{5}$$

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

where the state vector *xst*,*<sup>i</sup>* and the input vector *ust*,*<sup>i</sup>* are:

$$\mathbf{x}\_{\mathrm{st},i} = \begin{bmatrix} I\_{\mathrm{i}1} & V\_{\mathrm{F1}} & I\_{\mathrm{F1}} & I\_{\mathrm{i}2} & V\_{\mathrm{F2}} & I\_{\mathrm{F2}} & \dots & I\_{\mathrm{i},\mathrm{m}\_{i}} & V\_{\mathrm{Fi},\mathrm{m}\_{i}} & I\_{\mathrm{Fi},\mathrm{m}\_{i}} \end{bmatrix}, \quad \mathbf{u}\_{\mathrm{st},j} = \begin{bmatrix} \mathsf{U}\_{\mathrm{i}1} & \mathsf{U}\_{\mathrm{i}2} & \dots & \mathsf{U}\_{\mathrm{i},\mathrm{m}\_{i}} \end{bmatrix}^{\mathrm{T}}$$

Matrices *Ast*,*i*, *Bst*,*i*, and *Cst*,*<sup>i</sup>* are the state, input, and output matrices, respectively:

$$A\_{sI,i} = \begin{bmatrix} \frac{-R\_{L1}}{L\_{L1}} & -\frac{1}{L\_{L1}} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ \frac{2}{L\_{L1}} & 0 & \frac{-2}{L\_{L1}} & 0 & 0 & 0 & 0 & 0 & 0\\ 0 & \frac{1}{L\_{L1}} & -\frac{R\_{L1}}{L\_{L1}} & 0 & \frac{-1}{L\_{L1}} & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & \frac{-R\_{L2}}{L\_{L2}} & \frac{-1}{L\_{L2}} & 0 & 0 & 0 & 0\\ 0 & 0 & \frac{2}{L\_{L1} + L\_{L2}} & \frac{2}{L\_{L1} + L\_{L2}} & 0 & \frac{-2}{L\_{L1} + L\_{L2}} & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & \frac{1}{L\_{L2}} & \frac{-R\_{L2}}{L\_{L2}} & 0 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & \frac{-R\_{L1m}}{L\_{L1m}} & \frac{-1}{L\_{L1m}} & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & \frac{2}{L\_{L1m} + 1} & 0 & \frac{-2}{L\_{L1m} - 1} & 0\\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{L\_{L1m} - 1} & -\frac{R\_{L1m}}{L\_{L1m}} \end{bmatrix}$$

$$B\_{st,i} = \begin{bmatrix} \frac{1}{L\_{11}} & 0 & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ 0 & \frac{1}{L\_{12}} & 0\\ 0 & 0 & 0\\ 0 & 0 & 0\\ & & \ddots\\ 0 & 0 & \frac{1}{L\_{1L\_i}}\\ 0 & 0 & 0\\ 0 & 0 & 0 \end{bmatrix}, \quad \mathbf{C}\_{st,i} = \left[ \begin{array}{cccc} 0 & 0 & 0 & 0 & 0 & 0 & \dots & 0 & 0 & 1 \end{array} \right]$$

and:

$$\mathcal{W}\_{st,i} = \begin{bmatrix} 0 & 0 & 0 & 0 & 0 & 0 & \dots & 0 & 0 & \frac{-1}{L\_{Li,m\_i}} \end{bmatrix}^T$$

#### *2.2. Overall Admittance Model*

The admittance model of the complete WPP was obtained by the radial connection of the *n* strings to the PCC (represented by voltage *VF*). The state space dynamics of the overall system are:

$$
\dot{\mathbf{x}}\_o = A\_o \cdot \mathbf{x}\_o + B\_o \cdot \mathbf{u}\_o + \mathbf{W}\_o \cdot \mathbf{V}\_F \tag{6}
$$

$$\mathbf{y}\_{\mathcal{O}} = \mathbb{C}\_{\mathcal{O}} \cdot \mathbf{x}\_{\mathcal{O}} \tag{7}$$

where *Ao*, *Bo*, and *Co* are the state, input, and output of the overall system obtained by combining the string state space equations:

$$\begin{aligned} A\_o &= \begin{bmatrix} A\_{st1} & 0 & 0 \\ 0 & A\_{st2} & 0 \\ & & \ddots \\ & & & 0 \end{bmatrix}, \quad B\_o = \begin{bmatrix} B\_{st1} & 0 & 0 \\ 0 & B\_{st2} & 0 \\ & & \ddots \\ 0 & 0 & B\_{st,n} \end{bmatrix} \\\\ \mathbf{C}\_o &= \begin{bmatrix} \mathbf{C}\_{st1} & \mathbf{C}\_{st2} & \dots & \mathbf{C}\_{st,n} \end{bmatrix} \\\\ & \mathbf{W}\_o = \begin{bmatrix} & \mathbf{W}\_{st1} & \mathbf{W}\_{st2} & \dots & \mathbf{W}\_{st,n} \end{bmatrix}^T \end{aligned}$$

and:

#### **3. Aggregation Technique**

The structure of the aggregated equivalent WPP is shown in Figure 2. At this point, the aggregation technique based on multi-objective optimization was carried out by following these steps:


#### *3.1. Input Reduction*

The proposed admittance model had a total of *n* voltage inputs that should be reduced to a single voltage input. Using the admittance model of the overall system, the superposition principle can be applied in order to obtain a single-input single-output model. The equivalent single-input voltage source was obtained by introducing the vector *σ* in Equations (6) and (7), leading to Equations (8) and (9):

$$
\dot{\mathbf{x}} = A \cdot \mathbf{x} + B \cdot \boldsymbol{\sigma} \cdot \mathbf{u} \tag{8}
$$

$$y = \mathbb{C} \cdot \mathfrak{x} + D \cdot \sigma \cdot u \tag{9}$$

*σ* is a vector with values of one for connected WTGs or zero for disconnected WTG:

$$
\sigma = \left( \sigma\_{1,1}, \sigma\_{1,2}, \dots, \sigma\_{i,j'}, \dots, \sigma\_{n,m\_i} \right)^T
$$

The obtained SISO system had a single equivalent voltage input (*u*) and the WPP PCC current as an output (*IF* = *IF*1,*m*<sup>1</sup> + *IF*2,*m*<sup>2</sup> + ... + *IFn*,*mn* ) and had a total of 3 · *n* states, with *n* equal to the number of WTGs in the WPP.

#### *3.2. Multi-Objective Optimization Problem Statement*

The target aggregated model is shown in Figure 2, with dynamics:

$$
\dot{\mathbf{x}}\_{eq} = A\_{eq} \cdot \mathbf{x}\_{eq} + B\_{eq} \cdot \boldsymbol{u}\_{eq} \tag{10}
$$

$$y\_{c\eta} = \mathbb{C}\_{c\eta} \cdot \mathbf{x}\_{c\eta} + D\_{c\eta} \cdot \mathbf{u}\_{c\eta} \tag{11}$$

where the inputs of the aggregated WTG model *ueq* were the aggregated voltage source (*Ueq*) and the voltage at the point of common coupling (*VF*). The model output was the WPP current (*IFeq*), and the state variables were *Ieq*, *VFeq*, and *IFeq*, as shown in Figure 2.

Matrices *Aeq*, *Beq*, *Ceq*, and *Deq* are given below:

$$\begin{aligned} A\_{eq} &= \begin{pmatrix} -\frac{R\_{Lq}}{L\_{Lq}} & \frac{-1}{L\_{Lq}} & 0\\\\ \frac{1}{\frac{C\_{Lq}}{2}} & 0 & \frac{-1}{\frac{C\_{Lq}}{2}}\\\\ 0 & \frac{1}{L\_{Lq}} & -\frac{R\_{Lq}}{L\_{Lq}} \end{pmatrix} & B\_{eq} &= \begin{pmatrix} \frac{1}{L\_{Lq}} & 0\\\\ 0 & 0\\\\ 0 & -\frac{1}{L\_{Lq}} \end{pmatrix} \\\\ \mathbf{C}\_{eq} &= \left( \begin{array}{ccc} 0 & 0 & 1 \end{array} \right) & D\_{eq} &= \left( \begin{array}{ccc} 0 & 0 \end{array} \right) \end{aligned}$$

The reduced order model in Equations (10) and (11) had five parameters to be identified (*LTeq*, *RTeq*, *LLeq*, *RLeq*, and *CLeq*). In order to formulate the multi-optimization problem properly, the following considerations were made:


At this point, three objectives are proposed to identify the remaining three parameters for the aggregated model, namely:


Finally, we shall consider that the WTG model parameters to identify always shall be positive and within a defined range (*RTeq*,*min* < *RTeq* < *RTeq*,*max*, *RLeq*,*min* < *RLeq* < *RLeq*,*max*, and *LLeq*,*min* < *LLeq* < *LLeq*,*max*).

From the assumptions above, the following multi-objective optimization function was defined for the aggregated model:

$$f(R\_{Tcq}, L\_{Lcq}, R\_{Lcq}) = \begin{pmatrix} \omega\_{peak}(R\_{Teq}, L\_{Lcq}, R\_{Lcq}) \\ K\_{\omega \simeq 0}(R\_{Teq}, L\_{Lcq}, R\_{Lcq}) \\ K\_{\omega \simeq \omega 0}(R\_{Teq}, L\_{Lcq}, R\_{Lcq}, \omega\_0) \end{pmatrix} \tag{12}$$

The goal function (*f* ∗) was obtained from the system described in Equations (8) and (9). Therefore, the aim of the optimization problem was to minimize the maximum of:

$$\mathcal{W}\_i(f(\mathcal{R}\_{\text{Teq}}, \mathcal{L}\_{\text{Leq}}, \mathcal{R}\_{\text{Leq}}) - f^\*) \tag{13}$$

where *Wi* is a diagonal matrix containing the weights to scale each function of *f*(*RTeq*, *LLeq*, *RLeq*). Additionally, the minimization problem shall satisfy the constraints:

$$
\begin{pmatrix} R\_{Teq,min} \\ L\_{Leq,min} \\ R\_{Leq,min} \end{pmatrix} < \begin{pmatrix} R\_{Teq} \\ L\_{Leq} \\ R\_{Leq} \end{pmatrix} < \begin{pmatrix} R\_{Teq,max} \\ L\_{Leq,max} \\ R\_{Leq,max} \end{pmatrix}
$$

The multi-objective optimization problem was solved by the goal attainment method described in [19–22]. The criteria used to choose the initial parameter and weight estimates is included in Section 5.1.

#### **4. HVDC Diode Rectifier-Connected Wind Power Plant**

The system under study is shown in Figure 3. It consisted of three 400-MW WPPs connected to the on-shore network through a diode rectifier-based HVDC link [1,14].

Each WPP consisted of 50 Type-4 WTGs rated at 8 MW each. Each WTG generator had a full-scale back-to-back converter, a PWM filter, and a WTG transformer that were connected to the 66-kV off-shore AC grid.

Each one of the three WPPs (*k* = 1, 2, 3) was composed of 50 WTGs, distributed in *n* = 7 strings (*i* = 1, 2, ..*n*) with *mi* wind turbines per string (*j* = 1, 2, ..*mi*).

Finally, the rectifier station consisted of three diode rectifier platforms, DC-side series connected. The rectifier stations were connected by means of the HVDC cable to the onshore VSC station, as shown in Figure 3. Each DR station consisted of two 12-pulse DR bridges, together with the corresponding transformers and AC filters [14].

The implemented control is shown in Figure 4. The control consisted of a centralized controller for the total delivered active power *POWPP*:

$$P\_{OWPP} = \sum\_{k=1}^{3} \sum\_{i=1}^{n} \sum\_{j=1}^{m\_i} P\_{WTG,ijk}$$

and a distributed controller based on a standard *P*/*ω*, *Q*/*V* droop controller.

**Figure 3.** Off-shore WPP with three WPPs of 400 MW each connected to the on-shore grid through a diode-based HVDC link.

**Figure 4.** Left: centralized power control; right: distributed droop controls.

The centralized controller was based on a PI controller operating at a 20-ms sampling time. Using the total optimum power reference from the WTGs (*P*∗ *OWPP*) and the actual power being generated by the WPP (*POWPP*), the controller calculates the reference voltage (*U*<sup>∗</sup> <sup>0</sup> ) to be used by each individual WTG droop controller. *ω*∗ *<sup>F</sup>* is the grid frequency reference. A communication delay of 40 ms was considered for the WPP active power controller.The WPP active power controller had a bandwidth of 4 Hz.

Each WTG implemented a distributed *P*/*ω* and *Q*/*V* droop controller. The droop controller calculates the frequency *ωFijk* and the voltage *Uijk* of the WTGs grid side converters (GSC) using the measured and reference active power *Pijk* and the reactive power *Qijk*, respectively. Local *Pijk* and *Qijk* measurements were filtered with a 10-rad/s low pass filter. More advanced controllers can also be implemented for *P* and *Q* sharing amongst wind turbines [23].

The AC-cable parameters, as well as all those of the system in Figure 3 are listed in the Appendix A.

#### **5. Wind Power Plant Aggregation**

This section includes the application of the proposed technique to the described case and its comparison with two well-established aggregation strategies.

#### *5.1. Aggregation Based on Multi-Objective Optimization*

The state space admittance model of each one of the wind power plants was obtained following the procedure explained in Section 2, including input reduction with *σ* equal to a vector of ones. The state space model thus obtained was of order 150.

Before carrying out the multi-objective optimization, it is important to ensure that the full system can be reduced to a system with only three states (according to Equations (10) and (11)). To this avail, the contribution of each state has been analyzed obtaining the Hankel singular values of the system. Figure 5 shows that the dynamics of the complete system are mainly dominated by three states; therefore, a third order reduced system can effectively capture most of the dynamic behavior of the full-order system.

**Figure 5.** Hankel singular values (state contributions; the first 20 states are shown).

To carry out the multi-objective optimization problem, the weights *Wi* and the restrictions on *RTeq*, *LLeq*, and *RLeq* need to be defined. To set the weights *Wi* properly, it is important to notice that the first objective used in the optimization function (12) is related to a frequency, while the other two are related to gains. Therefore, the weight *Wi* related to the first optimization objective is 1000-times higher than the other two weights. The selected weights *Wi* are shown in the Appendix A.

The parameters *RTeq*, *LLeq*, and *RLeq* were restricted to be always positive. Moreover, their limits can be further refined. Restriction on the equivalent resistance can be obtained from the objective *K*∗ *<sup>ω</sup>*=0. Analyzing the circuit of Figure 2, we know that the sum of *RTeq* + *RLeq* shall be equal to 1/*K*<sup>∗</sup> *<sup>ω</sup>*=0. Therefore, the maximum value of the addition of both resistors was 1/*K*∗ *<sup>ω</sup>*=0.

The restrictions on *LLeq* were calculated by considering the simplified expression for the cable's first resonant frequency:

$$
\omega\_{\rm res} = \sqrt{\frac{1}{LC}} \to L\_{\rm ini} = \frac{1}{\mathcal{C}\_{\rm eq} \omega\_{\rm res}^2} \tag{14}
$$

where *Ceq* is the total array cable capacitance and *ωres* is the first cable resonance frequency obtained from the detailed state space admittance model. Therefore, the restrictions for *LLeq* were set to [0.01*Lini*, 100*Lini*]. The selected parameter ranges are shown in the Appendix A.

With all the restrictions set, the multi-objective optimization problem was solved, and the obtained parameters are shown in Table 1. Figure 6 shows the admittance of the obtained aggregated model as a function of frequency.

**Table 1.** Equivalent parameters obtained from different aggregation methods.


**Figure 6.** WPP grid admittance as a function of frequency: blue: 150 order system; red: multi-objective optimization; yellow: voltage drop; purple: power loses.

#### *5.2. Aggregation Based on Voltage Drop*

The voltage drop aggregation method presented in [12] was proposed in order to simplify frequency response analysis for WTG interaction with the grid. The following equations are a summary of the aggregation method in [12]. When WTGs are connected in a string, the total equivalent string impedance is:

$$Z\_{ik} = \frac{1}{m\_i} \sum\_{j=1}^{m\_i} j Z\_{ijk} \tag{15}$$

where *Zik* is the equivalent R-L impedance corresponding to the *i*th string of the *k*th WPP and *Zijk* is the cable R-L impedance between WTGs *j* and *j* + 1. Equation (16) was used for the parallel connection of several strings.

$$Z\_k = \frac{1}{\sum\_{i=1}^n \frac{1}{Z\_{ik}}} \tag{16}$$

where *Zk* is the equivalent impedance of WPP *k*. The aggregated capacitance of the system was obtained as follows:

$$\mathbf{C}\_{k} = \sum\_{j=1}^{m\_i} \sum\_{i=1}^{n} \mathbf{C}\_{ijk} \tag{17}$$

where *Cijk* is the cable capacitance between WTGs *j* and *j* + 1. Additionally, the WTG transformer impedance and PWM filter were scaled by using their per unit values and taking into consideration how many WTGs were in operation.

Following these steps, an aggregated model of the system was calculated. The obtained values are shown in Table 1.

#### *5.3. Aggregation Based on Power Losses*

WTG aggregation based on power losses was used for load-flow studies and also in electromechanical (RMS) stability. Moreover, this aggregation technique has been applied successfully in many real-life wind farm projects [7].

The power loss aggregation technique proposed in [9] was based on the following equations:

$$Z\_{ik} = \frac{1}{m\_i^2} \sum\_{j=1}^{m\_i} j^2 Z\_{ijk} \tag{18}$$

$$Z\_k = \frac{\sum\_{i=1}^n m\_i^2 Z\_{ik}}{[\sum\_{i=1}^n m\_i]^2} \tag{19}$$

$$\mathbf{C}\_{k} = \sum\_{j=1}^{m\_i} \sum\_{i=1}^{n} \mathbf{C}\_{ijk} \tag{20}$$

where *Zik* is the equivalent R-L impedance corresponding to the *i* th string of the *k*th WPP, *Zijk* is the R-L cable impedance between WTGs *j* and *j* + 1, *Cijk* is the cable capacitance between WTGs *j* and *j* + 1, and *Ck* is the total capacitance of WPP *k*.

In the same way as the voltage drop-based aggregation, WTG transformer impedance and PWM filter were scaled from their pu values, taking into consideration how many turbines were in operation.

The values obtained for this aggregation method applied to the considered case study are shown in Table 1.

#### *5.4. Analytical Frequency Response of Different Aggregation Techniques*

Table 1 shows the resulting parameters of the aggregated models obtained with the proposed aggregation technique and with two alternative state-of-the-art techniques.

Figure 6 shows the WPP grid admittance as a function of the frequency of the complete 150-state SISO model, as well as that of each one of the aggregated models (multi-objective optimization, voltage drop, and power losses' aggregation). All aggregation techniques showed a relative good agreement with the detailed model admittance for a wide range of frequencies.

However, the proposed technique showed a much closer match of low frequency admittance values and, more importantly, an excellent match of both main resonant peak frequency and amplitude. It is worth noting that both voltage drop- and power loss-based aggregation showed main resonant peak frequencies that were far off from the actual resonant peak. The largest error was obtained from the power loss error, which led to a resonant peak frequency of more than 500 Hz higher than its actual value.

Additionally, different WPP topologies have been studied, in order to verify that the proposed aggregation technique is valid for different WPP configurations. The studies covered different numbers of connected strings and also different numbers of connected WTGs per string. In all cases, the multi-objective aggregation technique always showed a better admittance frequency response match than the other two considered aggregation techniques.

Clearly, the results in Figure 6 represent small signal operation; however, the proposed aggregation technique had a very good agreement with the actual grid admittance for a very wide range of frequencies, including the main resonant peak. Therefore, the proposed aggregation technique is expected to show dynamic characteristics very close to the detailed 150-state system. The dynamic response of the proposed and state-of-the-art techniques is covered in Section 6.

#### **6. Results**

The frequency results in Figure 6 suggest that the proposed aggregation method should be able to represent the full system dynamics with less error than the other two considered aggregation techniques. Therefore, the EMTsimulation of the considered test cases was carried out in order to verify the dynamic behavior of each aggregation method.

Four scenarios have been considered, one with full detail models of each cluster (considering the full 150 WTGs of 8 MW each), and considering an aggregated equivalent for each one of the three 400-MW clusters (with multi-objective optimization, voltage drop, and power losses' aggregation methods).

The PSCAD-EMT simulations for the system shown in Figure 3 have been carried out in all four scenarios considered. However, for the sake of clarity, only the full simulation results corresponding to the detailed and multi-objective aggregated model are shown in the simulation section. The results corresponding to the other two aggregation methods are shown as errors with respect to the detailed model simulations.

The validation of the aggregated models has been carried out by comparing the active and reactive power and voltage and current (*PFk*, *QFk*, *VFk*, *IFk*) responses at the points of common coupling of each one of the wind power plants (buses PCC1,2,3 in Figure 3).

Three test cases have been considered:


Finally, a study of the computational load of each alternative has been carried out.

#### *6.1. Simulation Results*

#### 6.1.1. Case 1: Fast Changes in Active Power Reference

The three wind power plants were initially operated in islanding mode (DR not conducting) with a voltage reference of 0.87 pu. From this state, the active power reference was ramped up from 0–1 pu in 150 ms for all three WPPs at the same time. Once all the WPPs reached their steady state operation, several step transients on active power reference were applied from 1 down to 0.5 pu, by sequentially reducing the active power reference *P*∗ *Fk* of each individual wind power plant.

Figure 7 shows the behavior of the active power *PFk*. Initially, the active power reference *P*<sup>∗</sup> *Fk* rose to 1 pu at *t* = 0.2 s, but power production did not increase until currents began to flow through the DR at *t* = 0.25 s. The reason for this is that the AC-side DR voltage had to be higher than 0.9 pu for the DR station to start conducting. Then, at *t* = 0.8 s, the active power reference *P*∗ *Fk* was decreased to 0.5 pu in each cluster every 0.7 s. Such rapid active power transients might not be realistic, but have been chosen as an extreme case for validation purposes. Figure 7 clearly shows that both aggregated and detailed models led to very similar active power dynamic results.

**Figure 7.** Case 1—Active power *PFk* at PCC-k of each diode rectifier platform for detailed and aggregated models.

Figure 8 shows the reactive power responses for each WPP *QFk*. The reactive power produced by the AC filters was compensated by the WTGs while DRs were not delivering power. Conversely, when the DRs started conducting, the capacitor and AC filter banks' reactive power compensated that absorbed by the DRs, and hence, the reactive power *QFk* delivered by the wind farms was reduced to a value very close to zero. Reactive power (*QFk*) responses were also very similar for both aggregated and detailed models.

**Figure 8.** Case 1—Reactive power *QFk* at PCC-k of each diode rectifier platform for the detailed and aggregated models.

Figure 9 shows the behavior of the voltage *VFk* at the PCC of each WPP. It shows clearly how active power flow by DRs depends of the voltage *VFk*. Figure 10 shows the magnitudes of current *IFk* at the PCC of each WPP. The current behavior was similar to that of the active power *PFk*. Figures 9 and 10 clearly show that the voltage and current dynamics from the aggregated and detailed models agree to a great extent.

**Figure 9.** Case 1—Voltage amplitude at the PCC of each WPP (*VFk*) for detailed and aggregated models.

**Figure 10.** Case 1—Current amplitude at the PCC of each WPP (*IFk*) for detailed and aggregated models.

Three aggregation methods have been used to simulate the fast active power transients considered in Case 1. The results are shown in Figure 11, which shows the active and reactive power simulation errors between aggregated and detailed models at PCC-1, i.e., (*PF*<sup>1</sup>*Detailed* − *PF*<sup>1</sup>*Agg*) and (*QF*<sup>1</sup>*Detailed* − *QF*<sup>1</sup>*Agg*). Only the traces corresponding to PCC-1 are shown, as the other PCCs showed a very similar behavior.

The simulation error obtained using the proposed multi-objective optimization aggregation technique (in blue in Figure 11) was clearly lower than those with the alternative aggregation techniques.

The maximum active power simulation error with the proposed technique was about 0.3%, whereas the maximum simulation error for the reactive power simulation was approximately 0.2%. For the power loss technique, the active and reactive power maximum simulation errors were 0.9% and 0.5%, respectively, whereas for the voltage drop technique for the maximum simulation, the corresponding errors were 1.2% and 0.7%. Therefore, the proposed technique clearly showed a more accurate simulation during transients.

Moreover, for the complete case, the error variance using the multi-objective optimization aggregation technique was 0.0423, whereas the error variance with the voltage droop and power losses techniques was 0.4094 and 0.3120, respectively. Therefore, the error variance for the proposed technique was up to 9.5-times better than the other aggregation techniques.

**Figure 11.** Case 1—Simulation errors for different aggregation techniques: top: *PF*<sup>1</sup>*Detailed* − *PF*<sup>1</sup>*Agg*; bottom: *QF*<sup>1</sup>*Detailed* − *QF*<sup>1</sup>*Agg*.

#### 6.1.2. Case 2: Wind Power Plant Disconnection

The second test case consisted of the disconnection of WPP3, by opening the breaker at *VF*<sup>3</sup> in Figure 3 at t = 0.05 s. Initially, all WPPs were generating 0.5 pu rated power. Figures 12 and 13 show the behavior of the WPP voltage and current magnitude at the PCC of each wind power plant.

**Figure 12.** Case 2—Comparison of detailed and proposed aggregated simulations. WPP voltage *VFk*.

**Figure 13.** Case 2—Comparison of detailed and proposed aggregated simulations. WPP current *IFk*.

After disconnection, WPP3's voltage remained at its 0.8 pu reference value, as all considered WTGs were grid-forming (Figure 12). After the transient, the voltage of the WPPs that remained connected settled to a slightly smaller voltage, as the power transmitted through the diode rectifiers was now reduced by one third [5].

From Figures 12 and 13, it is clear that voltage and currents obtained from the detailed and from the proposed aggregated models showed an excellent agreement during the transient.

The same WPP3 disconnection case shown in Figures 12 and 13 has been repeated considering two alternative aggregation techniques (power loss and voltage drop). The results are shown in Figures 14 and 15, which show the difference between the detailed simulation and each one of the three considered aggregation techniques. Simulation errors are expressed as per unit with respect to rated values.

**Figure 14.** Case 2—Simulation errors for different aggregation techniques: top: *VF*<sup>1</sup>*Detailed* − *VF*<sup>1</sup>*Agg*; bottom: *IF*<sup>1</sup>*Detailed* − *IF*<sup>1</sup>*Agg*.

**Figure 15.** Case 2—Simulation errors for different aggregation techniques: top: *VF*<sup>3</sup>*Detailed* − *VF*<sup>3</sup>*Agg*; bottom: *IF*<sup>3</sup>*Detailed* − *IF*<sup>3</sup>*Agg*.

Figure 14 shows the errors on the simulated voltages and currents for the detailed and aggregated models, (*VF*<sup>1</sup>*Detailed* − *VF*<sup>1</sup>*Agg*) and (*IF*<sup>1</sup>*Detailed* − *IF*<sup>1</sup>*Agg*) respectively. These voltages and currents correspond to one of the two WPP that was not disconnected (WPP1). The proposed aggregation technique (in blue) clearly showed a better performance than existing aggregation techniques, although all aggregation techniques showed relatively small *VF*<sup>1</sup> and *IF*<sup>1</sup> simulation errors.

However, the voltage and current transients were much larger in the disconnected WPP3, as previously shown in Figures 12 and 13. The top graph of Figure 15 shows that the *VF*<sup>3</sup> simulation error for the proposed aggregation method was much smaller than that of the existing aggregation methods (which even reached simulation errors larger than 10%).

On the other hand, the *IF*<sup>3</sup> simulation errors were very small for all aggregation techniques, with only a small deviation just at the instant of WPP3 disconnection (Figure 15 bottom graph).

#### 6.1.3. Case 3: Disconnection and Re-Connection of Diode Rectifier AC Filters

The third case consisted of the disconnection and re-connection of the AC capacitor and filter banks of DRU3(Figure 3). The capacitor and filter banks of each DRU were rated at 0.4 pu (160 MVAr) and were connected in several steps. However, this transient considered the connection and disconnection of the full AC filter banks, in order to validate the dynamic performance of the aggregated models.

The voltages and current output of each WPP during this test case are shown in Figures 16 and 17. The DRU3 AC filter banks were disconnected at t = 0.1 s and re-connected at t = 0.8 s. Clearly, the harmonic contents of both voltage and current increased when the filter banks were disconnected. Moreover, voltage and current distortion was larger for WPP3, as it was electrically closer to the DRU station with no AC filters connected. WPP1 showed a relatively small voltage and current ripple when the filters were disconnected, whereas WPP2 showed intermediate harmonic contents.

**Figure 16.** Case 3—Comparison of detailed and proposed aggregated simulations. WPP voltage *VFk*.

**Figure 17.** Case 3—Comparison of detailed and proposed aggregated simulations. WPP voltage *IFk*.

In any case, there was a very large agreement between detailed and aggregated simulations during the complete transient, except for the large current peak during AC filter bank reconnection in *IF*<sup>3</sup> at t = 0.8 s (Figure 17), where the aggregated model overestimated the *IF*<sup>3</sup> current peak.

The comparison of the proposed aggregation technique with existing techniques is shown in Figures 18 and 19. These figures show the simulation errors for the three aggregation techniques, during the DRU3 AC filter disconnection and re-connection transient considered in this section.

Figure 18 shows the error on WPP3 active and reactive powers (*PF*<sup>3</sup> and *QF*3) when simulated with each one of the considered aggregation techniques. All aggregation models adequately simulated the steady state active and reactive power delivered by WPP3. However, it was clear that, during the transients, the proposed aggregation technique (shown in blue) showed a smaller simulation error than existing techniques.

Figure 19 (top graph) shows that the proposed aggregation technique also showed a smaller *VF*<sup>3</sup> error than other techniques when the filters were disconnected (i.e., between t = 0.1 and t = 0.8 s). Clearly, the proposed technique was relatively better at simulating *VF*<sup>3</sup> voltage harmonics.

On the other hand, Figure 19 (bottom graph) shows that current *IF*<sup>3</sup> simulation errors between t = 0.1 and t = 0.8 s were very similar for the three considered techniques, around 0.5% (albeit that the proposed technique error was slightly better).

It is worth noting that all aggregation techniques showed relatively large simulation errors for *VF*<sup>3</sup> and *IF*<sup>3</sup> when the AC filter banks were re-connected at t = 0.8 s (Figure 19). However, the proposed technique performed clearly better than voltage drop and power loss aggregation techniques.

**Figure 18.** Case 3—Simulation errors for different aggregation techniques: top: *PF*<sup>3</sup>*Detailed* − *PF*<sup>3</sup>*Agg*; bottom: *QF*<sup>3</sup>*Detailed* − *QF*<sup>3</sup>*Agg*.

**Figure 19.** Case 3—Simulation errors for different aggregation techniques: top: *VF*<sup>3</sup>*Detailed* − *VF*<sup>3</sup>*Agg*; bottom: *IF*<sup>3</sup>*Detailed* − *IF*<sup>3</sup>*Agg*.

#### *6.2. Simulation Performance*

The simulations of all considered cases have been carried out with a 3.3-GHz Intel Core i7 PC with 16 GB DDR3 memory, SSD hard drive, and MS-Windows 7.

Considering a simulation run time of 3 s, Table 2 shows the complexity and simulation times of both detailed and aggregated models. The detailed simulation took 2 h, whereas the aggregated simulation was more than 350-times faster.

Therefore, the use of the proposed aggregated models for large Type-4 wind power plants was validated considering the similarity of the obtained results and the improvement on simulation times.


#### **7. Discussion and Conclusions**

This paper has presented an aggregation approach based on model order reduction using a multi-objective optimization technique. The proposed method consisted of first obtaining full admittance model of the WPP (in our case using state space techniques), then the inputs to the system were reduced to a single input applying the superposition principle, and finally, multi-objective optimization has been used to reduce the order of the system. The proposed aggregation technique considered voltage source WTGs' line side converters; therefore, it is also applicable to current-controlled voltage source converters.

Using the multi-objective optimization aggregation technique, the achieved aggregated admittance model had the same DC gain, main resonant frequency, and the same gain at the operating frequency as the detailed admittance model.

The proposed aggregation technique was validated considering the PSCAD/EMTDC simulation of a 1.2-GW HVDC DR-connected system, consisting of three 400-MW WPP of 50 Type-4 grid-forming WTGs each. Hence, the detailed system consisted of 150 individual wind turbines and 2450 nodes.

The proposed aggregation technique and two commonly-used aggregation strategies have been used to reduce the WPPs to only three aggregated WPPs (each one equivalent to 400 MW).

The comparison has been carried out considering three different test cases, namely fast active power reference changes, disconnection of a wind power plant, and disconnection and re-connection of DRU AC filter banks. These cases cover a wide range of dynamic and transient conditions.

The response against active power changes achieved by using the proposed aggregation method showed a worst case error three-times and an error variance 9.5-times better than aggregation techniques based on voltage drop or on power losses.

During WPP disconnection, all aggregation techniques showed relatively good PCC voltage and current simulation accuracy for the two wind power plants that remained connected. However, for the disconnected WPP, the proposed method showed voltage simulation accuracy 10-times better than standard methods, during the disconnection transient.

During disconnection and re-connection of one of the DRU AC filter banks, the proposed technique also showed better simulation accuracy regarding PCC voltage, active, and reactive power.

Finally, it has been shown that the aggregated simulations were 350-times faster than their detailed counterpart.

Therefore, when only the behavior at the points of common coupling is of interest, the proposed aggregated model provides important simulation time savings while delivering accurate results and preserving the main resonant characteristics of the full array system.

**Author Contributions:** The present work was developed with the following contributions: conceptualization, methodology, software, validation, formal analysis, research, writing, original draft preparation, and data curation: J.M.-T., S.A.-V., S.B.-P., and R.B.-G.; writing, review and editing, and supervision: S.A.-V., S.B.-P., and R.B.-G.

**Funding:** The authors would like to thank the support of the Spanish Ministry of Science, Innovation and Universities and EU FEDER funds under Grants DPI2014-53245-R and DPI2017-84503-R. This project has received funding from the European Union's Horizon 2020 research and innovation program under Grant Agreement No. 691714.

**Acknowledgments:** The authors would like to thank the editorial board, as well as the anonymous reviewers for their valuable comments that have greatly improved the quality of this document.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; nor in the decision to publish the results.

#### **Appendix A**


#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
