Wind farms providing secondary frequency regulation: evaluating the performance of model-based receding horizon control

This paper is an extended version of our paper presented at the 2016 TORQUE conference (Shapiro et al., 2016). We investigate the use of wind farms to provide secondary frequency regulation for a power grid using a model-based receding horizon control framework. In order to enable real-time implementation, the control actions are computed based on a time-varying one-dimensional wake model. This model describes wake advection and wake interactions, both of which play an important role in wind farm power production. In order to test the control strategy, it is implemented in a large-eddy simulation (LES) model of an 84-turbine wind farm using the actuator disk turbine representation. Rotor-averaged velocity measurements at each turbine are used to provide feedback for error correction. The importance of including the dynamics of wake advection in the underlying wake model is tested by comparing the performance of this dynamic-model control approach to a comparable static-model control approach that relies on a modified Jensen model. We compare the performance of both control approaches using two types of regulation signals, “RegA” and “RegD”, which are used by PJM, an independent system operator in the eastern United States. The poor performance of the static-model control relative to the dynamic-model control demonstrates that modeling the dynamics of wake advection is key to providing the proposed type of model-based coordinated control of large wind farms. We further explore the performance of the dynamic-model control via composite performance scores used by PJM to qualify plants for regulation services or markets. Our results demonstrate that the dynamic-model-controlled wind farm consistently performs well, passing the qualification threshold for all fast-acting RegD signals. For the RegA signal, which changes over slower timescales, the dynamic-model control leads to average performance that surpasses the qualification threshold, but further work is needed to enable this controlled wind farm to achieve qualifying performance for all regulation signals.

1 Introduction

Recent market trends are rapidly changing the composition of power grid energy sources, replacing conventional dispatchable power sources with non-dispatchable, variable resources, such as wind energy. These changes are putting pressure on the power system by reducing the number of resources available to provide a wide range of grid services traditionally provided by conventional power plants (Aho et al., 2012). A particularly important example is grid frequency regulation, which is closely tied to short-term imbalances in active power generation and load over timescales ranging from milliseconds to tens of minutes (Rebours et al., 2007). In order to deal with this challenge, a number of independent system operators (ISOs) are beginning to consider requiring wind plants to provide frequency regulation services and expanding frequency regulation markets to include wind plants (Aho et al., 2012; Díaz-González et al., 2014).

Secondary frequency regulation, in which participating generators track a power signal sent by an ISO over tens of minutes, is an area of growing interest. Recent work (Aho et al., 2013; Jeong et al., 2014) has shown that stand-alone wind turbines can effectively provide secondary frequency regulation, but recent fluid dynamics simulations (Fleming et al., 2016) have shown that interactions between wakes can lead to poor tracking performance when these single-turbine control strategies are applied to an array of turbines (Aho et al., 2013; Jeong et al., 2014). This poor performance is not unexpected because aerodynamic interactions between turbines occur on timescales commensurate with those of the secondary frequency regulation signals. Such considerations have led to recent emphasis on approaches that consider the collective behavior of the farm (Annoni et al., 2016; Gebraad et al., 2015), but few of these studies provide real-time implementable algorithms that can respond to changing wind farm power output levels.

Our recent work (Shapiro et al., 2017a) sought to overcome these challenges by developing a time-varying extension of the classic Jensen wake model (Katić et al., 1986) that accounts for the dynamics of wake advection through the farm. This new wake model was incorporated into a predictive model-based receding horizon control framework to coordinate an array of wind turbines to provide secondary frequency regulation by modulating the thrust coefficients of individual turbines. This approach used predictions from the underlying model to iteratively solve an online optimization problem representing the power tracking goal. Feedback from measurements of the velocity at each turbine was used to correct modeling errors. This approach showed promising results when tested in a large-eddy simulation (LES) model of a wind farm in which turbines were represented using actuator disk models (Shapiro et al., 2017a). In these simulations, we used set-point reductions of only 50 % of the maximum regulation provided, but were able to track a sample regulation signal with the wind farm test system used. In previous studies (Aho et al., 2013; Jeong et al., 2014), in which the control was designed at the individual turbine level, successful power tracking required set-point reductions exactly equal to the maximum change in power production requested by the ISO. The ability to lower the set-point reduction represents an important advantage over single-turbine approaches, as the amount of set-point reduction corresponds directly to the amount of power that wind farms are sacrificing in the bulk energy market to provide regulation. In fact, previous studies have shown that set-point reductions equal to the one-sided regulation signal variation may not be economically prudent (Rose and Apt, 2014).

The feasibility of providing secondary frequency regulation with wind farms was demonstrated by our initial results (Shapiro et al., 2017a). In this work we further evaluate the performance of this approach and consider the effect of reducing the control design and wake model complexity. In particular, we evaluate the importance of explicitly modeling the dynamics of wake advection by comparing the performance of the dynamic-model approach to a similar static-model approach, i.e., we replace the dynamic wake model with a wake model that does not include wake advection (Katić et al., 1986). In order to make appropriate comparisons, the static-model controller solves an online optimization problem with feedback similar to that solved in the dynamic-model controller.

Both controllers are implemented in ”virtual wind farms” comprised of actuator disk turbine models in LES. We evaluate the two approaches with regulation test signals from PJM, an ISO in the United States Eastern Interconnection (PJM, 2012, 2015). PJM has two types of secondary frequency regulation signals that are based on the area control error (ACE) signal, a combined measure of the power imbalance and deviation of the frequency from its nominal operating value. The “RegA” signal is a low-pass filter of the ACE that is generally followed using traditional regulating resources, such as fossil fuel plants. The “RegD” signal is a high-pass filter of the ACE that can be followed by more quickly responding resources, such as energy storage devices. Our results show that the static wake model leads to poor tracking performance, which indicates that the complexity of this particular control design cannot be reduced by ignoring the dynamics of wake advection. We then evaluate the performance of the dynamic-model controller using PJM's performance evaluation criteria (PJM, 2012, 2015) to determine whether the controlled wind farm system can meet PJM's threshold for qualification in the two regulation markets. These computations allow us to evaluate whether wind farms with this control strategy are better suited to provide traditional or fast-acting regulation.

The remainder of this paper is organized as follows. The static and dynamic wake models are described in Sect. 2, and the respective model-based controller designs are outlined in Sect. 3. Sections 4 and 5 describe the virtual wind farm test system and the simulation cases. The two controllers are compared in Sect. 6. The performance of dynamic-model control is further explored in Sect. 7 using PJM's performance criteria. Finally, we present conclusions and discuss directions for future work in Sect. 8.

