This work introduces the area localized coupled (ALC) model, which extends the applicability of approaches that couple classical wake superposition models and atmospheric boundary layer models to wind farms with arbitrary layouts. Coupling wake and top–down boundary layer models is particularly challenging since the latter requires averaging over planform areas associated with turbine-specific regions of the flow that need to be specified. The ALC model uses Voronoi tessellation to define this local area around each turbine. A top–down description of a developing internal boundary layer is then applied over Voronoi cells upstream of each turbine to estimate the local mean velocity profile. Coupling between the velocity at hub-height based on this localized top–down model and a wake model is achieved by enforcing a minimum least-square-error in mean velocity in each cell. The wake model in the present implementation takes into account variations in wind farm inflow velocity and represents the wake profile behind each turbine as a super-Gaussian function that smoothly transitions between a top-hat shape in the region immediately following the turbine to a Gaussian profile downstream. Detailed comparisons to large-eddy simulation (LES) data from two different wind farms demonstrate the efficacy of the model in accurately predicting both wind farm power output and local turbine hub-height velocity for different wind farm geometries. These validations using data generated from two different LES codes demonstrate the model's versatility with respect to capturing results from different simulation setups and wind farm configurations.

## I. INTRODUCTION

Analytical modeling of wind farms has traditionally fallen into two main categories: (1) wake models that focus on the velocity deficits (wakes) of turbines and interactions between the turbine wakes and (2) top–down representations that examine the overall effect of the farm on the atmospheric boundary layer (ABL). Wake models describe the magnitude, profile, and evolution of the velocity deficits arising from the mean-flow kinetic energy and momentum extraction by individual turbines (the wake). This velocity deficit profile and wake growth law is then combined with a superposition law, which is used to capture wake interactions between upstream turbines in order to predict the mean velocity seen by each turbine in the array. There is a vast amount of literature proposing different velocity deficit representations that range from analytical functions^{1–5} to experimental^{6,7} or data-driven parameterizations.^{8} A number of different superposition laws for the deficits^{9–12} have also been proposed, with the most common ones being quadratic or linear superposition. For example, the widely used Jensen–Katic model^{1,13} assumes a top-hat velocity profile, linear wake growth, and a square superposition of wake velocity deficits.

On the other hand, top–down representations do not model detailed array velocity distributions but instead consider large scale interactions of the wind farm with the ABL. Many of these models are based on the notion that in the “fully developed” region of a large wind farm, the boundary layer is no longer growing and the wind farm is in balance with the ABL. This concept was introduced by Templin in 1974 (Ref. 14). A complete model was proposed by Frandsen,^{15} who specified the last rows (latter portions) of the farm as the fully developed region, where the farm causes a modification of the mean velocity profile in the ABL. A number of subsequent models build upon the Frandsen top–down model. Such models employ a more detailed description of the mean velocity profile in the turbine array boundary layer,^{16} while the theory of developing boundary layers^{17} enabled the model to be applied in both fully developed and developing regions of the wind farm. Top–down models have also been extended to include various atmospheric stratification conditions^{18,19} and used in larger atmospheric simulations.^{20}

In general, wake models are applied to predict the power of individual turbines and wake behavior in the farm, while top–down models provide an average view that gives information regarding the physics of the ABL flow and the limiting effects that begin to affect large wind farms (deep array effects). Coupled models that seek to leverage the relative strengths of both approaches have been shown to improve the accuracy of the hub-height velocity estimates and total wind farm power output predictions over either type of model in isolation. The Frandsen model^{21} represents a precursor to this approach in which it linked the top–down model in the fully developed region to a wake model at the front of a farm through an intermediate region. While this model recognized the benefits of both wake and top–down models for different wind farm flow regimes, the two models were not integrated and the information exchange between the regions was limited to the transition from one to the other over the intermediate region.

The coupled wake boundary layer (CWBL) model proposed by Stevens *et al.*^{22,23} demonstrated that two-way coupling of a wake and top–down model led to further improvements in power output predictions. This model and subsequent extensions showed excellent agreement with large-eddy simulations (LES) of operational wind farms.^{23} As with the Fransden^{21} model, this CWBL required *a priori* identification of developing and fully developed regions of the farm. In addition, the CWBL model^{23} assumed a regular rectangular array of wind turbines since regular streamwise (*s _{x}*) and spanwise (

*s*) spacings needed to be prescribed, precluding its direct application to wind farms with arbitrary turbine layouts.

_{y}Many wind farm sites have topographical features that dictate irregular wind turbine arrangements, which complicates the *a priori* identification of fully developed regions and motivates the development of models that can be applied to a wide-range of turbine configurations. In the work of Shapiro *et al.*,^{11} an earlier form of the area localized coupled (ALC) model was proposed, which built on the ideas of the CWBL model. This model took important steps toward full generalization to arbitrary layouts by localizing the top–down portion of the model to each turbine instead of averaging over an entire region of the wind farm. A local top–down model for the developing boundary layer^{17} then provides a description of the growth of internal boundary layers (IBLs) for each turbine. This method enables the model to be applied in both the developing and fully developed region of the farm and eliminates the need for an empirical specification of the evolution of the wake growth parameter as a function of distance into the farm. The Shapiro *et al.*^{11} model also improved upon the CWBL implementation through enhancements of the wake model including a super-Gaussian velocity deficit (wake) profile that first approximates a top-hat profile near the turbine and then transitions smoothly toward a Gaussian profile^{2,11,24} further downstream. Linear superposition was used to account for wake interactions.^{11}

In this paper, we further advance the goals of increasing model applicability to a wider range of wind farm geometries and inflow conditions. In particular, the ALC model introduced here redefines the identification of the freestream turbines and the planform thrust coefficient based on the localized turbine areas by combining a local formula with information about the flow upstream. This ALC model also directly incorporates non-uniform inflow velocity across the freestream turbines as a further improvement over the preliminary method of accounting for these differences through the addition of an estimation algorithm.^{11} Directly accounting for these variations in the inlet velocity across the farm (transverse direction) as a function of the wind inlet direction further improves the velocity deficit estimates at both freestream and downstream turbines.

