Journal topic
Wind Energ. Sci., 3, 11–24, 2018
https://doi.org/10.5194/wes-3-11-2018

Special issue: The Science of Making Torque from Wind (TORQUE) 2016

Wind Energ. Sci., 3, 11–24, 2018
https://doi.org/10.5194/wes-3-11-2018

Research article 22 Jan 2018

Research article | 22 Jan 2018

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

Wind farms providing secondary frequency regulation: evaluating the performance of model-based receding horizon control
Carl R. Shapiro1, Johan Meyers2, Charles Meneveau1, and Dennice F. Gayme1 Carl R. Shapiro et al.
• 1Department of Mechanical Engineering, Johns Hopkins University, 3400 N Charles St, Baltimore, Maryland, 21218, USA
• 2Department of Mechanical Engineering, KU Leuven, Celestijnenlaan 300A, 3001 Leuven, Belgium

Correspondence: Dennice F. Gayme (dennice@jhu.edu)

Abstract

This paper is an extended version of our paper presented at the 2016 TORQUE conference . 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 . 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 . 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 .

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 has shown that stand-alone wind turbines can effectively provide secondary frequency regulation, but recent fluid dynamics simulations 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 . 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 , but few of these studies provide real-time implementable algorithms that can respond to changing wind farm power output levels.

Our recent work sought to overcome these challenges by developing a time-varying extension of the classic Jensen wake model 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 . 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 , 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 .

The feasibility of providing secondary frequency regulation with wind farms was demonstrated by our initial results . 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 . 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 (PJM2012, 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 (PJM2012, 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 . 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=sn. Every turbine is assumed to have the same rotor diameter D.

Figure 1Diagram 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)}& {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)}& \mathit{\delta }{u}_{nm}\left(x\right)=\frac{\mathrm{2}{U}_{\mathrm{\infty }}{a}_{nm}}{{\left[{D}_{\mathrm{w}}\left(x-{s}_{n}\right)/D\right]}^{\mathrm{2}}},\end{array}$

where anm 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)}& {u}_{\mathrm{\infty }nm}={U}_{\mathrm{\infty }}-\sqrt{\sum _{\left(j,k\right)\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)}& {\stackrel{\mathrm{^}}{P}}_{nm}=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{C}_{\mathrm{P}nm}{u}_{\mathrm{\infty }nm}^{\mathrm{3}},\end{array}$

where CPnm 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 kn.

Furthermore, we express the turbine power production using the local thrust coefficient ${C}_{\mathrm{T}n}^{\prime }$ and modeled velocity at the turbine rotor ${\stackrel{\mathrm{^}}{u}}_{n}$. Simple momentum theory can be used to show that

$\begin{array}{}\text{(5)}& {a}_{n}=\frac{{C}_{\mathrm{T}n}^{\prime }}{\mathrm{4}+{C}_{\mathrm{T}n}^{\prime }},\phantom{\rule{0.25em}{0ex}}{\stackrel{\mathrm{^}}{u}}_{n}=\left(\mathrm{1}-{a}_{n}\right){u}_{\mathrm{\infty }n},\phantom{\rule{0.25em}{0ex}}\mathrm{and}\phantom{\rule{0.25em}{0ex}}{C}_{\mathrm{T}n}^{\prime }=\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}}/\left(\mathrm{1}-{a}_{n}{\right)}^{\mathrm{3}}$, from which we conclude ${C}_{\mathrm{T}n}^{\prime }$=${C}_{\mathrm{P}n}^{\prime }$. Replacing the induction factor an in Eq. (2), the modeled upstream velocity un in Eq. (3), and the power coefficient CPn in Eq. (4) with these equations, the power production can be rewritten as

