On the power and control of a misaligned rotor – beyond the cosine law

We present a new model to estimate the performance of a wind turbine operating in misaligned conditions. The model is based on the classic momentum and lifting-line theories, considering a misaligned rotor as a lifting wing of finite span, and accounts for the combined effects of both yaw and uptilt angles.

Improving on the classical empirical cosine law in widespread use, the new model reveals the dependency of power not only on the misalignment angle, but also on some rotor design parameters and – crucially – on the way a rotor is governed when it is yawed out of the wind. We show how the model can be readily integrated with arbitrary control laws below, above, and around the rated wind speed. Additionally, the model also shows that a sheared inflow is responsible for the observed lack of symmetry for positive and negative misalignment angles. Notwithstanding its simplicity and insignificant computational cost, the new proposed approach is in excellent agreement with large eddy simulations (LESs) and wind tunnel experiments.

Building on the new model, we derive the optimal control strategy for maximizing power on a misaligned rotor. Additionally, we maximize the total power of a cluster of two turbines by wake steering, improving on the solution based on the cosine law.

1 Introduction

Wind farm control by wake steering consists of deflecting the wake away from downstream rotors to boost the total power of a plant (Meyers et al., 2022). The effectiveness of this control strategy has been proven numerically (Jiménez et al., 2010), experimentally in the wind tunnel (Campagnolo et al., 2016), and in field tests (Fleming et al., 2019; Doekemeijer et al., 2021). At the core of the power-boosting ability of wake steering is a trade-off: on the one hand, there is an enhanced momentum of the inflow at a downstream turbine when a wake is shifted laterally away from it; on the other hand, some power is lost at the upstream misaligned rotor, because it does not point into the wind anymore. In general the trade-off budget is positive, in the sense that the power that is gained downstream is larger than the power upstream. The problem is, however, highly complex: downstream, power capture is determined by the interaction of the impinged rotor with the wake that, in turn, is influenced by the ambient conditions and those of the wake-shedding turbine; upstream, power losses depend on the inflow characteristics but also on the rotor and on the way it is governed. Understanding and controlling this delicate balance between upstream and downstream behavior is clearly of paramount importance for improving the power capture of wind farms by wake steering. Great progress has been made in recent years to understand, model, and control wakes (see, for example, the review in Meyers et al., 2022), i.e., on the downstream aspect of the problem. However, the upstream aspect remains much less explored and understood. How much power does a yawed turbine really lose? And what inflow, rotor, and rotor-control parameters influence its behavior? It is a major aim of this paper to try and answer these questions.

The aerodynamic power P of a wind turbine is customarily written as $P=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }A{u}_{\mathrm{\infty }}^{\mathrm{3}}{C}_{\mathrm{P}}$, where ρ is the air density, A is the rotor swept area, u_∞ is the ambient free-stream wind speed, and C_P is the power coefficient. When a turbine is misaligned with respect to the wind vector by an angle γ, the rotor-orthogonal velocity component becomes u_∞cos γ. Accordingly, one would expect the yaw-induced power loss to be ${\mathit{\eta }}_{\mathrm{P}}=P/P_{\mathrm{0}}={\cos }^{\mathrm{3}}\mathit{\gamma }$, where P₀ is the aerodynamic power produced for γ=0. Unfortunately, this is only a naive interpretation of the true behavior of a misaligned rotor, and its predictions are not confirmed by experimental and numerical observations (Liew et al., 2020). To reflect this fact, a pragmatic solution has been adopted by most of the literature, where power losses due to misalignment are assumed to obey the simple law ${\mathit{\eta }}_{\mathrm{P}}\approx {\cos }^{{\mathrm{p}}_{\mathrm{p}}}\mathit{\gamma }$, where p_p is a tunable parameter.

Unsurprisingly, since such a model is not based on actual physics, a large spread of values for p_p has been reported in the literature. In wind tunnel experiments with scaled models, Campagnolo et al. (2020) measured p_p=2.1, Krogstad and Adaramola (2012) and Bartl et al. (2018) reported p_p≈3, and Medici (2005) found a value p_p=2. Numerically, Fleming et al. (2015) measured p_p=1.88 on the NREL 5 MW wind turbine (Jonkman et al., 2009), whereas Draper et al. (2018) obtained values between 1.3 and 2.5 for scaled wind turbine models operating in waked inflow conditions. The power production in misaligned conditions has also been measured in multiple field tests. For example, Fleming et al. (2017) reported a value of 1.41 for an Envision 4 MW turbine, and Dahlberg and Montgomerie (2005) published a range of values for p_p between 1.9 and 5.1 at an offshore plant. More recently, Hulsman et al. (2022) observed $\mathrm{2}<p_{\mathrm{p}}<\mathrm{2.5}$ at an onshore wind farm in the north of Germany.

The large scatter characterizing the p_p coefficients reported in the literature is a relevant source of uncertainty, creating a significant hindrance to the development of wind farm control strategies and suggesting that some relevant phenomena are not captured by the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law. In hindsight, this is to be expected, because this simple model fails to explicitly represent how the power coefficient C_P changes when a turbine is misaligned, and it somehow absorbs this effect into the tunable exponent. Some indications that there is more to this problem than a simple power cosine law have already been reported by various authors. Based on experiments and numerical simulations, Campagnolo et al. (2023), Cossu (2021a, b), and Heck et al. (2023) suggested that the power of a misaligned rotor strongly depends on its loading, in the form of the thrust coefficient C_T; clearly, in turn this has a strong effect on the behavior of the wake (Cossu, 2021a, b). Other variables that have been shown to play a role on power losses are related to the inflow. Recently, Draper et al. (2018) and Liew et al. (2020) have observed that power losses in misaligned conditions differ depending on whether a rotor is waked or not. Howland et al. (2020) observed a significant influence of shear and veer, while Simley et al. (2021) measured a strong dependency on inflow speed. The behavior of power losses has also been shown to depend on the direction of yaw misalignment and not only on its magnitude as implied by the power cosine law. This asymmetric behavior of yaw misalignment has been observed by Fleming et al. (2015), Schottler et al. (2017), Fleming et al. (2018), and Campagnolo et al. (2020), among others. However, an agreement on which misalignment direction yields more or less power has not been reached yet.

In this paper, we present a new analytical model for misaligned wind turbine rotors. The proposed approach combines the classic momentum and lifting-line theories, considering a misaligned rotor as a lifting wing of finite span, in close parallel to the analysis conducted for helicopter rotors in forward flight (Johnson, 1995). While existing approaches do not explain the lack of symmetry with respect to yaw direction, the present model includes the effects of wind shear, which is shown to be the culprit for the observed break of symmetry with respect to the misalignment direction. For improved accuracy, the model also includes the effects of the uptilt angle, as it contributes to the overall misalignment of the rotor with respect to the wind vector. Thanks to this feature, the proposed methodology is also readily applicable to vertical wake-steering control, which could be implemented with floating wind turbines (Nanos et al., 2022) or downwind teetering rotors. The resulting model equations are integrated over the blade span and averaged over one rotor revolution, leading to a semi-analytical formulation of negligible computational cost that can be readily coupled with engineering wake models such as FLORIS (NREL, 2023b) or PyWake (Pedersen et al., 2019). However, the model governing equations could also be integrated numerically and embedded into blade element momentum (BEM) codes (Hansen, 2015), such as the AeroDyn package (NREL, 2023a) implemented in OpenFAST (NREL, 2023c).

Very recently, Heck et al. (2023) published a misaligned rotor model based on similar arguments. However, their approach does not include the effects of shear and therefore fails to capture the asymmetric behavior of yaw direction. More importantly, their formulation uses a modified thrust coefficient $C_{\mathrm{T}}^{\prime }$, which is assumed to remain constant between aligned and misaligned conditions. This hypothesis is indeed verified when the turbine operates in the partial-load region. Departing from this approach, the method proposed here is based on a completely general dependency of the thrust coefficient on the misalignment angle and therefore can readily accommodate arbitrary regulation strategies in the partial, full, and intermediate regulation regions, including thrust clipping and derating (Campagnolo et al., 2023). Moreover, the model of Heck et al. (2023) cannot predict power losses higher than cos ³γ, which have been reported in the literature. A detailed comparison of the new proposed model and the one of Heck et al. (2023) is developed in the following.

The proposed semi-analytical model shows that the behavior of a misaligned rotor does not follow the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law, contradicting this empirical formula in widespread use. Additionally, the new model clarifies the behavior of power capture with respect to some rotor design parameters and – even more importantly – with respect to the way a rotor is governed when it is misaligned. This is an effect that has been neglected so far, but that – as already noted by Howland et al. (2020) – most probably explains the large scatter observed by various authors. Building on the unique ability of the proposed method to handle arbitrary control policies, we derive the optimal strategy for maximizing power capture when pointing a rotor away from the wind. Finally, we implement the semi-analytical model in FLORIS, and we optimize the power of a cluster of two turbines. We obtain set points that differ from those that can be computed with the empirical ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law and that lead to a slight improvement of the cluster power.

The new models exhibits an excellent match with high-fidelity LESs obtained with a TUM-modified version of NREL's large eddy simulator actuator line model (LES-ALM) SOWFA (Fleming et al., 2014; Wang et al., 2019, 2018). Additionally, the model is further validated with wind tunnel data from experiments conducted with the TUM G1 scaled wind turbines (Bottasso and Campagnolo, 2022; Campagnolo et al., 2020).

The paper is organized as follows: Sect. 2 presents the new formulation, and Sect. 3 explains its implementation in an engineering wake model, including the integration with arbitrary control strategies. Next, Sect. 4 considers its validation with respect to simulated and experimental data, while Sect. 5 analyzes the effects of the new model on wake steering. Finally, Sect. 6 draws conclusions and offers an outlook towards future work.

2 Misalignment model

2.1 Frames of reference

Three reference frames are necessary to completely characterize a misaligned rotor interacting with the wind, as shown in Fig. 1: a ground-fixed reference frame and a nacelle-fixed reference frame, which together describe the relative orientation of the rotor with respect to the ground, and a wake-deflection intrinsic frame, which describes the relative orientation of the rotor with respect to the incoming wind vector.

Figure 1 Reference frames used in the derivation of the model. Ground-fixed wind-aligned reference frame (in red, subscript “g”) (a); intermediate frame, obtained by the tilt rotation δ about the horizontal axis y_g (in orange, subscript “t”) (b); nacelle-fixed reference frame, obtained from t by a yaw rotation γ about the vertical axis z_g (in green, subscript “n”) (c); plane Π formed by the wind vector u_∞ and the rotor axis x_n; the plane contains the x_d and y_d unit vectors of the wake-deflection intrinsic frame (in light blue, subscript “d”) (d); rotor plane Ψ, formed by the y_d and z_d unit vectors, with rotor azimuthal angle ψ and radial position $r/R$ (e). Rotations are positive according to the right-hand rule; notice that the value of the uptilt angle of the rotor in panel (b) is therefore negative.

The ground-fixed wind-aligned frame of reference is indicated with a subscript “g” and is defined by the right-handed triad of unit vectors ${\mathcal{F}}_{\mathrm{g}}=\mathit{\{}{\mathit{x}}_{\mathrm{g}},{\mathit{y}}_{\mathrm{g}},{\mathit{z}}_{\mathrm{g}}\mathit{\}}$. z_g points vertically down towards the ground; x_g is parallel to the terrain pointing downstream and is contained in the plane formed by the wind vector u_∞ and z_g; finally, y_g completes a right-handed triad. In the following, for simplicity we consider the wind vector to be parallel to the terrain, i.e., u_∞∥x_g, although this is not strictly necessary.

The nacelle-fixed frame of reference is indicated with a subscript “n” and is defined by the triad of unit vectors ${\mathcal{F}}_{\mathrm{n}}=\mathit{\{}{\mathit{x}}_{\mathrm{n}},{\mathit{y}}_{\mathrm{n}},{\mathit{z}}_{\mathrm{n}}\mathit{\}}$. ℱ_n is obtained from ℱ_g by two successive rotations: a first rotation by the tilt angle δ about the horizontal axis y_g, followed by a second rotation by the yaw angle γ about the vertical axis z_g. Both rotations are positive about their respective axes according to the right-hand rule (notice that, according to this definition, the typical uptilt of an upwind turbine results in a negative value for δ).

However, the interaction of the rotor with the flow depends only on their mutual orientation and not on how they are oriented with respect to the ground, which is a fundamental principle of fluid mechanics known as Galilean relativity. Therefore, a third frame is necessary, which is termed here wake-deflection intrinsic frame and is indicated with a subscript “d”. The frame is formed by a right-handed triad of unit vectors ${\mathcal{F}}_{\mathrm{d}}=\mathit{\{}{\mathit{x}}_{\mathrm{d}},{\mathit{y}}_{\mathrm{d}},{\mathit{z}}_{\mathrm{d}}\mathit{\}}$. Vector x_d is parallel to the rotor axis, i.e., x_d=x_n, while vectors y_d and z_d are contained in the rotor disk plane Ψ. Together, the rotor axis x_d and the wind velocity vector u_∞ define the Π plane. The angle in the Π plane between these two vectors is the true misalignment angle μ:

$$\begin{array}{}\text{(1)}& \cos \mathit{\mu }={\displaystyle \frac{{\mathit{u}}_{\mathrm{\infty }}}{u_{\mathrm{\infty }}}}\cdot {\mathit{x}}_{\mathrm{n}},\end{array}$$

where $u_{\mathrm{\infty }}=\left|{\mathit{u}}_{\mathrm{\infty }}\right|$ is the scalar ambient wind speed. The unit vector z_d is orthogonal to the Π plane, i.e., ${\mathit{z}}_{\mathrm{d}}={\mathit{x}}_{\mathrm{d}}\times {\mathit{u}}_{\mathrm{\infty }}/\left(u_{\mathrm{\infty }}\sin \mathit{\mu }\right)$, while unit vector y_d is finally chosen to form a right-handed triad. Using the coordinate transformations in Appendix A, it can be readily shown that cos μ=cos δcos γ; i.e., the total misalignment is caused by both the tilt and yaw angles, with the former typically being neglected in most wake models. Notice that, given its definition, the misalignment angle μ is always positive, because z_d flips from one side of the Π plane to the other, depending on the relative orientation of the wind velocity and rotor axis vectors. When the wind comes from the right looking upstream in the Π plane, z_d points downwards (see Fig. 1d), whereas it points upwards when the wind comes from the left.

