Mid-fidelity simulations and comparisons of five techniques for axial induction control of a wind turbine

As wind turbines are more frequently placed in arrays, the need to understand and mitigate problems arising from their wakes has increased. When downstream turbines are in the wakes of upstream ones, the downstream turbines produce less power, require more maintenance, and have shorter lifetimes. One wake mitigation technique is known as axial induction control (AIC) and it involves derating (operating suboptimally) upstream turbines such that more energy remains in their wakes for downstream turbines to harvest. While there has been considerable research on this technique, much of it has suffered from 5 a misunderstanding of the most important parameters in optimizing AIC. As such, the research has been largely inconclusive. Herein, we seek to rectify several perceived shortcomings of previous work by using mid-fidelity simulations to compare five different techniques for AIC at three different derate percentages against a baseline case and examining the recovery of the wake. We find that only the case with the lowest derate, 10%, and using maximum thrust exceeds the baseline when estimating the combined power of the simulated turbine and a virtual turbine five diameters downstream and that it produced 10% more 10 power. Furthermore, these results help to validate previous work that concluded that the excess energy that is in the wake of a derated turbine will be at the edges of the wake unless the wake can sufficiently recover before the next downstream turbine. Finally, all together this suggests that the precise combination of derate percentage and the method used to derate turbines (i.e., the precise combination of pitch and torque controls), as well as the spacing and arrangement of turbines, must all be considered when optimizing AIC, and that substantial power gains may be possible. 15


Introduction
As wind energy progresses, attention has shifted from optimizing individual turbines to optimizing wind farms to mitigate the negative effects of upstream turbines' wakes. In addition to turbine siting, there has been very promising research into coordinating turbine controls within a wind farm to maximize its total power. These methods have been recently and well-reviewed in Houck (2021) and Kheirabadi and Nagamune (2019). Axial induction control (AIC) is one method that has been heavily re-20 searched, and it involves derating upstream turbines for the benefit of downstream turbines. By reducing the power production of upstream turbines, the effects of their wakes are mitigated and downstream turbines produce more power as a result. The net result of this is an increase in total power of the multi-turbine system.
Previous research into AIC has produced mixed results and has been largely inconclusive as to its effectiveness. Some of the ambiguity of the results of previous studies can be attributed to one or more of three particular criticisms: the use of low-fidelity 5.5 methods at computational speeds that are much faster than LES or direct numerical simulations (DNS) and can include the effects of a turbine's pitch and TSR, and so its induction and thrust as well. To further investigate the multiple pathways to a given power target when derating, we evaluated the differences in the wake when using five different strategies to achieve derates of 10%, 20%, and 40% of the baseline power production: constant pitch, constant TSR, maximum rotation rate, minimum thrust, and maximum thrust. When considering the simulated turbine and a virtual downstream turbine, our results indicate that maximizing the thrust of the derated turbine yields the highest combined power of the two, though the placement of the downstream turbine is an important consideration in choosing a derate method.
The remainder of the article is organized as follows: In Section 2, details of the FVWM code, the turbine used, the simulation set up, and careful consideration of the code's limitations are provided. In Section 3, we provide results of the 16 simulations that were performed, and in Section 4 we evaluate the relevant differences in the wakes created by different derates and derate 70 methods. Finally, conclusions and recommendations for future work are given in Section 5.

