Increasing the Power Production of VerticalAxis WindTurbine Farms Using Synergistic Clustering
Abstract
Verticalaxis wind turbines (VAWTs) are being reconsidered as a complementary technology to the more widely used horizontalaxis wind turbines (HAWTs) due to their unique suitability for offshore deployments. In addition, field experiments have confirmed that verticalaxis wind turbines can interact synergistically to enhance the total power production when placed in close proximity. Here, we use an actuator line model in a largeeddy simulation to test novel VAWT farm configurations that exploit these synergistic interactions. We first design clusters with three turbines each that preserve the omnidirectionality of verticalaxis wind turbines, and optimize the distance between the clustered turbines. We then configure farms based on clusters, rather than individual turbines. The simulations confirm that verticalaxis wind turbines have a positive influence on each other when packed in welldesigned clusters: such configurations increase the power generation of a single turbine by about 10 percent. In addition, the cluster designs allow for closer turbine spacing resulting in about three times the number of turbines for a given land area compared to conventional configurations. Therefore, both the turbine and windfarm efficiencies are improved, leading to a significant increase in the density of power production per unit land area.
Introduction
Despite the concerted effort to improve energy efficiency and decouple economic growth from energy consumption, the U.S. Energy Information Administration projects that global total energy consumption will grow by about 45% between 2015 and 2040 (U.S. Energy Information Administration 2013). Mitigating the concomitant large increase in greenhouse gas emissions necessitates exploring alternative loweremission energy sources, particularly since the majority of the current fossilbased energy resources are finite and have other adverse side effects on the environment. Wind energy is expected to be one of the primary sources of clean, renewable energy that would allow a rapid transition away from fossilfuelbased energy. In the USA, for example, wind power is projected to contribute around 20% of electrical energy by the year 2030 (Marquis et al. 2011). As a result, increasingly larger wind farms are being deployed, and the continued spread and expansion of these farms poses a challenge since the required land area will increase. A major goal of current research is thus to increase the windfarm power density, i.e. how much energy can be produced per unit land area used.
In a wind farm, turbines should be far enough apart to allow wind speeds to recover, through lateral or vertical momentum entrainment, after deceleration by the upwind generator (Cortina et al. 2016). Spacing the turbines also reduces the fatigue load generated by turbulence from the upstream turbines and thus increases turbine lifetime (Chamorro and PortéAgel 2009). The large majority of existing farms use horizontalaxis wind turbines (HAWTs); the behaviour of horizontalaxis wind turbines in large wind farms, and the required spacing between them, have been extensively studied (Wu and PortéAgel 2012; Troldborg and Sørensen 2014). Calaf et al. (2010) investigated the vertical transport of momentum and kinetic energy in a fullydeveloped HAWTarray boundary layer (defined as the internal boundary layer developing above a wind farm). They showed that, for large wind farms, regeneration of the kinetic energy is mainly from downward vertical fluxes across the plane delineating the top of the farm, unlike farms with a limited number of windturbine rows where the streamwise advection of kinetic energy dominates. The concept of a windturbinearray boundary layer is particularly useful for wind farms where streamwise farm length is an order of magnitude larger than the depth of the atmospheric boundary layer (ABL) since the influence of such farms extends to the top of the ABL. Meyers and Meneveau (2010) used an actuatordisk model and largeeddy simulations (LES) to model large HAWT wind farms and investigate their interaction with the ABL. They have shown that a staggered wind farm can extract 5% more power than an aligned configuration and in a followup study (Meyers and Meneveau 2012), investigated the optimization of turbine spacing in fullydeveloped wind farms. They showed that the ratio of land cost to turbine cost in the financial optimization analysis (maximizing power per unit cost) influences the deduced optimal spacing. Meyers and Meneveau (2012) indicate that the optimal turbine spacing is higher than that currently being used in HAWT wind farms. Recently, it has also been shown that the highest achievable mean windfarm power is strongly dependent on the alignment of the turbine arrays relative to the mean wind direction, and the optimal alignment angle is significantly smaller than that in a perfectlystaggered farm (Stevens et al. 2014). For windfarm sites with a dominant wind direction, these findings can be implemented to improve windfarm performance.
All of the above, and other previous work, have focused on wind farms consisting of horizontalaxis wind turbines (Chamorro and PortéAgel 2010; Lu and PortéAgel 2011; Yuting 2011; Meyers and Meneveau 2012). However, recently Dabiri (2011) has suggested the possibility of an order of magnitude increase in power densities for wind farms when verticalaxis wind turbines (VAWTs) are used. Due to their axis of rotation, VAWT wakes and the flow in a VAWT farm are distinctly different from their HAWT counterparts. This potential increase in power density can be achieved by configuring VAWT farms with a closer spacing to better exploit the flow patterns created by upstream turbines. Dabiri (2011) and collaborators (Kinzel et al. 2012) performed experiments on various counterrotating configurations of 9m tall verticalaxis wind turbines and demonstrated that, unlike the typical performance reduction of horizontalaxis wind turbines with close spacing, there is an increase in VAWT performance when adjacent turbines are arranged to interact synergistically. However, high experimental costs and time requirements prevent the extension of these field investigations to large farm scales or the assessment of a large number of configurations. The previous findings thus only pertain to a limited number of turbines where the mean kinetic energy is primarily replenished by streamwise advection and crossstream turbulent transport, rather than by vertical transport as in large farms. Our aim here is to bridge this research gap and assess the feasibility of increasing power density in large VAWT farms using a synergistic clustering of turbines. Building on Hezaveh et al. (2016), where a LES model for verticalaxis wind turbines was extensively validated and the flow recovery in the wake of a single turbine investigated, here we simulate the interactions of multiple verticalaxis wind turbines in small clusters, and subsequently use these clusters to design large VAWT farms.
Numerical Model
In order to investigate verticalaxis wind turbines in the ABL, we used the LES model with a VAWT actuatorline model (ALM–LES) presented and validated in Hezaveh et al. (2016), and present here a brief overview only. In this model, which is a variant of a LES model that has been previously used and validated for flow around horizontalaxis wind turbines (Chamorro and PortéAgel 2009; Calaf et al. 2010, 2014; Lu and PortéAgel 2011) and other complex flows (BouZeid et al. 2007; Huang and BouZeid 2013; Li et al. 2016), the following continuity and Navier–Stokes equations are solved at each timestep for the large resolved scales assuming an incompressible flow with a mean in vertical hydrostatic equilibrium
$$ frac{{partial tilde{u}_{i} }}{{partial x_{i} }} = 0, $$
(1)
$$ frac{{partial tilde{u}_{i} }}{partial t} + tilde{u}_{j} left( {frac{{partial tilde{u}_{i} }}{{partial x_{j} }} – frac{{partial tilde{u}_{j} }}{{partial x_{i} }}} right) = – frac{1}{rho }frac{{partial tilde{p}^{*} }}{{partial x_{i} }} – frac{{partial tau_{ij} }}{{partial x_{j} }} + F_{i} + F_{i}^{t} , $$
(2)
where ( tilde{u}_{i} ) is the resolved velocity vector with the tilde denoting a filtered quantity, (u, v, w) are its streamwise, crossstream and vertical components, respectively. This instantaneous velocity is decomposed into a mean Ui and a resolved perturbation u′i; xi is the position vector with components (x, y, z) in the streamwise, crossstream and vertical directions respectively, ( tilde{p}^{*} ) is a modified pressure that includes the resolved and subgridscale turbulent kinetic energies, ρ is the air density, Fi is the mean pressure gradient driving the flow, τij is the deviatoric subgridscale stress tensor; and ( F_{i}^{t} ) represents the aerodynamic forces of the turbine blades on the airflow. Note the omission of the Coriolis force, which is assumed to have no significant impact at such small distances (about 10 m) from the Earth’s surface. At each timestep, ( F_{i}^{t} ) is computed using the actuatorline model as detailed in Hezaveh et al. (2016). The horizontal boundary conditions are numerically periodic, but nonperiodic flows can be simulated using an inlet sponge region as shown later. At the top boundary, zero vertical velocity and zero shear stress are imposed. The bottom boundary has zero vertical velocity, while the surface shear stress is imposed using an equilibrium logarithmiclaw wall model with a wall roughness length z0 = 10−6zi, where zi is the depth of the computational domain used to normalize all length scales in the model (zi = 25 m). The details of the wall and subgridscale models are provided in BouZeid et al. (2005). The model details are summarized in Fig. 1: an angle of attack (α, the angle between the blade chord and the flow velocity relative to the blade ( overrightarrow {V}_{text{rel}} )) is first computed by knowing the location of each blade represented as a vertical line in the actuatorline model, the upstream undisturbed flow velocity (( overrightarrow {U}_{infty } )), and the rotational speed of the turbine (ω). This then allows us to obtain the lift and drag force coefficients, CL and CD respectively, to be calculated from experimental data, bladeresolving Reynoldsaveraged simulations, or tabulated airfoil data after applying a dynamic stall correction, as in Hezaveh et al. (2016). CL and CD are then used to compute the normal and tangential force coefficients, CN and CT respectively,
$$ C_{{N}} = left {C_{{L}} } rightcos alpha + left {C_{{D}} } right sin alpha , $$
(3)
$$ C_{{T}} = left {C_{{L}} } rightsin alpha – left {C_{{D}} } right cos alpha , $$
(4)
which are then used to compute the corresponding forces
$$ {text{d}}F_{N} left( theta right) = frac{1}{2}rho c V_{text{rel}}^{2} C_{{N}} {text{d}}z, $$
(5)
$$ {text{d}}F_{T} left( theta right) = frac{1}{2}rho ;c V_{text{rel}}^{2} C_{{rm T}} {text{d}}z, $$
(6)
where c is the blade chord length, θ is the azimuthal angle of the blade at a given time, and dz is the vertical length of blade occupying each LES grid cell. Finally, the forces in the LES Cartesian coordinate system are obtained as,
$$ {text{d}}F_{x} = {text{d}}F_{N} left( theta right)cos theta + {text{d}}F_{T} left( theta right)sin theta , $$
(7)
$$ {text{d}}F_{y} = {text{d}}F_{N} left( theta right){ sin } theta – {text{d}}F_{T} left( theta right){ cos }theta . $$
(8)
In order to evaluate CL and CD herein, actual dynamic curves for lift and drag forces determined experimentally were implemented into the numerical model. Since we aim to mimic the turbines used in the field experiments for verticalaxis wind turbines reported by Dabiri (2011) and Kinzel et al. (2012), we simulate the same DU06W200 airfoil. Claessens (2006) has conducted several numerical and windtunnel experiments at various Reynolds numbers and tabulated the lift and drag curves for different airfoil blade types. Figure 2 depicts their measured values for the DU06W200 airfoil, which include the dynamicstall effect (directly accounted for in the windtunnel measurements). Based on the sign of dα/dt, the instantaneous rate of change of the angle of attack, different paths emerge in these curves in Fig. 2: CLm+ versus CLm− and CDm+ versus CDm−. We note that the DU06W200 blade type is redesigned from the NACA 0018 airfoil and is 2% thicker, and cambered rather than symmetric (Claessens 2006). This blade type and the experimental results reported in Fig. 2 are used throughout.
Adapted from Hezaveh et al. (2016)
Schematic twodimensional crosssection (top view) of the VAWT blade path, the forces on the blades, and representative LES grid cells. θ is the azimuthal angle denoting the angular location of the blade relative to the incoming flow direction (defined positive in the same direction as the turbine rotation); it continuously varies in time for each blade as the turbine rotates. The depicted relative scales of the blade chord length to the LES grid cell illustrate that they are comparable, but their exact ratio (might be larger or smaller than 1) varies in the different simulations and the figure is not to scale.
Dynamic CL and CD for the DU 06W200 blade type as measured by (Claessens 2006). The + subscript is for dα/dt > 0 and the − subscript is for dα/dt α/dt is the time rate of change of the angle of attack
Results and Discussions
As mentioned before, Dabiri (2011) and Kinzel et al. (2015) report on several field measurements with different configurations of 9m tall verticalaxis wind turbines using the same blade type and rotor configuration we adopt here. They have shown that, by using counter rotating verticalaxis wind turbines and special configurations, the turbines can exploit the flow deflection from upwind adjacent turbines and there is a potential of a one orderofmagnitude increase in power density. To complement the previous validation of this ALM–LES model performed against laboratory experiments (Hezaveh et al. 2016), and to ensure that the simulations accurately represent the flow in between multiple turbines and therefore within and in the wake of turbine clusters, we compare our LES results to the fieldmeasured data described in Kinzel et al. (2015) for two adjacent counterrotating turbines. This is the first validation of our model against data from realsized VAWT field measurements and, to the best of our knowledge, the first validation of any verticalaxis wind turbine ALM–LES against field data.
The details of the experimental setup are presented in Table 1, and the schematic configuration is shown in Fig. 3. The 1200W turbines are a modified version of a commercially available model from Windspire Energy Inc. (Dabiri 2011) and they were placed 1.6D apart (D is the rotor diameter). The velocity profiles were measured at 16 points with streamwise coordinates (relative to the line joining the centre of the turbines) of x = − 15D, − 1.5D, 2D, and 8D and elevations above ground of z = 3, 5, 7 and 9 m. All velocity components are normalized using a measured 10m wind speed from a meteorological tower in the vicinity of the experiments (Araya et al. 2014; Kinzel et al. 2015), which took place in Antelope Valley, north of Los Angeles, California. Further details about the measurements are provided in Appendix A—also see Dabiri (2011) for further information.
Table 1 Turbine characteristics from Dabiri (2011) and Kinzel et al. (2015)
Schematic of a the two turbines in the threedimensional flow domain, and b top view of the computational domain
The computational domain has Lx × Ly × Lz = 31.2 m × 15.6 m × 25 m, respectively spanned by 128 × 64 × 192 grid nodes. This resolution yields about 5 × 5 horizontal grid points spanning each turbine rotor (five in each direction). The distance between the domain inflow and turbines was set equal to the distance between the furthest upstream measurement point and the turbines in the experiment, that is 15D. In order to match the inflow conditions such as turbulence intensity and mean upstream wind speed profile in the LES to the observed field data, a precursor periodic simulation was run to generate the inflow. The roughness length and friction velocity of this precursor simulation were calibrated (with adopted values of 0.001 m and 0.5 m s−1, respectively) to yield the experimentallyobserved logarithmic velocity profile measured 15D upstream of the turbines. The inflow and validation simulations were conducted with neutral stability and, as detailed later, field experimental periods were selected during nearneutral conditions. Then, y–z slices of instantaneous velocity and pressure were saved at each timestep and fed to the simulation with the turbines as upwindinflow boundary conditions.
The results are shown in Fig. 4, and it is clear that the ALM–LES model is capable of closely reproducing the wake generated by the interactions of the two counterrotating verticalaxis wind turbines (the blades move towards the back when facing the other turbine so that the flow acceleration in between the two rotors is maximized). We should emphasize that it is essential to provide the LES with accurate inflow (left panel of Fig. 4, from the precursor simulation) for the experimental profiles near and behind the turbines (right three panels of Fig. 4) to be reproduced accurately. These results confirm that the ALM–LES produces realistic wakes even where turbines are interacting, and hence the model can be used to investigate large wind farms and VAWT clusters with confidence. It should also be mentioned that the ALM–LES model is capable of realistically capturing wake meandering, but this meandering does not appear in the figures herein since we only show mean velocities. Furthermore, the omission of the Coriolis force does not influence the result given the high Rossby number in the atmospheric surface layer at such low elevations. While the Coriolis force induces Ekman turning, for the omnidirectional verticalaxis wind turbines the effect on performance is smaller than for horizontalaxis wind turbines.
Profiles of incoming and wake wind speed for the ALM–LES model versus field experimental data
Clustering verticalaxis wind turbines in small arrangements have been shown to have several advantageous implications for power generation (Dabiri 2011). The global performance of the turbines is enhanced since the downstream turbines benefit from the flowdeflection effect and the resulting higher flow speed induced by upstream turbines. However, depending on the wind direction relative to the alignment of the turbine arrays in the farm, compact clustering might also have negative effects when one turbine is mainly in the wake/shadow of an upstream rotor. For example, if two turbines are clustered together, the range of wind directions for which one of the turbines is in the shadow (partially or fully) of the other is 2β, where β = tan−1(2D/2L) (Fig. 5, left), L being the turbine spacing (centre to centre) in a cluster. We note that this is a purely geometric consideration that does not account for the expansion of the wake. On the other hand, when the flow is approximately perpendicular to the centretocentre axis, the higher induced speed inbetween the two turbines is not being exploited.
Wind directions in which verticalaxis wind turbines are in the wake of an upstream rotor for two, three and four turbines (γ ≈ β)
By introducing one additional turbine, the range of wind directions where two turbines can directly shadow each other is increased to 6β (Fig. 5, middle). However, the third turbine can benefit from the higher wind speed induced inbetween the two upstream turbines or the two downstream rotors can benefit from the transverse flow deflection of the upstream turbine (depending on wind direction). This has the potential to result in power production from these three turbines that is greater than the power from three distant noninteracting ones (this improvement depends on L/D, as shown below). By increasing the number of turbines in the cluster beyond three, the flowrelated benefits decrease and the range of wind directions where the turbines shadow each other increases to n (n − 1) β, where n is the number of turbines in the cluster (e.g., Fig. 5, right). In Fig. 6, the variation of β with L/D and n for various clusters is shown; physically, β represents the total range of wind directions where shadowing influences the turbines. A value of β = π implies, for example, that one turbine is partially or fully shadowed for 50% of possible wind directions, or alternatively that two turbines are shadowed for 25% of possible wind directions. As such, β is the cumulative (partial or full) shadowing of all turbines from all possible wind directions and it can therefore exceed 2π. By increasing the value of L/D of a cluster, β is reduced, while on the other hand, increasing n results in higher β. For n > 3, the β/2π value can become larger than 1, which indicates that there is no wind direction for which the turbines are not casting at least partial shadows on each other. However, one notes that, for L/D > 5, the differences between the β values for n = 2 and n = 3 are minor. Moreover, a clustering with higher n has the important benefit of using a smaller land area. Therefore, the most efficient design for a cluster when there is no dominant wind direction at the site seems to be a triangle (n = 3) since it has a limited β, while at the same time allowing for compact clustering and synergistic interaction between the turbines. A value of n = 4 almost doubles the shadowing angle β, with no increase in the winddirection range for which synergistic interactions occur. Therefore, we focus on triangular clusters hereafter.
Variation of β, the cumulative shadowing angle, with the L/D ratio and the number of turbines in a cluster n
In order to investigate the characteristics of the proposed triangular cluster design, we conduct a suite of largeeddy simulations in a computational domain containing three of the same turbines defined in Table 1. The basic domain size is Lx × Ly × Lz = 72 m × 48 m × 25 m = 60 × 40 × 20.8D and is spanned by Nx × Ny × Nz = 288 × 192 × 192 grid nodes. At such resolution (dx = dy = 0.25 m), the turbine diameter is covered by about 5 × 5 horizontal grid points, which is similar to the validation runs. These remain constant for the analyses used in this subsection (except for the domainsize sensitivity analysis detailed later). The power coefficient CP of a single isolated turbine, as modelled by the LES, is 0.36. Since the wake deficit increases with D and decreases with the distance between the turbines L, L/D emerges as an important dimensionless number to consider; that is, in addition to its impact on the cumulative shadowing angle β illustrated in Fig. 6 and discussed above. In order to investigate the optimal distance, various L/D ratios ranging from 2 to 8 were simulated.
Before conducting these simulations however, the computational setup needed verification; therefore, for a fixed value of L/D = 6, an analysis of the sensitivity of the results to the domain size was performed (such that the domain size and number of grid points increase proportionally, and thus the grid resolution is unchanged). Two parameters were investigated for sensitivity to the domain size: the clusteraveraged power coefficient CP and the wake velocity deficit at 15D and 20D downstream of the cluster. The wake velocity deficit is averaged in time and over a y–z rectangle that is aligned in the xdirection with the turbine crosssection projected area. As can be seen from Fig. 7, if one compares the 80 × 27, 80 × 40, and 80 × 54 runs, changes in domain width Ly can be significant when Ly becomes very small (27D). For such narrow domains, the crossflow area blocked by the turbines becomes large and prevents correct sideways deflection of the streamlines. Since our domain is periodic in y, a small Ly allows the clusters to interact with “virtual” adjacent ones. The figure suggests that a minimal Ly ≥40D should be used since increasing the transverse domain size beyond that, to Ly =54D, results in insignificant changes in the average CP or in wake recovery.
Domain size sensitivity analysis. The adopted size is 60D × 40D. CP = Turbine power/(0.5ρAU 3∞ ), where A is the rotor area
Changes in domain length Lx have little impact on the average CP (compare the 54 × 54, 60 × 54, and 80 × 54 runs). However, the velocity deficit is clearly influenced by Lx (due to downstream boundarycondition effects). Based on the sensitivity analysis results in Fig. 7, a minimal Lx = 40D was deemed necessary (and sufficient) to avoid an impact of the domain length on the average CP and the wake velocity deficit (compare the 60 × 40 to the 60 × 54, or the 80 × 40 to the 80 × 54 runs). Therefore, a domain size of Lx = 60D by Ly =40D is adopted for the singlecluster simulations hereafter. All these simulations were conducted using an imposed laminar logarithmicprofile inflow with a surface roughness length z0 = 0.001 m and a friction velocity = 0.5 m s−1. To assess the influence of the turbulence levels in the inflow, a simulation was conducted using inflow planes from a precursor periodic turbulent run. As can be seen in Fig. 7, using a turbulent inflow reduces the deficit values at 15D and 20D downstream of clusters significantly, which is expected since the increase in turbulence intensity increases momentum entrainment into the wake and speeds up its recovery.
With the basic domain size set, simulations with triangular clusters for different L/D ratios were conducted using, at first, the unique wind direction of 60° depicted in Fig. 7. Based on simulation results for different cases (see Fig. 8), it is obvious that increasing L/D improves the performance of the first (upwind) turbine due to the decrease in upstream blockage effect from turbines 2 and 3. Due to the rotation direction of the first turbine (shown in Fig. 7), which deflects the flow towards the third turbine, the third turbine has a slightly higher CP value compared to the second turbine. On the other hand, Fig. 8 illustrates that the performance of the second and third turbines first improves as L/D increases from 2 to 3, then reaches a plateau until L/D = 5, and finally decreases again. When L/D > 5, these turbines are less able to utilize the higher wind speed induced from the flow deflection by the upstream rotor. The clusteraveraged CP value (related to the upstream wind speed U∞) thus peaks at an intermediate L/D value. The three cases with the highest average CP value, corresponding to L/D values of 3, 4 and 5, were hence selected for further analysis.
CP for each turbine in the cluster, and the average for the whole cluster, as a function of L/D
These analyses consisted of simulations where all parameters remain the same for a given L/D, but with a different incoming flow orientation. We aim to investigate the omnidirectionality of the proposed VAWT clusters, as well as to find the most efficient VAWT spacing averaged over all wind directions. Figure 9a shows the average CP value versus incoming wind direction; the case with L/D = 5 has the highest CP averaged over all turbines for all wind directions. This is confirmed in Fig. 9b that depicts the influence of L/D on CP averaged over all wind directions and all turbines, and normalized by the CP value of a single isolated turbine. The setup with L/D = 5 has improved performance because the angle at which the verticalaxis wind turbines cast shadows on downstream turbines (β) is reduced due increasing distance between turbines, as well as because wake recovery is improved when shadowing occurs at the larger recovery distance. Finally, a key observation from Fig. 9b is that the average CP is about 10% higher than that for a single isolated turbine when L/D = 5; this confirms our premise that the synergistic interaction between closelyspaced turbines can indeed result in a higher overall power generation when exploited adequately.
a Triangular clusteraverage CP versus wind direction ζ. b Clusteraveraged CP, also averaged over all wind directions, and normalized by the CP of a single isolated turbine (angled brackets denote averaging)
With an effective cluster design selected, we now turn our attention to the design of farms based on these clusters. An important parameter in designing and optimizing wind farms is the distance required for wind speed and power recovery downstream of turbines (Hezaveh et al. 2016). This applies for arrays consisting of individual turbines, as well as of clusters (unless there is a dominant known wind direction, which is not an assumption made here). The windspeed deficit (1 − U(x, y, z)/U∞(z)) was averaged over the y–z planes encompassing the whole cluster (projected flownormal area) at varying x distances from the hub using data from the same simulations described in the previous subsection. We also investigated various values of L/D to confirm that our choice of L/D = 5, made based on the power output of an isolated cluster, does not produce longer wakes than for other L/D values. The results reported in Fig. 10 indicate that increasing the distance between the turbines in each triangle cluster significantly reduces the distance needed for the wind speed to recover to 75% of its upstream value U∞. It is clear from the figure that the recovery distance to 75% speed is reduced from 25D for L/D = 3 to 15D for L/D = 5. The choice of the 75% recovery speed is somewhat arbitrary, and other thresholds can be selected, of course. However, the comparative analysis of the recovery distances would reach the same conclusions regarding the optimal L/D to adopt, irrespective of the exact recovery threshold.
Comparison between averaged velocity deficit for various L/D and different wind directions ζ
Recovery is an important criterion for designing a wind farm, which further confirms the selection of L/D = 5. In a wind farm, it is important that downstream turbines are placed at distances were the available flow has recovered to significant levels of its upstream undisturbed speed (e.g. to over 75%, although higher levels are advantageous) so that the powergeneration capacity of these turbines is not underutilized. Furthermore, as shown in Fig. 10, by increasing the distance between the turbines in the clusters, the effect of incoming wind direction on the recovery distance is reduced. The L/D = 3 recovery is sensitive to the change in incoming wind direction ζ; the recovery distance to 75% of upstream flow speed occurs anywhere between 18D and 28D as the wind direction changes. The recovery for L/D = 5 on the other hand is much less sensitive to wind direction and thus yields more omnidirectional farms. The results also indicate that when designing farms based on clusters with L/D = 5, the velocity recovery for a separation of 10D between the clusters is about 70–75%, while a separation of 20D allows a recovery to over 80% of upstream flow speed. Both of these separations are tested in the fullfarm simulations below. While other separations can be examined, the outcome of the testing of our hypothesis regarding the potential benefits from synergistic interaction between verticalaxis wind turbines remains the same.
Now we tackle the main question: can synergistic interactions between verticalaxis wind turbines increase windfarm power density? Practically, we need to investigate whether farms with synergistic clusters have improved performance (produce more power per unit land or per unit invested cost) relative to two prototypical layouts of wind farm, aligned and staggered regular arrays. Based on the size of the selected turbine and the results obtained above, eleven farm configurations were simulated. One configuration is illustrated in Fig. 11, while four more layouts can be visualized in the streamwise velocity plots of Fig. 12; the turbines are the same as those detailed in previous sections. The simulations used are all periodic (representing an infinite farm), with Nx × Ny × Nz = 320 × 160 × 336 nodes and Lx × Ly × Lz = 96 m × 48 m × 32 m. The resolution yields 4 × 4 horizontal grid nodes per rotor diameter, which is comparable to the validation tests presented earlier. The vertical height of the simulation domain is selected based on a sensitivity analysis performed for the 10D staggeredspacing wind farm. Domains with vertical heights of 32, 45 and 54 m were chosen and the total power coefficients (we use two definitions, CP and C *P , that are described below) of these three domains are compared. Due to the small blockage ratio of the wind turbines (projected area of turbines normal to the flow over the y–z area of the computational domain), which is 2). Therefore, the 32m high domain was deemed sufficient and chosen accordingly; this agrees with Sarlak et al. (2016) concerning appropriate blockage ratios for windenergy studies.
Schematic of the windfarm configuration with the VAWT triangular clusters, with turbinetoturbine separation distances of L = 5D and intercluster spacing of 20D for an aligned cluster configuration
Streamwise velocity magnitude in wind farms with 10D horizontal spacings, mean flow from left to right: a regular aligned, b regular staggered, c cluster staggered with 0° wind direction, and d cluster staggered with 60° wind direction
Table 2 Vertical height of the computational domain versus farmaveraged power coefficient for a wind farm with 10D staggered configuration
Six of the simulated cases are prototypical wind farms with staggered or regular array configurations and with separation distances of 5D, 10D or 20D; these configurations result in 128, 32 and 8 turbines in the computation domain, respectively. Four additional experiments were conducted using staggered clusters (with L/D = 5 separation within each cluster) with wind direction of 0° or 60° and with intercluster distances of 20D or 10D for each direction, corresponding to 24 and 96 turbines in the domain, respectively. It should be mentioned that since the LES model uses periodic boundary conditions in both x and y directions, these simulations correspond to infinite wind farms with an infinite number of turbines. Therefore, the number of turbines in the computational domain does not influence the results when normalized per turbine. Finally, one experiment is conducted using 20D spacing and aligned clusters with a 60° wind direction relative to the cluster. To visualize the differences in the flow patterns in these designs, the average streamwise velocity component for a few selected configurations is shown in Fig. 12. The lower values in the staggered configurations reflects the higher power extraction resulting from a larger number of turbines, and therefore metrics that allow consistent comparison of these configurations are needed. Finding such metrics is not straightforward as discussed in Meyers and Meneveau (2010) and Goit and Meyers (2015).
The most direct metric and the easiest to compute is the average windfarm CP value that uses as a reference speed the average streamwise velocity component in the whole windfarm volume containing the verticalaxis wind turbines (i.e. over a volume spanning a full x–y plane and the z domain from the bottom to the top of the blades). The comparison of this CP value for different layouts is shown in Fig. 13. As anticipated the staggered cases have higher CP values compared to the aligned ones, for both the clustered and regular designs. Of more interest and relevance is that the clustered designs consistently produce higher power than the prototypical design for any spacing. As indicated previously, the CP value for an isolated turbine is 0.36 and the clusterstaggered designs with an intercluster spacing of 20D surpass this value over the whole wind farm for both wind directions. This is due to the gain in average CP that clusters allow, and the large intercluster spacing that minimizes the effects of being in the wake of the upstream cluster. By reducing the distance between clusters to 10D, the average CP value decreases but remains significantly higher than for the corresponding regular wind farms.
Average windfarm CP values for various configurations (staggered or aligned, clustered or regular) and for various separation distances
Another important result is that the staggered configurations, even at small separation distances, consistently perform better than the aligned ones. The streamwise separation in the staggered 10D case, for example, is the same as the separation in the 20D aligned case, and yet the staggered 10D layout yields a higher CP. One physical reason for this improved performance is that in the staggered cluster farms, in addition to the synergistic interactions within each cluster, the clusters themselves probably interact favourably. One can observe, for example, in Fig. 12c, d that two adjacent clusters produce flow acceleration in between clusters, allowing the next staggered row to benefit from this flow deviation. This is exactly similar to the acceleration within a cluster but now occurs in between clusters, suggesting a fractal attribute to these synergistic interactions (although with only two fractal generations here).
The results in Fig. 13, however, exclude an important difference between these periodic simulations. Due to higher drag forces exerted on the ABL in cases with higher densities (5D spacing) or cases with more efficient farm layouts, the required pressure gradient imposed in the simulations to yield a steadystate mean flow will also be higher. In the LES, at each timestep, the drag exerted on the ABL by the whole wind farm and the ground surface is computed and the needed mean streamwise pressure gradient to balance this drag is imposed. This gradient eventually reaches a steady state as the mean flow equilibrates. As a result, the different cases have unique pressure gradients over the wind farm, and the supplied power input (estimated as the product of the pressure gradient force and the streamwise velocity magnitude) into the domain is not consistent across all cases. As found in Meyers and Meneveau (2010), Goit and Meyers (2015), several approaches can be used to overcome this potential inconsistency. Since our simulations omit the Coriolis force, the best approach is to normalize the power extracted by the power supplied in the windfarm volume to the simulations. With no Coriolis force and under steadystate conditions, the pressure gradient force has to balance the total (turbine + ground surface) drag. One can thus characterize these two equal and opposite forces with a squared friction velocity related to the total domain drag uτH (defined similarly to Calaf et al. 2010; Goit and Meyers 2015). The total power input is thus proportional to UTu 2τH , where UT is the streamwise velocity component averaged over the windfarm volume (average of the domain containing the blades as defined before). In large wind farms, this is the rate of mean kinetic energy input into the domain that can be extracted by the turbines. The velocity upwind of a specific row (used to define CP before) is an outcome of this input rather than the main source of energy as in very small farms. Thus, the kinetic energy that can be extracted in large wind farms scales with the pressure drop and with the farmvolumeaveraged streamwise velocity component UT, and since these farms influence the atmospheric pressure field as well as flow inside them significantly, they influence the power available to them. Therefore, in order to be able to compare the various configurations without this pressuredrop discrepancy, the following C *P relation normalized per unit power input, is introduced,
$$ C_{{P}}^{*} = frac{P_{T} }{frac{1}{2}rho Aleft( {u_{tau H} } right)^{2} U_{T} }, $$
(9)
$$ u_{tau H}^{2} = Delta P_{{drop}} frac{{L_{x} }}{{L_{z} }}, $$
(10)
where PT is the average power per turbine, uτH is the root square of the total mean drag (on ground + turbines), which scales with the total pressure drop ΔPdrop across the farm, and A is the rotor area of a single turbine. The comparison of this new performance metric for the various configurations is presented in Fig. 14. Since uτH is roughly about 10% of the wind speed, C *P is about 100CP and should not be interpreted in the same way as the classic power coefficient. Even after normalizing the total power generated in these layouts by the power input for each case, the cluster cases maintain the highest C *P value, implying that these cases are able to extract more energy from the applied pressure gradient in the field compared to regular wind farms. The relative differences in the performance of wind farms are expected to be closer to the differences depicted in Fig. 13 for smaller farms (CP is strictly applicable only for a single row), and closer to those in Fig. 14 for larger farms.
Average windfarm C *P values normalized per unit applied power input to the turbine domain for various configurations (staggered or aligned, clustered or regular) and for various separation distances
A comparison of the power density per unit land area used for the various configurations has also been performed, confirming that clustered designs increase the power density, and validating our hypothesis. However, the results have the caveat that the power density is invariably higher for smaller spacings, even when the turbines in the farm are not being used efficiently (low CP). Therefore, power density itself cannot be used as a metric for optimizing farm layout. In order to have a more realistic and practical metric, the total capital cost per unit power generation Ttotal is computed. Since the power generation for each farm is proportional to the sum of the C *P values of all the individual turbines in a given lot area of fixed size, we use this sum, denoted as C *Σ P , for normalization instead of the actual power. The capital costs consist mainly of the cost of the land and the turbines. The following normalized energy cost function can thus be computed (similar to Meyers and Meneveau 2012)
$$ frac{{T_{{total}} }}{{C_{{varSigma {{P}}}}^{*} }} = left( {frac{{T_{{land}} ,A_{{L}} + varGamma_{{A}} ,A_{{L}} ;T_{{turbine}} }}{{C_{{varSigma {{P}}}}^{*} }}} right) = left( {frac{{T_{{land}} }}{{T_{{turbine}} }} + varGamma_{{A}} } right)frac{{A_{{L}} ,T_{{turbine}} }}{{C_{{varSigma {{P}}}}^{*} }} , $$
(11)
where ΓA is the wind turbine density per unit area and AL is the total lot land area. Tland is the cost of land per unit area and Tturbine is the cost of a single turbine. Using different landcost to turbinecost ratios, and the cost for a typical individual turbine similar to the one simulated (≈ $US 10,000) (Dabiri 2011) in Eq. 11, the normalized energy costs were computed and plotted in Fig. 15. Using this comparison metric also indicates that the triangularcluster staggered layout has the lowest capital cost per projected unit power generated, and is therefore the optimum design among those investigated.
Total capital cost per “unit power” generated for the various cases
A similar analysis has been made using total CP and the results also indicate that wind farms with cluster designs are the most optimal amongst those investigated here. Again, we reiterate that the comparison with CP is more relevant for very small farms, while if one uses C *P , the results are more representative of large farms.
Conclusions
We have presented a novel concept for optimizing the layout of large verticalaxis windturbine farms, taking advantage of synergistic interactions between closelyspaced turbines that were previously shown to yield higher power for a limited number of turbines. Using an actuatorlinemodel representation of the turbines, embedded in a largeeddy simulation, the modelled wake generated by two counterrotating turbines is first successfully validated against observations from field experiments. To take advantage of the high wind speed created by the flow deflection of verticalaxis wind turbines when placed in close vicinity, we propose a triangular cluster design consisting of three verticalaxis wind turbines, which form the basis for larger wind farms. The triangular design is the one that best exploits flow acceleration, with a minimal increase in wake shadowing.
The influence of interturbine spacing relative to their diameter, L/D, was then investigated to optimize a single cluster in terms of the total generated power, the omnidirectionality of its performance, and the needed downstream wakerecovery distance. Changing the turbine spacing, the cases with L/D values of 3, 4 and 5 were shown to generate the highest clusteraveraged power. Further tests were then performed with these three spacings only, and the case with L/D = 5 emerged as the one with the highest clusteraveraged power over all wind directions: the generated power for this case is about 10% higher than that produced by three isolated turbines. Furthermore, L/D = 5 results in the lowest variation of the generated power with wind direction, and the downstream wakerecovery distance is the shortest (since the cluster is more “porous”). Therefore, this cluster design confirmed the use of synergistic verticalaxis windturbine interactions to increase power production, and would generate a higher power density (power generated per unit land used) due to the proximity of the rotors. It was hence adopted for configuring large VAWT wind farms.
Farms that use this advanced cluster design, and a sufficient distance for wake recovery between clusters of 10D and 20D, were then compared to prototypical aligned and staggered configurations for infinitelylarge wind farms, with different turbine horizontal spacings of 5D, 10D and 20D. For the very large wind farms chosen, the results show that the average windfarm power coefficient, using two distinct normalizations, is much higher for the staggeredtriangle clusters than for the wind farms with regular configurations. Using these average power coefficient results and a simple capital cost function for the whole wind farm, while varying the landtoturbine cost ratio, we also showed that the windfarm design with staggeredtriangle clusters is the optimal design (amongst the ones considered here) in terms of cost per unit power produced.
These results strongly indicate that VAWT farms can and should be configured using different approaches than those used for horizontalaxis rotors (although the potential benefits of HAWT clustering could also be investigated). A significant increase in power and decrease in capital costs can be achieved using the ability of verticalaxis wind turbines to positively boost the power production of nearby turbines if properly configured. A further important aspect of the results is that, in addition to turbine interactions within a cluster, the clusters themselves interact synergistically, further boosting power production. It should also be mentioned that one of our criteria in optimizing the clusters and farms was omnidirectionality. We sought to propose configurations with performances that are not strongly dependent on wind direction since this is also a major advantage of individual verticalaxis wind turbines. However, if this criterion is relaxed, for example in places where there is a dominant wind direction, the optimal cluster designs can be very different and can use this synergistic interaction between clusters as well, with potentially higher power densities.
Finally, one factor that plays an important role in modulating power output from a large wind farm is atmospheric stability. It has been shown that the diurnal cycle and a range of ABL stability strongly influence the performance of large HAWT wind farms (Lu and PortéAgel 2011; Abkar et al. 2016), and the same is expected for the VAWT farms used herein. However, we focused on the basic neutral case only, and ALM–LES model investigations of the effect of atmospheric stability on VAWT windfarm operation are left for future studies.
References