Figure 2 shows a visualization of the wakes developing behind a wind turbine rotor for two different pairs of tilt and yaw values: δ=0°, $\mathit{\gamma }=-\mathrm{30}$° and $\mathit{\delta }=-\mathrm{28.43}$°, $\mathit{\gamma }=-\mathrm{10}$°. Both pairs correspond to a same true misalignment angle μ=30°. The figure confirms that the wake is invariant for an observer on the ℱ_d frame. This is particularly evident in the images of the longitudinal speed on the Π plane (marked with a solid black border), which are identical in the two cases. Clearly, for large values of tilt the interaction of the wake with the ground or with a sheared inflow would break the Π-frame invariance.

Figure 2 Visualization of the wakes developing behind a wind turbine operating in steady inflow conditions for different pairs of tilt and yaw values, all corresponding to a same total (true) misalignment μ=30°. The deflection of the wake occurs in the Π plane, marked with a black border. Top row (a): δ=0°, $\mathit{\gamma }=-\mathrm{30}$°; bottom row (b): $\mathit{\delta }=-\mathrm{28.43}$°, $\mathit{\gamma }=-\mathrm{10}$°. Left column: isosurfaces of Q criterion; right column: image of the longitudinal flow speed $u/u_{\mathrm{\infty }}$ on the Π plane. Distances are expressed in rotor diameters D. Interactive 3D versions of the figures are available at the following links: https://tinyurl.com/btcl-fig-2-a (last access: 17 July 2024) and https://tinyurl.com/btcl-fig-2-b (last access: 17 July 2024) (b).

As a result, in the following the wake analysis is developed in the Π plane instead of the horizontal one as customarily done. Transformation matrices that map vector components from one frame to the other are reported in Appendix A.

2.2 Sheared inflow

Considering a linear vertical shear of the inflow, the ambient wind speed is written as

$$\begin{array}{}\text{(2)}& u_{\mathrm{\infty }}\left(z_{\mathrm{g}}\right)=u_{\mathrm{\infty },\mathrm{hub}}\left(\mathrm{1}-k{\displaystyle \frac{z_{\mathrm{g}}}{R}}\right).\end{array}$$

Here u_∞,hub is the ambient wind speed at hub height; k is the vertical linear shear coefficient; z_g is the vertical coordinate in the ground frame of reference, centered at the hub; and R is the rotor radius. The choice of a linear shear distribution was made just to simplify the derivations, and other choices are clearly possible, for example, to model the more common power law or the presence of low-level jets. Additionally, it would be interesting to also include the effects of a horizontal shear, to account for waked conditions, and of veer. These further model improvements are, however, deferred to a continuation of this study.

By applying the coordinate transformation of Appendix A, the ambient wind speed of Eq. (2) can be written in terms of the radial r and azimuthal ψ coordinates on the rotor plane, yielding

$$\begin{array}{}\text{(3)}& \begin{array}{rl}{u}_{\mathrm{\infty }}(r,\mathit{\psi },\mathit{\delta },\mathit{\gamma })& =u_{\mathrm{\infty },\mathrm{hub}}(\mathrm{1}-k{\displaystyle \frac{r}{R}}{\displaystyle \frac{\cos \mathit{\delta }}{\sin \mathit{\mu }}}(\sin \mathit{\gamma }\cos \mathit{\psi }\\ & -\cos \mathit{\gamma }\sin \mathit{\delta }\sin \mathit{\psi }\left)\right).\end{array}\end{array}$$

Here ψ is positive about x_d according to the right-hand rule (i.e., clockwise looking downstream), and it is measured starting from the z_d unit vector (which flips from one side of the Π plane to other depending on whether the wind blows from the right or left looking upstream, as explained in Sect. 2.1; see also Fig. 1).

2.3 Force and velocity components at a blade section

With reference to Fig. 3, the tangential F_t and normal F_n components of the aerodynamic force at a blade section are

$$\begin{array}{}\text{(4a)}& {\displaystyle }& {\displaystyle }{F}_{\mathrm{t}}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{u}^{\mathrm{2}}c\left(C_{\mathrm{L}}\sin \mathit{\phi }-C_{\mathrm{D}}\cos \mathit{\phi }\right),\text{(4b)}& {\displaystyle }& {\displaystyle }{F}_{\mathrm{n}}={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{u}^{\mathrm{2}}c\left(C_{\mathrm{L}}\cos \mathit{\phi }+C_{\mathrm{D}}\sin \mathit{\phi }\right),\end{array}$$

where $\mathit{\phi }={\tan }^{-\mathrm{1}}{u}_{\mathrm{n}}/u_{\mathrm{t}}$ is the inflow angle, $u=\sqrt{u_{\mathrm{t}}^{\mathrm{2}}+u_{\mathrm{n}}^{\mathrm{2}}}$ is the total flow speed at the blade section, c is the sectional chord length, and finally C_L and C_D are the lift and drag coefficients, respectively. Using the coordinate transformations of Appendix A, the tangential u_t and normal u_n velocity components are written as

$$\begin{array}{}\text{(5a)}& {\displaystyle }& {\displaystyle }{u}_{\mathrm{t}}=\mathrm{\Omega }r+u_{\mathrm{\infty }}\sin \mathit{\mu }\cos \mathit{\psi },\text{(5b)}& {\displaystyle }& {\displaystyle }{u}_{\mathrm{n}}=u_{\mathrm{\infty }}\cos \mathit{\mu }(\mathrm{1}-a),\end{array}$$

where Ω is the angular speed of the rotor and a is the axial induction factor, which expresses how much the rotor-orthogonal component of the free-stream speed u_∞ is slowed down at the rotor disk.

Figure 3 Blade cross section, with triangle of velocities (in blue), lift and drag (in light blue), and resulting aerodynamic force components (in red).

2.4 Induction model

It is well known that a non-uniform description of the induction is necessary in order to accurately capture the azimuthal variation of loads on a rotor operating in non-axial conditions (Johnson, 1995). However, it appears that this is not necessary when computing integral rotor quantities such as power, torque, and thrust, as in the present case. To show this, the induction is modeled here with an expansion limited to one-per-revolution (1P) harmonics, i.e.,

$$\begin{array}{}\text{(6)}& a=a_{\mathrm{0}}\left(\mathrm{1}+{\mathit{\kappa }}_{\mathrm{1}\mathrm{s}}{\displaystyle \frac{r}{R}}\sin \mathit{\psi }+{\mathit{\kappa }}_{\mathrm{1}\mathrm{c}}{\displaystyle \frac{r}{R}}\cos {\mathit{\psi }}_{\mathrm{g}}\right),\end{array}$$

where a₀ is the constant-over-the-rotor (0P) induction, while a₀κ_1s and a₀κ_1c are the 1P sine and cosine harmonic amplitudes, respectively.

Following the classical approach used for helicopter rotors in forward flight (Johnson, 1995), the sine term accounts for the tilting of the induction plane caused by the misalignment μ of the rotor with the incoming wind. As such, it is written in terms of ψ, which is measured starting from the z_d unit vector, and therefore it expresses a rotation of the induction plane about the axis normal to the wake-deflection intrinsic frame Π. The coefficient κ_1s can be modeled according to Coleman et al. (1945) and Pitt and Peters (1981), resulting in the expression

$$\begin{array}{}\text{(7)}& {\mathit{\kappa }}_{\mathrm{1}\mathrm{s}}=-{\displaystyle \frac{\mathrm{15}\mathit{\pi }}{\mathrm{32}}}\tan \left({\displaystyle \frac{\mathit{\chi }}{\mathrm{2}}}\right),\end{array}$$

where the initial wake skew angle is $\mathit{\chi }=\mathit{\mu }+\sin \mathit{\mu }{C}_{\mathrm{T}}/\mathrm{2}$ (Jiménez et al., 2010), and $C_{\mathrm{T}}=\mathrm{2}T/\left(\mathit{\rho }A{u}_{\mathrm{\infty },\mathrm{hub}}^{\mathrm{2}}\right)$ is the thrust coefficient. Notice that the definition of the skew angle differs from the one given by Eq. (20) of Jiménez et al. (2010), because of the different definition of the thrust coefficient used in that publication.

The cosine term is introduced to account for the effects on the induction caused by vertical shear. As such, it is written as a function of the azimuthal angle ψ_g, which is measured from the (vertical) z_g unit vector, and therefore it expresses a rotation of the induction plane about the (horizontal) unit vector y_g. Using Eq. (A4b), it is readily found that ψ_g=ψ cos μ. Following Meyer Forsting et al. (2018), the cosine term is proportional to both the shear k and the thrust C_T coefficients, i.e.,

$$\begin{array}{}\text{(8)}& {\mathit{\kappa }}_{\mathrm{1}\mathrm{c}}={\mathit{\kappa }}_{\mathrm{1}\mathrm{c}}^{*}k{C}_{\mathrm{T}}.\end{array}$$

The cosine term significantly complicates the analytical derivations of power, torque, and thrust, which must now be expressed in terms of Bessel functions (Abramowitz et al., 1988) because of the term cos (ψcos μ). Before attempting the modeling of the proportionality coefficient ${\mathit{\kappa }}_{\mathrm{1}\mathrm{c}}^{*}$, this term was numerically optimized to best fit the numerical simulations and experimental measurements, as explained later in Sect. 3.

The inclusion of the sine and cosine induction terms has only an extremely modest effect on the quality of the results. In fact the match of C_P improves by 0.35 % when the sine term is included and by 0.60 % when both terms are used, as more precisely shown in Appendix B. Because of their modest effects, these terms are dropped from the following discussion, to simplify the resulting expressions, and they were not used in the results reported later in this article. However, these terms are retained in the software implementation of the model (Tamaro et al., 2024a) and can be switched on if desired by the user.

It should also be noted that the present model neglects the effects of the tangential induction, which in fact does not appear in the tangential velocity component expressed by Eq. (5a). This is justified by the fact that the rotor swirl is concentrated close to the hub, and it is small for a large extent of the blade span (Burton et al., 2011), where most of the thrust and power are generated.

More in general, there are several other effects that are present in a rotor and that are not modeled here, such as radial drag, tip and root losses, blade sweep, prebend and cone, and others. All these effects can be taken into account in detailed BEM models (Hansen, 2015; Burton et al., 2011) but would significantly complicate the present simplified analytical method. However, notwithstanding these limitations, the results of Sect. 4 show a remarkable ability of the proposed approach in predicting the trends of power and thrust as functions of various operating and inflow conditions. Additionally, the ability of the model in predicting actual power and thrust values (instead of trends) can be improved by the use of loss functions, as explained in Sect. 3.2.

2.5 Streamtube model

An expression for the axial induction can be derived using the concept of a streamtube (Hansen, 2015), as shown in Fig. 4 with reference to the present case of a misaligned rotor. Four stations are considered along the stream tube: inlet i, outlet o, section r⁻ located immediately in front of the rotor, and section r⁺ located immediately behind the rotor.

Figure 4 Schematic view of a streamtube around a misaligned wind turbine. Cross-sectional stations: inlet i, outlet o, r⁻ immediately in front of the rotor, and r⁺ immediately behind the rotor.

The principle of impulse and momentum applied to the streamtube is written as

$$\begin{array}{}\text{(9)}& T\cos \mathit{\mu }=\dot{m}\left(u_{\mathrm{\infty }}-u_{\mathrm{o}}\right),\end{array}$$

where T is the thrust force; u_o is the longitudinal flow speed at the streamtube outlet; and $\dot{m}$ is the mass flux,

$$\begin{array}{}\text{(10)}& \dot{m}=\mathit{\rho }A{u}_{\mathrm{n}}.\end{array}$$

By using the thrust coefficient, Eq. (9) yields the non-dimensional longitudinal flow speed at the streamtube outlet:

$$\begin{array}{}\text{(11)}& {\displaystyle \frac{u_{\mathrm{o}}}{u_{\mathrm{\infty },\mathrm{hub}}}}=\mathrm{1}-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \frac{C_{\mathrm{T}}}{\mathrm{1}-a_{\mathrm{0}}}}.\end{array}$$

Next, Bernoulli's energy conservation theorem is applied between the streamtube inlet and the section immediately upstream of the rotor (stations i and r⁻ in Fig. 4) and between the section immediately downstream of the rotor and the streamtube outlet (stations r⁺ and o in the same figure):