Simulation Setup
This study used the Code for Axial and Cross-flow TUrbine Simulation (CACTUS), a mid-fidelity FVWM developed by Sandia National Labs (SNL), to simulate a turbine in a three-dimensional domain. The turbine's blades are represented by 75 discretized blade elements with lookup tables of two-dimensional airfoil polar data. The Kutta-Joukowski theorem relates the lift data to the bound vorticity at each blade element and the Helmholz theorem recovers the spanwise and trailing vorticities.
A potential flow model uses straight line vortex elements to represent the wake and the induced velocities are calculated via the Biot-Savart law. CACTUS does not attempt to approximate viscous processes as is common in other FVWMs (Ananthan et al., 2002;Leishman et al., 2002;Marten et al., 2020Marten et al., , 2016, so the vortex filaments are advected with the freestream flow 80 without dissipation or diffusion, though regularization is implemented with a cut-off radius for the vortex cores. More complete descriptions of CACTUS are given in Murray and Barone, 2011;Wosnik et al., 2016).
Extensive test cases were run in CACTUS to verify the best combination of many settings to ensure as accurate results as the tool allows. Performance coefficients of power, torque, and thrust are within 1% of their converged values within only approximately 20 revolutions for the majority of simulations, but wake statistics require more revolutions. All simulations were 85 run for 45 revolutions and the velocity data from the last 15 revolutions is saved for analysis based on convergence tests of many of the wake statistics used in section 3. The number of time steps per revolution was also given careful consideration. Bhagwat and Leishman (2000) showed in their FVWM code that increments of five degrees or less were required to achieve second order accuracy. Test cases with CACTUS indicated that 50 time steps per revolution, or 7.2 • , achieved convergence of wake statistics, but 72 time steps per revolution, or 5 • , was used to ensure greater accuracy and to help with the stability of the 90 solution by reducing the change in blade loads between timesteps.
CACTUS allows for user selection of either no vortex core model, a constant one, or a linear one. It was observed that having no core model accelerated wake recovery while the constant and linear core models produced wakes that looked fairly similar to each other. The linear core model was chosen because it is a closer approximation to the viscous core of a real vortex and because it showed a slightly slower, and so more conservative, wake recovery rate. In all of these models, regularization is 95 performed using a cut-off radius to avoid the singularity of the core at zero radius. This was left at the default of 10 −7 . Defaults were also used for the initial vortex core sizes, which are set as the maximum chord for the bound vorticity, the maximum blade element span for the trailing vorticity, while the core size for the spanwise vorticity is calculated from a reference distance between spanwise wake lines given the temporal discretization level (Murray et al.). These initial core sizes are in line with recommendations from Xu et al. (2019) who suggest keeping them above half the blade chord to ensure the power coefficient 100 is accurately predicted. Time step filtering was used to smooth bound vorticity over three timesteps in all simulations as it was found necessary for stability during simulations with high T SR. The simulated turbine is the National Rotor Testbed (NRT) rotor designed by SNL to fly on a modified V27 turbine while replicating the scaled wake of a 1.5 MW turbine . This turbine is part of the Scaled Wind Farm Test facility (SWiFT) operated by SNL at the Reese Technology Center in Lubbock, TX (Berg et al., 2013). The rotor has a diameter of 105 27 m with a 4 • shaft up-tilt and a hub height of 32 m. The CACTUS software includes several Matlab scripts for creating the necessary geometry files by inputting blade data and setting the desired pitch. Each blade was discretized into 42 uniformly distributed elements using nine different sets of polars that are provided to CACTUS as lookup tables over a range of Reynolds numbers. Polars for the lookup tables were generated in QBlade (Marten et al., 2013), which integrates XFOIL for the calculation of two-dimensional airfoil data and allows for extrapolation of polar data to 360 degrees with the Montgomerie method 110 (Montgomerie, 2004). The Leishman-Beddoes dynamic stall model is implemented in CACTUS to account for any transient stall effects (Leishman and Beddoes, 1989). It should be noted that the turbine's nacelle, including the hub, are not simulated and result in an unrealistic jet in the middle of the near-wake. A tower model adds an empirically defined wake deficit to the freestream flow to model the presence of the tower.
For purposes of the code's calculations, the turbine is in an infinite domain. The velocity field was saved over a domain 115 stretching two diameters (D) upstream of the turbine, 6D downstream, and 1.5D in either spanwise direction and upward as seen in Fig. 2. After 24D downstream of the turbine, the vortex elements are ignored in the induced velocity calculations to facilitate faster computation. Based on a grid sensitivity study, the velocity data are computed and saved on a grid with spacing of 0.05D in all directions.
The inflow is the same in all cases. CACTUS does not provide methods for a turbulent inflow. Simulating a sheared inflow 120 using a passive inflow model results in non-physical lifting of the turbine wake (Branlard et al., 2015). The inflow is uniform, laminar, and steady at 6.8 m/s, the density of air is 1.08 kg/m 3 , and its dynamic viscosity is 1.8 × 10 −5 (used for calculation of the Reynolds number and selection of the appropriate airfoil polar), which are the annual mean wind speed, density, and viscosity of air at the SWiFT site (Kelley and Ennis, 2016).