2 Wake models

The two wake models employed here are based on static and dynamic adaptations of the classic Jensen wake model (Katić et al., 1986). In this presentation of the Jensen model, we consider regularly arranged wind farms with N rows and M columns of turbines, where each column is aligned with the prevailing wind direction, as shown in Fig. 1. The streamwise coordinate is denoted as x, and the nth turbine row is located at x = s_n. Every turbine is assumed to have the same rotor diameter D.

Figure 1 Diagram of a regular wind farm (left panel) with N = 5 rows and M = 3 columns of turbines and the corresponding row-averaged wake model representation (right panel). Each column is aligned with the streamwise coordinate x, and each row is aligned with the spanwise coordinate y.

The standard Jensen model assumes each turbine generates a wake region that expands radially at a linear rate k with increasing downstream distance from the turbine. This leads to the following definition of the wake diameter:

$$\begin{array}{}\text{(1)}& {\displaystyle }{\displaystyle }{D}_{\mathrm{w}}\left(x\right)=D+\mathrm{2}kx,\end{array}$$

where x is the streamwise distance from the turbine rotor plane. Conservation of mass leads to the following velocity deficit in the wake of the mth turbine in the nth row:

$$\begin{array}{}\text{(2)}& {\displaystyle }{\displaystyle }\mathit{\delta }{u}_{nm}\left(x\right)={\displaystyle \frac{\mathrm{2}{U}_{\mathrm{\infty }}{a}_{nm}}{{\left[D_{\mathrm{w}}\left(x-s_n\right)/D\right]}^{\mathrm{2}}}},\end{array}$$

where a_nm is the induction factor and U_∞ is the velocity upstream of the wind farm. This representation yields top-hat profiles of velocity deficits in each cross-stream plane. The velocity field experienced by the each turbine is found by superimposing the squared velocity deficits:

$$\begin{array}{}\text{(3)}& {\displaystyle }{\displaystyle }{u}_{\mathrm{\infty }nm}=U_{\mathrm{\infty }}-\sqrt{\sum _{(j,k)\in {\mathcal{S}}_{nm}}\mathit{\delta }{u}_{jk}^{\mathrm{2}}\left(s_n-s_j\right)},\end{array}$$

where 𝒮_nm is defined as the set of turbines whose wakes lie within the swept area of the turbine rotor of the mth turbine in row n. The definition of these sets means that Eq. (3) reduces to u_∞1m = U_∞ for the first row of turbines. The power production of each turbine is subsequently found using

$$\begin{array}{}\text{(4)}& {\displaystyle }{\displaystyle }{\widehat{P}}_{nm}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{\displaystyle \frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}}{C}_{\mathrm{P}nm}{u}_{\mathrm{\infty }nm}^{\mathrm{3}},\end{array}$$

where C_Pnm is the power coefficient of the turbine in row n and column m.

For ease of implementation, each wake model used in this paper makes the following modifications to the standard Jensen model. First, we consider each row of turbines collectively (Shapiro et al., 2017a), which means that each modeled value is homogeneous in the spanwise direction and we neglect the spanwise merging of wakes. To reflect this modification, the column index m used in Eqs. (2)–(3) is dropped in subsequent equations. Second, to account for entrance effects in the farm and compensate for the neglected spanwise wake interactions, we allow each wind turbine row to have a separate wake expansion rate k_n.

Furthermore, we express the turbine power production using the local thrust coefficient $C_{\mathrm{T}n}^{\prime }$ and modeled velocity at the turbine rotor ${\widehat{u}}_n$. Simple momentum theory can be used to show that (Goit and Meyers, 2015; Meyers and Meneveau, 2010; Shapiro et al., 2017a)

$$\begin{array}{}\text{(5)}& {\displaystyle }{\displaystyle }{a}_n={\displaystyle \frac{C_{\mathrm{T}n}^{\prime }}{\mathrm{4}+C_{\mathrm{T}n}^{\prime }}},{\widehat{u}}_n=\left(\mathrm{1}-a_n\right)u_{\mathrm{\infty }n},\mathrm{and}{C}_{\mathrm{T}n}^{\prime }={\displaystyle \frac{C_{\mathrm{T}n}}{{\left(\mathrm{1}-a_n\right)}^{\mathrm{2}}}}.\end{array}$$

Similarly, one can show that $C_{\mathrm{P}n}^{\prime }$ = $C_{\mathrm{P}}/(\mathrm{1}-a_n{)}^{\mathrm{3}}$, from which we conclude $C_{\mathrm{T}n}^{\prime }$ = $C_{\mathrm{P}n}^{\prime }$. Replacing the induction factor a_n in Eq. (2), the modeled upstream velocity u_∞n in Eq. (3), and the power coefficient C_Pn in Eq. (4) with these equations, the power production can be rewritten as

$$\begin{array}{}\text{(6)}& {\displaystyle }{\displaystyle }{\widehat{P}}_n=M{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{\displaystyle \frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}}{C}_{\mathrm{T}n}^{\prime }{\widehat{u}}_n^{\mathrm{3}}.\end{array}$$

These idealized conditions assume that the electrical power generated by the turbine is proportional to the power extracted from the flow and the control actions do not significantly affect the aerodynamic efficiency of the blades (Goit and Meyers, 2015). Aerodynamic losses could also be taken into account by reducing the local power coefficient by a constant factor α (Stevens and Meneveau, 2014), i.e., $C_{\mathrm{P}}^{\prime }$ ≈ $\mathit{\alpha }{C}_{\mathrm{T}}^{\prime }$. For example, a wind turbine operating at a thrust coefficient of C_T = 0.75 and C_P = 0.45 would use α = 0.8. The following subsections describe the static and dynamic wake models in more detail.

2.1 Static wake model

The static wake model used in this work is the Jensen model with the modifications described above. In order to use this static wake model in a model-based controller for the farm power production, the model must be augmented to account for time-varying changes in the local thrust coefficient $C_{\mathrm{T}n}^{\prime }\left(t\right)$. Including time dependency in the thrust coefficient and replacing the induction factor in Eq. (2) with the expression in Eq. (5) gives the following expression for the velocity deficit for the nth turbine

$$\begin{array}{}\text{(7)}& {\displaystyle }{\displaystyle }\mathit{\delta }{u}_n(x,t)={\displaystyle \frac{C_{\mathrm{T}n}^{\prime }\left(t\right)}{\mathrm{4}+C_{\mathrm{T}n}^{\prime }\left(t\right)}}{\displaystyle \frac{\mathrm{2}{U}_{\mathrm{\infty }}}{{\left[\mathrm{1}+\mathrm{2}{k}_n\left(x-s_n\right)/D\right]}^{\mathrm{2}}}}.\end{array}$$