$\begin{array}{}\text{(6)}& {\stackrel{\mathrm{^}}{P}}_{n}=M\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{C}_{\mathrm{T}n}^{\prime }{\stackrel{\mathrm{^}}{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 Meyers2015). Aerodynamic losses could also be taken into account by reducing the local power coefficient by a constant factor α , i.e., ${C}_{\mathrm{P}}^{\prime }$$\mathit{\alpha }{C}_{\mathrm{T}}^{\prime }$. For example, a wind turbine operating at a thrust coefficient of CT= 0.75 and CP= 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)}& \mathit{\delta }{u}_{n}\left(x,t\right)=\frac{{C}_{\mathrm{T}n}^{\prime }\left(t\right)}{\mathrm{4}+{C}_{\mathrm{T}n}^{\prime }\left(t\right)}\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)}& {\stackrel{\mathrm{^}}{u}}_{n}\left(t\right)=\left(\mathrm{1}-\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 ${\mathbit{W}}_{\mathrm{s}}\left({\mathbit{C}}_{\mathrm{T}}^{\prime }$, qs)=0, where qs= [δu, $\stackrel{\mathrm{^}}{\mathbit{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 . The velocity deficit is governed by

$\begin{array}{}\text{(9)}& \frac{\partial \mathit{\delta }{u}_{n}}{\partial t}+{U}_{\mathrm{\infty }}\frac{\partial \mathit{\delta }{u}_{n}}{\partial x}=-{w}_{n}\left(x\right)\mathit{\delta }{u}_{n}\left(x,t\right)+{f}_{n}\left(x,t\right),\end{array}$

where wn(x) is the wake decay function and fn(xt) 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)}& {w}_{n}\left(x\right)=\mathrm{2}\frac{{U}_{\mathrm{\infty }}}{{d}_{n}\left(x\right)}\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 dn(x)=Dwn(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 }\right)$ . In order to retrieve this expected velocity reduction, the forcing function is specified as

$\begin{array}{}\text{(11)}& {f}_{n}\left(x,t\right)=\frac{\mathrm{2}{U}_{\mathrm{\infty }}^{\mathrm{2}}}{{d}_{n}^{\mathrm{2}}\left(x\right)}\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(xsn) is a smoothing function that integrates to unity, centered at the streamwise location of the turbine x=sn. A Gaussian function with characteristic width Δ,

$\begin{array}{}\text{(12)}& G\left(x-{s}_{n}\right)=\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 , the dimensionless diameter of the wake generated by turbine row n is dn(x)= 1 + 2kn(xsn)∕D, where kn is an empirical wake expansion coefficient. We make two modifications to this equation. First, the linear expansion is assumed to begin at x=sn+ 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)}& {d}_{n}\left(x\right)=\mathrm{1}+{k}_{n}\mathrm{ln}\left[\mathrm{1}+\mathrm{exp}\left(\frac{x-{s}_{n}-\mathrm{2}\mathrm{\Delta }}{D/\mathrm{2}}\right)\right].\end{array}$

As in the static model, squared deficits are superposed to calculate the estimated streamwise velocity ${\stackrel{\mathrm{^}}{u}}_{n}$ at the turbine:

$\begin{array}{}\text{(14)}& {\stackrel{\mathrm{^}}{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}}\left(x,t\right)\right)}^{\mathrm{1}/\mathrm{2}}G\left(x-{s}_{n}\right)\mathrm{d}x.\end{array}$