Abkar M, Sharifi A, PortéAgel F (2016) Wake flow variability in a wind farm throughout a diurnal cycle. J Turbul 17:420–441. https://doi.org/10.1080/14685248.2015.1127379

Araya DB, Craig AE, Kinzel M, Dabiri JO (2014) Loworder modeling of wind farm aerodynamics using leaky Rankine bodies. J Renew Sustain Energy 6:63118. https://doi.org/10.1063/1.4905127

BouZeid E, Meneveau C, Parlange M (2005) A scaledependent Lagrangian dynamic model for large eddy simulation of complex turbulent flows. Phys Fluids 17:25105. https://doi.org/10.1063/1.1839152

BouZeid E, Parlange MB, Meneveau C (2007) On the parameterization of surface roughness at regional scales. J Atmos Sci 64:216–227. https://doi.org/10.1175/JAS3826.1

Calaf M, Meneveau C, Meyers J (2010) Large eddy simulation study of fully developed windturbine array boundary layers. Phys Fluids 22:1–16. https://doi.org/10.1063/1.862466

Calaf M, Higgins C, Parlange MB (2014) Large wind farms and the scalar flux over an heterogeneously rough land surface. BoundaryLayer Meteorol 153:471–495. https://doi.org/10.1007/s1054601499596

Chamorro LP, PortéAgel F (2009) A windtunnel investigation of windturbine wakes: boundarylayer turbulence effects. BoundaryLayer Meteorol 132:129–149. https://doi.org/10.1007/s1054600993808