With this approach, thrust coefficient changes instantaneously affect the velocity deficit everywhere, i.e., the wakes implicitly have an infinitely fast advection speed. Finally, the velocity at the turbines of the nth row can be found by explicitly writing out the set of upstream turbines in Eq. (3) affecting the velocity at the nth turbine and using the equation for the rotor-averaged velocity in Eq. (5):

$$\begin{array}{}\text{(8)}& {\displaystyle }{\displaystyle }{\widehat{u}}_n\left(t\right)=\left(\mathrm{1}-{\displaystyle \frac{C_{\mathrm{T}n}^{\prime }\left(t\right)}{\mathrm{4}+C_{\mathrm{T}n}^{\prime }\left(t\right)}}\right)\left(U_{\mathrm{\infty }}-\sqrt{\sum _{m=\mathrm{1}}^{n-\mathrm{1}}\mathit{\delta }{u}_m^{\mathrm{2}}\left(s_n-s_m\right)}\right).\end{array}$$

Equations (7) and (8) are therefore the static wake model equations ${\mathit{W}}_{\mathrm{s}}({\mathit{C}}_{\mathrm{T}}^{\prime }$, q_s) = 0, where q_s = [δu, $\widehat{\mathit{u}}$] denotes the model states and boldface indicates vectors.

2.2 Dynamic wake model

The dynamic wake model also borrows from the classic Jensen model, but instead allows the wake velocity deficit to move downstream at a finite velocity. The resulting one-dimensional time-varying wake model, which assumes that the wake travels with the inlet velocity U_∞, was previously proposed and validated against LES of wind farms at startup (Shapiro et al., 2017a). The velocity deficit is governed by

$$\begin{array}{}\text{(9)}& {\displaystyle }{\displaystyle \frac{\partial \mathit{\delta }{u}_n}{\partial t}}+U_{\mathrm{\infty }}{\displaystyle \frac{\partial \mathit{\delta }{u}_n}{\partial x}}=-w_n\left(x\right)\mathit{\delta }{u}_n(x,t)+f_n(x,t),\end{array}$$

where w_n(x) is the wake decay function and f_n(x, t) is a forcing function used to account for the effect of the turbine on the flow field. The wake decay function

$$\begin{array}{}\text{(10)}& {\displaystyle }{\displaystyle }{w}_n\left(x\right)=\mathrm{2}{\displaystyle \frac{U_{\mathrm{\infty }}}{d_n\left(x\right)}}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}x}}{d}_n\left(x\right)\end{array}$$

is determined by assuming that the wake diameter normalized by the rotor diameter d_n(x) = D_wn(x)∕D at a fixed location x is constant in time. Momentum theory shows that as the air flows through the turbine rotor, the velocity reduces to U_∞ − 2$U_{\mathrm{\infty }}{C}_{\mathrm{T}n}^{\prime }/$(4 + $C_{\mathrm{T}n}^{\prime })$ (Shapiro et al., 2017a). In order to retrieve this expected velocity reduction, the forcing function is specified as

$$\begin{array}{}\text{(11)}& {\displaystyle }{\displaystyle }{f}_n(x,t)={\displaystyle \frac{\mathrm{2}{U}_{\mathrm{\infty }}^{\mathrm{2}}}{d_n^{\mathrm{2}}\left(x\right)}}{\displaystyle \frac{C_{\mathrm{T}n}^{\prime }\left(t\right)}{\mathrm{4}+C_{\mathrm{T}n}^{\prime }\left(t\right)}}G\left(x-s_n\right),\end{array}$$

where G(x − s_n) is a smoothing function that integrates to unity, centered at the streamwise location of the turbine x = s_n. A Gaussian function with characteristic width Δ,

$$\begin{array}{}\text{(12)}& {\displaystyle }{\displaystyle }G\left(x-s_n\right)={\displaystyle \frac{\mathrm{1}}{\mathrm{\Delta }\sqrt{\mathrm{2}\mathit{\pi }}}}{e}^{-\frac{{\left(x-s_n\right)}^{\mathrm{2}}}{\mathrm{2}{\mathrm{\Delta }}^{\mathrm{2}}}},\end{array}$$

maintains smoothness in the velocity deficit fields.

In the Jensen model (Katić et al., 1986), the dimensionless diameter of the wake generated by turbine row n is d_n(x) = 1 + 2k_n(x − s_n)∕D, where k_n is an empirical wake expansion coefficient. We make two modifications to this equation. First, the linear expansion is assumed to begin at x = s_n + 2Δ to prevent the wake expansion from occurring within the induction zone imposed by the Gaussian forcing. The second modification addresses the fact that the equation for the standard Jensen dimensionless wake diameter is ill-posed upstream of the turbine, where it can vanish or become negative. Therefore, we use the following modified function that smoothly approximates the linear expansion in the far wake while preventing the wake diameter from becoming less than unity close to the turbine:

$$\begin{array}{}\text{(13)}& {\displaystyle }{\displaystyle }{d}_n\left(x\right)=\mathrm{1}+k_n\ln \left[\mathrm{1}+\exp \left({\displaystyle \frac{x-s_n-\mathrm{2}\mathrm{\Delta }}{D/\mathrm{2}}}\right)\right].\end{array}$$

As in the static model, squared deficits (Katić et al., 1986) are superposed to calculate the estimated streamwise velocity ${\widehat{u}}_n$ at the turbine:

$$\begin{array}{}\text{(14)}& {\displaystyle }{\displaystyle }{\widehat{u}}_n\left(t\right)=U_{\mathrm{\infty }}-\underset{\mathrm{0}}{\overset{L}{\int }}{\left(\sum _{m=\mathrm{1}}^N\mathit{\delta }{u}_m^{\mathrm{2}}(x,t)\right)}^{\mathrm{1}/\mathrm{2}}G\left(x-s_n\right)\mathrm{d}x.\end{array}$$

Finally, the total estimated power ${\widehat{P}}_n$ of the M turbines in row n is found using Eq. (6). The dynamic wake model equations, Eqs. (9)–(14), are written as ${\mathit{W}}_{\mathrm{d}}({\mathit{C}}_{\mathrm{T}}^{\prime }$, q_d) = 0, where q_d = [δu, $\widehat{\mathit{u}}]$ denotes the model states.

3 Controlled wind farm designs