Finally, the total estimated power ${\stackrel{\mathrm{^}}{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 ${\mathbit{W}}_{\mathrm{d}}\left({\mathbit{C}}_{\mathrm{T}}^{\prime }$, qd) =0, where qd= [δu, $\stackrel{\mathrm{^}}{\mathbit{u}}\right]$ 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 Pref(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 un(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 2Block diagram of the controlled wind farm system for both controllers. Each controller computes a thrust coefficient command signal ${\mathbit{C}}_{\mathrm{T}n}^{\prime }\left(t\right)$ using the reference signal Pref(t) and an error correction term ϵn(t). The error correction is computed using the measured velocity un(t) and the predicted velocity ${\stackrel{\mathrm{^}}{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)}& & \underset{{\mathbit{C}}_{\mathrm{T}}^{\prime },\mathbit{q}}{minimize}\phantom{\rule{0.25em}{0ex}}\mathcal{J}\left({\mathbit{C}}_{\mathrm{T}}^{\prime },\mathbit{q}\right)+\mathcal{R}\left({\mathbit{C}}_{\mathrm{T}}^{\prime }\right)\text{(16)}& & \mathrm{subject}\phantom{\rule{0.25em}{0ex}}\mathrm{to}\phantom{\rule{0.25em}{0ex}}\mathbit{W}\left({\mathbit{C}}_{\mathrm{T}}^{\prime },\mathbit{q}\right)=\mathrm{0},\end{array}$

where $\mathbit{W}\left({\mathbit{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)}& {\mathcal{J}}_{\mathrm{s}}\left({\mathbit{C}}_{\mathrm{T}}^{\prime },{\mathbit{q}}_{\mathrm{s}}\right)=\frac{\mathrm{1}}{{\mathcal{P}}^{\mathrm{2}}}{\left(\sum _{n=\mathrm{1}}^{N}{\stackrel{\mathrm{^}}{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)}& {\mathcal{R}}_{\mathrm{s}}\left({\mathbit{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(\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−5 in this study.

The dynamic-model controller , 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 TA before re-solving the minimization problem. More details about this procedure can be found in and . With this predictive framework, the reference tracking goal is represented by the cost functional

$\begin{array}{}\text{(19)}& {\mathcal{J}}_{\mathrm{d}}\left({\mathbit{C}}_{\mathrm{T}}^{\prime },{\mathbit{q}}_{\mathrm{d}}\right)=\frac{\mathrm{1}}{{\mathcal{P}}^{\mathrm{2}}T}\underset{\mathrm{0}}{\overset{T}{\int }}{\left(\sum _{n=\mathrm{1}}^{N}{\stackrel{\mathrm{^}}{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}{\mathcal{R}}_{\mathrm{d}}\left({\mathbit{C}}_{\mathrm{T}}^{\prime }\right)& =\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)}& & +\mathit{\gamma }T\sum _{n=\mathrm{1}}^{N}\underset{\mathrm{0}}{\overset{T}{\int }}{\left(\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({\mathbit{C}}_{\mathrm{T}}^{\prime }\right)$=𝒥(q, ${\mathbit{C}}_{\mathrm{T}}^{\prime }\right)$ , instead of the cost functional defined in Eq. (15). Minimizations are performed using the gradient-based nonlinear Polak–Ribière conjugate gradient method (Press2007) combined with the Moré–Thuente line search method 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 . The full procedure is detailed in . 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 TA 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 TA= 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 . 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)}& {P}_{n}=\sum _{m=\mathrm{1}}^{M}{P}_{nm},\phantom{\rule{0.25em}{0ex}}\mathrm{and}\phantom{\rule{0.25em}{0ex}}{u}_{n}=\frac{\mathrm{1}}{M}{\left(\sum _{m=\mathrm{1}}^{M}{u}_{nm}^{\mathrm{3}}\right)}^{\mathrm{1}/\mathrm{3}},\end{array}$

where unm 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 Pn=$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)}& {\stackrel{\mathrm{^}}{P}}_{n}=M\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\frac{\mathit{\pi }{D}^{\mathrm{2}}}{\mathrm{4}}{C}_{\mathrm{T}n}^{\prime }{\left({\stackrel{\mathrm{^}}{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)}& {\mathit{ϵ}}_{n}^{p}={u}_{n}^{p-\mathrm{1}}-{\stackrel{\mathrm{^}}{u}}_{n}^{p-\mathrm{1}}.\end{array}$

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

$\begin{array}{}\text{(24)}& {\mathit{ϵ}}_{n}\left(t\right)=\left({u}_{n}\left({t}_{\mathrm{c}}\right)-{\stackrel{\mathrm{^}}{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 .

These simulations are performed using Johns Hopkins University's LES code LESGO , 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 . Inlet conditions for the wind farm are generated using the concurrent-precursor method . 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 unm and the thrust coefficient ${C}_{\mathrm{T}nm}^{\prime }$. The force is modeled as a drag force Fnm=$-\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 Pnm=Fnmunm. An instantaneous color contour plot of the streamwise velocity field from one of these simulations is shown in Fig. 3.

Figure 3Instantaneous 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−1, respectively. The inlet conditions are generated using a concurrent precursor simulation shown at the beginning of the plotted domain.

Figure 4Probability density functions of 𝒮1𝒮3 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 (PJM2012, 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 Pref(t)= [1 α+ 0.08r(t)]Pbase, where Pbase 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 Pbase.

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}& {\mathcal{S}}_{\mathrm{1}}=\frac{\mathrm{1}}{T}\underset{\mathrm{0}}{\overset{T}{\int }}r\left(t\right)\mathrm{d}t\\ & {\mathcal{S}}_{\mathrm{2}}=\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)}& & {\mathcal{S}}_{\mathrm{3}}=\frac{\mathrm{1}}{T}\underset{\mathrm{0}}{\overset{T}{\int }}{\left(\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, 𝒮1 is the mean of r(t), 𝒮2 is the variance of r(t), and 𝒮3 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 𝒮1 of the RegA signals, represents a more difficult test for the controlled wind farm because more energy is requested than the average.

Table 1Test 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 2Characteristics 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 kn 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−1 and a standard deviation σu 1 m s−1 as measured at the first row of turbines during the T0= 5 min prior to initiation of the control:

$\begin{array}{ll}& \stackrel{\mathrm{‾}}{u}=\frac{\mathrm{1}}{{T}_{\mathrm{0}}M}\underset{-{T}_{\mathrm{0}}}{\overset{\mathrm{0}}{\int }}\sum _{m=\mathrm{1}}^{M}\frac{\mathrm{4}+{C}_{\mathrm{T},\mathrm{ref}}^{\prime }}{\mathrm{4}}{u}_{\mathrm{1}m}\left(t\right)\mathrm{d}t\\ \text{(26)}& & {\mathit{\sigma }}_{u}={\left[\frac{\mathrm{1}}{{T}_{\mathrm{0}}M}\underset{-{T}_{\mathrm{0}}}{\overset{\mathrm{0}}{\int }}\sum _{m=\mathrm{1}}^{M}{\left(\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 (IEC2005). The simulations assume region 2 operation 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−1, 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 CP= 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 CP= 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 kn are also calculated for each initial condition using measurements from the T0= 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 }\right){T}_{\mathrm{0}}^{-\mathrm{1}}{\int }_{-{T}_{\mathrm{0}}}^{\mathrm{0}}{u}_{\mathrm{1}}\left(t\right)\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 for the last row. This operating point is likely the optimal power point for the Jensen model with constant wake expansion coefficients.

Figure 5Comparison 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 6Comparison 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 7Comparison 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 8Power 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 (PJM2012, 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 𝒮1, 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 9Box 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. 59 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 . 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 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 . 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.

Data availability
Data availability.

Data from the simulations, including the disk-averaged velocity, local thrust coefficient, and disk-center velocity, are provided in the data repository .

Author contributions
Author contributions.

CS designed the study, performed the simulations, analyzed results, and wrote the first draft of the paper. JM, CM, and DG helped design the study, analyzed results, and contributed to writing the paper.

Competing interests
Competing interests.

J. Meyers is a member of the editorial board of the journal.

Special issue statement
Special issue statement.

This article is part of the special issue “The Science of Making Torque from Wind (TORQUE) 2016”. It is a result of the The Science of Making Torque from Wind (TORQUE 2016), Munich, Germany, 5–7 October 2016.

Acknowledgements
Acknowledgements.

Carl R. Shapiro, Charles Meneveau, and Dennice F. Gayme are supported by the National Science Foundation (grant nos. ECCS-1230788, CMMI 1635430, and OISE-1243482, the WINDINSPIRE project). Johan Meyers is supported by the European Research Council (ActiveWindFarms, grant no. 306471). This research project was conducted using computational resources at the Maryland Advanced Research Computing Center (MARCC).

Edited by: George Kariniotakis
Reviewed by: two anonymous referees

References

Aho, J., Buckspan, A., Laks, J., Fleming, P., Jeong, Y., Dunne, F., Churchfield, M., Pao, L., and Johnson, K.: A tutorial of wind turbine control for supporting grid frequency through active power control, in: American Control Conference, Montreal, Canada, 3120–3131, https://doi.org/10.1109/ACC.2012.6315180, 2012. a, b

Aho, J., Buckspan, A., Pao, L., and Fleming, P.: An active power control system for wind turbines capable of primary and secondary frequency control for supporting grid reliability, in: AIAA Aerospace Sciences Meeting, Grapevine, TX, 2013. a, b, c, d

Annoni, J., Howard, K., Seiler, P., and Guala, M.: An experimental investigation on the effect of individual turbine control on wind farm dynamics, Wind Energy, 19, 1453–1467, https://doi.org/10.1002/we.1930, 2016. a

Bewley, T. R., Moin, P., and Temam, R.: DNS-based predictive control of turbulence: An optimal benchmark for feedback algorithms, J. Fluid Mech., 447, 179–225, 2001. a, b

Borzì, A. and Schulz, V.: Computational Optimization of Systems Governed by Partial Differential Equations, SIAM, Philadelphia, PA, 2011. a

Bou-Zeid, E., Meneveau, C., and Parlange, M.: A scale-dependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows, Phys. Fluids, 17, 025105, https://doi.org/10.1063/1.1839152, 2005. a

Calaf, M., Meneveau, C., and Meyers, J.: Large eddy simulation study of fully developed wind-turbine array boundary layers, Phys. Fluids, 22, 015110, https://doi.org/10.1063/1.3291077, 2010. a, b, c

De Rijcke, S., Driesen, J., and Meyers, J.: Power smoothing in large wind farms using optimal control of rotating kinetic energy reserves, Wind Energy, 18, 1777–1791, https://doi.org/10.1002/we.1790, 2015. a

Díaz-González, F., Hau, M., Sumper, A., and Gomis-Bellmunt, O.: Participation of wind power plants in system frequency control: Review of grid code requirements and control methods, Renew. Sustain. Energy Rev., 34, 551–564, https://doi.org/10.1016/j.rser.2014.03.040, 2014. a

Fleming, P., Aho, J., Gebraad, P., Pao, L., and Zhang, Y.: Computational fluid dynamics simulation study of active power control in wind plants, in: American Control Conference, Boston, MA, 1413–1420, https://doi.org/10.1109/ACC.2016.7525115, 2016. a

Gebraad, P. M. O., Fleming, P. A., and van Wingerden, J. W.: Wind turbine wake estimation and control using FLORIDyn, a control-oriented dynamic wind plant model, in: American Control Conference, Chicago, IL, 1702–1708, https://doi.org/10.1109/ACC.2015.7170978, 2015. a

Goit, J. P. and Meyers, J.: Optimal control of energy extraction in wind-farm boundary layers, J. Fluid Mech., 768, 5–50, 2015. a, b, c, d, e, f

IEC: Wind turbines – Part 1: Design requirements, Tech. Rep. IEC61400-1, Geneva, Switzerland, 2005. a

Jeong, Y., Johnson, K., and Fleming, P.: Comparison and testing of power reserve control strategies for grid-connected wind turbines, Wind Energy, 17, 343–358, https://doi.org/10.1002/we.1578, 2014. a, b, c, d

Johnson, K., Pao, L. Y., Balas, M., and Fingersh, L.: Control of variable-speed wind turbines: standard and adaptive techniques for maximizing energy capture, IEEE Control Syst. Mag., 26, 70–81, https://doi.org/10.1109/MCS.2006.1636311, 2006. a

Katić, I., Højstrup, J., and Jensen, N.: A simple model for cluster efficiency, in: European Wind Energy Association Conference and Exhibition, Rome, Italy, 407–410, 1986. a, b, c, d, e

Meyers, J. and Meneveau, C.: Large eddy simulations of large wind-turbine arrays in the atmospheric boundary layer, in: 50th AIAA Aerospace Sciences Meeting, Orlando, FL, 2010. a

Moré, J. J. and Thuente, D. J.: Line search algorithms with guaranteed sufficient decrease, ACM Trans. Math. Softw., 20, 286–307, 1994. a

PJM: PJM Manual 11: Energy & Ancillary Services Market Operations, PJM, Audubon, PA, 2012. a, b, c, d

PJM: PJM Manual 12: Balancing Operations, PJM, Audubon, PA, 2015. a, b, c, d

PJM: Ancillary Services, http://www.pjm.com/markets-and-operations/ancillary-services.aspx, last access: 27 August 2016. a

Press, W. H.: Numerical recipes 3rd edition: The art of scientific computing, Cambridge University Press, Cambridge, UK, 2007. a

Rebours, Y., Kirschen, D., Trotignon, M., and Rossignol, S.: A Survey of Frequency and Voltage Control Ancillary Services – Part I: Technical Features, IEEE T. Power Syst., 22, 350–357, https://doi.org/10.1109/TPWRS.2006.888963, 2007.  a

Rose, S. and Apt, J.: The cost of curtailing wind turbines for secondary frequency regulation capacity, Energy Systems, 5, 407–422, https://doi.org/10.1007/s12667-013-0093-1, 2014. a

Shapiro, C. R., Meyers, J., Meneveau, C., and Gayme, D. F.: Wind farms providing secondary frequency regulation: Evaluating the performance of model-based receding horizon control, J. Phys. Conf. Ser., 753, 052012, https://doi.org/10.1088/1742-6596/753/5/052012, 2016. a

Shapiro, C. R., Bauweraerts, P., Meyers, J., Meneveau, C., and Gayme, D. F.: Model-based receding horizon control of wind farms for secondary frequency regulation, Wind Energy, 20, 1261–1275, https://doi.org/10.1002/we.2093, 2017a. a, b, c, d, e, f, g, h, i, j

Shapiro, C. R., Meyers, J., Meneveau, C., and Gayme, D. F.: Dynamic wake modeling and state estimation for improved model-based receding horizon control of wind farms, in: American Control Conference, Seattle, WA, https://doi.org/10.23919/ACC.2017.7963036, 2017b. a

Shapiro, C. R., Meyers, J., Meneveau, C., and Gayme, D. F.: Data associated with publication “Wind farms providing secondary frequency regulation: Evaluating the performance of model-based receding horizon control”, Johns Hopkins University Data Archive Dataverse, https://doi.org/10.7281/T1/TLEBVW, 2017. a

Stevens, R. J. A. M. and Meneveau, C.: Temporal structure of aggregate power fluctuations in large-eddy simulations of extended wind-farms, J. Renew. Sustain. Energy, 6, 043102, https://doi.org/10.1063/1.4885114, 2014. a

Stevens, R. J. A. M., Gayme, D. F., and Meneveau, C.: Large eddy simulation studies of the effects of alignment and wind farm length, J. Renew. Sustain. Energy, 6, 023105, https://doi.org/10.1063/1.4869568, 2014a. a, b

Stevens, R. J. A. M., Graham, J., and Meneveau, C.: A concurrent precursor inflow method for Large Eddy Simulations and applications to finite length wind farms, Renewable Energy, 68, 46–50, https://doi.org/10.1016/j.renene.2014.01.024, 2014b. a, b, c

van Wingerden, J.-W., Pao, L., Aho, J., and Fleming, P.: Active Power Control of Waked Wind Farms, in: International Federation of Automatic Control World Congress, Toulouse, France, 2017. a

VerHulst, C. and Meneveau, C.: Altering Kinetic Energy Entrainment in Large Eddy Simulations of Large Wind Farms Using Unconventional Wind Turbine Actuator Forcing, Energies, 8, 370–386, 2015. a