The proposed localized approach employed in the ALC and Shapiro *et al.*^{11} models has the added advantage of providing significantly more information about the flow field than conventional wake models, top–down models, and the original CWBL implementation.^{23} These prior approaches do not typically furnish information about the velocity field such as turbulence conditions (friction velocities above and below turbine layer) within each cell. They instead provide only mean velocities or use empirical correlations for turbulence intensities, which make it more difficult to quantify loading characteristics for the turbines.

In this work, the ALC model is applied to two wind farms, and its predictions are compared to data from two LES. The first wind farm is arranged in a circular configuration that provides an example of a non-rectangular array with streamwise and spanwise spacings between turbines varying throughout the farm for any given wind direction. We provide a detailed analysis of the velocity predictions for each turbine in this circular farm over a range of wind inflow directions to illustrate the use of the local turbine area in the computations and to evaluate the accuracy of the predictions of the model at both the turbine and farm level. A second wind farm consisting of two regions, one with regular and the other with random turbine placements, is also considered under a single wind direction. The results illustrate the benefits of the area localized approach in coupling physics-based models for both the turbine wakes and the ABL and providing good predictions of the local velocity field and power produced by individual turbines as well as the total wind farm production.

An overview of the ALC model is presented in Sec. II with brief descriptions of both the wake and the top–down models that constitute the building blocks of the ALC model. In Sec. III, we test the model by comparing its predictions to LES data from a circular wind farm. In Sec. IV, we compare the model results with LES data from a wind farm that contains two regions, one with a regular array and another with a random distribution of turbines. Finally, in Sec. V, we present conclusions .

## II. AREA LOCALIZED COUPLED (ALC) MODEL

The ALC model couples a wake model with a super-Gaussian wake profile and linear wake superposition, with a top–down model for a developing wind farm. Each turbine has its own instance of the top–down model that is localized to its own “area” and described by the local planform thrust coefficient, which represents the momentum extracted from the flow in a developing internal boundary layer upstream of each of the turbines by accounting for their axial flow resistance (the localized internal boundary layer concept distinguishes the ALC model from the Shapiro *et al.*^{11} approach). The coupling is imposed through a minimization of the difference between the area average velocities computed from the two models, see Fig. 1. This matching condition is used to determine the wake expansion coefficients for each local turbine cell. Sections II A–II C summarize the wake model and the top–down model and their coupling in the ALC model, see Fig. 1.

### A. Wake model

The aim of the wake model is to calculate the mean streamwise velocity distribution in the entire wind farm as a function of position. The wake model used here is developed in the work of Shapiro *et al.*^{11} When evaluated at turbine hub locations, the velocity can be used to predict the power generated by each turbine. This flow field is described as a superposition of individual wakes by the following equation:

where *x* is the streamwise coordinate (in the direction of the incoming freestream wind), $x\u2032=x\u2212sn,x$ is the coordinate relative to the position of the turbine located at $x=sn,x$, and $\u2009r\u2032=[(y\u2212yn)2+(z\u2212zn)2]1/2$ is the radial distance from the center position of the $nth$ turbine. $U\u221e(y,z)$ is the incoming freestream mean velocity which can vary across the inflow plane (*y* is the horizontal coordinate transverse to the incoming wind) and, in general, also as function of vertical coordinate *z*. For the wake model used in this paper for simplicity, we will neglect the variations of $U\u221e(y,z)$ in the vertical direction and set $U\u221e(y)=U\u221e(y,zh)$, where *z _{h}* is the hub-height. $\delta un(x)$ is the deficit velocity of the n

^{th}turbine at downstream position

*x*, and $Wn(x,r)$ is the shape and superposition function of wake

*n*. The velocity deficit $\delta un(x)$ is calculated as the steady state solution to a one-dimensional partial differential equation (PDE) that describes the more generally valid time-dependent behavior of the wake,