Simulation Cases
Five different derate strategies were simulated with three different derates for a total of 16 simulations including a Baseline case during which the turbine is operating at its rated power for the freestream wind speed. During this simulation, the turbine operates with a pitch of zero degrees and a T SR of nine. The remaining cases can be seen in Table 1 and are plotted in Fig.   3 on top of contours of C P and C t . It should be noted that, in theory both the constant pitch and constant TSR cases present 130 two options: an increase or decrease. For the constant pitch cases, the the T SR had to be decreased as an increase would have exceeded the turbine's maximum rotation rate. For the constant T SR cases, a higher pitch was chosen as it maintains a positive angle of attack on the airfoil for which performance is more accurately predicted when calculating polars, which should aid in the accuracy of the simulations.
The three derates were chosen to extend over a wide range of possibilities. In a deep array, it is possible (though it remains to 135 be determined what is optimal) that the most upstream turbines would require significant derates to maximize the array power.   Table 1 is plotted along with the Baseline on top of contours of CP with contours of Ct (white lines) and their values overlaid. The five strategies considered were chosen for the following reasons: Keeping pitch or T SR constant while derating yields an easier control option and could reduce wear on components if frequent adjustment is required. The idea of maximizing the rotation rate for a given derate comes from literature studying the breakdown of tip vortices (Houck and Cowen, 2019;Iungo et al., 2013;Ivanell et al., 2010;Lignarolo et al., 2015). These studies have found that when turbines operate at higher rotation 140 rates, it reduces the pitch (streamwise distance) between the tip vortex helices, and this causes them to interact and break down sooner. This breakdown accelerates wake mixing and decay. Finally, the reasoning behind minimizing or maximizing the thrust was reviewed above.

145
As mentioned, the induced velocity is calculated at every timestep and output on a user-defined grid. This three-component velocity field is then used for further wake analyses as well as the turbine's coefficients of power, thrust, and torque. Care must be taken when examining turbine wakes simulated by FVWMs as CACTUS, in particular, has primarily been validated against experimental data for power and blade loads as opposed to wake data. In the near-wake, the individual vortex filaments remain neatly organized in helices and have been shown to accurately represent near wake dynamics in CACTUS (Ennis et al.,150 2015) and other FVWM codes (Corniglion et al., 2020;Marten et al., 2020Marten et al., , 2016. In the far-wake, however, and in particular after tip vortex breakdown, the filaments become highly disorganized and entangled and create a sort of numerical turbulence. Physically, there is a feedback loop in which vortex-vortex interactions cancel the spanwise and vertical components of the induced velocity, which causes the vortices to align more in the streamwise direction with varying orientations such that the spanwise and vertical induced velocity components further cancel and result in apparent wake recovery. This leads to an 155 increasingly streamwise-dominant flow, though it may be impossible to distinguish between physical and numerical instabilities and dissipation (Bhagwat and Leishman, 2001). The calculations performed by CACTUS are explicitly inviscid, i.e., there is no diffusion or dissipation, so solutions become less physically realistic in proportion to the wake age or its distance downstream, but this is difficult to quantify without comparison to validation data. Wind turbine wake decay, however, is a largely inviscid process (large scale turbulence is essentially inviscid), and dissipation is actually highest in the near-wake (Chamorro et al., 2014;Houck and Cowen, 2019). Since the effects of viscosity occur very close to surfaces and within the relatively small vortex cores, the bulk of the flow can be simulated by a potential flow model such as CACTUS.
To the author's knowledge, there is only one study that compares CACTUS to a higher fidelity method in the far-wake of a wind turbine. Ennis et al. (2015) compared CACTUS to the Virtual Wind Simulator (VWiS), an actuator line model LES code (Yang et al., 2015) simulating the same turbine and inflow conditions in both. Their results indicated that the velocity deficit 165 profiles from CACTUS and VWiS match well up through 7D downstream and that the wakes demonstrate similar growth rates as defined by the radius that contains 95% of the mass deficit from the freestream. A standard error calculation on their timeaveraged streamwise velocities shows that they have an error of less than 10% between approximately 3D and 7D downstream and that it remains below 15% between 0D and 8D. CACTUS, however, had greater Reynolds shear stress due to its lack of dissipation and, as a result, experienced a more rapid momentum recovery of its wake.