$$\begin{array}{}\text{(12a)}& {\displaystyle }& {\displaystyle }{p}_{\mathrm{i}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{V}_{\mathrm{i}}^{\mathrm{2}}=p_{{\mathrm{r}}^{-}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{V}_{{\mathrm{r}}^{-}}^{\mathrm{2}},\text{(12b)}& {\displaystyle }& {\displaystyle }{p}_{{\mathrm{r}}^{+}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{V}_{{\mathrm{r}}^{+}}^{\mathrm{2}}=p_{\mathrm{o}}+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{V}_{\mathrm{o}}^{\mathrm{2}},\end{array}$$

where p is pressure and $p_{\mathrm{i}}=p_{\mathrm{o}}=p_{\mathrm{\infty }}$, while p_∞ is the ambient value. Additionally, $V_{{\mathrm{r}}^{-}}=V_{{\mathrm{r}}^{+}}$ for continuity, and furthermore $V_{\mathrm{o}}^{\mathrm{2}}=u_{\mathrm{o}}^{\mathrm{2}}+v_{\mathrm{o}}^{\mathrm{2}}$ at the outlet section, where v_o is the lateral (sidewash) speed component.

Following a customarily textbook assumption used for helicopters' rotors in forward flight (Johnson, 1995), only more recently adopted also for wind turbines by Shapiro et al. (2018), the misaligned rotor can be seen as a lifting wing of finite span (albeit of a small aspect ratio $\mathrm{AR}=D^{\mathrm{2}}/A=\mathrm{4}/\mathit{\pi }$) operating at an angle of attack μ¹. The chord C(y_d) of the wing in the streamwise direction has an elliptic distribution: $\left(C\right(y_{\mathrm{d}})/\mathrm{2}{)}^{\mathrm{2}}+y_{\mathrm{d}}^{\mathrm{2}}=R^{\mathrm{2}}$. According to Prandtl's lifting-line theory (Tietjens and Prandtl, 1957; Katz and Plotkin, 2001), the wing has consequently an elliptic lift distribution, which induces a spanwise-constant downwash (in this case, sidewash) $v_{\mathrm{o}}=\mathrm{\Gamma }/\mathrm{4}R$. $\mathrm{\Gamma }=\stackrel{\mathrm{‾}}{L}/\mathit{\rho }{u}_{\mathrm{\infty },\mathrm{hub}}$ is the circulation at the wing midsection y_d=0, and $\stackrel{\mathrm{‾}}{L}=LC\left(\mathrm{0}\right)/A$ is the lift per unit span at that same location. Since the wing lift is the rotor side force, i.e., L=Tsin μ, it follows that the non-dimensional sidewash at the streamtube outlet can be expressed as (Heck et al., 2023)

$$\begin{array}{}\text{(13)}& {\displaystyle \frac{v_{\mathrm{o}}}{u_{\mathrm{\infty },\mathrm{hub}}}}={\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{C}_{\mathrm{T}}\sin \mathit{\mu }.\end{array}$$

Combining the previous equations yields an expression for the 0P axial induction a₀ as a function of the misalignment μ and thrust coefficient C_T:

$$\begin{array}{}\text{(14)}& \mathrm{1}-a_{\mathrm{0}}={\displaystyle \frac{\mathrm{1}+\sqrt{\mathrm{1}-C_{\mathrm{T}}-\frac{\mathrm{1}}{\mathrm{16}}{C}_{\mathrm{T}}^{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{\mu }}}{\mathrm{2}\left(\mathrm{1}+\frac{\mathrm{1}}{\mathrm{16}}{C}_{\mathrm{T}}{\sin }^{\mathrm{2}}\mathit{\mu }\right)}}.\end{array}$$

2.6 Thrust force

Equation (14) furnishes an expression for the 0P axial induction as a function of the thrust coefficient. To close the problem, an expression for the thrust coefficient in terms of the operating conditions of the turbine is necessary. To this end, the thrust force T is expressed in terms of the normal sectional force F_n as

$$\begin{array}{}\text{(15)}& T={\displaystyle \frac{B}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathrm{R}}{\int }}{F}_{\mathrm{n}}\mathrm{d}\mathit{\psi }\mathrm{d}r,\end{array}$$

where B indicates the number of blades. Using Eq. (4b) under the assumption of a small inflow angle (i.e., sin φ≈φ and cos φ≈1) yields

$$\begin{array}{}\text{(16)}& T={\displaystyle \frac{B}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathrm{R}}{\int }}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{u}^{\mathrm{2}}c\left(C_{\mathrm{D}}\mathit{\phi }+C_{\mathrm{L}}\right)\mathrm{d}\mathit{\psi }\mathrm{d}r.\end{array}$$

The lift coefficient is written as $C_{\mathrm{L}}=C_{\mathrm{L},\mathit{\alpha }}\mathit{\alpha }=C_{\mathrm{L},\mathit{\alpha }}({\mathit{\alpha }}_{\mathrm{g}}-{\mathit{\alpha }}_{\mathrm{0}})$, where C_L,α is the lift slope, and $\mathit{\alpha }=\mathit{\phi }-\mathit{\theta }$ is the angle of attack measured with respect to the zero-lift direction, whereas α_g is the angle of attack measured with respect to a generic direction. Without any loss of generality, the use of α is preferred in the following to avoid carrying along in the derivations the unnecessary extra term α₀. Furthermore, $\mathit{\theta }={\mathit{\theta }}_{\mathrm{p}}+\mathit{\beta }$ is the local pitch angle (see Fig. 3), where θ_p is the blade pitch rotation at the pitch bearing, and β indicates the blade twist referred to the zero-lift direction.

Neglecting swirl induction, the inflow angle is $\tan \mathit{\phi }=(\mathrm{1}-a_{\mathrm{0}})/(\mathit{\lambda }r/R)$, where $\mathit{\lambda }=\mathrm{\Omega }R/u_{\mathrm{\infty },\mathrm{hub}}$ is the tip speed ratio. For a turbine operating close to optimal induction (i.e., $a_{\mathrm{0}}=\mathrm{1}/\mathrm{3}$) and a typical tip speed ratio of 8.5, the inflow angle at three-quarter span is less than 6^∘, justifying the small-angle assumption. This assumption clearly becomes less accurate for small inductions and tip speed ratios or close to the blade root, where only a modest contribution to the thrust is generated.

Using again the small inflow angle assumption, it follows that $\mathit{\phi }\approx u_{\mathrm{n}}/u_{\mathrm{t}}$ and u≈u_t, and the thrust T becomes

$$\begin{array}{}\text{(17)}& T={\displaystyle \frac{B}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathrm{R}}{\int }}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{u}_{\mathrm{t}}^{\mathrm{2}}c\left(C_{\mathrm{D}}{\displaystyle \frac{u_{\mathrm{n}}}{u_{\mathrm{t}}}}+C_{\mathrm{L},\mathit{\alpha }}\left({\displaystyle \frac{u_{\mathrm{n}}}{u_{\mathrm{t}}}}-\mathit{\theta }\right)\right)\mathrm{d}\mathit{\psi }\mathrm{d}r.\end{array}$$

Using Eqs. (3) and (5), solving the double integral, and expressing T through the thrust coefficient $C_{\mathrm{T}}=C_{{\mathrm{T}}_{\mathrm{1}}}+C_{{\mathrm{T}}_{\mathrm{2}}}$ finally gives

$$\begin{array}{}\text{(18a)}& {\displaystyle }\begin{array}{rl}& C_{{\mathrm{T}}_{\mathrm{1}}}={\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}\left(C_{\mathrm{D}}+C_{\mathrm{L},\mathit{\alpha }}\right)\cos \mathit{\mu }\left(\mathit{\lambda }-k\cos \mathit{\delta }\sin \mathit{\gamma }\right)\left(\mathrm{1}-a_{\mathrm{0}}\right),\end{array}\text{(18b)}& {\displaystyle }\begin{array}{rl}{C}_{{\mathrm{T}}_{\mathrm{2}}}& =-{\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}{C}_{\mathrm{L},\mathit{\alpha }}\mathit{\theta }({\sin }^{\mathrm{2}}\mathit{\mu }+{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}{\mathit{\lambda }}^{\mathrm{2}}-{\displaystyle \frac{k\cos \mathit{\delta }}{\mathrm{12}}}(\mathrm{8}\mathit{\lambda }\sin \mathit{\gamma }\\ & -k\cos \mathit{\delta }\left({\cos }^{\mathrm{2}}\mathit{\gamma }{\sin }^{\mathrm{2}}\mathit{\delta }+\mathrm{3}{\sin }^{\mathrm{2}}\mathit{\gamma }\right)\left)\right),\end{array}\end{array}$$

where $\mathit{\sigma }=BcR/A$ is the rotor solidity.

For null shear, i.e., k=0, this expression simplifies to

$$\begin{array}{}\text{(19)}& \begin{array}{rl}{C}_{\mathrm{T}}& ={\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}(\left(C_{\mathrm{D}}+C_{\mathrm{L},\mathit{\alpha }}\right)\cos \mathit{\mu }\left(\left(\mathrm{1}-a_{\mathrm{0}}\right)\mathit{\lambda }\right)\\ & -C_{\mathrm{L},\mathit{\alpha }}\mathit{\theta }\left({\sin }^{\mathrm{2}}\mathit{\mu }+{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}{\mathit{\lambda }}^{\mathrm{2}}\right)).\end{array}\end{array}$$

Notice that the terms in Eq. (18) depending on shear k also depend on the angles γ and δ, whereas Eq. (19) only depends on the total misalignment angle μ. This is because the wind shear is defined with respect to the ground frame, which is mapped into the nacelle frame by the γ and δ angles, whereas μ only depends on the relative orientation of the wind vector with the rotor axis, as explained in Sect. 2.1.

Furthermore, we note that – differently from the approach of Heck et al. (2023) – these expressions for thrust are applicable to any desired control policy, as they depend explicitly on the tip speed ratio λ and the pitch setting θ.

2.7 Power

The aerodynamic power P generated by a wind turbine is P=QΩ, where the aerodynamic torque Q is written as

$$\begin{array}{}\text{(20)}& Q={\displaystyle \frac{B}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathrm{R}}{\int }}{F}_{\mathrm{t}}r\mathrm{d}\mathit{\psi }\mathrm{d}r.\end{array}$$

Considering small angles, power can be written as

$$\begin{array}{}\text{(21)}& P={\displaystyle \frac{B}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\underset{\mathrm{0}}{\overset{\mathrm{R}}{\int }}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{u}_{\mathrm{t}}^{\mathrm{2}}c\left(-C_{\mathrm{D}}+C_{\mathrm{L},\mathit{\alpha }}\left({\displaystyle \frac{u_{\mathrm{n}}}{u_{\mathrm{t}}}}-\mathit{\theta }\right){\displaystyle \frac{u_{\mathrm{n}}}{u_{\mathrm{t}}}}\right)\mathrm{\Omega }r\mathrm{d}\mathit{\psi }\mathrm{d}r.\end{array}$$

Using Eqs. (3) and (5), expressing the angular velocity as $\mathrm{\Omega }=\mathit{\lambda }{u}_{\mathrm{\infty },\mathrm{hub}}/R$, and solving the double integral yields the power coefficient C_P:

$$\begin{array}{}\text{(22)}& {\displaystyle }\begin{array}{rl}{C}_{\mathrm{P}}& ={\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}\mathit{\lambda }(C_{\mathrm{L},\mathit{\alpha }}\left(\mathrm{1}-a_{\mathrm{0}}\right)\cos \mathit{\mu }\left(\left(\mathrm{1}-a_{\mathrm{0}}\right)\cos \mathit{\mu }-{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\mathit{\lambda }\mathit{\theta }\right)\\ & -{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{C}_{\mathrm{D}}\left({\mathit{\lambda }}^{\mathrm{2}}+{\sin }^{\mathrm{2}}\mathit{\mu }\right)+k\cos \mathit{\delta }\sin \mathit{\gamma }\\ & \left({\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\cos \mathit{\mu }{C}_{\mathrm{L},\mathit{\alpha }}\mathit{\theta }\left(\mathrm{1}-a_{\mathrm{0}}\right)-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\lambda }{C}_{\mathrm{D}}\right)\\ & +{\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{k}^{\mathrm{2}}{\cos }^{\mathrm{2}}\mathit{\delta }({\cos }^{\mathrm{2}}\mathit{\mu }{C}_{\mathrm{L},\mathit{\alpha }}{\left(\mathrm{1}-a_{\mathrm{0}}\right)}^{\mathrm{2}}\\ & -{\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}{C}_{\mathrm{D}}\left({\sin }^{\mathrm{2}}\mathit{\mu }+\mathrm{2}{\sin }^{\mathrm{2}}\mathit{\gamma }\right)\left)\right).\end{array}\end{array}$$

For null shear, i.e., k=0, the power coefficient simplifies to

$$\begin{array}{}\text{(23)}& \begin{array}{rl}{C}_{\mathrm{P}}& ={\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}\mathit{\lambda }(C_{\mathrm{L},\mathit{\alpha }}\left(\mathrm{1}-a_{\mathrm{0}}\right)\cos \mathit{\mu }\left(\left(\mathrm{1}-a_{\mathrm{0}}\right)\cos \mathit{\mu }-{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\mathit{\lambda }\mathit{\theta }\right)\\ & -{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{C}_{\mathrm{D}}\left({\mathit{\lambda }}^{\mathrm{2}}+{\sin }^{\mathrm{2}}\mathit{\mu }\right)).\end{array}\end{array}$$

Due to the explicit dependency of C_P on cos μ and of 1−a₀ on sin ²μ (see Eq. 14), it follows that – in an unsheared inflow – the aerodynamic power does not depend on the misalignment direction.

2.8 Dependency on misalignment direction

The power model reveals that vertical shear is the culprit for the observed lack of symmetry with respect to yaw misalignment. In fact, the term responsible for the asymmetry in Eq. (22) is

$$\begin{array}{}\text{(24)}& C_{\mathrm{P}}^{\mathrm{asymm}}=k\cos \mathit{\delta }\sin \mathit{\gamma }\left({\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}\cos \mathit{\mu }{C}_{\mathrm{L},\mathit{\alpha }}\mathit{\theta }\left(\mathrm{1}-a_{\mathrm{0}}\right)-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\lambda }{C}_{\mathrm{D}}\right).\end{array}$$

This term shows that yawing a rotor out of the wind in a sheared inflow will produce a non-symmetric behavior with respect to positive and negative yaw angles γ, i.e., $P(-\mathit{\gamma })\ne P(+\mathit{\gamma })$. In fact, the following can be noted.

Figure 5 Normalized thrust $C_{\mathrm{T}}/C_{\mathrm{T}}(\mathit{\gamma }=\mathrm{0}\mathit{°})$ (a) and power $C_{\mathrm{P}}/C_{\mathrm{P}}(\mathit{\gamma }=\mathrm{0}\mathit{°})$ (b) coefficients, plotted as functions of the misalignment angle γ, for different shear coefficients k and pitch angles θ_p. Pitch values: colors; k=0: solid lines; k=0.3: dotted lines. The plots were generated considering the following values: λ=8.5, $C_{\mathrm{D}}=\mathrm{5.2}\times {\mathrm{10}}^{-\mathrm{3}}$, $C_{\mathrm{L},\mathit{\alpha }}=\mathrm{4.76}$, β=3.35°, σ=4.16 %, $\mathit{\delta }=-\mathrm{5}$°, and R=65 m. An interactive version of the figure that allows one to plot the thrust and power coefficients for user-defined values of the model parameters is available as a Jupyter notebook at the link https://tinyurl.com/btcl-fig-5 (last access: 12 July 2024).

• For small θ, the two terms within the parenthesis of Eq. (24) are small. Consequently, the asymmetry is small. Furthermore, either positive or negative yaw angles could produce more power, depending on the balance of these two terms.

• On the other hand, for larger θ the thrust coefficient C_T decreases, in turn increasing the term (1−a₀). As a consequence, the first of the two terms of Eq. (24) prevails. The origin of this prevailing term can be traced to the contribution to power of the fraction of lift that depends on the pitch angle, which is proportional to C_L,αθu_tu_n (see Eq. 21). When the rotor is misaligned, u_t exhibits a cos ψ variation (see Eq. 5a) that is in phase with or in opposition to a similar cos ψ variation of u_n caused by shear (see Eqs. 5b and 3). The integral of the product of these two terms over one revolution is different from zero, and it has the sign of the misalignment, resulting in some extra power for positive yaw angles and some losses for negative ones, when the thrust is low.

This complex balance of effects is probably the cause for the lack of agreement in the literature on which misalignment direction yields more. As shown by the model, there is no simple answer, and the behavior depends on the rotor design and on how it is operated.

Very similar conclusions apply also to the thrust coefficient. According to Eq. (18), $C_{{\mathrm{T}}_{\mathrm{1}}}(-\mathit{\gamma })>C_{{\mathrm{T}}_{\mathrm{1}}}(+\mathit{\gamma })$, because $C_{{\mathrm{T}}_{\mathrm{1}}}$ depends on −ksin γ, whereas $C_{{\mathrm{T}}_{\mathrm{2}}}(-\mathit{\gamma })<C_{{\mathrm{T}}_{\mathrm{2}}}(+\mathit{\gamma })$ (when θ>0), because $C_{{\mathrm{T}}_{\mathrm{2}}}$ depends on +θksin γ. Therefore, one can expect a slightly higher thrust for negative yaw angles at low pitch settings and vice versa at the higher pitch values (with the effect being more pronounced at larger tip speed ratios).

To illustrate these findings in an example case, Fig. 5 shows the thrust and power coefficients as functions of the misalignment angle γ, for different shear coefficients k and blade pitch angles θ_p. All coefficients have been normalized by their respective value in aligned conditions.

The lack of symmetry of the rotor with respect to misalignment direction is in general rather small. In a typical field implementation of wake steering, various uncertainties – due to limits in the knowledge of the ambient conditions, actual yaw orientation of the rotor, asymmetric behavior of the onboard anemometry, etc. – and other model errors probably dominate the problem, making the asymmetric behavior of misalignment a negligible effect, especially for small pitch values and moderate tip speed ratios.

3 Implementation

3.1 Integration with a wake model

The analytical model derived in the previous pages can be readily implemented in engineering flow models (NREL, 2023b; Pedersen et al., 2019). An open-source implementation in FLORIS (NREL, 2023b) is available on GitHub (Tamaro et al., 2024a). As a summary, we report here for convenience the governing Eqs. (14) and (18):

$$\begin{array}{}\text{(25a)}& {\displaystyle }\begin{array}{rl}& \mathrm{1}-a_{\mathrm{0}}={\displaystyle \frac{\mathrm{1}+\sqrt{\mathrm{1}-C_{\mathrm{T}}-\frac{\mathrm{1}}{\mathrm{16}}{C}_{\mathrm{T}}^{\mathrm{2}}{\sin }^{\mathrm{2}}\mathit{\mu }}}{\mathrm{2}\left(\mathrm{1}+\frac{\mathrm{1}}{\mathrm{16}}{C}_{\mathrm{T}}{\sin }^{\mathrm{2}}\mathit{\mu }\right)}},\end{array}\text{(25b)}& {\displaystyle }\begin{array}{rl}{C}_{\mathrm{T}}=& {\displaystyle \frac{\mathit{\sigma }}{\mathrm{2}}}(\left(C_{\mathrm{D}}+C_{\mathrm{L},\mathit{\alpha }}\right)\cos \mathit{\mu }\left(\mathit{\lambda }-\cos \mathit{\delta }\sin \mathit{\gamma }k\right)\\ & \left(\mathrm{1}-a_{\mathrm{0}}\right)-C_{\mathrm{L},\mathit{\alpha }}\mathit{\theta }({\sin }^{\mathrm{2}}\mathit{\mu }+{\displaystyle \frac{\mathrm{2}}{\mathrm{3}}}{\mathit{\lambda }}^{\mathrm{2}}-{\displaystyle \frac{k\cos \mathit{\delta }}{\mathrm{12}}}\\ & \left(\mathrm{8}\mathit{\lambda }\sin \mathit{\gamma }-k\cos \mathit{\delta }\left({\cos }^{\mathrm{2}}\mathit{\gamma }{\sin }^{\mathrm{2}}\mathit{\delta }+\mathrm{3}{\sin }^{\mathrm{2}}\mathit{\gamma }\right)\right)\left)\right).\end{array}\end{array}$$

This represents a closed system of equations that can be solved for the axial induction a₀ and thrust coefficient C_T, given the tip speed ratio λ, the pitch setting ${\mathit{\theta }}_{\mathrm{p}}=\mathit{\theta }-\mathit{\beta }$, and the yaw misalignment γ. Crucially, the presence of λ and θ_p enable using any desired control policy when misaligning the turbine. Having obtained a₀ and C_T, the power coefficient $C_{\mathrm{P}}(\mathit{\gamma },{\mathit{\theta }}_{\mathrm{p}},\mathit{\lambda })$ is obtained by Eq. (22), and, finally, the rotor power is computed as

$$\begin{array}{}\text{(26)}& P\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }A{u}_{\mathrm{\infty },\mathrm{hub}}^{\mathrm{3}}{C}_{\mathrm{P}}\left(\mathit{\gamma },{\mathit{\theta }}_{\mathrm{p}},\mathit{\lambda }\right).\end{array}$$

3.2 Improved accuracy by the use of loss functions

The analytical derivation of the equations implies that the model lacks many of the features that are present in more sophisticated BEM implementations, such as radially and azimuthally non-uniform induction, swirl induction, tip and root losses, radial flow, spanwise varying geometric characteristics, prebend, and others. Clearly, this lack of accuracy could be resolved by numerically implementing the same model in a BEM code (Hansen, 2015). However, this way the use of the misaligned rotor model in combination with an engineering wake model would become much more complex and numerically expensive.

To address this problem, the following implementation is recommended, which was used in the validation and the examples reported in the following section. A power loss function η_P is computed by using Eq. (22) to yield

$$\begin{array}{}\text{(27)}& {\mathit{\eta }}_{\mathrm{P}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)={\displaystyle \frac{C_{\mathrm{P}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)}{C_{\mathrm{P}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathrm{0}\right)}}.\end{array}$$

Next, a refined power coefficient $C_{\mathrm{P}}^{\mathrm{rf}}$ is obtained as

$$\begin{array}{}\text{(28)}& C_{\mathrm{P}}^{\mathrm{rf}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)={\mathit{\eta }}_{\mathrm{P}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)C_{\mathrm{P}}^{\mathrm{hf}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathrm{0}\right),\end{array}$$

where $C_{\mathrm{P}}^{\mathrm{hf}}$ is the power coefficient computed in aligned conditions through a higher-fidelity model, for example, based on a sophisticated BEM implementation or even on experimental measurements, when available. In other words, the analytical model is used not to predict the actual power coefficient but only the fraction of it that is lost by misalignment. The actual power coefficient is obtained by applying the loss model to a more accurate power coefficient model in aligned conditions.

The same approach is adopted for thrust. First, a thrust coefficient loss factor is computed by using the proposed model as ${\mathit{\eta }}_{\mathrm{T}}(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma })=C_{\mathrm{T}}(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma })/C_{\mathrm{T}}(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathrm{0})$. Next, a refined estimate in misaligned conditions is obtained by computing $C_{\mathrm{T}}^{\mathrm{rf}}$ from a higher-fidelity aligned value $C_{\mathrm{T}}^{\mathrm{hf}}$, i.e.,

$$\begin{array}{}\text{(29)}& C_{\mathrm{T}}^{\mathrm{rf}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)={\mathit{\eta }}_{\mathrm{T}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)C_{\mathrm{T}}^{\mathrm{hf}}\left(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathrm{0}\right).\end{array}$$

In the following, we always adopt this approach for the power and thrust coefficients. However, to simplify the notation, we drop the superscript (⋅)^rf. Hence, for example, when we write $C_{\mathrm{P}}(\mathit{\lambda },{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma })$, we in reality imply that Eq. (28) is used; the same holds for C_T.

3.3 Implementation of arbitrary control strategies

When a wind turbine yaws out of the wind, the inflow seen by the rotor changes with respect to the aligned condition. The controller reacts to the changed inflow, modifying the set point, which in turn affects the power captured by the rotor and its loading. Hence, the problem is implicit, in the sense that the misalignment model has to be solved together with the controller. This general implicit approach should be contrasted with the explicit one proposed in Heck et al. (2023), which assumes a rotor performance parameter, $C_{\mathrm{T}}^{\prime }=\mathrm{2}T/\left(\mathit{\rho }A{u}_{\mathrm{n}}^{\mathrm{2}}\right)$, remains constant even in misaligned conditions. This section explains how arbitrary control laws can be integrated with the present more general model.

The control of a modern variable-speed wind turbine is typically based on the definition of two or three main operational regions.

In region II (also called the below-rated or partial-load regime), the turbine should maximize its power output. This is achieved by operating at the maximum power coefficient $C_{\mathrm{P}}^{*}({\mathit{\lambda }}^{*},{\mathit{\theta }}_{\mathrm{p}}^{*})$, which corresponds to the optimal tip speed ratio λ^* and pitch setting ${\mathit{\theta }}_{\mathrm{p}}^{*}$. As the tip speed ratio must remain constant at its value λ^* throughout this control region, the rotor speed increases linearly with wind speed, i.e., $\mathrm{\Omega }={\mathit{\lambda }}^{*}{u}_{\mathrm{\infty }}/R$. The aerodynamic torque Q_a is readily computed as Q_a=K(ρ)Ω², with $K\left(\mathit{\rho }\right)=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }{\mathrm{AR}}^{\mathrm{3}}{C}_{\mathrm{P}}^{*}/{\mathit{\lambda }}^{{*}^{\mathrm{3}}}$. Once the aerodynamic torque is known, the torque provided by the generator Q_g is obtained from the expression Q_a=η_mη_eQ_g, where η_m and η_e are the mechanical and electrical efficiencies, respectively.

When the ambient wind speed is above the rated value $u_{{\mathrm{\infty }}_{\mathrm{r}}}={\mathrm{\Omega }}_{\mathrm{r}}R/{\mathit{\lambda }}^{*}$, there is enough power carried by the wind for the turbine to produce its maximum (rated) output P_r. This is called region III (also termed the above-rated or full-load regime), and the turbine operates at the constant (rated) rotor speed Ω_r. Hence, the aerodynamic torque is constant, i.e., $Q_{\mathrm{a}}=P_{\mathrm{r}}/{\mathrm{\Omega }}_{\mathrm{r}}$, whereas blades are progressively pitched into the wind to reduce C_P as wind speed increases.

To implement these standard region II and III control strategies in the proposed model, the following power equation is introduced:

$$\begin{array}{}\text{(30)}& {\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }A{u}_{\mathrm{\infty }}^{\mathrm{3}}{C}_{\mathrm{P}}\left(\mathit{\lambda }\left(\mathrm{\Omega }\right),{\mathit{\theta }}_{\mathrm{p}},\mathit{\gamma }\right)=Q_{\mathrm{a}}\left(\mathrm{\Omega }\right)\mathrm{\Omega }.\end{array}$$

The equation has three unknowns: the rotor speed Ω, the blade pitch θ_p, and the aerodynamic torque Q_a. Given ambient conditions u_∞ and ρ, two additional conditions are necessary before the three unknowns can be computed. To this end, one can first assume that the machine operates in region II. Hence, Eq. (30) is solved by appending to it the following two constraints:

$$\begin{array}{}\text{(31a)}& {\displaystyle }& {\displaystyle }{\mathit{\theta }}_{\mathrm{p}}={\mathit{\theta }}_{\mathrm{p}}^{*},\text{(31b)}& {\displaystyle }& {\displaystyle }{Q}_{\mathrm{a}}=K\left(\mathit{\rho }\right){\mathrm{\Omega }}^{\mathrm{2}}.\end{array}$$

If the computed rotor speed exceeds the rated value, i.e., Ω>Ω_r, then it means that – for the given ambient conditions and misalignment angle – the turbine operates in region III and not region II. Hence, the solution is discarded, and Eq. (30) is solved again by appending this time the following two constraints:

$$\begin{array}{}\text{(32a)}& {\displaystyle }& {\displaystyle }\mathrm{\Omega }={\mathrm{\Omega }}_{\mathrm{r}},\text{(32b)}& {\displaystyle }& {\displaystyle }{Q}_{\mathrm{a}}=P_{\mathrm{r}}/{\mathrm{\Omega }}_{\mathrm{r}}.\end{array}$$

This same approach can be used for curtailment and derating strategies (Juangarcia et al., 2018).

Figure 6 Behavior of the IEA 3.4 MW reference turbine as it yaws out of the wind, transitioning between the above-rated region III and the below-rated region II. Thrust coefficient C_T (a), power coefficient C_P (b), blade pitch angle θ_p (c), and tip speed ratio λ (d). An interactive version of the figure is available as a Jupyter notebook at the link https://tinyurl.com/btcl-fig-6 (last access: 12 July 2024).

Figure 6 shows the application of this approach to the IEA 3.4 MW reference wind turbine, a typical onshore machine with contemporary design characteristics (Bortolotti et al., 2019). In this example the turbine is exposed to an inflow characterized by a wind speed u_∞=10.5 m s⁻¹, an air density ρ=1.22 kg m⁻³, and a linear vertical shear coefficient k=0.2. For these ambient conditions, the turbine operates in region III when it is aligned with the wind (γ=0°). As the turbine starts yawing out of the wind, it initially keeps operating in region III. Accordingly, the tip speed ratio λ (Fig. 6d) and the power coefficient C_P (Fig. 6b) remain constant, while the thrust coefficient increases (Fig. 6a) and the blades pitch back (Fig. 6c). However, at around $\left|\mathit{\gamma }\right|\approx \mathrm{15}$°, the turbine enters into region II, because the rotor-orthogonal component of the wind speed is no longer large enough to maintain the rated power output. As the misalignment keeps increasing, the pitch angle remains fixed at its optimal value (Fig. 6c), whereas the tip speed ratio drops on account of the slowing rotor speed (Fig. 6d), in accordance with the region II policy.

Often turbines present an additional intermediate operating regime, called region II$\mathrm{1}/\mathrm{2}$, which occupies a wind speed interval across the rated value $u_{\mathrm{\infty }}^{\mathrm{lb}}\le u_{{\mathrm{\infty }}_{\mathrm{r}}}\le u_{\mathrm{\infty }}^{\mathrm{ub}}$. In such cases, the turbine operates in region II when $u_{\mathrm{\infty }}<u_{\mathrm{\infty }}^{\mathrm{lb}}$, in region III when $u_{\mathrm{\infty }}>u_{\mathrm{\infty }}^{\mathrm{ub}}$, and in region II$\mathrm{1}/\mathrm{2}$ when $u_{\mathrm{\infty }}^{\mathrm{lb}}\le u_{\mathrm{\infty }}\le u_{\mathrm{\infty }}^{\mathrm{ub}}$.

Differently from regions II and III, where controllers only require knowledge of the rotor speed, the control policy in region II$\mathrm{1}/\mathrm{2}$ typically requires prescribing the desired pitch and torque settings as functions of wind speed. In other words, one has to provide the schedules θ_p=θ_p(u_∞) and Q_a=Q_a(u_∞) in the desired range $u_{\mathrm{\infty }}^{\mathrm{lb}}\le u_{\mathrm{\infty }}\le u_{\mathrm{\infty }}^{\mathrm{ub}}$. Two common examples of region II$\mathrm{1}/\mathrm{2}$ control policies are provided by load and noise alleviation techniques.

Load alleviation is often necessary because thrust reaches a sharp maximum at the rated wind speed, $T_{\mathrm{r}}=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }A{u}_{{\mathrm{\infty }}_{\mathrm{r}}}^{\mathrm{2}}{C}_{\mathrm{T}}({\mathit{\theta }}_{\mathrm{p}}^{*},{\mathit{\lambda }}^{*})$. To reduce the effects of this large load on the sizing of various turbine components, thrust clipping (or peak shaving) is used, where blades are pitched to feather according to a desired schedule θ_p=θ_p(u_∞). This has the effect of reducing the angle of attack and hence the thrust (Zalkind et al., 2022), at the cost of some reduced power. To minimize power losses for a given pitch schedule, the optimal power coefficient schedule can be computed as $C_{\mathrm{P}}\left(u_{\mathrm{\infty }}\right)={\max }_{\mathrm{\Omega }}{C}_{\mathrm{P}}(\mathrm{\Omega }R/u_{\mathrm{\infty }},{\mathit{\theta }}_{\mathrm{p}}(u_{\mathrm{\infty }}\left)\right)$ under the constraint Ω≤Ω_r, which also returns the rotor speed schedule Ω(u_∞). Consequently, the torque schedule becomes $Q_{\mathrm{a}}\left(u_{\mathrm{\infty }}\right)=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }A{u}_{\mathrm{\infty }}^{\mathrm{3}}{C}_{\mathrm{P}}\left(u_{\mathrm{\infty }}\right)/\mathrm{\Omega }\left(u_{\mathrm{\infty }}\right)$.