The model-based controllers implementing the static and dynamic wake models discussed above are designed to track the power reference signals P_ref(t) sent by an ISO by modulating the thrust coefficients of each turbine row $C_{\mathrm{T}n}^{\prime }\left(t\right)$. Thrust modulation control is used as a proxy for direct actuation of blade pitch angle and generator torque. Explicit actuation of these control variables is the subject of future work. In both cases, feedback is included by measuring the row-averaged, rotor-averaged wind speed u_n(t). The resulting feedback term ϵ_n(t) is fed into the controller and used to correct the predicted power output of the wake model. A block diagram of the controlled wind farm system is shown in Fig. 2.

Figure 2 Block diagram of the controlled wind farm system for both controllers. Each controller computes a thrust coefficient command signal ${\mathit{C}}_{\mathrm{T}n}^{\prime }\left(t\right)$ using the reference signal P_ref(t) and an error correction term ϵ_n(t). The error correction is computed using the measured velocity u_n(t) and the predicted velocity ${\widehat{u}}_n\left(t\right)$ from the underlying wake model.

3.1 Controller designs

Each controller calculates the local thrust coefficient trajectories by repeatedly solving a minimization problem of the form

$$\begin{array}{}\text{(15)}& {\displaystyle }& {\displaystyle }\underset{{\mathit{C}}_{\mathrm{T}}^{\prime },\mathit{q}}{minimize}\mathcal{J}\left({\mathit{C}}_{\mathrm{T}}^{\prime },\mathit{q}\right)+\mathcal{R}\left({\mathit{C}}_{\mathrm{T}}^{\prime }\right)\text{(16)}& {\displaystyle }& {\displaystyle }\mathrm{subject}\mathrm{to}\mathit{W}\left({\mathit{C}}_{\mathrm{T}}^{\prime },\mathit{q}\right)=\mathrm{0},\end{array}$$

where $\mathit{W}({\mathit{C}}_{\mathrm{T}}^{\prime }$, q) = 0 and q are placeholders for the static and dynamic wake model and states, which were previously indicated by the subscripts “s” and “d”, respectively. The cost functional 𝒥 represents the reference tracking goal, and the functional ℛ contains regularizations to maintain well-behaved control trajectories. These regularizations include a penalization for fast changes in the thrust coefficients to avoid excessive oscillations in the control and a penalization for deviations away from the pre-control reference power to prevent the thrust coefficients from moving outside of physical bounds.

Although both controllers solve an online optimization problem, the mechanics of the implementation are quite different. Since the equations for the static model have no dynamics, every instance in time is uncoupled. Therefore, the static-model controller can consider each point in time as a separate minimization problem, and the cost functionals can be written solely in terms of the current state. With this approach, the power tracking cost functional at time t is written as

$$\begin{array}{}\text{(17)}& {\displaystyle }{\displaystyle }{\mathcal{J}}_{\mathrm{s}}\left({\mathit{C}}_{\mathrm{T}}^{\prime },{\mathit{q}}_{\mathrm{s}}\right)={\displaystyle \frac{\mathrm{1}}{{\mathcal{P}}^{\mathrm{2}}}}{\left(\sum _{n=\mathrm{1}}^N{\widehat{P}}_n\left(t\right)-P_{\mathrm{ref}}\left(t\right)\right)}^{\mathrm{2}},\end{array}$$

and the regularization terms are

$$\begin{array}{}\text{(18)}& {\displaystyle }{\displaystyle }{\mathcal{R}}_{\mathrm{s}}\left({\mathit{C}}_{\mathrm{T}}^{\prime }\right)=\mathit{\eta }\sum _{n=\mathrm{1}}^N{\left(C_{\mathrm{T}n}^{\prime }\left(t\right)-C_{\mathrm{T},\mathrm{ref}}\right)}^{\mathrm{2}}+\mathit{\gamma }{T}^{\mathrm{2}}\sum _{n=\mathrm{1}}^N{\left({\displaystyle \frac{\mathrm{d}{C}_{\mathrm{T}n}^{\prime }}{\mathrm{d}t}}\right)}^{\mathrm{2}}.\end{array}$$

The constants 𝒫 and T in Eqs. (17)–(18) are used to make each term in the power tracking cost functional and regularization functional dimensionless and of comparable magnitude. Here we choose 𝒫 = $M\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{U}_{\mathrm{\infty }}^{\mathrm{3}}$ and time T as the time horizon of the reference signal considered. The constants η and γ are the respective weights of each regularization term, which are set to η = 0.005 and γ = 2.083 × 10⁻⁵ in this study.

The dynamic-model controller (Shapiro et al., 2017a), however, accounts for the time-dependent advection of turbine wakes. We therefore employ a model-based receding horizon framework, which is a predictive approach that uses the model to plan future control actions. The receding horizon method works by iteratively solving a finite-time minimization problem over a time horizon T. The solution is implemented for a shorter period T_A before re-solving the minimization problem. More details about this procedure can be found in Bewley et al. (2001) and Goit and Meyers (2015). With this predictive framework, the reference tracking goal is represented by the cost functional

$$\begin{array}{}\text{(19)}& {\displaystyle }{\displaystyle }{\mathcal{J}}_{\mathrm{d}}\left({\mathit{C}}_{\mathrm{T}}^{\prime },{\mathit{q}}_{\mathrm{d}}\right)={\displaystyle \frac{\mathrm{1}}{{\mathcal{P}}^{\mathrm{2}}T}}\underset{\mathrm{0}}{\overset{T}{\int }}{\left(\sum _{n=\mathrm{1}}^N{\widehat{P}}_n\left(t\right)-P_{\mathrm{ref}}\left(t\right)\right)}^{\mathrm{2}}\mathrm{d}t,\end{array}$$

and the regularization functional is defined as

$$\begin{array}{ll}{\displaystyle }{\mathcal{R}}_{\mathrm{d}}\left({\mathit{C}}_{\mathrm{T}}^{\prime }\right)& {\displaystyle }={\displaystyle \frac{\mathit{\eta }}{T}}\sum _{n=\mathrm{1}}^N\underset{\mathrm{0}}{\overset{T}{\int }}{\left(C_{\mathrm{T}n}^{\prime }\left(t\right)-C_{\mathrm{T},\mathrm{ref}}^{\prime }\right)}^{\mathrm{2}}\mathrm{d}t\\ \text{(20)}& {\displaystyle }& {\displaystyle }+\mathit{\gamma }T\sum _{n=\mathrm{1}}^N\underset{\mathrm{0}}{\overset{T}{\int }}{\left({\displaystyle \frac{\mathrm{d}{C}_{\mathrm{T}n}^{\prime }}{\mathrm{d}t}}\right)}^{\mathrm{2}}\mathrm{d}t.\end{array}$$