170
This last point suggests some caution in interpreting results herein as the turbine is operated in several distinct ways. Some of the test cases are likely to (indeed, intended to) create differences in wake recovery that are likely to mitigate or exacerbate the shortcomings of CACTUS depending on the levels of "turbulence" they generate. Despite the fact that CACTUS simulations are not as high fidelity as other methods, this is somewhat mitigated by the control and treatment approach. Any errors induced by the method relative to higher fidelity or experimental methods can likely be viewed as bias errors that are negligible when 175 restricting results to differences among the cases. Additionally, the drastically lower computational cost of CACTUS relative to higher fidelity LES allows for the broader exploratory study presented here, which can be followed by higher fidelity methods to verify any promising results found herein.

Verification of simulations 180
Before exploring the results related to the main objective, the data produced by CACTUS were used to verify it against theory in lieu of higher fidelity or experimental data. In particular, a control volume of energy was calculated as the primary concern in a simulation without viscosity is that energy in the domain will become imbalanced due to the lack of dissipation. The complete equation for the control volume of energy is given as 185 where E is energy; V is the control volume; U is the velocity vector; n is the outward facing normal of each face of the control volume, A; ν is viscosity; and S ij is the strain rate tensor. The transport, T , is where p is the fluid pressure, and ρ is the fluid density.Ẇ is the work done by the turbine and iṡ where C P is the coefficient of power defined as where P mech is the mechanical power produced by the turbine, U ∞ is the freestream velocity, and A T is the turbine rotor-swept area. While it would not typically apply to flow through a wind turbine, the assumptions of the simulation, namely that it is inviscid, incompressible, and steady (in the mean), allow us to apply Bernoulli's Principle to calculate the pressure field at each 195 time step, which is where U res is the resultant velocity at each point in the velocity field and is where u, v, and w are the instantaneous streamwise, transverse, and vertical velocities, respectively.

200
Ignoring viscous terms, the resulting energy balance is shown in Fig. 4. The first term of Eqn. 1 is the "Volume" and the second term results in a sum of energy flux for each surface of the control volume. The work done by the turbine is negligible compared to the other terms due its relatively small area compared to the control volume surfaces and its efficiency. It can be seen that, in an average sense, energy is fairly balanced in most simulations, with 10maxCt being the largest exception due to a larger Volume term. The larger Volume terms indicate a lack of complete convergence as one would expect the temporally converge to near-zero given enough time, but it is clear that the simulations are not all equally converged. Concern for this is mitigated by looking at Fig. 6, which shows the temporally averaged energy residual as a function of x. Here we see that the average residual never has a magnitude greater than 2% of the total energy. Within the confines of CACTUS, further convergence would most likely be achieved by simply simulating additional rotor revolutions.
Finally, the quality of the apparent turbulence generated by CACTUS can be evaluated by looking at the spectra of the 215 centerline resultant velocity as a function of x. Under realistic flow conditions, we would expect to see a large inertial subrange and small dissipation subrange in the near-wake with the former decreasing and latter increasing as the wake progresses downstream (Chamorro and Porté-Agel, 2009). Figure 7 shows the spectra at x/R = 0, 1, 3, 6, and 12 and shows slopes of −5/3 and −7/3. The −5/3 slope indicates that CACTUS captures an inertial subrange in the apparent turbulence that is generated. The −7/3 slope, however, does not conform to any known turbulence theory for this type of flow and is a direct 220 result of the non-physical lack of dissipation. The spectra at x/R = 0 and 1 have little to no inertial subrange and are even flat in many cases indicating that nothing like turbulence has developed from vortex-vortex interactions yet. By x/R = 3, most   cases begin to demonstrate an inertial subrange and by x/R = 12, they all exhibit a −7/3 slope. This is likely the result of the flow developing smaller and smaller velocity fluctuations without any energy actually being dissipated while large scale turbulence is continuously generated, but this would need to be verified with further investigation. Furthermore, these results 225 are within the wake, though any ambient turbulence is created by the induced velocities of vortex-vortex interactions of the wake. In a sense, the ambient turbulence develops along with the wake-generated turbulence.