The equation describes the rate of change of $\delta un$ when moving downstream at speed $U\u221e$ in response to two effects: (1) the transverse diffusion expressed as a wake expansion term (first term on the equation's right-hand-side) leading to a decrease in $\delta un$ and (2) the creation of a wake deficit at and near the turbine locations due to the axial induction (second term). Under the simplifying assumption of steady state flow conditions, the solution to Eq. (2) has the form

where $\delta u0,n$ is the initial velocity deficit at the turbine location, $dw,n(x)$ is the wake expansion function (the wake diameter divided by the turbine diameter), and Δ is the characteristic width of the Gaussian function used to smooth the forcing in Eq. (2) along the streamwise direction. The wake expansion function is assumed to be of the following form:

which at large *x*/*R* tends to the classical linear growth $dw,n(x)=kw,n(x/R)$ but smoothly merges to $dw,n(x)=1$ in the near-turbine region, thus preventing values that are not physical. The parameter $kw,n$ is the turbine specific wake expansion coefficient for turbine *n* that determines the rate at which the wake expands as it travels downstream. The expansion of each individual turbine wake depends on local flow conditions. For instance, one would expect it should increase with increasing local turbulence intensity. In the ALC model, the coefficient will be calculated with an expression of the form

where $u\u221e,n$ is the *y*-dependent upstream velocity at turbine n and $u*,n$ is related to the local friction velocity representing the state of turbulence. The precise form for the local friction velocity to be used in the model is derived from the top–down model as described in Sec. II B. The coefficient *α* is a model parameter that will be determined through the coupling of the wake and top–down models.

We incorporate effects of spanwise variations in the inflow velocity through the “upstream” velocity $u\u221e,n$ of each turbine. To calculate $u\u221e,n$, i.e., the velocity that would have existed at the disk of wind turbine *n* without the turbine there, we use the following expression:

where *R* is the radius of the turbine, $sn,x$ is, again, the position of the turbine in the streamwise direction, and $r\u2032=[(y\u2212yn)2+(z\u2212zn)2]1/2$ is the radial coordinate of point (*y*, *z*) for turbine *n* with rotor center at (*y _{n}*,

*z*). This quantity is computed as follows. We first define the turbines that are not in the wake of others for a given wind direction as “freestream turbines,” to which we ascribe a

_{n}*y*dependent freestream velocity $u\u221e,n=U\u221e(yn)$. Here, $U\u221e(yn)$ denotes the average velocity across the turbine disk based on the incoming flow field. The sum is then performed over the remaining (waked) turbines and includes all turbines upstream of turbine

*n*. The nonzero terms in the sum are determined by the wake shape function

*W*, which takes the form of a super-Gaussian,

_{n}where $r2=y2+z2$. This function smoothly transitions from a top-hat wake profile immediately following the turbine into a Gaussian wake profile downstream. The following functional form for *p*(*x*) was proposed to reproduce the shape and rate of transition from a top-hat to a Gaussian seen in experimental wake data:

The function *C*(*x*) can then be shown^{11} to have the following form assuming mass conservation (or linearized momentum conservation):

Finally, the initial velocity deficit is obtained from an inviscid model and has the traditional form based on actuator disk momentum theory,

where $CT,n$ is the coefficient of thrust for the turbine, which determines the amount of thrust the turbine extracts from the flow. $CT,n\u2032$ is the local coefficient of thrust and is related to the coefficient of thrust through the relation $CT,nU\u221e,n2=CT,n\u2032ud,n2$, where $ud,n$ is the average velocity over the turbine disk area, which is defined later in Eq. (23).

In the present work, we use the steady state solution to Eq. (2) and we do not take into account the time-dependence. However, if $\delta un(x,t)$ were to be determined from the time-dependent PDE, the ALC model can be readily extended to time-varying applications.

### B. Top–down model

As previously discussed, the top–down model^{16,21} regards the wind farm as a homogenized surface roughness in the atmospheric boundary layer. Assuming a fully developed (horizontally homogeneous) wind turbine array boundary layer (WTABL),^{16} the velocity profile can be shown to consist of two constant shear stress regions, one below the turbines and one above the turbines. Each layer has friction velocities, $u*,lo$ and $u*,hi$, and roughness heights, $z0,lo$ and $z0,hi$, respectively. The model connects the two regions through the wind turbine layer.^{16} A momentum balance results in the following relationship between the two friction velocities and the planar average velocity $u\xafh$ at the hub-height *z _{h}* of the turbines:

where *c _{ft}* is the planform thrust coefficient that represents the momentum extracted from the flow by the turbines' axial flow resistance per unit horizontal area. This coefficient is obtained from the ALC model coupling, as described in Sec. II C.

In the wind turbine layer, which is defined as $zh\u2212R\u2264z\u2264zh+R$, where *R* is the radius of the turbines, the turbulent flow is assumed to include the effects of an additional eddy viscosity $\nu w*$. The total eddy viscosity in this region is defined in terms of the friction velocities in the top and bottom halves of the turbine region. In the lower region, where $zh\u2212R\u2264z\u2264zh$, the eddy viscosity is represented as $(1+\nu w*)\kappa u*,lo$, and in the upper region, $zh\u2264z\u2264zh+R$, it is represented as $(1+\nu w*)\kappa u*,hi$. The extra eddy viscosity needs to be modeled and it was assumed^{16} to be given as $\nu w*=2812cft$. Using the known value of the lower surface roughness height $z0,lo$ and the momentum balance [Eq. (11)] and enforcing velocity continuity across the vertical profile, one can derive the roughness height due to the wind farm as,^{16}

where $\zeta =ln\u2009[(zh/z0,lo)(1\u2212R/zh)\beta ]$ and $\beta =\nu w*/(1+\nu w*)$. Using this result, we can also calculate the friction velocities in both of the layers

where $u*$ is the friction velocity of the incoming ABL and *δ* is the overall ABL height.

Note that the top–down model is based on the planform-averaged momentum equation and therefore requires averaging over extended horizontal areas. In order to apply the model in a more localized fashion, horizontal areas associated with each of the turbines must be defined. We define these areas using Voronoi tessellation (also used in the previous model^{11}) which naturally associates cells and their areas to each turbine *n*. In particular, each turbine is a node, and the vertices are defined as points that are equidistant from three separate nodes. Figure 2 shows the local areas resulting from our Voronoi tessellation procedure applied to several different configurations of turbines: a regular array, a staggered array, a circular configuration, and an entirely random distribution. The edge cells are defined by using ghost points projected outside of the array.^{11}

The Voronoi cells separate the wind farm into planform areas belonging to each turbine in a way that can be generalized to any wind farm layout. The ALC model exploits this local turbine area definition to compute friction velocities and roughness heights individually in each Voronoi cell obtaining information on flow conditions in the areas around each turbine. Information related to the turbulence of the local flow field is then used for each turbine and cell to determine appropriate wake expansion coefficients for the wake model.

In order to compute the local friction velocities, we use the notion of a developing IBL over a wind farm^{17,25} that depends on the streamwise position *x*. The IBL height is modeled according to

where $x\u2212x0$ represents the distance to the beginning of the IBL, i.e., the distance to the start of the farm directly upstream to turbine *n*, and *δ* is the height of the overall ABL. This implies that the friction velocities evolve as a function of *x* through the farm, until *δ _{ibl}* reaches the final boundary layer height

*δ*in the fully-developed region. Analysis of the developing wind farm internal boundary layer

^{17,26}shows that the friction velocity evolves similarly to that of a boundary layer flow over a surface with a smooth-to-rough transition, i.e., the friction velocity increases sharply at the transition (the front of the wind farm) but then gradually decreases again as the IBL grows and the mean velocity gradient decreases. The analysis

^{17}shows that the position-dependent friction velocity $u*,hi(x)$ can be obtained from Eq. (13) by replacing the overall

*δ*by $\delta ibl(x)$ for turbine

*n*. That is, to say, we will evaluate $u*,hi,n$ for each individual turbine by replacing

*δ*in Eq. (13) by $\delta ibl(sn,x)$ evaluated at the position of turbine

*n*.

When dealing with an irregular wind farm, special care must be taken to properly define the starting position of the internal boundary layer. Here, we choose to define the IBL start to be located at the freestream turbines, i.e., the turbines in the farm that only see the upstream incoming velocity and are not affected by the wakes of other turbines. The process to identify freestream turbines is described in Sec. II C. Figure 3 shows how we model this as a boundary layer initialization, defined for the case of a circular wind farm. In this case, the start (“trip line”) of the IBL is defined at the first set of turbines to see the freestream velocity. Let $xn,0$ be the *x*-location on the trip-line directly upstream of turbine *n*. The streamwise distance between any given turbine and the boundary layer trip line for turbine *n* is thus given by $sn,x\u2212xn,0$, which is shown with the blue arrow for a sample turbine in the farm.

In the localized top–down model, calculations can now be performed for each turbine in its own turbine cell, enabling the localization needed to represent both the spatial non-uniformities due to arbitrary turbine placements as well as nonuniform (*y*-dependent) velocity inflow. In the case of the freestream turbines, the freestream values of the flow are used. For the remaining turbines, we use the developing boundary layer framework. For a standard boundary layer, the friction velocity of a boundary layer flow over a rough surface with roughness $z0,lo$ and inflow velocity $U\u221e$ (mean velocity at hub-height *z* = *z _{h}* upstream of the wind farm) would be given by $u*=U\u221e\kappa /\u2009ln\u2009(zh/z0,lo)$. If the inflow is

*y*-dependent, and if the freestream inflow velocity corresponding to turbine

*n*is denoted as $U\u221e,n$ (it can be computed by averaging the inflow over the cell), we can write $u*,n=U\u221e,n\kappa /\u2009ln\u2009(zh/z0,lo)$. Replacing in Eq. (13) to obtain the actual friction velocity for turbine

*n*, we can write

Once the friction velocity is known, the top–down model provides a prediction for the mean velocity at hub height,^{16} and when applied individually to each cell, the mean velocity in the cell associated with turbine n is given by

Note that in order to evaluate $u\xafh,ntd$ we require $z0,hi,n$ which according to Eq. (12) depends on the planform thrust coefficient *c _{ft,n}* in the cell associated with turbine n. This value can differ from turbine to turbine since the total turbine forces and momentum exchanges affecting the development of the local internal boundary layer may differ across the wind farm. In order to determine $cft,n$ for each individual turbine

*n*, information from the wake model (Sec. II A) is required as described in Sec. II C where the coupling of the constituent models is presented.

### C. Coupling of wake and top–down models: ALC Model

The wake model and the top–down model are coupled by comparing their respective predictions of the average planar velocity values in each waked Voronoi cell. In the top–down model, these are the velocities $u\xafh,ntd$ calculated in each Voronoi cell according to Eq. (17). The output of the wake model is a velocity field $u(x,y,z)$ from Eq. (1). In order to compare the wake model to the top–down model predictions, a cell-averaged wake model velocity is defined by averaging the velocity field predicted by Eq. (1) at *z* = *z _{h}* over each cell for each turbine

*n*. We denote this cell averaged wake model velocity as $u\xafh,nwm$. The steps are illustrated in Fig. 4. The average planar velocities in each cell are illustrated in Fig. 4(c).

Ideally one would want the velocities predicted by the top–down and wake methods to yield the same cell averaged velocities, for each turbine. Free model parameters are adjusted to minimize the least-square error between both predictions from the two approaches. Specifically, recall that the parameter *α* required to specify the wake expansion parameter *k _{w,n}* in Eq. (5) was left unspecified. We now follow the basic idea of the CWBL approach

^{22}to determine the wake expansion coefficient using the condition that both approaches match. Since in the ALC model has many cells, we find the

*α*that minimizes the square difference between the average planar velocities from the two models over the cells of turbines that are in the wake of upstream turbines (called “waked cells”), according to

where *wake* defines the set of waked cells. The summation does not include freestream turbines because for those turbines the top–down and wake models would use the same value of *k _{w}* which is unaffected by wakes and expected to be different from those inside the farm. The freestream turbines are found by setting a threshold on the velocity deficit of 1%. Specifically, turbine

*n*is denoted as a freestream turbine if

where $U\u221e,n$ is the average freestream velocity over the disk for the *n*th turbine according to the inflow profile and $u\u221e,n$ is the inflow to the n^{th} turbine defined by Eq. (6) from the wake model, which takes into account the wakes of upstream turbines. The procedure is iterative and begins by treating only turbines that fall along the convex hull of the farm as freestream turbines. After a first iteration (see below), $u\u221e,n$ can be determined and additional freestream turbines [those complying with the condition in Eq. (19)] that may be located further inside the farm and yet be unaffected by upstream turbines can be identified.

For the waked turbines, we compute a turbine-specific wake expansion coefficient $kw,n$ in terms of the arithmetic mean of the high and low cell-specific friction velocities [found from Eqs. (16) and (14)], according to

The denominator ($u\xafh,ntd$) represents the advection velocity for which we use the average hub-height velocity in the cell from the top–down model.

For the freestream turbines, we use the same expression but since they are subjected to the incoming freestream advection and friction velocities, their expansion coefficient is evaluated according to Eq. (5), or simply

using the same *α* as for the entire wind farm.

Evaluation of the top–down model for each cell requires specification of the planform thrust coefficient *c _{ft,n}* in Eq. (12), which in turn determines the evolution of $\delta ibl(x)$ as well as the top–down mean velocity from Eq. (17). The planform thrust coefficient

*c*is defined as the total force per unit horizontal area. We argue that the relevant area affecting the turbulence at turbine

_{ft,n}*n*is the area of all the Voronoi cells upstream of turbine

*n*, as illustrated in the sketch in Fig. 5. This is the region over which the IBL evolves until it reaches turbine

*n*, and we therefore consider the vertical momentum flux averaged over this entire upstream area to characterize the evolution of the boundary layer reaching turbine

*n*.

The local planform area thrust coefficient is then given by

where $Apf,n$ is the planar area of the Voronoi cell for the *n*th turbine, and *line* refers to the cells that would be crossed if a line were drawn from the current cell to the front of the farm. This area is illustrated in Fig. 5, where the line is denoted in red and the cells included in the calculation are shaded. Moreover, $ud,i$ is the wake-model velocity field $u(x,y,z)$ [obtained from Eq. (1)] averaged over the rotor disk area *A _{d}* for turbine

*i*,

The coefficient $cft,n$ can be thought of as the ratio of the total force applied by the turbines in the cells upstream, including turbine *n* that affects the local boundary layer development, and the total horizontal momentum flux associated with the horizontally averaged mean velocity over the same planform area. The addition of the cells in front of the current cell is meant to represent the effect of the developing boundary layer over the farm and links the downstream cells to those upstream, since the turbines in the upstream cells influence those in their wakes. The coefficient $cft,n$ is then used to determine a local roughness height $z0,hi,n$ that is then used in the evaluation of $\delta ibl(xn)$ and $u*,hi,n$ according to Eqs. (15) and (16).

Figure 6 summarizes the coupling in the ALC model using an example of the planar average velocity view from each model. Using the developing boundary layer approach for the entire farm enables the model to find one *α* value that includes input from all the waked turbine cells. The top–down model provides the information required to determine the wake expansion coefficient in the wake model, while the wake model provides the planform thrust coefficient to be used in the top–down model.

After the determination of *α*, the power at each turbine is found using the velocity field produced by the wake model with the optimized wake expansion coefficients from the coupling. In order to compute the power from turbine *n*, we use

where *ρ* is the density of the air, $CP\u2032$ is the local power coefficient, and $ud,n$ is the average disk velocity of turbine n calculated using Eq. (23).

From the model description, it is apparent that the ALC model provides significantly more information about the flow than conventional wake models. The ALC model calculates friction velocities for each of the cells, which gives local flow information regarding turbulence conditions [assuming that turbulent root mean square (RMS) velocities are proportional to the friction velocity]. This local information can be used in combination with operating conditions to predict other quantities such as the unsteady loading characteristics for each individual turbine. Additionally, each cell is associated with an individual wake expansion coefficient even though only one *α* is selected for the entire farm. Having individual *k _{w}* values for wakes emanating from each turbine can represent, for example, situations in which the expansion rate for turbines close to the entrance is lower than the wake expansion rate further downstream, where turbulence levels are enhanced.

Sections III and IV demonstrate the application of the ALC model to two wind farms that are not arranged in a rectangular lattice pattern. In Sec. III, we consider a circular farm, while in Sec. IV we examine a wind farm that includes random turbine spacings to demonstrate the ability of the model to capture a range of farm geometries.

## III. MODEL VALIDATION FOR CIRCULAR WIND FARM

In this section, we compare the ALC model predictions to LES results generated using the National Renewable Energy Laboratory (NREL) SOWFA code^{27,28} for the circular wind farm configuration shown in Fig. 7. Comparisons include both total wind farm power and power from individual wind turbines for 12 wind inflow directions, 30° apart, spanning a full 360°. Figure 7 indicates the orientation of the wind directions; the $0\xb0$ direction correspond to a wind inflow direction from the north. Angles are then measured in the clockwise direction (i.e., at $90\xb0$ the flow is going from right to left, and for 270° the flow is going left to right). Given this orientation, the wind farm is symmetric around the east–west axis but has more complicated wake interactions than a rectangular lattice as the streamwise and spanwise spacings between turbines are not uniform for any given wind direction.

### A. LES and ALC model setup

The LES of the farm comprises a set of 38 NREL 5-MW reference turbines^{29} represented by actuator disk models (ADM). Velocity wake profiles from LES using ADM have been shown to be in good agreement with those using actuator line models (ALM).^{30} Simulations employing ADM can be accomplished using coarser LES grid-spacings as compared to those using ALM, and thus the choice of ADM enables us to run many more cases (inflow angles) and run simulations for longer times (aiming for statistical convergence). The simulation domain comprises ($Nx\xd7Ny\xd7Nz=500\xd7500\xd7100$) grid-points. The inflow has a mean velocity of 8 m/s and was generated by a precursor simulation. The simulation uses a roughness height $z0=0.15$ m. The details of this simulation are available in the work of Thomas *et al.*^{31}

Figure 8 shows the average power as a function of wind direction obtained from LES at 30° intervals beginning at 10°, denoted by the blue markers. The model data are also on the figure, represented by the red line, and will be discussed in more detail later. The LES data points were found by averaging the total power output of the LES farm over a time interval spanning approximately 0.75 h of real time farm operation, which represents approximately five flow-through times through the farm.

To provide an accurate comparison, the parameters used to evaluate the ALC model are selected to closely match those of the LES. We therefore use the values for the NREL 5-MW reference turbine,^{29} which have a hub height of *z _{h}* = 90 m and a diameter of

*D*= 126 m. The inflow velocity distribution $U\u221e(y)$ is found using average inflow velocity data from the LES. All of the LES have very similar inflow profiles, so an average over the 12 simulations was used as the ALC inflow profile for all wind directions. For the lower roughness height in the model, we use the roughness value from the LES: $z0,lo=0.15$ m. The maximum boundary layer height was set to

*δ*= 750 m since the LES had a temperature inversion layer at this height, which caps the growth of the boundary layer.

In addition to the turbine parameters and flow conditions, we match the local thrust and power coefficients ($CT\u2032$ and $CP\u2032$, respectively) of the turbines in the LES for the ALC model evaluation. In general, the thrust and power coefficients are determined by the turbine type and its operating conditions. In our case, we have the power produced as well as the velocity at the disk for the wind directions $10\xb0\u2212190\xb0$ (in 30° increments: seven directions) for each turbine from the LES. We extract the coefficients by rearranging Eq. (24),

Using this expression, we average over the time series to find an average local coefficient of power for each turbine. Figure 9 shows the local coefficients of power calculated for all seven directions, plotted by turbine number. The red triangles represent unwaked (or freestream) turbines in that orientation and the blue circles represent the waked turbines. There is significant spread in the local coefficients of power for the turbines. However, all of the freestream turbines are located at the bottom of the plot, around $CP\u2032\u22481.4$. Since all of the spread occurred in the waked turbines, where turbulent fluctuations may skew calculations of the local power coefficient, the waked turbines were judged to be a less accurate measure for this quantity. Therefore, we used the average value from the unwaked turbines, $CP\u2032=1.387$, for the local power coefficient in the ALC model.

Once the local power coefficient was calculated, blade element momentum theory is used to find the relationship between $CP\u2032$ and $CT\u2032$,

where the prefactor $\gamma =0.9032$ was determined using the specifications of the NREL 5 MW turbine [we are grateful to Dr. Carl Shapiro (personal communication) for performing this calculation]. In the ALC model, $CT\u2032$ is used to calculate the planform thrust coefficient and the local forcing of velocity deficit in the wake model, while $CP\u2032$ is used to calculate the power from the velocity computed at the disk from Eq. (24).

### B. Results

The ALC model was run using the parameters outlined above for the cases where the wind approaches from $0\xb0\u2013360\xb0$ in increments of 5°. Figure 8 compares the ALC model power output predictions (red line) to the LES results (blue dots). The farm is symmetric about the east–west axis, but due to the nonuniform inflow used in the simulation, the model results are not quite symmetric. We can see that in most cases the ALC model results match those of the LES, i.e., they fall reasonably close to the blue symbols. A noticeable outlier is for the inflow direction of 10°. In this case, the model over-predicts the power; however, the results for a few degrees to the side provide a better match. This highlights how small changes in angle can significantly impact the power output of the wind farm. Such changes are especially apparent, for example, over the range 90°–110°.

Figure 10 shows the free parameter *α* obtained from the error minimization over the waked portion of the farm, as a function of wind direction. The fact that the results are on the order of unity (around 2) confirms the validity of the scaling assumptions underlying Eq. (20). However, a non-trivial dependence on angle can be observed, which is a result of the complicated relationships involved in the model, the turbine geometric layout, and the transverse variations in inflow velocity.

We next examine the post-optimization hub-height velocity fields obtained from the two constituent models of the ALC model, and compare predicted turbine power predictions to the LES data, the Jensen model, and the model in Shapiro *et al.*^{11} Figure 11 shows (a) the average velocity field from the LES and (b) the velocity field calculated by the ALC model when the incoming wind direction is 190°. We can see that allowing nonuniform inflow enables the ALC model to reproduce the variations of higher and lower velocities seen across the LES, which continue back through the wind farm. The ALC model also captures wake superpositions, as can be seen by examining the wakes at the back of the farm. The LES shows some flow acceleration in between the wakes, which the ALC model does not capture.

Continuing the comparison, the top row of Fig. 12 compares the planar averaged velocities in each Voronoi cell over the farm, calculated using the method shown earlier in Fig. 4 for (a) the LES, (b) the resulting ALC wake model, and (c) the resulting ALC top–down model velocity fields for the same 190° wind inlet direction. We can see that the wake model, which provides the output of the ALC model, and the LES average velocity fields compare well with each other. The freestream turbines agree well, which is a result of the use of a nonuniform inflow velocity profile. Further back in the wind farm, the LES average velocity plot has a pattern of cells where some of the cells have a slower average velocity and thus more wake effects. The wake model is able to capture the overall patterns of cell velocities well although there are some cells with lower velocity in the LES. Note that the wake model and top–down model average cell velocities do not agree exactly since we minimize the difference over the entire farm to find a global *α* value rather than computing individual values of *α* on a cell-by-cell basis. Using only one variable in the minimization rather than cell-specific *α* values enables the models to inform each other with significantly reduced the computation time over a multiple *α* approach.

The differences between the LES and the wake and top–down portion of the models can be quantified by computing the root mean squared (RMS) error between the cell velocity values from the LES and from the model. The RMS error was calculated according to

where $uh,nLES$ denotes the average cell velocity from LES and $u\xafh,nm$ represents the average cell velocity from the ALC model [in the figure, (b) wake model and (c) top–down model]. These values are reported in the caption of Fig. 12.

Next, we compare the ALC approach to some prior analytical models. In the second row of Fig. 12, we compare scatter plots of the power produced by each turbine in the farm computed with the Jensen model, the model in Shapiro *et al.*,^{11} and the ALC model. Here, the LES power is on the x-axis and the model power is on the y-axis, so points along 45° line (marked as a black dashed line) represent $1:1$ agreement. The red triangles denote the freestream (unwaked) turbines, while the blue circles are the waked turbines. The left plot (d) shows the scatterplot using the Jensen model with a wake expansion coefficient derived from the average friction velocity of the inflow, which was *k* = 0.0653. The middle plot (e) uses the model in Shapiro *et al.*^{11} with Voronoi cells and a uniform inflow. In this case, we can see that the uniform inflow prevents the freestream turbines from capturing the behavior of the LES, causing all of the unwaked turbines to give almost the same value. In the right plot (f), we show results from the ALC model, where the freestream turbines are much closer to the 45° line as a result of the nonuniform inflow profile used. The ALC model predictions for the waked turbines are also closer to the LES results. The improvement in this direction is especially prominent in turbines with higher power, whose prediction is aided by the inclusion of the nonuniform inflow velocity profile. We also calculated the RMS error between each model and the individual power from LES, evaluated according Eq. (27) for power rather than velocity, these values are given in the caption of Fig. 12.

Figure 13 provides a detailed comparison of the LES field to the wake and top–down models comprising the ALC model for four different incoming wind directions: $70\xb0,\u2009130\xb0,\u2009190\xb0$, and $250\xb0$. Each direction corresponds to a row in the figure, while each column depicts a different quantity for that direction. The first column shows the time averaged LES velocity field with the Voronoi cells drawn in to show the area over which the planar velocity values are calculated. The first column also has arrows to indicate the direction of the incoming wind. The second column shows the cell average LES planar velocities calculated from the time averaged LES velocity fields. The third column shows the average planar velocities from the wake model of the ALC model, and the last column shows the same quantity from the top–down model. In all the cases, the wake model is able to predict the pattern of faster and slower velocity in the cells in the waked region, the back half of the wind farm, of the LES farm quite well. The pattern of the cells is important in seeing how the turbines and wakes interact with each other. The accurate prediction is consistent with the power output comparison in Fig. 8, where the ALC model matches the LES for the 130° and 250° directions quite well. However, in the 190° direction, the model over-predicts the power, which is clear from the slightly higher velocities predicted from both the top–down and wake models for that wind direction.

We can further examine wake interactions between turbines by examining the wake expansion coefficient determined as part of the ALC model. Even though only one parameter is matched across the whole farm (*α*), the model gives varying wake expansion coefficients across the farm because each cell has a different freestream velocity and friction velocity $u*,hi,n$, with the freestream friction velocity $u*,n$ used for the freestream turbines. Figure 14 shows how the wake expansion coefficient differs across the farm for the same four wind directions shown in Fig. 13, this time with the wind direction changing across the row. A lighter value indicates a higher wake expansion coefficient, and thus a wider wake, more expansion and velocity recovery, and more turbulence in the flow. A darker value indicates a lower value and a stronger wake effect and velocity deficit on the following turbines. Here, we can see that the wake expansion coefficient is higher overall for the first direction (a) 70°, but is lower for the later three directions:(b) 130°, (c) 190°, and (d) 250°. The significantly higher wake expansion coefficient for the case of 70° for turbines inside the wind farm can be attributed to the fact that for this wind direction, several inside turbines are freestream since they do not fall in wakes. The freestream turbines have smaller expansion coefficients in their wakes, i.e., the wakes dissipate more slowly initially. However, once these relatively stronger wakes interact with downstream turbines inside the farm, the large friction velocity and low mean velocity associated with the cells with strong wakes translate into large expansion coefficients for the turbines further downstream. These interactions explain the suddenly large *k _{w}* values for the 70° case.

In all cases, the freestream turbines have a significantly lower wake expansion coefficient than the waked turbines. To capture the larger velocity deficit that occurs behind freestream turbines, we use the freestream friction velocity for the freestream wake expansion coefficient since they are the first to interact with the incoming flow. The waked turbines wake expansion coefficients are determined based on Eq. (20), which depends on the local friction velocity and inflow in each cell. The variations in these values are significantly smaller than the difference between the unwaked and waked wake expansion coefficients. Another visible trend is that the wake expansion coefficient is typically the highest in the back center of the farm, as referenced from the freestream wind direction. This trend reproduces the realistic effect that an increase in turbulence further back in the farm is due to cumulative wake interactions.

In addition to the planar velocity comparisons presented above, we can also examine a more quantitative view of the model performance using scatter plots of the individual turbine velocities for each wind direction. Figure 15 shows the scatter plots for all 12 directions with LES data. As in Fig. 12, the freestream turbines are represented by the red triangles, the waked turbines are represented by blue circles and the $1:1$ correspondence line with LES data are represented by the black dashed lines. As can be expected, the ALC model provides relatively accurate power predictions for the freestream turbines compared to the LES results, while the results for the waked turbines have more spread overall. Some of the discrepancies can be explained by the spread in the local power coefficients calculated from the waked turbines, which is evident in Fig. 9, shown earlier. Since we chose the power coefficient based on the freestream turbines, the waked turbines could deviate from this constant value in the LES somewhat, possibly due to wake interactions, causing more spread in the results. Overall, however, the ALC model reproduces the LES trends quite well for any of the angles considered, with an overall RMS error of 0.197 MW.

## IV. MODEL VALIDATION USING A MIXED REGULAR AND NONUNIFORM WIND FARM

In this section, we compare the model to LES of a wind farm layout that mixes a regular array region with a random region. More specifically, the wind farm starts with an array consisting of four staggered rows, and then has a total of fourteen additional turbines placed in a random fashion behind the staggered-array turbines. The LES is performed with the JHU LESGO code^{32} with an actuator disk model using a correction factor accounting for finite grid resolution.^{33} The code has been validated on several wind energy applications.^{34–36} As the inflow, the simulation used a field generated by the concurrent-precursor approach.^{37} The turbines have a diameter and hub height of $D=zh=100$ m, and a ground surface roughness height of $z0,lo=0.1$ m is used to prescribe the bottom boundary condition of the LES. The main simulation domain uses $(Nx\xd7Ny\xd7NZ=384\xd7256\xd7128$) grid-points. The Lagrangian scale dependent dynamic subgrid-scale model^{38} is used to determine the eddy-viscosity without adjustable parameters. The turbine (local) coefficient of thrust is kept constant at $CT\u2032=1.33$ throughout the simulation. The same parameters are used in the ALC model. The LES is run for a single inflow direction of 270° for a time period of approximately ten flow-through times for the farm, which translates to roughly 1.75 h in real time. The data also provide a well-converged averaged inflow condition for the ALC model application.

Figure 16 shows the comparison of the time averaged streamwise velocity field at hub-height given by (a) the LES and (b) the ALC model. The inflow used in the ALC model was taken from the time-average LES velocity field at *x* = 0 m in front of the turbines. The dimensionless wake expansion parameter found through the error minimization was $\alpha =1.68$, which is slightly lower than the values found in the circular wind farm study, but still of order unity. We can see that the ALC model captures the flow well, and represents the wakes, particularly in the rear of the farm, quite accurately. The wakes in the ALC model appear sharper, and the velocity is slightly lower in the wake as compared to the LES results, but their relative speed of decay is similar in both cases. The comparison also shows the importance of the nonuniform inflow condition in accurately modeling the overall flow especially near the inlet. The lower panels of Fig. 16 show the average velocities in each of the Voronoi cells obtained from (c) the time and cell-averaged field from the LES and the (d) the mean velocity obtained from the ALC model (the wake model part). We can see that the cell velocities match well, showing the same general trends in both cases.

We can also analyze the wake expansion coefficients calculated for the hybrid wind farm, shown in Fig. 17. In this farm, the wake expansion coefficients start out the lowest at the front of the farm, then increase to their highest values in the rows immediately after the first rows. The values then gradually decrease from this point as a function of downstream distance, which follows trends expected for the friction velocity development in a smooth-to-rough transition. The boundary layer flow downstream of the transition, in the third and fourth rows at around 1–1.5 km in the streamwise direction, displays larger mean shear near the surface (here at the turbine height) and thus larger turbulence levels and wake expansion coefficients. Further downstream, as the boundary layer approaches equilibrium conditions where the turbines are randomly placed, the shear decreases and so do the turbulence levels causing a slight decrease in the wake expansion coefficients.

Finally, as a more quantitative comparison between the ALC model results and the LES results, Fig. 18 shows a scatterplot comparing the power output from each turbine. As before, the freestream turbines (first two rows at *x* = 0.5 and 0.85 km) are denoted by the red triangles, while the waked turbines further downstream are denoted by the blue circles. We can see that the points follow the $1:1$ correspondence line rather closely, implying very good agreement of the ALC predictions with the LES results. There is substantially less scatter compared to the circular wind farm case considered in Sec. III, with this case having an RMS error of 0.062 MW. This could be due to the longer time period used for averaging which for the mixed regular-nonuniform wind farm was 1.75 h (approximately 11 wind farm flow through times), while it was only 45 min (approximately five wind farm flow through times) for the circular wind farm. We recall that the ALC model is designed to predict the mean flow, which requires significant time averaging in LES. In summary, we can conclude that the ALC model is able to capture the mean flow and power production of the mixed regular-nonuniform wind farm quite well.

## V. CONCLUSIONS

The area localized coupled (ALC) model combines a wake model and a top–down model in a localized fashion to create a more generally applicable approach that predicts several useful quantities about wind farms. The model uses Voronoi cells to divide the wind farm up into areas that belong to each turbine. Since the calculations can be applied to each cell, the model can be applied to regular as well as irregularly arranged wind farms. The wake expansion coefficient is obtained as a function of the local turbulence properties as described by the ratio of the friction velocity to the mean advection velocity. The dependence of the wake expansion coefficient on the local friction velocity determined from the growth of a local internal boundary layer enables a localized result involving different position-dependent wake expansion coefficients across the wind farm. A global parameter *α* is selected to minimize the difference in the cell-averaged mean velocity predictions between the top–down and wake models. In addition, the model can be implemented using a nonuniform inflow for both the wake and top–down model constitutive parts. This feature improves the performance of the model for all of the turbines, specifically for turbines near the inflow that are directly exposed to the position-dependent inflow velocity.

The model was validated using LES data for two different wind farms. The first was a circular wind farm. In this case, the model was able to accurately predict the cells of higher and lower velocity across the wind farm when compared to the time average LES data. It also reproduced the velocity patterns in the cells further back in the farm, where the wake interactions are most important. The ALC model also captured general trends in the variation of total wind farm power with wind direction as obtained from LES data, and predicted the power of individual turbines well.

The ALC model was also applied to a mixed regular-nonuniform wind farm for which LES data were also generated. The wind farm started with regular staggered rows and then continued with a set of randomly placed turbines downstream. The ALC model predicted velocity field matched the time-average velocity from the LES very well. The model was able to predict the power of individual turbines quite accurately, as well.

The newly proposed ALC model provides more information than typical wake models since physically relevant information about position-dependent friction velocity, internal boundary layer, and roughness height becomes available. Since the Voronoi tessellation is generally applicable to any turbine arrangement, there are no limitations to apply the ALC to any wind farm geometry, regular, random or mixed. It can account for position-dependent inflow velocity distributions. Additionally, the model has the potential to be extended to dynamic, time dependent situations if, instead of using the steady state solution of the partial differential equation for the velocity deficit, one uses the full time-dependent solution. It should be noted that the top–down model assumes a quasi-steady description of the boundary layer and so further improvements to include time-dependence in the top–down model may still be required. At any rate, such possible extensions make the ALC model an attractive and versatile option to be used in wind farm design and control. Future work will include the validation of this model on operational wind farms and possible extensions of this framework to a dynamic setting.

## ACKNOWLEDGMENTS

G.M.S., C.M., and D.F.G. gratefully acknowledge partial funding support from the National Science Foundation (Grant Nos. DGE-1746891, CBET-1949778, and CMMI-1635430) and the Maryland Advanced Research Computing Center for computing resources.

This work was co-authored by the National Renewable Energy Laboratory, operated by Alliance for Sustainable Energy, LLC, for the U.S. Department of Energy (DOE) under Contract No. DE-AC36–08GO28308. Funding provided by the U.S. Department of Energy Office of Energy Efficiency and Renewable Energy Wind Energy Technologies Office. The views expressed in the article do not necessarily represent the views of the DOE or the U.S. Government. The U.S. Government retains and the publisher, by accepting the article for publication, acknowledges that the U.S. Government retains a nonexclusive, paid-up, irrevocable, worldwide license to publish or reproduce the published form of this work, or allow others to do so, for U.S. Government purposes.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding author upon reasonable request.