Chamorro LP, PortéAgel F (2010) Effects of thermal stability and incoming boundarylayer flow characteristics on windturbine wakes: a windtunnel study. BoundaryLayer Meteorol 136:515–533. https://doi.org/10.1007/s1054601095121

Claessens MC (2006) The design and testing of airfoils for application in small vertical axis wind turbines. Masters Thesis, McMaster University, Hamilton, Ontario, Canada

Cortina G, Calaf M, Cal RB (2016) Distribution of mean kinetic energy around an isolated wind turbine and a characteristic wind turbine of a very large wind farm. Phys Rev Fluids 1:74402. https://doi.org/10.1103/PhysRevFluids.1.074402

Dabiri JO (2011) Potential orderofmagnitude enhancement of wind farm power density via counterrotating verticalaxis wind turbine arrays. J Renew Sustain Energy 3:43104. https://doi.org/10.1063/1.3608170

Goit JP, Meyers J (2015) Optimal control of energy extraction in windfarm boundary layers. J Fluid Mech 768:5–50. https://doi.org/10.1017/jfm.2015.70

Hezaveh SH, BouZeid E, Lohry MW, Martinelli L (2016) Simulation and wake analysis of a single vertical axis wind turbine. Wind Energy 20:713–730. https://doi.org/10.1002/we.2056