Qualitative results
Moving on to results focused on the main objective, the wakes generated by each case are shown in Fig. 8. There are many ways to visualize wind turbine wakes (see (Quon et al., 2020)) and these results are only meant to be qualitative. The contours used 230 to generate these wake visualizations were determined by the following steps: Contours of the temporally averaged streamwise velocity deficit calculated as 1 − U /U ∞ were created in each y − z plane. These were then filtered to remove, first, open contours, then contours with areas smaller than the rotor-swept area, and then contours with a centroid coordinate that was more than 0.7R from the center of the rotor-swept area. From the remaining contours in each plane, the contour level closest to zero was chosen corresponding to a more recovered (smaller deficit) wake. The red dot at the middle of each contour is the 235 centroid of that contour and provides an indication of wake center. Looking down the columns of Fig. 8, we see that, for a (a) Baseline given derate technique, wake recovery happens more slowly for higher derates, though the difference can be subtle. This may be expected as smaller velocity deficits would create shear of a smaller magnitude and the shear-induced mixing is largely responsible for wake decay (Castillo and Newman, 2013;Chamorro et al., 2012). Already it appears that there might be an exception to this in the maxCt cases as the 40% derate appears to have the shortest near-wake of those cases. Looking across 240 the rows, we can compare the different techniques for a given derate and here we see much larger differences. The maxCt cases stand out as those wakes appear to decay more quickly while the cstP cases barely decay at all.
Contours of the velocity deficit in the x − z plane shown in Fig. 9 provide additional detail regarding the evolution of the wake in each case. In general, it appears that higher derates (lower power setpoints) produce less mixing and a slower wake  limited, which slows wake recovery. This is also seen to have the effect of stretching the region of high velocity (low velocity deficit) that appears where the rotor hub and nacelle would be. The maxCt cases have a behavior opposite of the others in that the highest derate (40%) appears to recover faster than the lower derates. In these cases, there is a substantially larger deficit about the edges of the near-wake with a very sharp gradient to increased velocities in the ambient flow. The high thrust in 250 these cases, which is higher for the higher derates, effectively creates a blockage behind the rotor leading to an acceleration of ambient flow around the near-wake. This area of high gradient produces rapidly growing shear (as shown below) and facilitates rapid mixing and wake recovery.
The time-averaged vertical vorticity is shown in Fig. 10. Tip vortices play an important role in wake decay as their highly stable nature tends to shield the wake from significant mixing until they destablize and decay at which point increased mixing 255 between the wake and the ambient flow occurs (Brown et al., 2021;Houck and Cowen, 2019;Ivanell et al., 2010;Lignarolo et al., 2015). Simulating the evolution of the tip vortices is also one of the strengths of an FVWM. In addition to the vertical vorticity, a metric was created to track the strength of the vorticity along the lower side of the rotor in each subplot. First, at   Figure 11. The smoothed wake width determined from contours shown in Fig. 8 as twice the average distance from the centroid to the contour. Most cases show a very similar pattern of wake growth, but the maxCt and Baseline cases grow and fluctuate much more. each x, we find the maximum vorticity between 50% and 150% of the rotor radius and that is marked with a small red dot.
The red star represents the first x location downstream of the rotor where the local maximum is less than or equal to 50% of 260 the global maximum for that case. The precise definition is arbitrary, but it helps to indicate the rate of decreasing vorticity strength. How far downstream 50% of the maximum vorticity occurs tracks with the qualitative observations made from the wake visualizations in Fig. 8 and we see that there is no significant wake decay until the vorticity significantly decays. In particular, we see that vorticity in the maxCt cases decays the fastest and that the maxRR method appears to have the predicted effect.

265
A final qualitative metric is shown in Fig. 11, which shows a calculation of wake width. This width was calculated from the same routine as the wake visualizations in Fig. 8 with the width taken as twice the average distance from a contour to its centroid. These widths as a function of x are then smoothed for clarity in plotting. The same x locations marking 50% of the maximum vorticity are also plotted for each case to show how it corresponds to the growth of the wake. The Baseline and maxCt cases all grow much faster at first with subsequent reductions suggesting wake decay. It should be noted that the contour 270 selecting routine described above is a bit arbitrary and imperfect, especially when the flow becomes more chaotic, so larger fluctuations should be considered with caution. The rest of the cases stack up very neatly in order of derate with lower derates (higher power) producing wider wakes due to the larger velocity deficits and resulting shear layers. Figure 12 shows the results of all simulated cases in terms of the achieved C P , C t , and C q , the coefficients of power, thrust, 275 and torque, respectively. C P was defined previously, and C t and C q are defined as     Figure 13. The rotor-equivalent wind speed provides an indication of velocity recovery, which largely correlates to the derate, though the maxCt and Baseline cases evolve very differently from the rest.