One effective way of constraining noise emissions is to limit the rotor speed (Leloudas et al., 2007; Bottasso et al., 2012) to a maximum noise-acceptable value Ω_n. When Ω_n<Ω_r, the increase in rotor speed as a function of wind speed that characterizes region II is stopped before the machine reaches the rated power. In this case, region II$\mathrm{1}/\mathrm{2}$ is entered when $u_{\mathrm{\infty }}>u_{\mathrm{\infty }}^{\mathrm{lb}}={\mathrm{\Omega }}_{\mathrm{n}}R/{\mathit{\lambda }}^{*}$, and, above this wind speed, the rotor operates at the constant speed Ω_n. To minimize power losses, the blade pitch setting can be computed for each wind speed u_∞ as ${\mathit{\theta }}_{\mathrm{p}}\left(u_{\mathrm{\infty }}\right)={\mathrm{argmax}}_{{\mathit{\theta }}_{\mathrm{p}}}{C}_{\mathrm{P}}\left(\mathit{\lambda }\right(u_{\mathrm{\infty }}),{\mathit{\theta }}_{\mathrm{p}})$, where $\mathit{\lambda }\left(u_{\mathrm{\infty }}\right)={\mathrm{\Omega }}_{\mathrm{n}}R/u_{\mathrm{\infty }}$. The corresponding aerodynamic torque schedule is readily obtained as $Q_{\mathrm{a}}\left(u_{\mathrm{\infty }}\right)=\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }\mathrm{AR}{C}_{\mathrm{P}}\left(\mathit{\lambda }\right(u_{\mathrm{\infty }}),{\mathit{\theta }}_{\mathrm{p}}(u_{\mathrm{\infty }}\left)\right)/\mathit{\lambda }\left(u_{\mathrm{\infty }}\right)$. The end of region II$\mathrm{1}/\mathrm{2}$ is reached when, for sufficiently high u_∞, the turbine reaches the rated power, finally entering into region III.