Consistent with the distinction between the non-predictive and predictive nature of the static-model and dynamic-model controllers, the functionals for the static wake model are not integrated forward in time. In other words, the static wake model is not a receding horizon method because the modeled system does not include dynamics.

All minimizations are solved using the modified unconstrained reduced cost functional $\stackrel{\mathrm{̃}}{\mathcal{J}}\left({\mathit{C}}_{\mathrm{T}}^{\prime }\right)$ = 𝒥(q, ${\mathit{C}}_{\mathrm{T}}^{\prime })$ (Bewley et al., 2001; Goit and Meyers, 2015), instead of the cost functional defined in Eq. (15). Minimizations are performed using the gradient-based nonlinear Polak–Ribière conjugate gradient method (Press, 2007) combined with the Moré–Thuente line search method (Moré and Thuente, 1994) and are terminated after 100 iterations. For the static wake model, gradients are obtained using finite differencing, which can be implemented efficiently because there are only N control variables per minimization. For the dynamic wake model, gradients are obtained using one backward simulation of the adjoint equations of the wake model using the formal Lagrangian method (Borzì and Schulz, 2011; Goit and Meyers, 2015). The full procedure is detailed in Shapiro et al. (2017a). This approach was chosen because it is computationally efficient for systems with large state spaces, such as the discretized partial differential equation system described by the dynamic wake model.

In this work, we use horizon and advancement times of T = 40 min and T_A ≈ 10 s, respectively. With these values, the optimization takes approximately 1 min on a single processor, which is roughly 6 times as long as the advancement time of T_A = 10 s. However, several modifications can reduce the optimization time significantly. For example, a previous implementation reduced the optimization time to a fraction of the advancement time by employing a quasi-Newton minimization method, reducing the horizon and advancement times, and limiting the number of minimization iterations (Shapiro et al., 2017b). As a result, this approach is feasible for real-time control.

3.2 Measurement feedback

As shown in Fig. 2, the controller employs closed-loop feedback for velocity measurements at each turbine to correct modeling errors and assumptions. The row-averaged power and row- and rotor-averaged wind velocities are defined as

$$\begin{array}{}\text{(21)}& {\displaystyle }{\displaystyle }{P}_n=\sum _{m=\mathrm{1}}^M{P}_{nm},\mathrm{and}{u}_n={\displaystyle \frac{\mathrm{1}}{M}}{\left(\sum _{m=\mathrm{1}}^M{u}_{nm}^{\mathrm{3}}\right)}^{\mathrm{1}/\mathrm{3}},\end{array}$$

where u_nm is the velocity measured at the turbine in the nth row and mth column of the wind farm. The definition of the row-average velocity at the turbine disk is necessary to ensure that P_n = $M\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{C}_{\mathrm{T}n}^{\prime }{u}_n^{\mathrm{3}}$. These measurements are used to calculate an error term ϵ_n and provide feedback by replacing Eq. (6) with

$$\begin{array}{}\text{(22)}& {\displaystyle }{\displaystyle }{\widehat{P}}_n=M{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{\displaystyle \frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}}{C}_{\mathrm{T}n}^{\prime }{\left({\widehat{u}}_n+{\mathit{ϵ}}_n\right)}^{\mathrm{3}}.\end{array}$$

For the static wake model, the error term for turbine row n at time step p is calculated using the difference between the measured and estimated velocity from the previous time step p − 1

$$\begin{array}{}\text{(23)}& {\displaystyle }{\displaystyle }{\mathit{ϵ}}_n^p=u_n^{p-\mathrm{1}}-{\widehat{u}}_n^{p-\mathrm{1}}.\end{array}$$

For the dynamic wake model, the error correction at the receding horizon iteration starting at time t_c is

$$\begin{array}{}\text{(24)}& {\displaystyle }{\displaystyle }{\mathit{ϵ}}_n\left(t\right)=\left(u_n\left(t_{\mathrm{c}}\right)-{\widehat{u}}_n\left(t_{\mathrm{c}}\right)\right)e^{-\left(t-t_{\mathrm{c}}\right)/\mathit{\tau }}.\end{array}$$

The exponential decay accounts for the reduced future accuracy of the error term in the receding horizon prediction and is set to τ = 120 s in this study.

4 Virtual wind farm test system

A LES model of a wind farm with wind turbines represented using actuator disk models is used to test the two control approaches. The wind farm is composed of N = 7 rows of M = 12 aligned columns of turbines. Each turbine has a 100 m rotor diameter D and a 100 m hub height. The spacing between turbines is 7D in the streamwise direction and 5D in the spanwise direction. Prior to the initiation of the control actions, all of the turbines are operated at a constant reference local thrust coefficient of $C_{\mathrm{T},\mathrm{ref}}^{\prime }$ = 1.33, which is assumed to be representative of wind turbines operating in region 2 (Calaf et al., 2010; Stevens et al., 2014a).

These simulations are performed using Johns Hopkins University's LES code LESGO (Calaf et al., 2010; Stevens et al., 2014b; VerHulst and Meneveau, 2015), which uses pseudo-spectral discretization in the horizontal directions with periodic boundary conditions. The code also employs second-order Adams–Bashforth time integration, second-order finite differencing in the vertical direction, and the dynamic scale-dependent Lagrangian Smagorinsky subgrid stress model (Bou-Zeid et al., 2005). Inlet conditions for the wind farm are generated using the concurrent-precursor method (Stevens et al., 2014b). The force exerted by and the power extraction of the mth turbine of the nth row are both a function of the filtered rotor-averaged velocity u_nm (Calaf et al., 2010) and the thrust coefficient $C_{\mathrm{T}nm}^{\prime }$. The force is modeled as a drag force F_nm = $-\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{C}_{\mathrm{T}nm}^{\prime }{u}_{nm}^{\mathrm{2}}$, and the power extraction is P_nm = −F_nmu_nm. An instantaneous color contour plot of the streamwise velocity field from one of these simulations is shown in Fig. 3.

Figure 3 Instantaneous streamwise velocity contours for a large-eddy simulation with actuator disk turbine models, which are indicated by black dashes. Each turbine has a rotor diameter of D = 100 m and hub height of 100 m. The mean and maximum inlet velocities are approximately 9.5 and 12 m s⁻¹, respectively. The inlet conditions are generated using a concurrent precursor simulation (Stevens et al., 2014b) shown at the beginning of the plotted domain.