Quantitative results
and where T is the rotor thrust, Q is the rotor torque, and R is the rotor radius. In Fig. 12a, we see that many cases did not perfectly 280 achieve the desired derate but that they are all very close. Figure 12b shows first that all methods except maxCt result in reduced thrusts relative to the Baseline and that, in most cases, a greater derate corresponds to less thrust. Notably, the maxCt cases are the exception and their thrusts increase with greater derate. Finally, Fig. 12c shows the torque coefficient. Here, it is the constant pitch cases that stand out. As mentioned previously, because the Baseline operation is near the maximum T SR for the conditions, the T SR was reduced to derate the turbine while keeping the pitch constant. Since the turbine is operating at a 285 lower speed, the torque is higher to achieve the desired power.
We can also evaluate the results through more quantitative metrics such as the rotor-equivalent wind speed, or REW S, shown in Fig. 13. The REW S is calculated for the same rotor location as the simulated turbine at each x as where U is the temporally averaged streamwise velocity, and the sum is taken over discrete ∆y and ∆z that fall on or within 290 A T at each x. This provides a sense of the velocity recovery rate for each case. Here again, we see that the Baseline and maxCt cases stand out producing a larger drop in REW S in the near-wake because of their larger velocity deficits, but a much faster recovery thereafter. The maxCt cases are the only cases to fully recover to the freestream velocity. It must also be noted that the locations marking 50% of the maximum vorticity correspond approximately to the beginning of the recovery phase of the wake after the initial decay in REW S. In the near-wake, the other cases are organized, first, by derate with lower derates maintaining 295 a higher REW S on account of having smaller deficits, and, second, by technique with minCt, cstP, cstT, and maxRR having  Figure 14. The excess energy, calculated as the spatially averaged difference between each case and the Baseline, shows how the evolution of energy in the domain differs among the cases.
decreasing REW S in that order within a derate. In the far-wake, some of the lower derates begin to exceed the REW S of higher derates, most notably the maxRR cases. This is the first suggestion of a break-even between producing a large velocity deficit and facilitating wake recovery by shear production.
Another way to view wake recovery is in terms of energy, and more specifically in terms of the excess energy relative to 300 the Baseline case. Figure 14 shows the difference between the mean kinetic energy (MKE) of each case at each x from the corresponding MKE of the Baseline spatially averaged in the y − z plane and normalized by the spatially average MKE of the Baseline. It is essentially a percent difference in energy between each case and the Baseline as a function of x. The MKE is where V and W are the temporally averaged transverse and vertical velocities, respectively. It should be noted that this metric 305 would be more meaningful if it were calculated across an area representing some physical parity such as a wake envelope or return to freestream velocity (similar to a boundary layer height) rather than spatially averaging across an entire plane. This gets at some of the inherent difficulty in defining a wind turbine wake, but this metric maintains consistency in that all cases use the same domain, so there is a parity in that. It is likely no surprise by now that the maxCt cases stand apart from the rest. Their large induction zones produce a deficit of excess energy upstream of the turbine, which is further evidence of their 310 blockage effect, and, opposite to all other cases, there is an almost immediate increase in excess energy downstream of the turbine. Though the REW S calculation would suggest a reduced MKE, recall that it is only calculated within a rotor-swept area while this excess energy calculation includes the entire y −z plane. Recalling Fig. 9, it is the acceleration of the freestream around the wake that causes the excess energy observed here. This effect, however, is short-lived, and the excess energy in the maxCt cases decays again before exhibiting another small increase that can be attributed to the increased mixing. All other 315 20 https://doi.org/10.5194/wes-2021-122 Preprint.  Figure 15. The absolute value of the Reynolds shear stress along the hub height centerline shows the stark difference in shear development between the maxCt and Baseline cases and all the other cases as well as the relation to vorticity breakdown.
cases are organized as they were in Fig. 13 and exhibit excess energy mostly by virtue of being derated relative to the Baseline.
Interestingly, the other cases have a decreasing excess energy throughout the entirety of their wake recoveries and all have less energy relative to the Baseline by 5D, which suggests reduced mixing relative to the Baseline.
The mixing between wake and ambient flow is well-represented by the Reynolds shear stress, which is plotted in Fig. 15 as a function of x. Again, in lieu of a well-defined wake envelope, the shear is spatially averaged over each y − z plane. Most 320 cases exhibit only mild shear within the domain length, but the maxCt cases produce substantially more shear than all cases including the Baseline. Notice again the effect of increased thrust in accelerating shear production as the higher derates in the maxCt cases produce more shear. In Fig. 16, contours of the difference in Reynolds shear stress between each case and the Baseline in the y − z plane at 5D downstream reveal the spatial distribution of shear. Here we see that, not only do the maxCt cases produce greater shear, it has also penetrated the wake and does not remain along the shear layer typically along 325 the edges of the wake envelope. This indicates that the maxCt cases have achieved net entrainment of energy into their wakes by 5D downstream while all other cases lag behind. Furthermore, we see that while 40maxCt had substantially larger spatially averaged shear, it does not appear to be as concentrated within the rotor disc area as it is for 10maxCt and 20maxCt. We can also notice the 10maxRR appears to be close behind the maxCt cases in this respect, perhaps due in part to its accelerated vorticity breakdown.