In order to implement a given control policy for region II$\mathrm{1}/\mathrm{2}$ with the proposed model, we assume that the controller will implement the desired schedules θ_p(u_∞) and Q_a(u_∞) (whether computed as explained above or according to different criteria) by reacting to the rotor-orthogonal wind speed component u_∞cos γ. In the absence of specific details on the implementation, this is a reasonable assumption, as region II$\mathrm{1}/\mathrm{2}$ controllers are typically based on rotor-effective wind speed estimates (Bottasso et al., 2012; Zalkind et al., 2022) that, in a misaligned condition, will sense u_∞cos γ and not u_∞. Therefore, when $u_{\mathrm{\infty }}^{\mathrm{lb}}\le u_{\mathrm{\infty }}\cos \mathit{\gamma }\le u_{\mathrm{\infty }}^{\mathrm{ub}}$, Eq. (30) is solved by appending to it the following two constraints:

$$\begin{array}{}\text{(33a)}& {\displaystyle }& {\displaystyle }{\mathit{\theta }}_{\mathrm{p}}={\mathit{\theta }}_{\mathrm{p}}\left(u_{\mathrm{\infty }}\cos \mathit{\gamma }\right),\text{(33b)}& {\displaystyle }& {\displaystyle }{Q}_{\mathrm{a}}=Q_{\mathrm{a}}\left(u_{\mathrm{\infty }}\cos \mathit{\gamma }\right).\end{array}$$

In summary, using Eq. (30) in combination with the constraint Eqs. (31), (32), or (33) yields the set point achieved by the turbine for given ambient conditions and a given misalignment, no matter what region it corresponds to and what control strategy is implemented by the controller.

3.4 Model calibration

Through Eq. (25b) and (26), the analytical model depends on C_D, C_L,α, and β. These are average parameters, which represent in the model the equivalent effect caused by corresponding quantities that in reality exhibit a spanwise variability.

When numerical or experimental measurements are available, the parameters C_D, C_L,α, and β can be calibrated to minimize the error produced by the model in the prediction of the power loss factor η_P and of the thrust coefficient C_T. Notice that C_T is preferred to η_T for this scope, because it was found that the informational content of η_T is very similar to the one of η_P, reducing the quality of the tuning. Calibration is here performed by numerically solving the following minimization problem:

$$\begin{array}{}\text{(34)}& \begin{array}{rl}{\min }_{C_{\mathrm{D}},C_{\mathrm{L},\mathit{\alpha }},\mathit{\beta }}& \sqrt{{\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N{\left({\mathit{\eta }}_{\mathrm{P},i}^{\mathrm{obs}}-{\mathit{\eta }}_{\mathrm{P},i}^{\mathrm{mod}}\right)}^{\mathrm{2}}}\\ & +\sqrt{{\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N{\left(C_{\mathrm{T},i}^{\mathrm{obs}}-C_{\mathrm{T},i}^{\mathrm{mod}}\right)}^{\mathrm{2}}},\end{array}\end{array}$$

where $(\cdot {)}_i^{\mathrm{obs}}$ represents N numerical or experimental observations, and $(\cdot {)}_i^{\mathrm{mod}}$ represents the corresponding model predictions. For each data set, tuning was performed solving $N/\mathrm{2}$ times the problem expressed by Eq. (34) using a gradient-based optimization, each time with a different random 50 % subset of the available data, and finally averaging the resulting parameters.

3.5 Simplified choice of model parameters in the absence of calibration data

When calibration is not possible, the equivalent model parameters C_D, C_L,α, and β must be estimated from the corresponding actual spanwise distributions $C_{\mathrm{D}}(r/R)$, $C_{\mathrm{L},\mathit{\alpha }}(r/R)$, and $\mathit{\beta }(r/R)$.

Examining the expression for thrust given by Eq. (16), neglecting drag, it appears that the lift force has roughly a spanwise triangular distribution. In fact, inspecting Eq. (5), u_t is proportional to r, whereas u_n does not depend on r. Additionally, the leading term of the Taylor series of the optimal twist distribution is $\mathrm{1}/r$ (Burton et al., 2011). Similarly, inspecting the expression for power given by Eq. (21), again neglecting the contribution of drag, it appears that – for the same reasons – also the spanwise power capture has a triangular distribution. This suggests the need to evaluate the spanwise integrals at the centroid of the triangle, which is located at $r/R=\mathrm{2}/\mathrm{3}$. Adopting this approach, the model parameters are then set to the following values:

$$\begin{array}{}\text{(35a)}& {\displaystyle }& {\displaystyle }{C}_{\mathrm{D}}=f_{\mathrm{d}}{C}_{\mathrm{D}}(\mathrm{2}/\mathrm{3}),\text{(35b)}& {\displaystyle }& {\displaystyle }{C}_{\mathrm{L},\mathit{\alpha }}=f_{\mathrm{l}}{C}_{\mathrm{L},\mathit{\alpha }}(\mathrm{2}/\mathrm{3}),\text{(35c)}& {\displaystyle }& {\displaystyle }\mathit{\beta }=\mathit{\beta }(\mathrm{2}/\mathrm{3}).\end{array}$$

Coefficient f_d is a correction factor for drag, while f_l is a knockdown factor for lift, which accounts for the finite span of the blades. Based on comparisons with calibrated values (see Sect. 4.2), we recommend a drag correction factor f_d=1 for moderate yaw (up to 20°) and pitch values and a smaller value of 0.45 if the model has to be used also for large yaw and pitch settings. This smaller value is probably due to the approximation of a small inflow angle used in the model, which is partially corrected by a smaller drag coefficient. For lift, we recommend the value $f_{\mathrm{l}}=\mathrm{2}/\mathrm{3}$.

The performance of this simplified choice of model parameters is demonstrated later in Sect. 4.2.

4 Model validation

4.1 Validation with respect to LES-ALM simulations

LES-ALM simulations are used for testing the accuracy of the model in representing misaligned conditions, similarly to what has been done by other authors (Gebraad et al., 2016; Liew et al., 2020; Nanos et al., 2022). The effects of the rotor on the flow are modeled with the filtered ALM of Troldborg et al. (2007) and Martínez-Tossas and Meneveau (2019) by projecting forces computed along the lifting lines onto the LES grid. The Cartesian mesh consists of approximately 3.5 million cells and uses four refinement levels. The smallest cells measure 1 m and are used to discretize the area surrounding the rotor.

Simulations were conducted for the IEA 3.4 MW reference wind turbine, whose complete technical specifications are reported in Bortolotti et al. (2019). Here we only note that the turbine has a 5° uptilt angle; i.e., $\mathit{\delta }=-\mathrm{5}$°. The parameters of the proposed model were calibrated as explained in Sect. 3 based on the LES-ALM simulations described in Sect. 4.1.2 but not the ones of Sect. 4.1.1, obtaining the values $C_{\mathrm{D}}=\mathrm{0.0052}\pm \mathrm{0.0001}$, $C_{\mathrm{L},\mathit{\alpha }}=\mathrm{4.759}\pm \mathrm{0.007}$ rad⁻¹, and $\mathit{\beta }=-\mathrm{3.345}\pm \mathrm{0.007}$°, for a 95 % confidence level.

Figure 7 Thrust loss factor η_T (a), power loss factor η_P (b), blade pitch θ_p (c), and tip speed ratio λ (d), plotted as functions of the misalignment angle γ, for wind speeds below, above and around the rated one. Proposed model: lines; LES-ALM: solid markers; Heck et al. (2023) (only for the 8.5 m s⁻¹ case): dashed orange lines.

4.1.1 Simulations in control regions II and III

First we demonstrate the integration with a standard controller, including operations in region II, III, and the transition between the two as the turbine is progressively yawed out of the wind. For the LES-ALM simulations, set points were computed using a controller in the loop, based on an implementation similar to the one of Bortolotti et al. (2019). For the proposed model, set points were obtained from Eq. (30) in combination with the constraint Eqs. (31) and (32).

We consider laminar inflows with four wind speeds: one below-rated speed of 8.5 m s⁻¹ and three above-rated speeds of 10.5, 11, and 13 m s⁻¹. Figure 7 reports the results in term of the thrust loss factor η_T (Fig 7a), power loss factor η_P (Fig. 7b), blade pitch θ_p (Fig. 7c), and tip speed ratio λ (Fig. 7d), all plotted as functions of the misalignment angle γ. The solid markers are the results of the LES-ALM simulations, while the lines represent predictions of the proposed model.

For the lowest wind speed of 8.5 m s⁻¹, the turbine always operates in region II. This is the only case where the method of Heck et al. (2023) is strictly applicable. In fact, their method does not contain a generic thrust model, but rather it is formulated in terms of the modified thrust coefficient $C_{\mathrm{T}}^{\prime }$, which is constant when a turbine yaws out of the wind in region II. The results of this alternative model are shown with a dashed orange line in the figure. The benchmark LESs feature a negative uptilt $\mathit{\delta }=-\mathrm{5}$°, which is not modeled in the approach of Heck et al. (2023). To avoid cluttering the results with this additional effect, here and in the following examples the total true misalignment angle μ (instead of γ) is provided as input to the model of Heck et al. (2023). The figure shows that both methods are in excellent agreement with the CFD results.

For the highest wind speed of 13 m s⁻¹, the turbine operates in region III for all misalignment angles. On the other hand, for a wind speed of 11 m s⁻¹ the machine enters region II around γ=27°, and for a wind speed of 10.5 m s⁻¹ the machine enters region II at about γ=16°.

In general, there is a very good agreement of the model with the higher-fidelity CFD results, not only in terms of loss factors, but also on the calculation of the set points.

4.1.2 Simulations with fixed tip speed ratio and pitch setting

Next, we present a second set of results obtained by varying the misalignment while keeping the tip speed ratio and pitch constant, for different inflows. These conditions are meant to provide a more general view of the performance of the proposed method in a variety of conditions, although – as explained in Sect. 3.3 – tip speed ratio and pitch in general do not both remain constant when a turbine yaws out of the wind.

The four operational scenarios of Table 1 are considered. The flow is laminar and steady in all scenarios. Cases 1 and 2 have no shear and different tip speed ratios, whereas cases 3 and 4 are sheared and have the same λ.

Table 1 Operational scenarios for the LES-ALM simulations.

Figures 8 and 9 report the power loss factor η_P in the range of yaw misalignment angles $-\mathrm{30}\mathit{°}<\mathit{\gamma }<\mathrm{30}\mathit{°}$ for different pitch settings, each corresponding to a different thrust coefficient C_T,0 in aligned conditions. Figure 8 corresponds to scenarios 1 and 2 of Table 1, i.e., no shear, while Fig. 9 reports the solution for the sheared cases 3 and 4 of that same table.

Figure 8 Power loss factor η_P vs. misalignment angle γ in the unsheared scenarios 1 (λ=8) (a–d) and 2 (λ=9.5) (e–h). Each subplot corresponds to a different value of the thrust coefficient in wind-aligned conditions: $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.74}$ (a); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.54}$ (b); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.44}$ (c); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (d); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.86}$ (e); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.60}$ (f); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.47}$ (g); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.35}$ (h).