Figure 4 Probability density functions of 𝒮₁–𝒮₃ defined in Eq. (25) for RegA (black) and RegD (red) during 2015. The three selected historical hours are shown in the PDFs as vertical dashed lines.

5 Test cases

The performance of the controlled wind farms is evaluated using a series of PJM's published RegA and RegD test signals as well as historical RegA and RegD signals from three hours in 2015 (PJM, 2012, 2015, 2016). For each regulation signal, we use three initial conditions for the wind farm and two levels of power set-point reduction, which we also refer to as derates. For each controller test, the reference signal is defined as P_ref(t) = [1 − α + 0.08r(t)]P_base, where P_base is the 5 min average power prior to initiation of the control, α is the derate amount, and r(t) ∈ [−1, 1] is the regulation signal from the ISO. As a result, the reference signal varies by ±8 % of the baseline power P_base.

The combination of test signals and initial conditions leads to 48 unique test cases, each of which is given a unique identifier that is a combination of identifiers for each of the variable types shown in Table 1. “Signal” refers to the regulation signal type (RegA or RegD), “Derate” refers to the derate amount (4 or 6 %), “Initial condition” refers to the initial condition of the controlled plant simulation, and “Period” refers to the regulation signal period, which is either the PJM test signals or one of the selected hours in 2015. For example, the test case “RegA.D4.IC1.TS” refers to the case with the RegA test signal, 4 % derate, and the first initial condition.

5.1 Historical PJM regulation signals

The number of historical hours used to test the controlled wind farm is constrained by the computational cost of running the model wind farm LES. As a result, it is impractical to select enough hours to sample the entire range of possible regulation signals provided by PJM. To prevent systematic bias, the three hours were selected without considering the characteristics of the regulation signals during those periods.

In order to evaluate whether the selected signals are representative cases, we compare them to the range of all possible regulation signals provided by PJM in 2015 using three statistics:

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }{\mathcal{S}}_{\mathrm{1}}={\displaystyle \frac{\mathrm{1}}{T}}\underset{\mathrm{0}}{\overset{T}{\int }}r\left(t\right)\mathrm{d}t\\ {\displaystyle }& {\displaystyle }{\mathcal{S}}_{\mathrm{2}}={\displaystyle \frac{\mathrm{1}}{T}}\underset{\mathrm{0}}{\overset{T}{\int }}{r}^{\mathrm{2}}\left(t\right)\mathrm{d}t-{\mathcal{S}}_{\mathrm{1}}^{\mathrm{2}}\\ \text{(25)}& {\displaystyle }& {\displaystyle }{\mathcal{S}}_{\mathrm{3}}={\displaystyle \frac{\mathrm{1}}{T}}\underset{\mathrm{0}}{\overset{T}{\int }}{\left({\displaystyle \frac{\mathrm{d}r}{\mathrm{d}t}}\right)}^{\mathrm{2}}\mathrm{d}t,\end{array}$$

where r(t) is the regulation signal, T = 60 min, 𝒮₁ is the mean of r(t), 𝒮₂ is the variance of r(t), and 𝒮₃ is the variance of $\frac{\mathrm{d}r}{\mathrm{d}t}$. For each of these statistical measures, the probability density function (PDF) is calculated using all possible hourly signals provided by PJM in 2015 and is shown in Fig. 4. These PDFs demonstrate the differences between the RegA and RegD signals. The RegA signals have a larger mean and variance than the RegD signals, but the variance of $\frac{\mathrm{d}r}{\mathrm{d}t}$ is smaller. The values of these statistics for the three selected hours is compared to the PDFs over the entire year in Fig. 4. These figures show that the selected historical signals represent a reasonable cross section of the possible PJM regulation signals. The only exception, the high percentile ranking in 𝒮₁ of the RegA signals, represents a more difficult test for the controlled wind farm because more energy is requested than the average.

Table 1 Test case identifiers describing the signal type, derate amount, initial condition of the wind farm, and regulation signal period. For example, the test case “RegD.D6.IC1.H2” refers to the case with the second RegD historical signal, 6 % derate, and the first initial condition.

Table 2 Characteristics of wind farm initial conditions, including mean inlet velocity $\stackrel{\mathrm{‾}}{u}$, standard deviation of inlet velocity σ_u, and turbulence intensity TI. The corresponding wake model inlet velocity U_∞ and wake expansion coefficients k_n are also shown.

5.2 Wind farm initial conditions

We set the initial conditions of the controlled wind farm simulations to correspond to uncontrolled simulations with a local thrust coefficient of $C_{\mathrm{T},\mathrm{ref}}^{\prime }$ = 1.33, as previously discussed. The inflow characteristics for the three initial conditions of interest are provided in Table 2. The inflow velocities of the initial conditions have a mean $\stackrel{\mathrm{‾}}{u}$ ≈ 9.5 m s⁻¹ and a standard deviation σ_u ≈ 1 m s⁻¹ as measured at the first row of turbines during the T₀ = 5 min prior to initiation of the control:

$$\begin{array}{ll}{\displaystyle }& {\displaystyle }\stackrel{\mathrm{‾}}{u}={\displaystyle \frac{\mathrm{1}}{T_{\mathrm{0}}M}}\underset{-T_{\mathrm{0}}}{\overset{\mathrm{0}}{\int }}\sum _{m=\mathrm{1}}^M{\displaystyle \frac{\mathrm{4}+C_{\mathrm{T},\mathrm{ref}}^{\prime }}{\mathrm{4}}}{u}_{\mathrm{1}m}\left(t\right)\mathrm{d}t\\ \text{(26)}& {\displaystyle }& {\displaystyle }{\mathit{\sigma }}_u={\left[{\displaystyle \frac{\mathrm{1}}{T_{\mathrm{0}}M}}\underset{-T_{\mathrm{0}}}{\overset{\mathrm{0}}{\int }}\sum _{m=\mathrm{1}}^M{\left({\displaystyle \frac{\mathrm{4}+C_{\mathrm{T},\mathrm{ref}}^{\prime }}{\mathrm{4}}}{u}_{\mathrm{1}m}\left(t\right)-\stackrel{\mathrm{‾}}{u}\right)}^{\mathrm{2}}\mathrm{d}t\right]}^{\mathrm{1}/\mathrm{2}}.\end{array}$$