Virtual rotor analysis
To better say what the effects of the different derates and techniques would have on a downstream turbine, a virtual rotor is created. Figures 17 and 18 show the combined power of the simulated turbine and one in its wake normalized by the same using the Baseline. Since CACTUS is restricted to simulating a single turbine, the combined power of two turbines was estimated by calculating the power a turbine would produce at a given location based on the REW S calculated at that location. Specifically, 335 for a rotor-swept area equal to the simulated turbine's, the REW S was calculated and then used to calculate the power of a virtual turbine at any location in the domain as where 4 m/s is the rotor's cut-in speed and C P is taken as 0.46, the rotor's theoretical coefficient of power at the range of wind speeds found in the wake. Considering Fig. 17, we see that, due to their relatively faster recovery, the maxCt cases 340 produce more total power than other cases. In particular, 10maxCt produces more total power than the Baseline from about 2.5D and onward reaching a peak of almost 1.2 at 3D. From 4-6D, 10maxCt is approximately constant at about 1.1 times the  Figure 17. The total power of each case and the power of another turbine operating optimally aligned and in its wake as a function of x shows that higher derates offer little to no advantage and that the maxCt cases perform best of all techniques.
power produced in the Baseline. Curiously, the 20% and 40% derates with maxCt end up at about the same total power at 6D downstream despite taking very diffent paths, especially in the near-wake. Indeed, 40maxCt technically produces the highest total power of all cases, but it is not sustained. This suggests a balance between maximizing thrust and sacrificing power to 345 a greater derate. That is, the added thrust achievable with a 40% derate produced wake dyanmics that made up for the initial power lost compared to a 20% derate with a lower thrust. While it does not surpass the Baseline, 10maxRR is best among the remaining cases and could, perhaps, exceed the Baseline given more distance downstream.
Finally, Fig. 18 shows the same total power calculation, but now as a function of z at 5D downstream. In this view, 10maxCt is still the best choice at practically all z's and best when aligned with the upstream turbine. This result may appear unexpected 350 given previous studies that showed that aligned turbines may not be well positioned to benefit from AIC (Annoni et al., 2015), but it is directly in line with those predictions. Because the maxCt cases produce wakes that are recovered enough for the excess energy from derating to have penetrated the wake, the best position for a downstream turbine in these cases is aligned with the upstream turbine. This is also fortuitous since it is the position for which wake mitigation is most needed. In contrast, all other cases suggest a peak total power for a downstream turbine that is offset transversely from the upstream turbine, which 355 is exactly what was shown by Annoni et al. (2015).  Figure 18. The total power of each case and another turbine operating optimally in its wake at different spanwise locations and effective wind directions at five diameters downstream shows that the highest combined power is achieved with a downstream turbine aligned with an upstream turbine operating with maximized thrust at a 10% derate.
where B is the number of blades and C t and C q are evaluated as functions of rotor radius (Manwell et al., 2009).