Figure 9 Power loss factor η_P vs. misalignment angle γ in the sheared scenarios 3 (k=0.06) (a–d) and 4 (k=0.19) (e–h). Each subplot corresponds to a different value of the thrust coefficient in wind-aligned conditions: $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.77}$ (a); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.56}$ (b); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.45}$ (c); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (d); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.77}$ (e); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.56}$ (f); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.45}$ (g); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (h).

In the figures, LES-ALM results (shown with circular black markers) are compared with the proposed approach (shown with solid blue lines) and the method proposed by Heck et al. (2023) (shown with dashed orange lines). For the latter, the modified thrust coefficient $C_{\mathrm{T}}^{\prime }$ was obtained directly from each LES at the corresponding γ value.

Overall, there is a very good match between the predictions of the proposed model and numerical simulations. For the null shear cases of Fig. 8, results are reported only for positive yaw angles, as power is symmetric. On the other hand, power is not symmetric for the sheared inflow cases of Fig. 9, which shows clear evidence of the complex behavior described in Sect. 2.8. At high C_T (low pitch), curves are very nearly symmetric with respect to γ. However, as thrust is decreased (and pitch increased), power capture is larger at positive γ values compared to negative γ values.

The model of Heck et al. (2023) performs similarly well at high and moderate rotor loading, when $C_{\mathrm{T}}^{\prime }$ is roughly constant. However, as the C_T is reduced, the model tends towards the solution cos ³γ, and therefore its accuracy is compromised. Moreover, the model fails to predict the shear-induced asymmetry (see Fig. 9).

As predicted by the proposed model, the power asymmetry increases with shear (see Table 1 and the explanation given in Sect. 2.8). To facilitate the visualization of this effect, Fig. 10 shows the difference $\mathrm{\Delta }{\mathit{\eta }}_{\mathrm{P},\mathit{\gamma }=\pm \mathrm{30}\mathit{°}}$ between the two values of η_P at $\mathit{\gamma }=\pm \mathrm{30}$° as a function of shear for varying thrust coefficients. The asymmetry also exhibits a noticeable dependency on the thrust coefficient, with larger asymmetries being observed for lower values of C_T, as predicted in Sect. 2.8 by examining Eq. (22).

Figure 10 Difference $\mathrm{\Delta }{\mathit{\eta }}_{\mathrm{P},\mathit{\gamma }=\pm \mathrm{30}\mathit{°}}$ between power loss factors η_P evaluated at misalignments γ=30° and at $\mathit{\gamma }=-\mathrm{30}$°, as a function of vertical linear shear coefficient k, for varying thrust coefficient C_T. Proposed model: solid lines; LES-ALM simulations: markers.

Figures 11 and 12 report the thrust loss factor η_T as a function of yaw misalignment for the same four scenarios and different thrust settings. Here again model predictions are indicated with lines and LES-ALM results with markers. There is a consistently good match for all scenarios and for all yaw and pitch values. The lack of symmetry is again consistent with the model, similarly to the case of power discussed above. Figure 12a–d show a higher thrust for positive yaw angles at low thrust coefficients (high pitch values) because of the high tip speed ratio of scenario 3, indicating that term $C_{{\mathrm{T}}_{\mathrm{2}}}$ prevails over $C_{{\mathrm{T}}_{\mathrm{1}}}$. The opposite happens in Fig. 12e–h due to the lower λ of scenario 4.

Figure 11 Thrust loss factor η_T vs. misalignment angle γ in scenario 1 with λ=8 (a–d) and scenario 2 with λ=9.5 (e–h). Each subplot corresponds to a different value of the thrust coefficient in wind-aligned conditions: $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.74}$ (a); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.54}$ (b); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.44}$ (c); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (d); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.86}$ (e); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.60}$ (f); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.47}$ (g); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.35}$ (h).

Figure 12 Thrust loss factor η_T vs. misalignment angle γ in scenarios 3 with k=0.06 (a–d) and scenario 4 with k=0.19 (e–h). Each subplot corresponds to a different value of the thrust coefficient in wind-aligned conditions: $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.77}$ (a); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.56}$ (b); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.45}$ (c); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (d); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.77}$ (e); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.56}$ (f); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.45}$ (g); $C_{\mathrm{T},\mathrm{0}}=\mathrm{0.36}$ (h).

Table 2 Characteristics of the wind tunnel campaigns.

Overall, it appears that the performance of the rotor is strongly dependent on the thrust coefficient and tip speed ratio and hence on the way it is controlled when it yaws out of the wind. Therefore, the standard power law ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ may oversimplify the complex aerodynamics that are typical of this problem. On the other hand, notwithstanding its simplicity, the proposed model is in very good agreement with sophisticated CFD simulations, and it is capable of describing even relatively minor effects of the complex behavior of a misaligned wind turbine rotor in a sheared inflow.

4.2 Validation of the simplified choice of model parameters

The simplified choice of model parameters described in Sect. 3.5 is based on Eq. (35), which includes the correction factors f_d and f_l. To verify the existence of typical values for these factors, we considered four different wind turbines: IEA 3.4 MW (Bortolotti et al., 2019); NREL 5 MW (Jonkman et al., 2009); G178, which is a modified version of the DTU 10 MW (Bak et al., 2013; Wang et al., 2021); and the small-scale G1 turbine (Bottasso and Campagnolo, 2022; Campagnolo et al., 2020). For the three full-scale machines, the model parameters were first calibrated using the LES-ALM simulation results of Sect. 4.1.2, whereas for the G1 model the calibration was performed using wind tunnel measurements (see later Sect. 4.3).

The parameters calibrated this way were then compared to the ones based on the simplified approach of Eq. (35), leading to the recommended values reported in Sect. 3.5. The $C_{\mathrm{D}}(r/R)$ and $C_{\mathrm{L},\mathit{\alpha }}(r/R)$ coefficients were obtained by averaging over the interval of angles of attack 2° below the negative and positive stall limits. In all cases, the calibrated value of the twist corresponded remarkably well with the actual twist at $\mathrm{2}/\mathrm{3}$ span, i.e., $\mathit{\beta }(\mathrm{2}/\mathrm{3})$.

The simplified choice of the model parameters was then applied to the NREL 5 MW and G178 10 MW wind turbines. Simulations were performed with a steady inflow, at different misalignments, and for two different blade pitch settings. The proposed model was calculated with the parameters based on Eq. (35), using the default correction coefficients f_d=0.45 and $f_{\mathrm{l}}=\mathrm{2}/\mathrm{3}$ (in other words, without using LES-ALM calibrated values, replicating what one could do in the absence of suitable tuning data).

The results in terms of η_T and η_P for the two turbines are reported in Fig. 13 and compared with LES-ALM simulations. For both turbines there is an excellent match between model predictions and CFD results. This seems to indicate that the even a simplified choice of the model parameters is sufficient for a good performance of the model.

Figure 13 Thrust η_T (a) and power η_P (b) loss factors as functions of the misalignment angle for two pitch settings. Results obtained with the NREL 5 MW (Jonkman et al., 2009) and G178 10 MW (Bak et al., 2013; Wang et al., 2021) wind turbines, using the simplified calculation of the model parameters based on Eq. (35). Proposed model: lines; LES-ALM simulations: markers.

4.3 Validation with respect to wind tunnel measurements

Next, the model is compared to data recorded during wind tunnel experimental campaigns performed with a G1 wind turbine (Campagnolo et al., 2016). This scaled machine has a diameter of 1.1 m, a rated rotor speed of 850 rpm, and null tilt. The design of the G1 is described in Bottasso and Campagnolo (2022), and its rotor aerodynamic and wake characteristics have been reported in Wang et al. (2021) and references therein.

Tests were performed in a boundary layer wind tunnel (Bottasso et al., 2014) with three different inflows: the first one, termed Low-TI, has no shear and a very low turbulence intensity (approximatively equal to 1 %); the other two, termed Mod-TI and High-TI, have TIs of about 6 % and 13 % at hub height and vertical linear shears in the rotor region equal to k=0.11 and k=0.15, respectively. Figure 14a reports the vertical profiles of the longitudinal wind speed component u measured by means of CTA probes (Bottasso et al., 2014), normalized by the wind speed u_pitot measured by a pitot tube placed at hub height. Figure 14b shows the vertical profiles of the turbulence intensity, as measured with the same instrumentation.

Figure 14 Wind tunnel inflows. Vertical wind speed profiles, with corresponding best-fitted linear shears (dotted lines) (a), and vertical profiles of turbulence intensity (b).

The experimental characterization of power losses in misaligned conditions was performed based on the three campaigns of Table 2, totaling 119 observations of a duration of 2 min each. The measured average θ_p and λ are reported in Appendix C for campaigns 1 and 2 and in Fig. 19 for campaign 3.

The first campaign was conducted in region II (i.e., the below-rated partial-load regime), using the classical variable-rotor-speed maximum-power-tracking strategy. Tests were conducted in all three inflow conditions – Low-TI, Mod-TI, and High-TI – with hub-height wind speeds of 5.86, 5.69, and 5.40 m s⁻¹, respectively.

The second campaign was also conducted in region II, but in this case the turbine was derated in the range $P_{\mathrm{d}}\in [\mathrm{50},\mathrm{100}]$ % while adopting two different strategies: iso-λ, where the tip speed ratio is held constant (Campagnolo et al., 2023), and min-C_T, where the thrust coefficient is minimized (Juangarcia et al., 2018). Tests were conducted only in the Mod-TI inflow, with a hub-height wind speed of 5.62 m s⁻¹,

The third campaign was conducted in region III (i.e., the above-rated full-load regime) in the Low-TI and Mod-TI inflows, with hub-height wind speeds of 6.97 and 6.11 m s⁻¹, respectively.

Various sources of error affect the experimental observations. These include measurements of the wind speed u_pitot upstream of the model (obtained by a pitot tube placed at hub height 3 D in front of the turbine), of the air density ρ, of the rotor speed Ω, of the shaft torque Q, of the bending moment at tower base (which is used to estimate thrust), of the blade pitch angle, and of the nacelle orientation with respect to the wind tunnel (i.e., of the misalignment angle). The error in u_pitot is related to the uncertainty associated with the measurements of flow density and dynamic pressure. This latter quantity is measured with a MKS Baratron-type 226A transducer (MKS Instruments Inc., 2022) with full span equal to 1 Torr, characterized by an accuracy of ±0.4 Pa. Density is instead derived from measurements of air pressure, temperature, and humidity, and it is affected by an error equal to ±0.01 kg m³ (Wang et al., 2020). Torque is measured with a load cell installed on the rotor shaft, and it is affected by an uncertainty of ±0.005 Nm. The rotor speed measurement, provided by an optical incremental encoder, is instead affected by an error equal to ±1.5 rpm. The measurement uncertainty on power P=QΩ is derived by adding in quadrature the uncertainties on Q and Ω. Thrust T is obtained by correcting the measurements of the bending moments at tower base by the effects induced by the drag of the tower, nacelle, and hub spinner (Wang et al., 2020). The calibration of the load cell at tower base revealed an uncertainty in the thrust of ±0.14 N. Blade pitch and nacelle orientation are measured by optical encoders, affected by uncertainties of ±0.2°. In turn, all these effects are used to quantify uncertainties in the tip speed ratio λ, as well as yaw-induced power and thrust losses η_P and η_T, again by adding errors in quadrature. In the following, the resulting uncertainties are reported for a 95 % confidence level.

Uncertainties in some experimental measurements affect also the predictions of the proposed model. The uncertainties of the four model inputs – tip speed ratio, blade pitch, rotor speed, and yaw misalignment – were propagated forward throughout the model by Latin hypercube sampling with 10 000 sample points, using the UQLab software (Marelli and Sudret, 2014).

Tuning of the model parameters C_D, C_L,α, and β was performed with Eq. (34), using a total of 94 observations from experimental campaigns 1 and 2.

Given the small size of the G1 wind turbine, the Reynolds number at its blade sections is particularly low. Although special low-Reynolds airfoils are used in the design of the G1 blades (Bottasso and Campagnolo, 2022), their aerodynamic characteristics are particularly sensitive to the operating conditions of the turbine (Wang et al., 2020). In fact, the Reynolds number has a significant effect on the drag and on the zero-lift direction, which in turn affects the parameter β, whereas the effect on the lift slope C_L,α is negligible (Wang et al., 2020). Accordingly, the model parameters C_D and β are assumed to depend on the rotational speed Ω, since the relative speed at the airfoils is close to the tangential speed (u≈u_t; see Fig. 3). The values of C_D and β at $\mathrm{\Omega }=[\mathrm{850},\mathrm{625},\mathrm{400}]$ rpm were assumed as unknowns, and a piecewise linear interpolation was used at other intermediate values of the rotor speed.

Figure 15 reports the tuned C_D and β parameters, with whiskers representing the corresponding 95 % confidence intervals. As expected, drag decreases for increasing rotor speed, i.e., for increasing sectional Reynolds number. The twist β also exhibits the same trend, since the zero-lift direction rotates nose up as the Reynolds number increases (Wang et al., 2020). The tuned parameter C_L,α is equal to 4.5033±0.0459.