The turbulence intensity as measured at the center of each of the turbine rotors is approximately 13 %, which corresponds to low to medium International Electrotechnical Commission turbulence levels (IEC, 2005). The simulations assume region 2 operation (Johnson et al., 2006) with idealized aerodynamic characteristics of $C_{\mathrm{P}}^{\prime }$ = $C_{\mathrm{T}}^{\prime }$. In order to avoid any significant interaction with the rated regime, we assume wind turbines with a rated wind speed of at least 12 m s⁻¹, which corresponds to the 99th percentile of the LES velocity field. Wind turbines with a diameter of D = 100 m and a power coefficient of C_P = 0.5625, which corresponds to $C_{\mathrm{P}}^{\prime }$ = $C_{\mathrm{T}}^{\prime }$, therefore have a rated power of approximately 4.5 MW and an average total farm power of approximately 100 MW. Under nonideal aerodynamic conditions (Stevens et al., 2014a), a power coefficient of C_P = 0.45 would yield a rated wind turbine power of 3.6 MW.

The required parameters of the static and dynamic wake models, inlet velocity U_∞, and wake expansion coefficients k_n are also calculated for each initial condition using measurements from the T₀ = 5 min prior to initialization of the control. The inlet velocity is set using the relationship U_∞ = $\frac{\mathrm{1}}{\mathrm{4}}$(4 + $C_{\mathrm{T},\mathrm{ref}}^{\prime }\left)T_{\mathrm{0}}^{-\mathrm{1}}{\int }_{-T_{\mathrm{0}}}^{\mathrm{0}}{u}_{\mathrm{1}}\right(t)\mathrm{d}t$, and the wake expansion coefficients are found using a least-squares fit between the measured power and the power predicted by the static model. Note that the inlet velocity for the model is defined using the average power and therefore the average inflow velocity is not equal to the inlet velocity for the models, i.e., $\stackrel{\mathrm{‾}}{u}$ ≠ U_∞. The resulting parameters are also shown in Table 2.

6 Comparison of control methods

The power tracking performance and control trajectories of the controlled wind farm, represented by the LES described in Sect. 4, are shown in Figs. 5 and 6. The left and right panels of these figures show the performance of the static- and dynamic-model controllers, respectively. Figure 5 shows the response of the controlled farms to the RegA test signals, and Fig. 6 shows the response of the controllers to the RegD test signals. The dynamic-model control demonstrates good overall tracking performance, although it has some trouble tracking the reference signal during the last 5–10 min of the RegA.D4.IC1.TS and RegA.D4.IC3.TS cases. Conversely, the static-model control demonstrates poor overall tracking performance, although it is able to track the signal for certain down regulation events, e.g., around minute 20 in all cases in Fig. 5.

The static-model control method appears to switch between two distinct operating points, depending on the characteristics of the regulation signal. Downregulation trajectories are often successfully tracked by increasing the thrust coefficient of the first row of turbines to values above $C_{\mathrm{T}}^{\prime }$ = 2. This change in operating conditions increases the magnitude of the velocity deficits throughout the farm, thereby reducing the overall wind speed and total power production. When there is a period of upregulation approaching or the wind farm is slightly underproducing, the controller quickly reduces the upstream thrust coefficients and moves to the Betz optimal thrust coefficient $C_{\mathrm{T}}^{\prime }$ = 2 (Goit and Meyers, 2015) for the last row. This operating point is likely the optimal power point for the Jensen model with constant wake expansion coefficients.

Figure 5 Comparison of static-model (a, c, e) and dynamic-model (b, d, f) control methods for RegA test signals with 4 % derates. All three initial conditions are shown from top to bottom. The top plot in each panel shows the controlled LES wind farm model power production (black) compared to the reference signal (red). The bottom plot in each panel shows the local thrust coefficients calculated with control methods by row: row 1 (dark blue), row 2 (orange), row 3 (yellow), row 4 (purple), row 5 (green), row 6 (light blue), and row 7 (brown).

Figure 6 Comparison of static-model (a, c, e) and dynamic-model (b, d, f) control methods for RegD test signals with 4 % derates. All three initial conditions are shown from top to bottom. The top plot in each panel shows the controlled LES wind farm model power production (black) compared to the reference signal (red). The bottom plot in each panel shows the local thrust coefficients calculated with control methods by row: row 1 (dark blue), row 2 (orange), row 3 (yellow), row 4 (purple), row 5 (green), row 6 (light blue), and row 7 (brown).

Figure 7 Comparison of static-model (a) and dynamic-model (b) control methods for the RegA.D4.IC2.TS simulation case. Each panel shows the controlled LES wind farm model power production, rotor-averaged velocity, and thrust coefficients by row: row 1 (dark blue), row 2 (orange), row 3 (yellow), row 4 (purple), row 5 (green), row 6 (light blue), and row 7 (brown).

Figure 8 Power tracking performance of dynamic-model-controlled wind farms comparing simulated farm power from a controlled LES wind farm model (black), an uncontrolled LES wind farm model (gray), and power reference signals (red) for 4 % derates and initial condition 3.

The performance of the static-model control provides an interesting demonstration of the importance of including time dependency in the wake model used in this type of control scheme. In an attempt to track the changing reference signal, the controller switches quickly between the two operating points discussed above. The static Jensen model erroneously models these transitions between operating points as an instantaneous change of the wake velocity deficit throughout the farm. In reality, the air around the turbine will slowly respond to a sudden change in the thrust coefficient and the reduced wake deficit must travel through the farm before the effects of changing upstream thrust coefficients on downstream power production and wind speeds are realized. Detailed trajectories of the power and rotor-averaged velocity of each row in Fig. 7 show that the LES wind farm does not respond instantaneously to the change in operating point. Instead, power production slowly increases between minutes 5 and 15.

As a result of these modeling errors, the static-model controller produces large transient variations in power production when moving between operating points. When moving to the upregulation operating point identified by the controller, the power production of the farm plunges. In some cases, the total power production slowly recovers to the desired set point. Furthermore, all of the static-model control cases in Figs. 5 and 6 demonstrate significant overshoot in the power production during the first 30 s of the simulations as the thrust coefficients quickly move away from the pre-control level. The frequency of the changes in the control actions, which is determined by the advancement time of approximately 10 s, is faster than these dynamics of wake advection and is therefore unlikely to explain these effects. Instead, neglecting the wake advection time in the Jensen model best explains the poor performance of the static model controller.

The dynamic-model control uses strategies similar to those of the static-model controller, including increasing the thrust coefficient during downregulation periods and moving toward a Jensen model optimal power point for upregulation periods. However, by including the time-dependent effects of wake advection, the controller avoids large transient changes when changing between states. The underlying dynamic model can correctly predict the time-varying effect of changing upstream thrust coefficients on downstream power production. In the next section we further study the performance of this dynamic-model control approach.

7 Performance evaluation of dynamic-model control