Huang J, BouZeid E (2013) Turbulence and vertical fluxes in the stable atmospheric boundary layer. Part I: a largeeddy simulation study. J Atmos Sci 70:1513–1527. https://doi.org/10.1175/JASD120167.1

Kinzel M, Mulligan Q, Dabiri JO (2012) Energy exchange in an array of verticalaxis wind turbines. J Turbul 13:N38. https://doi.org/10.1080/14685248.2012.712698

Kinzel M, Araya DB, Dabiri JO (2015) Turbulence in vertical axis wind turbine canopies. Phys Fluids 27:115102. https://doi.org/10.1063/1.4935111

Li Q, BouZeid E, Anderson W, Grimmond S, Hultmark M (2016) Quality and reliability of LES of convective scalar transfer at high Reynolds numbers. Int J Heat Mass Transf 102:959–970. https://doi.org/10.1016/j.ijheatmasstransfer.2016.06.093

Lu H, PortéAgel F (2011) Largeeddy simulation of a very large wind farm in a stable atmospheric boundary layer. Phys Fluids 23:65101. https://doi.org/10.1063/1.3589857

Marquis M, Wilczak J, Ahlstrom M, Sharp J, Stern A, Smith JC, Calvert S (2011) Forecasting the wind to reach significant penetration levels of wind energy. Bull Am Meteorol Soc 92:1159–1171. https://doi.org/10.1175/2011BAMS3033.1