Figure 15 Tuned model parameters C_D and β as functions of rotational speed Ω.

For maximum-power tracking operation in region II (test campaign 1), Fig. 16 reports a comparison between model-predicted and measured power and thrust losses. The present model results are indicated with solid blue lines, while measurements are indicated by black circles. Whiskers indicate the respective 95 % confidence intervals. There is a very good match between experimental measurements and the present model, with the latter falling within the uncertainty range of the measurements in most cases. As predicted by the model, the sheared inflow conditions Mod-TI and High-TI exhibit the expected non-symmetric behavior with respect to positive and negative yaw angles.

Figure 16 Experimental campaign 1. Power loss factor η_P (a–c) and thrust loss factor η_T (d–f) vs. yaw misalignment γ, in region II operation for the three different inflow conditions. Whiskers indicate the 95 % confidence intervals.

In principle, the model of Heck et al. (2023) would be applicable to these tests in region II. However, their method does not directly consider the aerodynamic characteristics of the blades, as it expresses their behavior through the single parameter represented by $C_{\mathrm{T}}^{\prime }$. Therefore, it is blind to the variability of twist and drag with respect to the Reynolds number, which drops significantly as the misalignment angle increases. Since this strong Reynolds dependency is specific to the small scale of wind tunnel models, the results of the method of Heck et al. (2023) are not shown here, because its poor match with the measurements would be misleading, as these effects would not be present at full scale. The present method is not affected by this issue, because it uses a lifting-line approach and specifically includes the blade aerodynamic characteristics in the governing Eqs. (22) and (25b).

The loss factors are reported for derated operation (test campaign 2) in Fig. 17. Here again the match between experimental data and predictions of the present model is very good, with the latter being mostly within the uncertainty band of the measurements. Slightly larger deviations are observed for the min-C_T case at P_d=50 %. This can be explained by the fact that the machine operates at significantly low λ values, with consequent low rotational speeds. This results in particularly high angles of attack (Juangarcia et al., 2018) and very low chord-based Reynolds numbers. Both have significant impacts on the airfoil performance, which are likely not properly captured by the analytical model. Overall, it appears that the model is capable of capturing the reduction in the thrust coefficient as derating P_d increases, as well as the lack of symmetry with respect to the misalignment angle.

Figure 17 Experimental campaign 2. Power loss factor η_P (a–h) and thrust loss factor η_T (i–p) vs. yaw misalignment γ, in derated operation in region II for the Mod-TI inflow case. Whiskers indicate the 95 % confidence intervals.

The effects of thrust and shear are visualized in Fig. 18 in terms of the average and the difference of the power loss factors at $\mathit{\gamma }\pm {\mathrm{20}}^{\circ }$, respectively denoted ${\stackrel{\mathrm{‾}}{\mathit{\eta }}}_{\mathrm{P},\mathit{\gamma }=\pm \mathrm{20}\mathit{°}}$ and $\mathrm{\Delta }{\mathit{\eta }}_{P,\mathit{\gamma }=\pm {\mathrm{20}}^{\circ }}$. It appears that power losses tend to decrease with increasing thrust coefficients (Fig. 18a), whereas there is no significant dependency on shear (Fig. 18b). The power loss asymmetry grows with increasing shear (Fig. 18d). On the other hand, the asymmetry is roughly constant with respect to the thrust coefficient (Fig. 18c).

Figure 18 Power loss factor ${\stackrel{\mathrm{‾}}{\mathit{\eta }}}_{\mathrm{P},\mathit{\gamma }=\pm \mathrm{20}\mathit{°}}$ (averaged between $\mathit{\gamma }=\pm \mathrm{20}$°) as a function of C_T(γ=0°) (a) and of shear k (b). Difference $\mathrm{\Delta }{\mathit{\eta }}_{\mathrm{P},\mathit{\gamma }=\pm \mathrm{20}\mathit{°}}$ between power loss factors η_P evaluated at misalignments γ=20° and at $\mathit{\gamma }=-\mathrm{20}$°, as a function of C_T(γ=0°) (c) and of shear k (d). Proposed model: solid lines; experimental measurements: black circles. Whiskers indicate the 95 % confidence intervals.

In the third test campaign, the wind turbine is operated above rated conditions. Figure 19 reports the 2 min average tip speed ratio and pitch angles measured during the experiment and plotted as functions of the misalignment angle γ. For the Low-TI case, the turbine operates in region III for all misalignment angles. Recalling that the tip speed ratio is defined as $\mathit{\lambda }=\mathrm{\Omega }R/u_{\mathrm{\infty },\mathrm{hub}}$, since both the ambient wind speed u_∞,hub and rotor speed Ω are constant, when the turbine yaws away from the wind λ (indicated by red circles in Fig. 19a) remains constant, while the blades are pitched back (red circles in Fig. 19b) in order to keep the power output equal to the rated value. The same happens for the Mod-TI case. However, when $\mathit{\gamma }<-\mathrm{25}$°, the turbine exits region III and enters region II. Therefore, as blade pitch (purple squares in Fig. 19b) reaches the value for the maximum power coefficient, λ starts decreasing (purple squares in Fig. 19a).

Figure 19 Experimental campaign 3. Blade pitch angle θ_p (a) and tip speed ratio λ (b) in Low-TI and Mod-TI inflows.

Figure 20 shows the results for the power and thrust loss factors in this scenario. Once again the proposed model exhibits a very good match with the experiments, falling within the uncertainty bands in most cases.

Figure 20 Experimental campaign 3. Power loss factor η_P (a, b) and thrust loss factor η_T (c, d) vs. yaw misalignment γ, in above-rated-speed operation for the Low-TI and Mod-TI inflow cases. Model predictions: solid lines; experimental measurements: black circles. Whiskers indicate the 95 % confidence intervals.

These results confirm the ability of the method to correctly represent the effects of different control approaches, covering both regions II and III, including derating. This is crucially important because, as shown, control laws have a strong impact on the behavior of power and trust in misaligned conditions.

5 Optimal wake steering

The insight provided by the new model suggests two questions.

• What is the power-optimal way to yaw a single turbine out of the wind?

• And does the new model affect the way wake steering should be conducted?

We try to give some initial answers to these questions in the following two sections.

5.1 Optimal power capture of a single misaligned turbine

The new model was used to compute the optimal power of a wind turbine when it is misaligned with respect to the wind. The analysis was conducted for the same IEA 3.4 MW wind turbine used for the previous numerical validation of the model.

Figure 21 Comparison between standard and optimal control strategies for different yaw angles γ. Tip speed ratio λ (a); pitch angle θ_p (b); percent thrust difference between the optimal and standard strategies (c); percent power difference between the optimal and standard strategies (d).

The optimal control strategy was computed by numerical optimization using an adaptive Nelder–Mead algorithm (Gao and Han, 2012), and results are shown in Fig. 21. The figure also reports the standard region II control approach, which consists in holding the pitch angle fixed while the generator torque is varied proportionally to the square of the rotor speed, i.e., Q∼Ω², as explained in Sect. 3.3. Assuming as a first approximation that $P=\cos {\mathit{\gamma }}^{{\mathrm{p}}_{\mathrm{p}}}$ and considering that $Q=P/\mathrm{\Omega }$, it follows that $\mathrm{\Omega }\sim \cos {\mathit{\gamma }}^{{\mathrm{p}}_{\mathrm{p}}/\mathrm{3}}$. The figure reports the solution computed for a coefficient p_p=1.88, following Fleming et al. (2015).

As the yaw misalignment increases, the tip speed ratio drops for the standard control strategy, driven by the reduced rotor-orthogonal component of the wind. Since the pitch angle remains fixed, the reduced λ also leads to a decreased thrust coefficient. Although this might be beneficial for reducing loading on the yawed turbine, the resulting drop in power is significant. On the other hand, the optimal strategy governs the turbine to keep a much more constant tip speed ratio and thrust coefficient, while the blade pitches back a little. This results in some power boost, which is small for moderate angles but reaches above 3 % around ±30°. These findings are in line with results presented by Cossu (2021a) and Heck et al. (2023). While the higher C_T implies that the turbine is loaded more than in the standard case, it also has an effect on the wake that will be felt downstream, as explored in the next section.

5.2 Optimal power capture of two turbines

The previous section showed that a single turbine can extract more energy from the wind in misaligned conditions when its C_T is increased compared to a standard region II control approach. This so-called overinductive yaw control (Cossu, 2021a) increases the velocity deficit in the wake, but it also affects its recovery and enhances its deflection. It is therefore necessary to find the optimal trade-off among these complex effects when considering wake-steering wind farm control (Meyers et al., 2022).

FLORIS v3 (NREL, 2023b), modified with the present model, was used to optimize the power capture of a cluster of two IEA 3.4 MW wind turbines placed at a distance of 5 D. The wake was modeled with the Gauss–curl-hybrid model (King et al., 2021). The inflow is characterized by an ambient wind speed $u_{\mathrm{\infty },\mathrm{hub}}=\mathrm{9.7}$ m s⁻¹, a shear of 0.12, and a turbulence intensity of 6 %. A 60° range of wind directions Φ was considered in order to realize different degrees of overlap between the wake and the downstream rotor. The optimal wind farm control strategy was computed by numerically maximizing the cluster power with the same adaptive Nelder–Mead algorithm used for the single-turbine case of the previous section.

Results are shown in Fig. 22. The plots report in dotted green lines the results obtained with greedy control (i.e., each turbine maximizes its own power capture), in solid blue lines the solution obtained with wake-steering control based on the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law using p_p=1.88, and in dashed red lines with wake-steering control based on the present model. For the three control strategies, results were validated with LES-ALM simulations run for five different wind directions, namely $\mathrm{\Phi }=\mathit{\{}\mathrm{270}\pm \mathrm{5.74},\mathrm{270}\pm \mathrm{2.5},\mathrm{270}\mathit{\}}$°, corresponding to rotor overlaps of 50 %, 78.2 %, and 100 %, respectively. The LES-ALM results are indicated in the figure with markers, where the colors correspond to the control strategy.

Figure 22 Control of a cluster of two turbines in wake interference conditions. Green: greedy policy; blue: optimal wake-steering solution based on ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$; red: optimal wake-steering solution based on the proposed model. Lines: FLORIS engineering wake model; markers: LES-ALM CFD. Tip speed ratio λ (a); thrust coefficient C_T (b); absolute yaw misalignment $\left|\mathit{\gamma }\right|$ (c); percent power changes with respect to the greedy policy for the upstream wake-steering turbine (d); percent power changes with respect to the greedy policy for the downstream turbine (e); overall percent power changes for the cluster of two turbines (f).

Figures 22a–c show the front turbine tip speed ratio λ, pitch angle θ_p, and absolute misalignment angle $\left|\mathit{\gamma }\right|$, respectively, all plotted as functions of wind direction Φ. The solution for the present model is characterized by a fairly constant tip speed ratio that, in conjunction with some pitch-back of the blades at the highest misalignments, also results in a roughly constant thrust coefficient (not shown for brevity). This is in contrast to the solution based on the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law, where both the tip speed ratio and the thrust coefficient drop at the higher misalignments that correspond to the strongest wake overlap conditions. In addition, the present model also results in slightly larger misalignment angles, as shown by Fig. 22c.

The bottom three plots show the effects of the various control strategies on power as a function of wind direction Φ. Figure 22d and e report the power changes with respect to the greedy strategy for the front and back turbines, respectively. It appears that the upstream machine, due to a larger misalignment, looses more power than in the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ case. Conversely, it also appears that the second machine gains more power with the strategy based on the new model, thanks to the larger misalignment of the upstream turbine but also due to its larger thrust coefficient. Finally, Fig. 22f shows the overall gain at the cluster level. Results indicate a fairly consistent improvement, in excess of roughly 1 %, for almost the entire wake overlap range.

The LES-ALM results confirm the findings based on the FLORIS engineering wake model: more power losses for the front turbine and more gains for the downstream one, resulting in a positive net gain for the cluster.

6 Conclusions

We have presented a new model to estimate the power performance of a misaligned wind turbine rotor. The model is a modified version of the classical blade element momentum theory, where the rotor is considered a lifting wing of finite span operating at an angle of attack.

The new model reveals the following characteristics of the behavior of a misaligned rotor.

• Power does not depend on the misalignment angle according to the ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ law, a formula in widespread use in the literature.

• The true effective misalignment angle that drives wake behavior is a combination of both yaw and tilt. Therefore, a two-dimensional wake model should be described in the plane formed by the rotor axis and wind vectors, not on the horizontal plane as commonly assumed.

• Power depends on the true misalignment angle but – crucially – also on the way the rotor is governed as it is pointed out of the wind, a fact that probably explains the widely different performance observed by various authors. This fact also means that power losses due to misalignment can be mitigated by using a suitable control strategy.

• According to the model, the observed lack of symmetry between positive and negative misalignment angles is caused by the interaction with a sheared inflow. In these conditions, there is a complex interplay of various effects that may lead to different outcomes in terms of which yaw sign yields more or less power. In general, one can expect a small asymmetry at high thrust coefficients, while a more pronounced asymmetry emerges for low thrust and high tip speed ratios, where a higher power is generated for positive yaw angles. However, in general the behavior of power (but of thrust too, which also exhibits an asymmetric behavior) depends not only on the rotor design characteristics but also on the way it is governed, through the values of the pitch setting and of the tip speed ratio. Additionally, in the field other effects may be present (e.g., due to an asymmetric behavior of the onboard wind vane), which may add to the phenomena described by the model.

• A constant-over-the-rotor induction is sufficient to accurately describe the power and thrust behavior of a misaligned rotor in a sheared inflow. In fact, under classical small-angle assumptions, the tilting of the inflow due to misalignment and shear has only a negligible effect on the quality of the results.

The model was derived in a semi-analytical form, leading to a closed system of equations that can be directly integrated with engineering wake models, at an irrelevant computational cost. To improve its accuracy, we proposed a specific implementation that overcomes the intrinsic limits brought by the analytical solution of some model integrals. The proposed implementation corrects for the effects of misalignment a higher-fidelity power model obtained in aligned conditions, and it calibrates the model parameters based on measurements.