The time evolution of the total LES wind farm power is compared to the reference signals for initial condition 3 and a 4 % derate in Fig. 8, which shows all regulation signals (RegA or RegD) and regulation period combinations. The controlled wind farm power production is also compared to the uncontrolled case, in which the wind farm is kept at the constant pre-control thrust coefficient. These results demonstrate the good overall tracking performance of the controlled wind farm, except for a few specific periods of underperformance. Furthermore, the results demonstrate that the dynamic-model-based receding horizon control method is also able to reduce the natural turbulent fluctuations in the wind farm power production. Indeed, the root mean square (RMS) of the controlled wind farm power production about the reference signal is 1.06 MW, which is almost a quarter of the 3.93 MW RMS of uncontrolled power production about the baseline power.

Quantitative measures of the performance of each regulation signal type (RegA or RegD) for derate values of 4 and 6 % are shown in terms of PJM's performance scores in Fig. 9. In order to participate in PJM's regulation market, power plants must pass the regulation qualification test for the particular regulation signal being supplied. This test is carried out over a 40 min period, and the tracking capability is quantified using a composite performance score, which is the weighted sum of accuracy, delay, and precision scores (PJM, 2012, 2015). The accuracy score measures the ability of the signal to respond to a change in the ISO regulation signal. The delay score measures the delay in the plant's response to the regulation signal. The precision score measures the difference between the requested power and the plant's power output. A minimum composite score of 75 % is needed to qualify to participate in each of the two regulation services. Once qualified for a particular service, a plant is continuously evaluated; if its average score over the last 100 h drops below 40 %, then the plant is disqualified from providing the service and must retake the initial performance test to re-qualify.

The controlled wind farm performs better for the RegD signals, meeting the composite score threshold for qualification of 75 % in all cases. The performance of the controlled farm in tracking the RegA signals is also satisfactory for PJM participation, but the controlled farm would not have qualified in all tests. These lower composite scores may be explained by the large values in 𝒮₁, which represent the total energy requested in the signals, compared to other PJM signals. However, in cases in which the controlled wind farm had poor performance for the RegA signal with a 4 % derate, increasing the derate to 6 % markedly improved the overall performance.

Figure 9 Box plots of PJM performance scores for a dynamic-model-controlled wind farm for all regulation signal types (RegA or RegD) with derate values of 4 and 6 %. The qualification threshold of 75 % for the composite score (red dashed line) is shown in the lower right panel. The average controlled system performance exceeds 75 % for both signal types, but only the RegD cases pass all of the time.

The results shown in Figs. 5–9 provide important insights into the possible strengths and limitations of the proposed approach to wind farm control for frequency regulation. These results suggest that wind farms may be well suited to act as a quickly responding resource for grid regulation services. For example, the consistent passing of the composite performance score for the RegD signals indicates that dynamic-model-controlled wind farms are able to provide this service reliably.

The power tracking results in Fig. 8 demonstrate that the controller is able to track the upregulation portions of the RegA signals at the beginning of the control period, such as during the first 5–10 min of the first two historical signals. In several cases the controlled LES wind farm is able to produce more power than the uncontrolled case, such as after minute 20 of the “RegA.D4.IC3.H1” and “RegD.D4.IC3.H1” simulations. However, when upregulation is requested for prolonged periods or towards the end of the control interval, such as the last 10 min of the “RegA.D4.IC3.TS” and “RegA.D4.IC3.H3” cases, the controller does not perform as well. A possible explanation is that the available energy in the wind is slowly changing as the atmospheric boundary layer evolves, as demonstrated by the declining power production of the uncontrolled simulations during these time periods. Since estimates of available energy are readily available over short time horizons, more frequent market clearing may allow wind farms to more effectively provide regulation. Ultimately, future work is needed to determine whether this is a fundamental limitation of the wind farm dynamics or the control strategy.

8 Conclusions

In this study we further characterize the performance of wind farms, providing secondary frequency regulation using the dynamic-model control framework proposed in Shapiro et al. (2017a). This model-based receding horizon approach relies on a simple one-dimensional time-varying wake model to provide thrust coefficient trajectories for individual turbines within a wind farm. As in previous work, the control approach is tested using a “virtual wind farm” represented by LESs of an 84-turbine wind farm with turbines modeled as actuator disks.

First, we evaluate the relative importance of including the dynamics of wake advection in the control scheme by comparing the performance of the dynamic-model controller to a comparable static-model controller. Tests using regulation signals from PJM indicate that the dynamic-model control demonstrates good overall tracking performance, whereas static-model control failed to match the reference signal for all simulated cases. These results indicate that the complexity of including the dynamics of wake advection is indeed required in model-based coordinated wind farm controls.

The tracking performance of the dynamic-model control method is then further quantified using PJM's performance metrics. Tests for both regulation signal types, RegA and RegD, exceed the PJM threshold for regulation participation on average, but only the RegD signal exceeds the threshold in all cases. These results indicate that this model-based receding horizon controller design could allow wind farms to meet industry design standards and allow wind farms to fully participate in regulation markets, particularly in fast-acting regulation markets.

The potential for reducing the derate required to participate in these regulation markets was also explored. Participating in frequency regulation markets currently requires a trade-off between revenue losses in the bulk power market and revenues generated in the frequency market. Previous approaches (Aho et al., 2013; Jeong et al., 2014) required power set-point reductions of an amount equal to the regulation amount, which directly reduces bulk power revenue by this amount. For this study we took a more aggressive approach by reducing the power set point by only 75 %, and even only 50 %, of the maximum regulation provided. For both of these derates, the controller is able to track fast-acting RegD signals. The potential for reducing the required derate has important economic implications for wind farms participating in both energy and regulation markets, a situation that will become increasingly common as more ISOs require wind farms to contribute to this grid service.

Although the dynamic-model controller design showed promising results, more work is needed to push this approach towards the implementation phase. We used the local thrust coefficient as a surrogate for real turbine control variables, such as generator torque and blade pitch angle. Improvements to our representation of these variables through actuator line methods and the inclusion of drivetrain dynamics in the control method are needed. Including rotational inertia may allow for further reductions in the amount of derate because rotational kinetic energy can compensate for short-term power shortages (De Rijcke et al., 2015). Furthermore, we assumed that the ISO provided the regulation signal at the beginning of the control period; however, PJM provides this reference at a 2 s scan rate. This shortcoming could be addressed by adding estimated reference trajectories to the control design. Finally, a systematic study of the relative advantages of all of the emerging control designs for wind farms to provide secondary frequency regulation (van Wingerden et al., 2017) is needed to identify which strategies are appropriate under various market, geographic, and technical constraints.