Meyers J, Meneveau C (2010) Large eddy simulations of large windturbine arrays in the atmospheric boundary layer. In: 48th AIAA aerospace sciences meeting, Orlando, Florida, USA, 4–7 January 2010, pp 1–10

Meyers J, Meneveau C (2012) Optimal turbine spacing in fully developed wind farm boundary layers. Wind Energy 15:305–317. https://doi.org/10.1002/we.469

Sarlak H, Nishino T, MartínezTossas LA, Meneveau C, Sørensen JN (2016) Assessment of blockage effects on the wake characteristics and power of wind turbines. Renew Energy 93:340–352. https://doi.org/10.1016/j.renene.2016.01.101

Stevens RJAM, Gayme DF, Meneveau C (2014) Large eddy simulation studies of the effects of alignment and wind farm length. J Renew Sustain Energy 6:1–14. https://doi.org/10.1063/1.4869568

Troldborg N, Sørensen J (2014) A simple atmospheric boundary layer model applied to large eddy simulations of wind turbine wakes. Wind Energy 17:657–669. https://doi.org/10.1002/we.1608

U.S. Energy Information Administration (2013) International energy outlook 2013. Technical Report, U.S. Energy Information Administration. https://www.eia.gov/outlooks/ieo/pdf/0484(2013).pdf