365
For the virtual turbine, the OOP moment was estimated by using a control volume approach to derive the forces over the rotor-swept area. At 5D, the temporally averaged streamwise velocity was extracted from the 10 grid points corresponding to a blade pointed straight up, down, left, or right and these 40 points were averaged to provide an approximately azimuthally averaged streamwise velocity for a rotor area aligned with the simulated turbine and freestream flow. Because the pressure difference created by an actual rotor would be significantly larger than the pressure difference across the virtual rotor, the 370 pressure difference was not included. The OOP moment was then calculated as where A E (r) are the span-dependent element areas approximated as the product of the element span length and the chord length at the center of each element and r is the spanwise distance along the rotor.
Looking first at the simulated turbine (Figures 19a and 19b), the price for increasing thrust is an approximately 25% increase 375 in the OOP moment. All other cases would decrease the OOP moment since they also decrease thrust and moreso for higher derates. The IP moments correlate to the results for C q shown in Fig. 12c as expected. For the cstP cases to produce the given power without adjusting pitch, the T SR must be very low, which requires a compensation in torque. It can also be inferred from Fig. 12b and Fig. 19a that the tower fore-aft moment would increase by around 25% as it is proportional to the turbine thrust. The OOP moments of the virtual turbine (Fig. 19c) indicate that it would experience an approximately 50% increase in 380 moments during the maxCt cases. Since these blade moments are produced by forces in the streamwise direction, the virtual turbine would likely also experience a similar increase in tower fore-aft moment. There is obvious uncertainty in the method for computing the OOP moments of the virtual turbine, which is why they have been normalized by the Baseline rather than showing the magnitudes that were calculated. It is interesting as well that there is no clear trend between derate and moments for the virtual turbine. Finally, the likely fatigue loading can be compared among the cases by using the turbulence intensity 385 (T I) as a proxy as it is the cause of most fatigue (Kanev et al., 2018). This is done by averaging the streamwise turbulence intensity across a rotor-swept area at 5D downstream as shown in Fig. 19d    that there is a price to pay in added fatigue loading. The cost of an increase in combined power output of the two turbines considered may be an increase in loads in both.

390
Within the limitations of the simulations, these results suggest that operating an upstream turbine at a small derate and maximized thrust could increase the total power of it and a downstream turbine by over 10% of a Baseline, conventional, operation with streamwise turbine spacing of about 4D to at least 6D. While many previous studies have identified minimum thrust as the optimal method (Mirzaei et al., 2015;Pedersen and Larsen, 2020;Schaak, 2006;Vitulli et al., 2019), several have also found that overinduction (increasing thrust) produces higher total power than underinduction Meyers, 2016, 2018;395 Cossu, 2021;Dilip and Porté-Agel, 2017). The latter studies, however, use very different methods both in the details of the derating and the simulations so it is difficult to compare their results to these.
It remains to be seen if the positive effects of high thrust can be carried on farther downstream and/or across additional rows of turbines, or if this would require something like a cascade of derates with higher derates upstream decreasing to conventional operation at the last turbine downstream as has been previously suggested (Bitar and Seiler, 2013;Corten and Schaak, 2004).

400
Though the domain does not extend as far as typical transverse turbine spacing, Fig. 18 suggests that this method works best when the turbines are aligned, which is contrary to conventional siting in which turbines are offset with respect to the dominant wind direction. This is, however, in agreement with when wake mitigation techniques such as AIC are most relevant, i.e., when turbines are aligned with the wind and subject each other to their wakes. It should also be noted that only the maxCt cases have peak combined power outputs around z/R = 0 while all other cases appear to have peaks at least 1.3D away from z/R = 0, 405 which is in agreement with results from Annoni et al. (2015). Indeed, Fig. 16 reinforces the importance of sufficient wake mixing when turbines are aligned such that the excess energy generated by derating can penetrate the wake and be captured by the downstream turbine.
The mechanism by which this improvement works appears to be the effective blockage that is created by higher turbine thrust forces. The high thrust rapidly slows the flow passing through the rotor, and, due to conservation, the volume it occupies 410 expands. Likewise, this forces the ambient flow to accelerate around the near-wake setting up a high gradient and rapid increase in Reynolds shear stress. This in turn aids in mixing between the ambient flow and wake facilitating faster recovery. It also appears that this process is intimately connected to the decay of the tip vortices. Their decay is likely aided by the turbulence that is generated by the high gradient across the near-wake, which in turn allows for greater mixing into the wake. Finally, there is an apparent connection between the degree of derate and the increase in thrust. The 10maxCt case was the only one of the 415 maxCt cases whose wake recovery could compensate for the derate. The precise balance between these variables and whether or not thrust needs to be maximized or just increased should be explored in future studies.
The verification study suggests that cases with higher thrust may benefit from simulating more revolutions for better convergence, though the average residuals of the energy balance remained small. The virtual rotor method is also imperfect, in part because it does not capture the feedback between the turbine and its inflow through its induction zone. The biggest limitation