The model was validated in a broad range of cases, considering LES-ALM numerical simulations of various multi-MW machines as well as experimental observations on a scaled wind tunnel model, in different inflows (from unsheared laminar to sheared highly turbulent conditions), operating with controllers in the loop in regions II and III, as well as in derated conditions. In all cases, the model achieved a very satisfactory agreement with the numerical and experimental reference power and thrust values. Additionally, we demonstrated how the model can be integrated with given control laws, achieving an excellent match in the calculation of the set points. The model was also compared with a similar model recently developed by Heck et al. (2023), limitedly to the control region II where it is applicable, consistently improving on its predictions and exhibiting a wider applicability to arbitrary control strategies.

Using the proposed model, we maximized the power capture of a wind turbine for a range of misalignment angles, obtaining the optimal power strategy in terms of pitch setting and tip speed ratio. Results indicate that the maximum power extraction is obtained by keeping an almost constant tip speed ratio and by slightly reducing the blade pitch as the turbine yaws out of the wind. This also implies a roughly constant thrust coefficient, which will increase the loading on the yawed turbine, but will also have an effect on its wake.

Next, we applied the new model to the maximization of power by wake steering for a cluster of two turbines. The resulting control strategy was compared to the one obtained by the classical ${\cos }^{{\mathrm{p}}_{\mathrm{p}}}$ power loss model and validated by means of LES-ALM simulations for a few selected cases. Results indicate that the proposed model results in slightly greater power losses for the wake-steering turbine, which are more than compensated for by greater power gains for the wake-affected one, achieving a small but consistent gain in power at the cluster level for the full range of possible wake overlaps.

Future work should investigate the effects of the new model and its resulting control strategy in more complex conditions. Of particular interest is the analysis of the effects on loads, which might increase because of the eliminated drop in thrust coefficient as the turbine is yawed out of the wind.

Appendix A: Transformation matrices

A1 Transformation from ground to nacelle frame of reference

The nacelle-fixed frame of reference is obtained from the ground-fixed frame by a first rotation δ about the horizontal axis y_g, followed by a rotation γ about the vertical axis z_g. The components of a generic vector v are denoted ${\underset{\mathrm{¯}}{v}}_{\mathrm{g}}$ when measured in the ground frame ℱ_g and ${\underset{\mathrm{¯}}{v}}_{\mathrm{n}}$ when measured in the nacelle frame ℱ_n. Combining the two successive rotations, one obtains the transformation of components from one frame to the other as

$$\begin{array}{}\text{(A1)}& \begin{array}{rl}& {\underset{\mathrm{¯}}{v}}_{\mathrm{n}}=\left[\begin{array}{ccc}\cos \mathit{\delta }& \mathrm{0}& -\sin \mathit{\delta }\\ \mathrm{0}& \mathrm{1}& \mathrm{0}\\ \sin \mathit{\delta }& \mathrm{0}& \cos \mathit{\delta }\end{array}\right]\left[\begin{array}{ccc}\cos \mathit{\gamma }& -\sin \mathit{\gamma }& \mathrm{0}\\ \sin \mathit{\gamma }& \cos \mathit{\gamma }& \mathrm{0}\\ \mathrm{0}& \mathrm{0}& \mathrm{1}\end{array}\right]\\ & {\underset{\mathrm{¯}}{v}}_{\mathrm{g}}=\left[\begin{array}{ccc}\cos \mathit{\delta }\cos \mathit{\gamma }& -\cos \mathit{\delta }\sin \mathit{\gamma }& -\sin \mathit{\delta }\\ \sin \mathit{\gamma }& \cos \mathit{\gamma }& \mathrm{0}\\ \sin \mathit{\delta }\cos \mathit{\gamma }& -\sin \mathit{\delta }\sin \mathit{\gamma }& \cos \mathit{\delta }\end{array}\right]{\underset{\mathrm{¯}}{v}}_{\mathrm{g}}.\end{array}\end{array}$$

The inverse transformation is simply given by the matrix transpose.

Using Eq. (A1), the components of the ambient velocity vector u_∞ in the nacelle-attached frame are readily found to be ${\underset{\mathrm{¯}}{u}}_{{\mathrm{\infty }}_{\mathrm{n}}}=u_{\mathrm{\infty }}\mathit{\{}\cos \mathit{\delta }\cos \mathit{\gamma },\sin \mathit{\delta },\sin \mathit{\delta }\cos \mathit{\gamma }{\mathit{\}}}^T$, where the scalar wind speed is $u_{\mathrm{\infty }}=\left|{\mathit{u}}_{\mathrm{\infty }}\right|$, while ${\underset{\mathrm{¯}}{x}}_{{\mathrm{n}}_{\mathrm{n}}}=\mathit{\{}\mathrm{1},\mathrm{0},\mathrm{0}{\mathit{\}}}^T$. Hence, it follows that

$$\begin{array}{}\text{(A2)}& {\displaystyle \frac{{\mathit{u}}_{\mathrm{\infty }}}{u_{\mathrm{\infty }}}}\cdot {\mathit{x}}_{\mathrm{n}}=\cos \mathit{\mu }=\cos \mathit{\delta }\cos \mathit{\gamma }.\end{array}$$

A2 Transformation from wake-deflection to ground frame of reference

The nacelle and wake-deflection frames share the same unit vector x_n=x_d, which corresponds to the rotor axis of rotation. The z_d unit vector is orthogonal to the plane composed by u_∞ and x_d, and therefore it can be written as ${\mathit{z}}_{\mathrm{d}}=z{\mathit{x}}_{\mathrm{d}}\times {\mathit{u}}_{\mathrm{\infty }}$, where z is a normalization scalar such that $\mathit{z}\cdot \mathit{z}=\mathrm{1}$. Performing the cross product and the normalization, one finds ${\underset{\mathrm{¯}}{z}}_{{\mathrm{d}}_{\mathrm{n}}}=z\mathit{\{}\mathrm{0},-\sin \mathit{\delta }\cos \mathit{\gamma },\sin \mathit{\gamma }{\mathit{\}}}^T$ and $z=\mathrm{1}/\sqrt{{\cos }^{\mathrm{2}}\mathit{\gamma }{\sin }^{\mathrm{2}}\mathit{\delta }+{\sin }^{\mathrm{2}}\mathit{\gamma }}=\mathrm{1}/\sin \mathit{\mu }$. A right-handed triad is completed by setting ${\mathit{y}}_{\mathrm{d}}=-{\mathit{x}}_{\mathrm{d}}\times {\mathit{z}}_{\mathrm{d}}$, which yields ${\underset{\mathrm{¯}}{y}}_{{\mathrm{d}}_{\mathrm{n}}}=z\mathit{\{}\mathrm{0},\sin \mathit{\gamma },\sin \mathit{\delta }\cos \mathit{\gamma }{\mathit{\}}}^T$. The transformation matrix between the wake-deflection and nacelle-fixed components is therefore readily obtained as

$$\begin{array}{}\text{(A3)}& {\underset{\mathrm{¯}}{v}}_{\mathrm{n}}=\left[\begin{array}{ccc}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& z\sin \mathit{\gamma }& -z\sin \mathit{\delta }\cos \mathit{\gamma }\\ \mathrm{0}& z\sin \mathit{\delta }\cos \mathit{\gamma }& z\sin \mathit{\gamma }\end{array}\right]{\underset{\mathrm{¯}}{v}}_{\mathrm{d}}.\end{array}$$

Finally, the transformation between wake-deflection and ground-fixed components follows by using Eqs. (A3) and (A1), which yields

$$\begin{array}{}\text{(A4a)}& {\displaystyle }\begin{array}{rl}{\underset{\mathrm{¯}}{v}}_{\mathrm{g}}=& \left[\begin{array}{ccc}\cos \mathit{\delta }\cos \mathit{\gamma }& \sin \mathit{\gamma }& \sin \mathit{\delta }\cos \mathit{\gamma }\\ -\cos \mathit{\delta }\sin \mathit{\gamma }& \cos \mathit{\gamma }& -\sin \mathit{\delta }\sin \mathit{\gamma }\\ -\sin \mathit{\delta }& \mathrm{0}& \cos \mathit{\delta }\end{array}\right]\\ & \left[\begin{array}{ccc}\mathrm{1}& \mathrm{0}& \mathrm{0}\\ \mathrm{0}& z\sin \mathit{\gamma }& -z\cos \mathit{\gamma }\sin \mathit{\delta }\\ \mathrm{0}& z\cos \mathit{\gamma }\sin \mathit{\delta }& z\sin \mathit{\gamma }\end{array}\right],\end{array}\text{(A4b)}& {\displaystyle }\begin{array}{rl}& {\underset{\mathrm{¯}}{v}}_{\mathrm{d}}=\left[\begin{array}{ccc}\cos \mathit{\delta }\cos \mathit{\gamma }& \mathrm{1}/z& \mathrm{0}\\ -\cos \mathit{\delta }\sin \mathit{\gamma }& z{\cos }^{\mathrm{2}}\mathit{\delta }\sin \mathit{\gamma }\cos \mathit{\gamma }& -z\sin \mathit{\delta }\\ -\sin \mathit{\delta }& z\sin \mathit{\delta }\cos \mathit{\delta }\cos \mathit{\gamma }& z\cos \mathit{\delta }\sin \mathit{\gamma }\end{array}\right]{\underset{\mathrm{¯}}{v}}_{\mathrm{d}}.\end{array}\end{array}$$

The inverse transformation is simply given by the matrix transpose.

Using Eq. (A4b), the longitudinal (given by Eq. 11) and lateral (sidewash, given by Eq. 13) flow velocity components at the streamtube outlet can be transformed into the corresponding longitudinal and lateral components in the ground frame:

$$\begin{array}{}\text{(A5a)}& {\displaystyle }& {\displaystyle \frac{v_{\mathrm{o},\mathrm{g}}}{u_{\mathrm{\infty }}}}={\displaystyle \frac{C_{\mathrm{T}}}{\mathrm{4}}}\cos \mathit{\delta }\sin \mathit{\gamma },\text{(A5b)}& {\displaystyle }& {\displaystyle \frac{w_{\mathrm{o},\mathrm{g}}}{u_{\mathrm{\infty }}}}={\displaystyle \frac{C_{\mathrm{T}}}{\mathrm{4}}}\sin \mathit{\delta }.\end{array}$$

Appendix B: Effect of a non-uniform axial induction on the rotor disk

As mentioned in Sect. 2.4, misalignment and shear cause a non-uniform distribution of the induction over the rotor disk. Following the classical approach used for helicopter rotors in forward flight (Johnson, 1995), the simplest model of non-uniform axial induction is the one expressed by Eq. (6), based on the 1P harmonics κ_1s and κ_1c. For the present application, however, it appears that the inclusion of these terms is not necessary. To show this, we consider the κ_1s term, which is triggered by the misalignment μ and results in the largest induction in the most downstream portion of the rotor disk. Figure B1 presents the loss factors η_T and η_P predicted by the model with (dotted red line) and without (solid blue line) the 1P sine term κ_1s. Differences appear to be negligible, especially for η_T. Neglecting κ_1s leads to a slight decrease in η_P as the misalignment increases, reaching a maximum difference of 0.21 % at $\mathit{\gamma }=-\mathrm{30}$°.

Figure B1 Thrust loss factor η_T (a) and power loss factor η_P (b) computed with the proposed model, with and without the 1P sine harmonic κ_1,s in Eq. (6). Results are computed for the IEA 3.4 MW reference wind turbine, subjected to a vertical shear k=0.2, operating at a tip speed ratio λ=8.5, and with a blade pitch angle θ=1°.

Appendix C: Experimental data set

Figures C1 and C2 report the blade pitch θ_p and tip speed ratio λ measured in experimental campaigns 1 and 2, respectively. The same quantities for experimental campaign 3 are reported in Fig. 19.

Figure C1 Experimental campaign 1. Blade pitch angle θ_p (a) and tip speed ratio λ (b) in Low-TI, Mod-TI, and High-TI inflows.

Figure C2 Experimental campaign 2. Blade pitch angle θ_p (a) and tip speed ratio λ (b) in the Mod-TI inflow.

Appendix D: Nomenclature

A	Rotor disk area
a	Axial induction
B	Number of blades
CD	Drag coefficient
CL	Lift coefficient
CL,α	Lift slope
CP	Power coefficient
CT	Thrust coefficient
$C_{\mathrm{T}}^{\prime }$	Modified thrust coefficient of Heck et al. (2023)
D	Rotor diameter
Fn	Normal force
Ft	Tangential force
K	Coefficient relating aerodynamic torque
	and squared rotor speed in control region II
k	Linear vertical wind shear coefficient
$\dot{m}$	Mass flux
p	Pressure
P	Power
Pd	Power demand (derating)
Q	Rotor torque
R	Rotor radius
r	Spanwise coordinate
T	Thrust force
u∞	Free-stream wind speed

u∞,hub	Free-stream wind speed at hub height
u	Longitudinal velocity component
un	Rotor-orthogonal velocity component
ut	Rotor-tangential velocity component
v	Lateral velocity component
x	Cartesian coordinate
y	Cartesian coordinate
z	Cartesian coordinate
α	Angle of attack
β	Blade twist angle
δ	Rotor tilt angle
ηP	Power loss factor
ηT	Thrust loss factor
γ	Rotor yaw angle
λ	Tip speed ratio
μ	Rotor total (true) misalignment angle
Ω	Rotor angular speed
φ	Inflow angle
Φ	Wind direction
ψ	Rotor azimuth angle
ρ	Air density
θ	Local pitch angle
θp	Blade pitch rotation at the pitch bearing
${\underset{\mathrm{¯}}{v}}_{\mathrm{f}}$	Components of vector v in frame f
(⋅)d	Quantity evaluated in the wake-deflection
	intrinsic frame of reference
(⋅)g	Quantity evaluated in the ground frame of
	reference

(⋅)n	Quantity evaluated in the nacelle frame of
	reference
(⋅)1c	Once-per-revolution cosine harmonic
(⋅)1s	Once-per-revolution sine harmonic
0P	Zeroth (constant) harmonic
1P	One-per-revolution harmonic
ALM	Actuator line method
BEM	Blade element momentum
CFD	Computational fluid dynamics
CTA	Constant temperature anemometry
FLORIS	FLOw Redirection and Induction in Steady
	State
LES	Large eddy simulation
rpm	Revolutions per minute
SOWFA	Simulator fOr Wind Farm Applications
TI	Turbulence intensity