Wu YT, PortéAgel F (2012) Simulation of turbulent flow inside and above wind farms: model validation and layout effects. BoundaryLayer Meteorol 146:181–205. https://doi.org/10.1007/s105460129757y

Yuting WPF (2011) Largeeddy simulation of windturbine wakes: evaluation of turbine parametrisations. BoundaryLayer Meteorol 138:345–366. https://doi.org/10.1007/s105460109569x
Acknowledgements
This work was supported by the Siebel Energy Challenge and the Andlinger Centre for Energy and the Environment of Princeton University. The simulations were performed on the supercomputing clusters of the National Centre for Atmospheric Research through project P36861020 and UPRI0007, and of Princeton University.
Author information
 Department of Civil and Environmental Engineering, Princeton University, Princeton, NJ, 08540, USA
 Seyed Hossein Hezaveh
 & Elie BouZeid
 Department of Civil and Environmental Engineering, Department of Mechanical Engineering, Stanford University, Stanford, CA, 94305, USA
 John Dabiri
 Graduate Aerospace Laboratories, California Institute of Technology, Pasadena, CA, 91125, USA
 Matthias Kinzel
 Department of Mechanical Engineering, The University of Utah, Salt Lake City, UT, 84112, USA
 Gerard Cortina
 Department of Mechanical and Aerospace Engineering, Princeton University, Princeton, NJ, 08540, USA
 Luigi Martinelli
Correspondence to Elie BouZeid.
Appendix A: Field Experimental SetUp
Experiments were conducted at the Caltech Field Laboratory located in the Antelope valley of northern Los Angeles County in California, USA (Kinzel et al. 2015). The surroundings are flat desert terrain with no obstacles for at least 1.5 km horizontally in all directions. The turbines were in operation, and data were collected, mostly in the middle of the day when the ABL was not neutrally stable. Therefore, buoyancy generation contributes to the turbulence levels of the flow. For the validation, periods with nearneutral stability were identified and used. The turbines are rated at 1.2 kW (Windspire Energy, Inc., Reno, Nevada, USA), with a total height of 9.1 m, a rotor height of 6.1 m, and a diameter of 1.2 m. Their maximum rotation rate is 420 revolutions min−1 at a wind speed of 12 m s−1, yielding a tipspeed ratio λ = 2.3. The cutin and cutout speeds for this turbine are 3.8 and 12 m s−1, respectively.
The wind velocity was measured from a movable 10m high meteorological tower (Aluma Towers, Inc., Vero Beach, Florida, USA) with seven, verticallystaggered, threecomponent ultrasonic anemometers (Campbell Scientific CSAT3, Logan, Utah, USA). The sensors are equally spaced by 1 m between the top and the bottom of the turbine rotor, i.e., between heights of 3 and 9 m. The data collection frequency was 10 Hz and the tower was moved to the various horizontal locations where data were to be collected.
The velocity profiles were measured at several locations along the centreline of the turbine arrays. For the configuration shown herein, these positions were 15D and 1.5D upstream, as well as 2D and 8D downstream, from the front of the array. The measurements were taken for at least 150 h at each position. The dataset was filtered for times when the freestream wind speed was within the cutin and cutout wind speeds of the turbines and the wind direction was within ± 10° from the array centreline. The mean horizontal wind speed as measured by the reference sensor was 8.2 m s−1 during the times when the doublet configurations used here were tested, leading to an average Reynolds number of ≈ 106 based on rotor diameter. The prevailing wind direction was south–southwest, i.e., 223° for doublet configurations.
Rights and permissions
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.