# A CFD Study on High-Thrust Corrections for Blade Element Momentum Models

Federico Zilic de Arcos, Aidan Wimshurst, Richard H J Willden, Grégory

Pinon, Christopher R Vogel

► To cite this version:

Federico Zilic de Arcos, Aidan Wimshurst, Richard H J Willden, Grégory Pinon, Christopher R Vogel. A CFD Study on High-Thrust Corrections for Blade Element Momentum Models. Wind Energy, 2024, 10.1002/we.2937. hal-04674174

HAL Id: hal-04674174

https://hal.science/hal-04674174v1

Submitted on 21 Aug 2024

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

RESEARCH ARTICLE

OPEN ACCESS

# A CFD Study on High-Thrust Corrections for Blade Element Momentum Models

Federico Zilic de Arcos<fnref n="1" /> | Aidan Wimshurst² | Richard H. J. Willden³ | Grégory Pinon<fnref n="1" /> | Christopher R. Vogel³

**Correspondence:** Federico Zilic de Arcos (federico.zilic@univ-lehavre.fr)

**Received:** 13 July 2023 | **Revised:** 22 January 2024 | **Accepted:** 4 June 2024

**Funding:** This research project was supported in part by CONICYT PFCHA/BECAS CHILE DOCTORADO EN EL EXTRANJERO 2016/72170292. This project also received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement no. 101034329 and recipient of the WINNINGNormandy Program supported by the Normandy Region.

**Keywords:** blade element momentum | CFD | Mexico | Oxford benchmarking experiment | tidal rotors | turbines | turbulent wake model

## ABSTRACT

This paper presents a reanalysis of four axial-flow rotor simulation datasets to study the relationship between thrust and axial induction factor. We concentrate on high-thrust conditions and study variations in induction factor and loads across the span of the different rotor blades. The datasets consist of three different axial-flow rotors operating at different tip-speed ratios and, for one dataset, also at different blockage ratios. The reanalysis shows differences between the blade-resolved CFD results and a wide-spread empirical turbulent wake model (TWM) used within blade element momentum (BEM) turbine models. These differences result in BEM models underestimating thrust and especially power for axial-flow rotors operating in high-thrust regimes. The accuracy of BEM model predictions are improved substantially by correcting this empirical TWM, producing better agreement with blade-resolved CFD simulations for thrust and torque across most of the span of the blades of the three rotors. Additionally, the paper highlights deficiencies in tiploss modelling in common BEM implementations and highlights the impact of blockage on the relationship between thrust and axial induction factors.

## 1 | Introduction

The interest in engineering models for wind and tidal turbines arises from the need to perform extensive and complex analyses across a range of operating and environmental conditions. Examples include design optimisation studies, analysis of rotors operating under transient environmental loads, fluid-structure interaction studies where blade deformations play an important role in rotor thrust and performance, or rotors mounted on dynamic floating platforms. In all of these cases, the computational cost of blade-resolved computational fluid dynamics (CFD) simulations is often very high, typically limiting the number of modelled cases to just a few of relatively short duration.

In light of these computational limitations, engineering models, often based on simplifications of complex phenomena and empirical corrections, are a fundamental element of engineering practice. Blade element momentum (BEM) models are the default engineering models used for analysis of axial-flow rotors, being used widely in both academia and industry to calculate loads and power performance, with often thousands of BEM solutions considered during rotor design.

Most modern implementations of BEM models are based on the classical work of Glauert [1], a model that stems from the slipstream momentum theory developed by Rankine [2] and the blade element model proposed by Froude[3]. Glauert

<footnote n="1">¹Laboratoire Ondes et Milieux Complexes, Université Le Havre Normandie, Le Havre, France | ²Frazer-Nash Consultancy, Bristol, UK | ³Department of Engineering Science, University of Oxford, Oxford, UK</footnote>

<footnote>This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.</footnote>

<footnote>© 2024 The Author(s). *Wind Energy* published by John Wiley & Sons Ltd.</footnote>
presents the problem as a system of two-dimensional blade element and momentum equations, which are solved over a series of concentric and radially independent annuli across the swept area of a rotor represented by an actuator disc. For each annulus, the flow is described by axial and tangential momentum balances across the actuator, where the influence of the rotor is captured as an opposing force that changes the axial and tangential flow speeds, that is, it induces velocities that are represented by axial and tangential induction factors. The rotor blades are modelled as a collection of two-dimensional aerofoils whose effects are azimuthally averaged around the annulus. The local relative flow velocity then determines the angle of attack on the aerofoils and thus the forces and rotor performance [4].

The original derivation of the BEM method has been extended to account for fluid-dynamic effects that were not considered in its basic formulation, increasing the model's fidelity and predictive capabilities. Recent developments for axial-flow rotors include: spanwise modifications to the two-dimensional polar coefficients to account for rotational [5] and near-tip effects (e.g., Shen et al. [6], Wimshurst and Willden [7]); blockage [8]; yawed flows [9]; and the inclusion of other transient flow effects [10]. However, possibly two of the most notable advances were the early introduction of a correction for rotors with a finite number of blades (Prandtl [11], Betz [12]), and the correction of the momentum equation for high axial induction factors (Glauert [4]).

One-dimensional inviscid linear momentum theory can be used to model the reduction in flow speed at the turbine plane, $aV_∞$, that results from the resistive turbine thrust, $C_T$ (ratio of thrust per disc area to the dynamic pressure of the undisturbed flow). The analytical inviscid momentum model predicts a parabolic relationship between the thrust coefficient, $C_T$, and the axial induction factor, $a$, with a maximum rotor thrust coefficient of $C_T = 1.0$ at an axial induction factor $a = 1/2$. While good agreement in predictions of rotor performance is found between the analytical model and higher fidelity studies for lower induction factors [13], at higher induction factors the rotor wake becomes characterised by large vortices, flow recirculation and viscous mixing effects that are not captured in the inviscid model. This operational regime is termed the turbulent wake or turbulent windmill state [14]. Predictions of the analytical momentum model depart from experiments (e.g., Lock et al. [15], Parra et al. [16]), as shown in Figure 1, as well as from higher fidelity simulations (e.g., Sørensen and Kock [17], Madsen [18], Martínez-Tossas et al. [13]).

Glauert proposed a turbulent wake model (TWM) as an empirical correction to the momentum equation to estimate the rotor thrust coefficient at high axial induction factors based on the experimental data published by Lock et al. [15]. Following Glauert's work, other turbulent wake corrections have been proposed for BEM models, including Buhl [14], Burton et al. [19] and Spera [20], among others. Buhl's TWM, a parabolic function fitted at a critical (transition) axial induction factor $a_c = 0.400$, is almost identical to Glauert's except that it guarantees smoothness and continuity between the analytical and empirical thrust coefficient functions in the presence of tip-loss corrections. Burton et al.'s and Spera's models are linear

<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | Variation of thrust coefficient $C_T$ with axial induction factor $a$ for different momentum models with turbulent wake corrections, together with empirical data from Lock et al. [15] and Parra et al. [16]. Lock et al.'s empirical data are reproduced from Buhl [14].</figcaption></figure>

functions tangent to the analytical $C_T$ curve at $a_c = 0.326$ and $a_c = 0.200$, respectively. These models can be seen in Figure 1 alongside the data in Lock et al. reproduced from the work of Buhl [14], and the data in Parra et al. [16] extracted at different radial locations (25% to 92% of blade span) from experiments of the Model Experiments in Controlled Conditions (MEXICO) rotor [21].

Developing and using turbulent wake corrections consistent with the BEM method framework is challenging due to the inherent unsteadiness and three-dimensional nature of the flow in the turbulent wake regime, to which the scatter in Lock et al.'s data has been attributed [22]. Similarly, the scatter in the data of Parra et al. [16] at high induction factors was partially attributed, by the authors of that study, to the loading and wake regime. The difficulties of extracting flow data and associated uncertainty from experiments must also be acknowledged. These challenges arise not only due to the problem of sampling the flow around three-dimensional blades and processing the information in terms that are compatible with the BEM framework, such as axial induction factors, but also due to the dependence of the induced velocities on axial and radial coordinates. In particular, Parra et al. [16] showed that radial location has a substantial impact on the relationship between thrust coefficient and axial induction factor for the New MEXICO experiments, while Zilic de Arcos et al. [23] showed the importance of the methods and sampling locations used for interpolating flow velocities at the rotor plane.

Lock et al.'s data and the empirical momentum corrections derived therefrom relate the whole-rotor thrust coefficient to the rotor axial induction factor. As BEM models are solved across a series of independent annular streamtubes, an assumption is made that the same TWM relationship between whole-rotor thrust and axial induction factor can also be applied at each annulus. Moriarty and Hansen [24] note that this assumption is necessitated by the paucity of rotor data with radial dependency that would allow the development of alternative models. The increasingly widespread use of CFD has since enabled higher fidelity simulation of highly loaded rotors which
can provide additional data for TWMs, particularly at the values of axial induction factor where the TWM is likely to be important.

Sørensen and Kock [17] and Madsen [18] used actuator disc CFD (AD-CFD) simulations to study rotors operating at high axial induction factors. Similar to BEM, AD-CFD models the rotor as a disc, a finite surface embedded within a CFD domain, where a pressure discontinuity occurs. Unlike traditional BEM models, however, an AD-CFD model captures the interactions between different annular sections, radial flows and other three-dimensional flow phenomena that affect rotor and wake. The results of Sørensen and Kock [17] and Madsen [18] divert from the analytical $C_T(a)$ function, showing similar trends to those observed by Lock et al. [15] and Parra et al. [16]. In addition, they show changes in axial speeds across discs with uniform loading, suggesting changes in the relationship between thrust and axial induction factors across different radial locations.

Blade-resolved CFD simulations have also been used to study axial-flow rotors. In particular, Reynolds-averaged Navier-Stokes (RANS) CFD models have been used extensively and validated for different rotor designs and operational regimes, including high-thrust conditions. These models are a compromise between model accuracy and computational cost; while they are not capable of accurately capturing transient flow dynamics, they can reproduce mean forces and flow quantities in agreement with higher order turbulence closure models [25]. Examples that compare RANS CFD models of axial-flow rotors with experiments, showing good agreement in terms of integrated and spanwise forces, are the works of Sørensen et al. [26], Wimshurst and Willden [27], Länger-Möller et al. [28] and Willden et al. [29].

In this paper we re-evaluate BEM theory in conjunction with observations made from blade-resolved CFD simulations of different axial-flow rotors. The results gathered from CFD are used to assess empirical momentum corrections, and evidence is presented that allows improved engineering modelling based on BEM theory for highly loaded axial-flow rotors.

The study is based on the analysis of four blade-resolved CFD datasets that correspond to three different rotors and that were generated originally for different studies. The CFD results were reanalysed so that they can be compared with the momentum component of BEM models using an analysis technique referred to as line-average [30]. This process allows some of the shortcomings and causes of accuracy losses in BEM when modelling highly loaded rotors to be identified, as well as highlighting potential improvements for the model. Following this analysis, a modified TWM is proposed and tested, showing an improvement in the power and thrust predictions made with the BEM model across the blade spans of the three different rotors.

The final section of the paper addresses the impact of blockage on the turbulent wake state. Blockage is defined as the aero- or hydrodynamic interaction between a body (e.g., an aerofoil, a ship, or an axial-flow rotor) and its surroundings (e.g., the walls in a confined flow), and normally affects rotor operation by

increasing power and thrust due to flow confinement. The evidence in this section shows that blockage can have a significant impact on the relationship between thrust and axial induction factors.

## 2 | Methods

### 2.1 | Rotor Characteristics and Available Data

Three different rotors and four CFD datasets were analysed in this study. The MEXICO rotor is a three-bladed 4.5m diameter wind turbine designed using three different aerofoils (DU-91-W2-250, RISØ A1-21 and NACA 64-418) along the blade span with geometric transitions between them. The MEXICO rotor simulation dataset [7] used in this study considers a fixed blade pitch angle of -2.3 degrees, a blockage $B = 0.01$ (defined as the ratio of the rotor swept area to the channel cross-sectional area), and tip-speed ratios $\lambda \in \{4.17, 5.56, 6.67, 8.34, 10.00, 13.3\}$, with $\lambda = \omega R/V_\infty$, $\omega$ the rotational speed, $R$ the rotor radius and $V_\infty$ the flow speed. The simulations were performed and the modelling approach was validated originally by Wimshurst and Willden [7] to analyse and develop tip-loss correction models for actuator line simulations.

The Sch15B is a 20-m diameter three-bladed tidal rotor originally designed by Schluntz and Willden [31] to operate at a blockage ratio $B = 0.196$ at an optimum tip-speed ratio $\lambda = 5.5$. The rotor was designed using a single RISØ A1-24 aerofoil along the span, with aerofoil section quarter-chords aligned along the blade axis. The original root section defined by Schluntz and Willden [31] was modified by Wimshurst and Willden [32] and by Zilic de Arcos et al. [33] to allow for a less substantial transition between the foil sections and the cylindrical cross-section at the attachment point on the nacelle [34].

Two simulation datasets of the Sch15B rotor were used. The main dataset corresponds to the simulations of the rotor at a flow speed of $V_\infty = 4.5$ m/s, at seven different tip-speed ratios ($\lambda \in \{3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0\}$) with a blockage of 1% and with different blade deformations, as described by Zilic de Arcos et al. [35] and Zilic de Arcos et al. [36]. Only the cases with undeformed rotor blades were considered in this study. A secondary dataset corresponds to the rotor simulated at the same flow speed but for different isotropic blockage ratios ($B \in \{0.01, 0.05, 0.10, 0.20, 0.40\}$), with each case simulated at four different tip-speed ratios ($\lambda \in 4.0, 5.0, 6.0, 7.0$), as described by Zilic de Arcos et al. [37].

Finally, the rotor designed, used and tested in the Tidal Turbine Benchmarking Project and labelled as OxB is a 1.6m diameter, three-bladed tidal rotor [38] designed by Tucker Harvey et al. [39] using the method of Schluntz and Willden [31]. The rotor uses a NACA63-415 profile projected radially on each blade section. The experimental campaign was performed at QinetiQ's Haslar towing tank, which has a test section of 12.2m by 5.4m providing a low blockage of 3.05%. The rotor was installed with a tip clearance of 0.354m relative to the free surface. The data used in this study correspond to the CFD simulations conducted by the Laboratoire Ondes et Milieux Complexes to participate in the blind prediction exercise workshop [29, 39]. The
cases analysed match the turbulence intensity of $T_i = 3.1\%$, the constant inflow speed $V_\infty = 0.9207$ [m/s], and five tip-speed ratios ($\lambda \in 4.46, 5.37, 5.82, 6.37, 7.37$). The simulations assume an equivalent but isotropic blockage ratio of 3.05%, and neglect interactions with the free surface and support structure.

The solidity and geometric twist of the three rotors, $\sigma$ and $\beta$ respectively, are presented in Figure 2, with solidity defined as $\sigma = cN_b/2\pi r$, where c is the local chord, $N_b$ the number of blades, and r the local radius. Twist is defined positive in the sense that reduces the flow's angle of attack to the blade. A comparison of the rotor's general characteristics is presented in Table 1.

## 2.2 | Computational Fluid Dynamics Models

The three rotors were simulated for different studies using steady-state Reynolds-averaged Navier-Stokes (RANS) models, with the Multiple Reference Frame approach [40] for blade rotation and the $k-\omega$ SST model for turbulence closure [41]. The Sch15B simulations were performed using the commercial code ANSYS Fluent 19.0 with second-order upwind discretisation schemes and a coupled pressure-velocity algorithm. The MEXICO and OxB rotors were simulated with OpenFOAM (2.3.1 and v2006, respectively) with the SIMPLE algorithm

for pressure-velocity coupling and second-order discretisation schemes. The simulations exploited the azimuthal periodicity of the rotors in uniform flow by simulating only one-third of the domain with periodic boundaries, and simulation residuals were reduced by at least 5 orders of magnitude for both the velocity and turbulence scalars.

The three rotors were discretised using a block-structured approach using entirely hexahedral cells. The first cell heights were chosen to target $y^+$ in the range of $0 < y^+ < 5$ for the MEXICO and OxB simulations, while a $30 < y^+ < 300$ range was used for the Sch15B. Mesh convergence details, validation of the simulation strategies, and further details can be found in inWimshurst and Willden [7], Zilic de Arcos [34] and Willden et al. [29] for the MEXICO, Sch15B and OxB rotors, respectively.

To determine angles of attack and axial induction factors, the line-average flow sampling method was used. This is a flow-sampling technique for CFD computations used to determine the flow-speed at different spanwise locations along a rotor blade. The method was first introduced by Josh et al. [30] to determine the undisturbed flow velocity that defines the angle of attack on a foil section. This is achieved through the implicit removal of the circulatory velocity component $\bar{U}_\Gamma$, induced by the blade's bound circulation $\Gamma$, from the sampled flow field $\bar{U}_M$. The velocity of the flow measured at an arbitrary point in a CFD

<figure><img src="image_4.png" /><figcaption>**FIGURE 2** | Spanwise distributions of blade solidity (left) and twist (right) for the MEXICO (blue dash dot), Sch15B (red dash), and OxB (black dot) rotors. Note that the MEXICO twist distribution is shown without additional pitch setting angle.</figcaption></figure>

<table><caption>**TABLE 1** | Key design and operational parameters of the Sch15B, MEXICO, and OxB rotors, displaying type of turbine, number of blades $N_b$, rotor and hub diameters $\emptyset$, original design tip-speed ratio $\lambda_D$ and velocity $V_D$, simulated tip-speed ratios $\lambda$ and flow velocities $V_\infty$, chord-based Reynolds numbers $Re_{0.75}$ evaluated at $r/R = 0.75$, and range of simulated thrust coefficients $C_T$.</caption><thead><tr><th>Turbine</th><th>Type</th><th>$N_b$</th><th>Rotor Ø (m)</th><th>Hub Ø (m)</th><th>$\lambda_D$</th><th>$V_D$ (m/s)</th><th>$\lambda$</th><th>$V_\infty$ (m/s)</th><th>$Re_{0.75}$</th><th>$C_T$</th></tr></thead><tbody><tr><td>Sch15B</td><td>Tidal</td><td>3</td><td>20.00</td><td>3.00</td><td>5.5<sup>a</sup></td><td>2.0</td><td>3.00-9.00</td><td>4.5000</td><td>~16-31 × 10<sup>6</sup></td><td>0.62-1.07</td></tr><tr><td>MEXICO</td><td>Wind</td><td>3</td><td>4.50</td><td>0.54</td><td>6.7</td><td>15.0</td><td>4.17-13.30</td><td>7.520-24.000</td><td>~0.5 × 10<sup>6</sup></td><td>0.40-1.10</td></tr><tr><td>OxB</td><td>Tidal</td><td>3</td><td>1.60</td><td>0.20</td><td>6.0<sup>b</sup></td><td>1.0</td><td>4.46-7.37</td><td>0.9207</td><td>~0.1-0.4 × 10<sup>6</sup></td><td>0.72-1.01</td></tr></tbody></table>

<footnote><sup>a</sup>Design tip-speed ratio for a blockage of $B = 0.1960$.</footnote>

<footnote><sup>b</sup>Design tip-speed ratio for a blockage of $B = 0.0305$.</footnote>
simulation, $\overline{U}_M$, is assumed to be the superposition of the approach flow velocity, $\overline{U}$, and the local blade-induced circulatory velocity due to the foil's bound circulation, $\overline{U}_Γ$:

$$ \overline{U}_M = \overline{U} + \overline{U}_Γ. \quad (1) $$

For the case of a rotating blade, this approach velocity $\overline{U}$ defines the angle of attack to a foil section and includes contributions from the free stream axial-flow speed $V_∞$, the rotational speed $ωr$, as well as non-circulatory induced velocity components resulting from momentum effects.

The line-average method removes $\overline{U}_Γ$ by averaging $\overline{U}_M$ over a closed circular path centred about the blade's centre of pressure line. The averaging of the flow velocity at points located on opposite sides of the circular segment cancels the circulation-induced velocity $U_Γ$. For a three-dimensional blade, the sampling is conducted over a circular path on a surface at the constant radial coordinate of the corresponding blade section to maintain consistency with the rotational reference frame seen by each blade section. Schematics of the method are shown in Figure 3 for a foil section and for a three-dimensional blade.

The line-average sampling method was extensively tested and compared with other flow-sampling techniques by Zilic de Arcos et al. [23]. The method showed robustness in the determination of angles of attack and axial induction factors, good agreement with other methods that had been described in the literature, and low sensitivity to sampling distance. Thus, this method was selected for the analysis of the relationship between thrust and blade-induced velocity as well as angles of attack due to its practicality and robustness. A local sampling radius of 1.1 chords was used.

The axial velocities induced by the blade in the flow are described through the axial induction factor $a$, defined along the span of the blade as

$$ a(r) = 1 - U_x(r)/V_\infty \quad (2) $$

where $U_x$ is the axial component of $\overline{U}$, as determined using the line-average method.

The spanwise force distributions $\overline{F} = [F_x, F_θ, F_Γ]$, with $x$, $θ$ and $r$ the axial, tangential and radial components, were extracted from the CFD simulations by integrating pressure $p$ and shear stress $\overline{S}$ around each foil section $C$, with surface normals $\overline{n}$:

$$ \overline{F}(r) = \oint (p\overline{n} + \overline{S})dC \quad (3) $$

Unless otherwise stated, thrust coefficients in this paper refer to local rather than integrated values. The local thrust coefficient $C_t$ is determined from the CFD simulations, for a specific radial location, as follows:

$$ C_t(r) = \frac{F_x(r)N_b}{1/2\rho V_\infty^2 2\pi r} \quad (4) $$

with $ρ$ the flow density, and $F_x(r)$ expressed as force per unit span. Finally, and where explicitly stated, the integrated power and thrust coefficients are defined as

$$ C_P = \frac{P}{1/2\rho V_\infty^3 A} \quad (5) $$

$$ C_T = \frac{T}{1/2\rho V_\infty^2 A} \quad (6) $$

where $P$ and $T$ are the whole-rotor integrated power and thrust, respectively, and $A$ is the rotor swept area.

## 2.3 | Blade Element Momentum Model

We used the BEM implementation of Ning [42], proposed as an alternative formulation to Glauert's original BEM model, and which the reader is referred to the original publication for full details of derivation and implementation of the method.

Instead of manipulating the blade element and momentum equations to derive a system of two equations for the axial and

<figure><img src="image_5.png" /><figcaption>FIGURE 3 | Two-dimensional sketch showing the symmetry of the velocity induced by a foil's bound vorticity $U_Γ$ (left), and a three-dimensional realisation of the line-average method [30] applied at different spanwise sections of a rotor blade (right). Inset contours show flow speed.</figcaption></figure>
tangential induction factors, as in Glauert's original derivation, Ning proposed a parametrisation of the equations based only on the inflow angle $\phi$, the angle between the inflow velocity and the rotational plane. The derivation starts from the definition of the thrust coefficient, which is defined not as a whole-rotor coefficient but for each annulus, for both the blade element (Equation 7) and momentum (Equation 8) models:

$$C_t(r) = \left( \frac{1-a}{\sin \phi} \right)^2 C_x \sigma \quad (7)$$

$$C_t(r) = 4a(1 - a)T_L \quad (8)$$

where $a$ is the axial induction factor, $C_x$ the normal force coefficient, $\sigma$ the local blade solidity, and $T_L$ a tip-loss correction (Note that $\sigma$ and $T_L$ are referenced as $\sigma'$ and $F$, respectively, in the work of Ning [42.]), with all variables corresponding to the annular section at $r$. Equating Equations (7) and (8) and solving for $a$ leads to an expression for the axial induction factor expressed as a function of $\phi$:

$$a(r) = a(\phi) = \frac{\kappa(\phi)}{(1 + \kappa(\phi))} \quad (9)$$

with

$$\kappa(\phi) = \frac{\sigma C_x}{4T_L \sin^2 \phi} \quad (10)$$

Analogously, the tangential induction factor $a'$ is defined by equating the torque evaluated by the momentum and blade element theories:

$$a'(r) = a'(\phi) = \frac{\kappa'(\phi)}{(1 - \kappa'(\phi))} \quad (11)$$

$$\kappa'(\phi) = \frac{\sigma C_{\theta}}{4T_{L} \sin \phi \cos \phi} \quad (12)$$

where $C_{\theta}$ is the tangential force coefficient. The normal and tangential force coefficients are found by resolving the local aerofoil lift and drag coefficients into the axial and tangential directions, respectively.

The axial and tangential inductions, which are both functions of the inflow angle, are related by the following geometric relationship:

$$\tan(\phi) = \frac{1 - a}{(1 + a')\lambda_r} \quad (13)$$

where $\lambda_r = \omega r / V_{\infty}$ is the local-speed ratio. Equation (13) can be rearranged to form a residual function, from which a solution for $\phi$ may be obtained:

$$f(\phi) = \frac{\sin(\phi)}{(1 - a)} - \frac{\cos(\phi)}{(1 + a')\lambda_r} = 0 \quad (14)$$

Ning's implementation also includes Prandtl's tip-loss correction and Buhl's [14] empirical momentum model for high axial inductions when $a > 0.40$ [42], which are discussed further in the following sections.

The advantage of this implementation of the BEM model lies in Ning's mathematical proof of the existence of an angle $\phi$ that fulfils Equation (14) and in the solution of a single residual

equation rather than a system of equations for the axial and tangential induction factors, as in other derivations of the BEM model (e.g., Burton et al. [19]). Solving a single equation enables the use of fast and robust root-finding algorithms, guaranteeing convergence to a solution for each annular streamtube along the entire span of the blade. Further details of the model derivation, implementation and the mathematical proof for the existence of a solution $\phi$ to Equation (14) can be found in the original work [42].

## 2.4 | Empirical Momentum Model

The general expression for Buhl's momentum model consists of a function with two branches that replaces Equation (8) in the BEM model. The branches are separated at a critical induction factor $a_c$. When the induction factor is small, the momentum equation corresponds to the analytical model of Equation (8), whereas when the induction factor is above $a_c$, an empirical formula is used:

$$C_t(r) = \begin{cases} 4a(1 - a)T_L & \text{if } a < a_c, \\ b_0 + b_1a + b_2a^2 & \text{if } a \ge a_c \end{cases} \quad (15)$$

with $b_0, b_1$ and $b_2$ empirical constants. The continuity and smoothness constraints are imposed by enforcing the value and first derivative between the two models in Equation (15) at $a_c$, where the transition between functions occurs. These continuity and smoothness constraints lead to

$$b_2 = \frac{b_0}{a_c^2} - 4T_L \quad (16)$$

$$b_1 = 4T_L - 8T_L a_c - 2b_2 a_c \quad (17)$$

with $b_0$ a free variable that can be fitted to empirical data or determined otherwise. Buhl, instead of directly fitting the data, closed the system by further assuming that $C_t = 2$ at $a = 1$, and established a value of $a_c = 0.40$ for the transition between models, leading to $b_0 = 8/9$.

Despite the widespread acceptance of Buhl's TWM, this does not guarantee accurate predictions for highly loaded rotors. Wind turbines normally operate at relatively low induction factors (e.g., the DTU 10 MW rotor operates at an axial induction factor on the order of 0.30 at design tip-speed ratio [43], the AVATAR rotor at around 0.24 [43], and the MEXICO rotor at approximately 0.26 [23]), an operational regime in which this empirical TWM is unlikely to be triggered. However, empirical evidence suggests that the onset of the turbulent wake state should occur at a lower $a_c$, as will be shown later. Furthermore, as the tip-speed ratio (TSR) increases (and, thus, thrust), wind rotors are prone to large errors when modelled using BEM methods, as shown in the following sections. This is a problem that also affects tidal rotors, which are conceptually similar to wind turbines but are often designed to operate at higher thrust regimes.

## 2.5 | Lift and Drag Coefficients

The aerofoil lift and drag coefficients used in this study were calculated using two-dimensional RANS CFD
simulations with a $k-\omega$ SST turbulence closure. Block structured meshes were generated for the foils used in each rotor, and simulated with a target $y^+ < 1$ for all cases. Polars for the Sch15B and OxB rotors were simulated for different Reynolds numbers using the turbulent kinetic energy and length scale corresponding to the blade-resolved CFD and tank conditions. The chords of the two-dimensional foils were set based on the spanwise average chords of the three-dimensional rotors.

For the two tidal rotors, Reynolds-number interpolation of lift and drag data was used in the BEM computations due to the large variations in relative flow speed observed for the different TSRs. These variations principally occur as the rotational speed was used to modify the TSR. Contrarily, for the MEXICO rotor, each of the different aerofoils was simulated at a single Reynolds number, as described in Wimshurst and Willden [7]. Reynolds number interpolation was not considered necessary as the rotor experiments were performed with a fixed rotational speed and, thus, Reynolds number variations were limited in comparison to the other two rotors (see Table 1). The polar coefficients used in the study are provided in Appendix A.

<figure><img src="image_6.png" /><figcaption>**FIGURE 4** | Comparison of thrust variation with axial induction factor with Lock's experimental data (black crosses, reproduced from Buhl [14]) and available CFD data at four radial locations, $r/R = 0.6, 0.7, 0.9, 0.98$. The CFD data include simulations of the Sch15B (red triangles), the MEXICO (blue triangles), and the OxB turbines (green circles). The analytical and empirical model curves are presented without tip-loss correction to highlight the impact of tip proximity on the empirical results.</figcaption></figure>

# 3 | Results and Discussion

## 3.1 | Analysis and Correction of the Turbulent Wake Model

### 3.1.1 | General CFD Results

The CFD results for the MEXICO, Sch15B ($B = 0.01$), and OxB rotors are shown in Figure 4 for four different radial locations. The plots compare the simulation data with the experimental data in Lock et al., the analytical momentum equation, Buhl's TWM, and a modified TWM. The CFD data are presented at four different radial locations: from the mid-span region of the blade ($r/R = 0.60 - 0.70$) to close proximity with the tip ($r/R = 0.90 - 0.98$). The momentum models are presented without tip-loss corrections to highlight the impact of tip proximity in the CFD results.

The OxB simulations include a modest blockage $B = 0.0305$ that could affect the relationship between the thrust coefficient and axial induction factor. However, no significant influence was observed in the results presented when compared with the MEXICO and Sch15B results. This is further reinforced by the
results discussed in Section 3.3, where only small differences were observed in the momentum plots for the Sch15B rotor at $B = 0.01$ and $B = 0.05$ (Figure 10). Thus, the impact of blockage for the OxB rotor in the thrust-induction relationship is considered negligible for this part of the analysis.

The plots that correspond to the inboard sections in Figure 4 show that divergence between the analytical momentum function and the empirical data occurs at a critical induction factor $a_c \approx 0.20$, earlier than the transition $a_c = 0.40$ specified by Buhl. The relatively small difference between the CFD data and Buhl's momentum function below $a_c = 0.40$, however, suggests that effects on BEM modelling would be limited at these relatively low induction factors. The results at $r/R = 0.90$ and $r/R = 0.98$ show that the sectional thrust coefficients decrease as the tip of the rotor is approached. This occurs due to the tip-loss effects and would normally be accounted for, within most BEM implementations, using a tip-loss factor $T_L$.

The CFD data also show less scatter than the results of Lock et al. [15], possibly due to the more controlled flow conditions and flow sampling achievable for the simulated flow, and the nature of the simulations (i.e., steady-state RANS turbulence modelling approach and the consequential lack of discrete velocity fluctuations). Nevertheless, and although generalisation is limited by the available simulations, the CFD data align relatively well with each other for the three different rotors at different tip-speed ratios, and at all radial locations, despite the rotors having different diameters and design characteristics. This suggests some generality in the results and their potential applicability to rotors of larger dimensions, including MW-scale wind turbines.

### 3.1.2 | Modifications to the Turbulent Wake Model

A simple modification to Buhl's empirical model is proposed to fit the available data. This model was derived using axial induction factors and thrust coefficients from the Sch15B dataset without blockage, the MEXICO data, and, owing to the limited blockage effects observed, the OxB data.

The modified TWM is fitted to the data at $r/R = 0.60$ with $a_c = 0.17$. This region of the blade is considered to be relatively free from root- and tip-induced three-dimensional flow effects, and simulated data points are available across a broad range of axial induction factors. Thus, fitting the free variable of the empirical equation in Equation (15) to the aforementioned data using the least-squares method leads to

$$b_0 = 0.0705 \qquad (18)$$

The modified empirical model for the turbulent wake regime is plotted in Figure 4 alongside the data at the different radial locations and Buhl's original TWM.

As shown for $r/R = 0.6$ and $r/R = 0.7$, the modified TWM provides a better agreement with the CFD data than Buhl's TWM, addressing the discrepancy for axial induction factors between $0.20 \lesssim a \lesssim 0.60$. The agreement between the empirical model

without any tip corrections and the CFD data points, however, is reduced as the tip is approached.

A comparison between Buhl's original model and the recalibrated TWM is shown in Table 2.

### 3.1.3 | Tip Proximity Effects

Ning's BEM implementation employs Prandtl's tip-loss correction model as a vertical scaling factor for the empirical momentum correction (Equation 15). The $T_L$ factor varies across the blade span and takes values from 1 (where no corrections are required) to 0. This modelling approach is consistent with the trends observed from the CFD simulations, as shown in Figure 5. The left plot shows the $C_L$ variation with $a$ at different radial locations for the three rotors, alongside the modified TWM presented previously ($b_0 = 0.0705$) scaled by an empirically fitted tip-loss factor $T_L$, as shown in Equations (15)–(17). This empirical $T_L$ is determined by a best-fit to the available data for each blade section.

The left plot in Figure 5 highlights that the use of this empirical scaling factor $T_L$ maintains the level of agreement between the TWM and the CFD data points, as the tip is approached, at similar levels as those observed for the inboard sections, without an evident influence of rotor design or TSR. This suggests that this form of correction for the momentum equations in the context of BEM models, such as the one proposed by Ning [42], has broad applicability and is appropriate.

The right plot in Figure 5 shows the best-fitted empirical $T_L$ as a function of blade radius, for all available radial sections, complementing the data from the left plot. The figure compares this empirically derived factor with Prandtl's tip-loss correction. Prandtl's model is a function of $\lambda$; thus, it is shown as shaded regions calculated for the Sch15B and MEXICO rotors at a range of different tip-speed ratios. The plot highlights significant discrepancies between the best-fitted empirical $T_L$ and Prandtl's model for most TSRs. Particularly, while the empirically derived $T_L$ is not sensitive to changes in TSR, the correction for a finite number of blades that Prandtl derived moves further outboard as $\lambda$ is increased, thus increasing the differences to the empirical $T_L$ at outboard locations.

The results shown in Figure 5 suggest that tip-loss factors applied to the momentum model should not be dependent on the tip-speed ratio or rotor design. This highlights a potential source of discrepancies between blade-resolved CFD and BEM

<table><caption>**TABLE 2 |** Model summary and comparison between Buhl's and the modified TWM.</caption><thead><tr><th></th><th>Buhl's TWM</th><th>Modified TWM</th></tr></thead><tbody><tr><td>a<sub>c</sub></td><td>0.4000</td><td>0.1700</td></tr><tr><td>b<sub>0</sub></td><td>8/9</td><td>0.0705</td></tr><tr><td>b<sub>1</sub></td><td>4T<sub>L</sub> - 8T<sub>L</sub>a<sub>c</sub> - 2b<sub>2</sub>a<sub>c</sub></td><td></td></tr><tr><td>b<sub>2</sub></td><td>b<sub>0</sub>/a<sub>c</sub><sup>2</sup> - 4T<sub>L</sub></td><td></td></tr><tr><td>C<sub>t</sub>(r)</td><td>4a(1 - a)T<sub>L</sub> if a &lt; a<sub>c</sub>,<br>b<sub>0</sub> + b<sub>1</sub>a + b<sub>2</sub>a<sup>2</sup> if a &ge; a<sub>c</sub></td><td></td></tr></tbody></table>
<figure><img src="image_7.png" /><figcaption>**FIGURE 5** | Left: Best fit of the scaling factor of the momentum function $T_L(b_1 = -1.561, a_c = 0.17)$ to the CFD data at different radial locations. Right: comparison of Prandtl's tip-loss factor for a range of TSRs with an empirical $T_L$ factor fitted to the CFD data. Prandtl's tip-loss factors are calculated for the Sch15B rotor for $4.0 \le \lambda \le 8.0$ (orange shaded area) and for the MEXICO rotor at TSRs in $6.67 \le \lambda \le 13.0$ (blue shaded area). The arrow indicates the direction in which the Prandtl correction evolves as TSR increases.</figcaption></figure>

models near the tip and especially at high tip-speed ratios when Prandtl's tip-loss model is used. The consequences of this are further discussed in Section 3.2.

## 3.2 | Evaluation of the Modified Turbulent Wake Model

Modifying the TWM affects the BEM predictions for thrust and torque on the three rotors. Figure 6 shows a comparison between rotor thrust and power coefficients calculated with Ning's BEM model using either Buhl's TWM or the modified TWM, and the blade-resolved CFD results for verification. Both BEM models follow the original implementation that includes Prandtl's tip-loss correction applied as a scaling factor to the momentum equations.

A first inspection of the integrated power and thrust reveals an apparently and misleading good agreement between Buhl's BEM and the CFD results for $C_T$ for the three rotors at the mid-range of $\lambda$ (relative error in the order of ~0.3% to 5.0%) which is adversely affected by using the modified TWM (relative error increases to a range of ~3.2% to 13.0%). $C_p$ shows a relatively better agreement with the CFD for the modified TWM over the low- to mid-range of $\lambda$.

The integrated quantities do not necessarily show or assess the full predictive quality of the BEM models. This can be further considered through the spanwise force distributions plotted in Figure 7. In this figure, different tip-speed ratios are plotted for each rotor, comparing the two BEM models and the CFD results. The small differences in integrated $C_T$ between the standard BEM and the CFD simulations observed for the MEXICO and Sch15B rotors are a consequence of a substantial overestimation of the loads near the tip (i.e., due to deficiencies in tip-loss modelling) in conjunction with a consistent underestimation of the loads through the mid and inboard sections of the blade. The two opposing errors negate each other, leading to a coincidental and misleading agreement in $C_T$ between the original TWM BEM and the CFD results.

The modified TWM shows better agreement across the mid and inboard sections ($r/R \le 0.70$) and an overestimation of thrust

near the tip which increases with $\lambda$. The discrepancy near the tip is, again, associated with limitations in tip-loss modelling. As discussed previously, the onset of Prandtl's correction moves further outboard as $\lambda$ increases, reducing the correction applied in the outboard region at higher $\lambda$. This explains the increasing difference in integrated thrust between the CFD and the modified TWM results for these two rotors, despite improved predictions across the inboard regions of the blade.

Another source of discrepancy that affects both TWMs is the use of two-dimensional lift and drag coefficients. As the tip of a blade is approached, the pressure difference between both sides of the blade tends to equalise, with lift and drag coefficients increasing and decreasing, respectively, thus modifying the shape of the $C_L$ and $C_D$ functions beyond the changes in angles of attack [44]. These changes in polar coefficients, not captured by typical BEM models nor by the model used in this paper, are likely to drive a further drop in local power and thrust near the tip.

A substantial improvement in the predictions of the spanwise torque distributions is observed using the modified TWM. While Buhl's TWM under-predicts torque across the entire span, the modified TWM shows better agreement with the CFD results across the mid and inboard sections of the blade ($r/R \le 0.70$). Following the same trend as for thrust, differences are observed between the CFD and the modified TWM near the tip, which ultimately leads to an overestimation of integrated torque, and therefore power coefficient, at the higher tip-speed ratios.

To extend the previous analysis, the relative errors in integrated thrust and torque between the BEM and CFD results, at mid-span sections of the blade, are shown in Figure 8 for a range of TSRs for the three rotors. The region of the blade where $r/R \in [0.4, 0.7]$ was selected due to the limited influence of root- or tip-loss effects. The figure highlights the substantially improved predictions produced with the modified TWM at mid-span of the blades across most TSRs, with the original BEM showing a greater propensity to under-predict both thrust and torque.

The CFD results for the OxB rotor shown in Figure 7, do not exhibit the same trend as the Sch15B and MEXICO rotor, despite
<figure><img src="image_8.png" /><figcaption>**FIGURE 6** | Rotor thrust and power coefficients calculated for the three rotors with blade-resolved CFD (red square), BEM with Buhl's turbulent wake empirical model [14, 42] (black cross), and BEM with the modified TWM (blue triangle). The shaded region highlights the difference between the results calculated with the original and modified TWMs.</figcaption></figure>

better agreement between the CFD and the modified TWM results for the whole-blade power and thrust. This is attributed to the mild effect of blockage in the simulated domain, where the blockage ratio was *B* = 0.031, which is expected to increase loads relative to fully unblocked analyses, as in the BEM solutions.

To provide an indication of whether blockage could be the cause for the difference in trend between the OxB and the other two rotors, we apply Bahaj's blockage correction [45] to the CFD results to estimate forces in a blockage-free scenario. This is a method extensively used to compare integrated thrust and power [37] for cases with different blockage. We used the method as a constant factor applied across the blade. The results, presented in Figure 9, show the same trends observed between the CFD and modified TWM for the other rotors: good agreement across the mid and in-board sections of the blade and an overestimation of thrust and torque towards the tip. Nevertheless, the use of Bahaj's correction as a constant factor across the blade has not been described in the literature to the best of our knowledge and needs to be tested in more detail before this indicative assessment can be substantiated.

## 3.3 | Effect of Blockage on the Turbulent Wake State

For axial-flow rotors, blockage generally results in increased thrust and power for a given undisturbed approach velocity $V_∞$

and rotational speed due to the imposed flow constraints that result in an increased mass flow through the rotor.

Blockage is a relevant factor in the design and operation of tidal rotors, rotor farms or fences, and in laboratory testing. However, most BEM models (e.g., Ning [42] or Burton et al. [19]) have been developed for wind energy applications and do not consider blockage effects. Vogel et al. [8] proposed a BEM implementation that incorporates blockage in its formulation through modifications to the momentum equations. This confined BEM incorporates Prandtl's tip-loss correction model and the Buhl [14] turbulent wake model. The authors showed an improvement in the predictive capabilities of their model compared with a standard BEM implementation when modelling a tidal rotor with blockage, although they recognise uncertainty regarding the turbulent wake state and the need for further studies in that area.

Analysis of the CFD results of Zilic de Arcos et al. [37], for the Sch15B rotor, provides evidence of the implications of blockage for BEM models. Figure 10 shows the relationship between thrust and axial induction factors for the Sch15B rotor operating at different tip-speed and blockage ratios, alongside Buhl's TWM, the inviscid momentum equation model, and the modified TWM developed in the previous section.

Figure 10 shows that blockage can have a very significant impact on the relationship between thrust and axial induction
<figure><img src="image_9.png" /><figcaption>**FIGURE 7** | Comparison of the spanwise thrust (top row) and torque (bottom row) distributions for different tip-speed ratios $\lambda$ calculated for the MEXICO (left), Sch15B (centre) and OxB (right) rotors with blade-resolved CFD (solid line), BEM with Buhl's TWM [14, 42] (dashed line), and the modified TWM (dot-dash line). The shaded region marks the sections that are less significantly affected by blade root- and tip-loss effects.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**FIGURE 8** | Relative error between the BEM models and CFD results for integrated thrust (left) and torque (right) on mid-span sections of the blade. The relative errors are calculated based on the loads integrated over the $r/R \in [0.40, 0.70]$ region to avoid undue influence from blade tip and root effects.</figcaption></figure>

factor for relatively high blockage ratios. While the effects are small up to $B = 0.05$, in agreement with the previous observations comparing the mildly blocked and unblocked OxB results,

the data still show an increase in thrust with increasing blockage even for small blockage, and a reduction in axial induction factor due to the increased mass flux through the rotor.
OxB Rotor - Blockage adjusted

<figure><img src="image_11.png" /><figcaption>**FIGURE 9** | Comparison of the spanwise thrust (left column) and torque (right column) distributions for different tip-speed ratios $\lambda$ calculated for the OxB rotor. The blade-resolved CFD results are adjusted for the effect of blockage using the method proposed by Bahaj et al. [45].</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>**FIGURE 10** | Impact of blockage on the relationship between thrust coefficient and axial induction factor for the Sch15B rotor at four different spanwise locations, $r/R = 0.6, 0.7, 0.9, 0.98$. The analytical and empirical model curves are presented without tip-loss correction.</figcaption></figure>

The limited data points are not sufficient to propose a modified TWM as a function of blockage especially at high induction factors. However, they suggest that the onset of the divergence

between the empirical data and the analytical thrust function occurs at lower $a_c$ values as blockage increases and, potentially, a change in the shape of the curve.
Although the CFD results show a change in the C<sub>t</sub>-a relationship as a result of changing blockage, we did not observe an impact on the foil's sectional lift and drag coefficients, and so changes in thrust and power with blockage are primarily a result of changes in rotor plane flow speed and its distribution rather than on modified behaviour of the foil sections themselves. This discussion falls outside the main scope of this paper and is thus presented in Appendix B.

## 4 | Conclusions

The analysis of the CFD datasets presented in this paper shows that turbulent wake models are necessary for accurate BEM computations even at relatively low axial induction factors. Buhl's widely used turbulent wake model does not agree well with the CFD-derived local thrust coefficients C<sub>t</sub> for a > 0.20 in regions of the blade that are relatively free from tip- and hub-losses.

Although divergence between Buhl's TWM and the CFD data starts to occur at a<sub>c</sub> ≈ 0.20, it becomes substantial at a<sub>c</sub> > 0.30. This is observed to affect BEM computations, resulting in an under-prediction of thrust and especially power at high tip-speed ratios for the three analysed rotors.

Buhl's TWM showed a relatively good prediction of the integrated thrust, with less good agreement for the integrated power. Torque was substantially underestimated by the BEM model using Buhl's TWM for all three rotors at high values of λ. The relatively good agreement of the integrated C<sub>T</sub> with the CFD results was found to be the result of counteracting influences from underprediction of thrust across mid-span and inboard regions and overprediction of loads towards the blade tip.

Prandtl's tip-loss formulation was originally derived to correct for the effects of a finite number of blades rather than a continuous disc and not necessarily to capture all the phenomena occurring close to blade tips. The CFD-derived data are consistent with requirement for a tip-loss model to act as a scaling factor on the momentum equations, as applied in the BEM model of Ning [42] used in this work. However, the fitting of an empirical scaling factor T<sub>L</sub> to the momentum model suggests that such correction should not depend on tip-speed ratio, as Prandtl's correction does. Note, however, that this refers exclusively to momentum corrections and not to performance corrections that account for changes in lift and drag coefficients as the tip is approached [6, 7, 44].

The modified TWM, an update of Buhl's model fitted to the CFD data at mid-board regions of the blade, significantly improved the results for both thrust and torque across most of the span of the blade compared with the original model. Note that the modified TWM was not derived as a best fit for the CFD power and thrust curves, but to provide a relationship between thrust coefficients and axial induction factors over an area of the blade where three-dimensional flow effects (tip- and hub-losses) would be limited.

While the integrated thrust and power were often overpredicted with the modified TWM, this was a consequence of the

aforementioned deficiencies in tip-loss modelling. Despite the overestimation, the CFD trends for integrated C<sub>p</sub> were much better followed with the modified TWM rather than with the original TWM, especially at higher tip-speed ratios.

The analysis of blockage effects showed a relatively limited impact of modest blockage ratios (B < 0.05) on the relationship between thrust and axial induction factors, and a substantial increase in C<sub>T</sub> at fixed axial induction factors as blockage becomes higher, affecting the shape of the C<sub>T</sub>-a relationship. Nevertheless, the observed impact on lift and drag coefficients (and, thus, the development of pressure changes across the chord and span of the blade) was negligible and suggests that only corrections to the momentum equations and, potentially, tip-loss models, would be required for adequate modelling of blockage effects within BEM models.

Finally, the discussion in this paper was presented without including any performance corrections to lift and drag coefficients (e.g., Shen et al. [6] or Wimshurst and Willden [7]), as well as only using the tip-loss correction of Prandtl. As extensively discussed throughout the paper, this modelling approach is not considered ideal. However, we limit the discussion to the momentum-thrust empirical relationship, instead of extending into an area where multiple corrections could coexist and, potentially, interact with each other. So far, this work has demonstrated the need for better models for the turbulent wake regime, and the good performance of the modified TWM in areas of the blade where three-dimensional flow effects are limited. Future work should address momentum and performance corrections to improve the prediction of tip-loss effects.

## Acknowledgements

This research project was supported in part by CONICYT PFCHA/BECAS CHILE DOCTORADO EN EL EXTRANJERO 2016/72170292. This project also received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 101034329, and recipient of the WINNINGNormandy Program supported by the Normandy Region. RHJW would like to acknowledge EPSRC who support his fellowship through grant number EP/R007322/1. CRV acknowledges the support of the UKRI through his Future Leaders Fellowship MR/ V02504X/1.

## Data Availability Statement

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

## References

1. H. Glauert, *Airplane Propellers* (Springer, 1935).

2. W. J. M. Rankine, “On the Mechanical Principles of the Action of Propellers,” *Transactions of the Institution of Naval Architects* 6 (1865).

3. W. Froude, “On the Elementary Relation Between Pitch, Slip, and Propulsive Efficiency,” *Transactions of the Royal Institution of Naval Architects* (1878).

4. H. Glauert, *The Elements of Aerofoil and Airscrew Theory* (Cambridge University Press, 1947).
5. G. Bangga, T. Lutz, E. Jost, and E. Krämer, "CFD Studies on Rotational Augmentation at the Inboard Sections of a 10 MW Wind Turbine Rotor," *Journal of Renewable and Sustainable Energy* 9, no. 2 (2017): 23304.

6. W. Z. Shen, R. Mikkelsen, J. N. Sørensen, and C. Bak, "Tip Loss Corrections for Wind Turbine Computations," *Wind Energy* 8, no. 4 (2005): 457–475.

7. A. Wimshurst and R. H. J. Willden, "Analysis of a Tip Correction Factor for Horizontal Axis Turbines," *Wind Energy* 20, no. 9 (2017): 1515–1528.

8. C. R. Vogel, R. H. J. Willden, and G. T. Houlsby, "Blade Element Momentum Theory for a Tidal Turbine," *Ocean Engineering* 169 (2018): 215–226.

9. M. Hansen, *Aerodynamics of Wind Turbines*, 2nd ed. (International Institute for Environment and Development, 2008).

10. R. R. Damiani and G. Hayman, "The Unsteady Aerodynamics Module for FAST8," National Renewable Energy Lab. (NREL), Golden, CO (United States), (2019).

11. L. Prandtl, "Applications of Modern Hydrodynamics to Aeronautics," *US Government Printing Office* (1921).

12. A. Betz, "Das maximum der theoretisch möglichen ausnützung des windes durch windmotoren," *Zeitschrift für das gesamte Turbinenwesen* 26 (1920).

13. L. A. Martinez-Tossas, E. Branlard, K. Shaler, et al., "Numerical Investigation of Wind Turbine Wakes Under High Thrust Coefficient," *Wind Energy* 25, no. 4 (2022): 605–617.

14. M. L. Buhl, "A New Empirical Relationship Between Thrust Coefficient and Induction Factor for the Turbulent Windmill State A New Empirical Relationship Between Thrust Coefficient and Induction Factor for the Turbulent Windmill State," NREL/TP-500-36834. National Renewable Energy Lab. (NREL), Golden, CO (United States), (2005).

15. C. N. H. Lock, H. Bateman, and H. C. H. Townend, *An Extension of the Vortex Theory of Airscrews With Applications to Airscrews of Small Pitch, Including Experimental Results* (HM Stationery Office, 1926).

16. E. A. Parra, K. Boorsma, J. G. Schepers, and H. Snel, "Momentum Considerations on the New MEXICO Experiment," *Journal of Physics: Conference Series*, Vol. 753. IOP Publishing, (2016): pp. 72001.

17. J. N. Sørensen and C. W. Kock, "A Model for Unsteady Rotor Aerodynamics," *Journal of Wind Engineering and Industrial Aerodynamics* 58, no. 3 (1995): 259–275.

18. H. A. Madsen, "A CFD Analysis of the Actuator Disc Flow Compared With Momentum Theory Results," 10th IEA Meeting on Aerodynamics. Technical University of Denmark. Department of Fluid Mechanics, (1997): pp. 109–124.

19. T. Burton, N. Jenkins, D. Sharpe, and E. Bossanyi, *Wind Energy Handbook* (John Wiley & Sons, 2011).

20. D. A. Spera, *Wind Turbine Technology* (ASME Press, 1994).

21. H. Snel, J. G. Schepers, and B. Montgomerie, "The MEXICO Project (Model Experiments in Controlled Conditions): The Database and First Results of Data Processing and Interpretation," *Journal of Physics: Conference Series* 75, no. 1 (2007): 12014.

22. C. Crawford, "Re-Examining the Precepts of the Blade Element Momentum Theory for Coning Rotors," *Wind Energy* 9, no. 5 (2006): 457–478.

23. F. Zilic de Arcos, C. Vogel, and R. Willden, "Extracting Angles of Attack From Blade-Resolved Rotor CFD Simulations," *Wind Energy* 23, no. 9 (2020): 1868–1885.

24. P. J. Moriarty and A. C. Hansen, "AeroDyn Theory Manual," Technical report. National Renewable Energy Lab., Golden, CO (US), (2005).

25. I. Afgan, J. McNaughton, S. Rolfo, D. D. Apsley, T. Stallard, and P. Stansby, "Turbulent Flow and Loading on a Tidal Stream Turbine by LES and RANS," *International Journal of Heat and Fluid Flow* 43 (2013): 96–108.

26. N. N. Sørensen, K. Boorsma, and G. Schepers, "CFD Computations of the Second Round of MEXICO Rotor Measurements," *Journal of physics: Conference series*, Vol. 753. IOP Publishing, (2016): pp. 22054.

27. A. Wimshurst and R. H. J. Willden, "Extracting Lift and Drag Polars From Blade-Resolved Computational Fluid Dynamics for Use in Actuator Line Modelling of Horizontal Axis Turbines," *Wind Energy* 20, no. 5 (2017): 815–833, https://doi.wiley.com/10.1002/we.2065.

28. A. Länger-Möller, J. Löwe, and R. Kessler, "Investigation of the NREL Phase VI Experiment With the Incompressible CFD Solver THETA," *Wind Energy* 20, no. 9 (2017): 1529–1549.

29. R. H. J. Willden, X. Chen, S. W. Tucker Harvey, et al., "Tidal Turbine Benchmarking Project: Stage I - Steady Flow Blind Predictions," Proceedings of the Fifteenth European Wave and Tidal Energy Conference. EWTEC Bilbao, (2023).

30. E. Jost, L. Klein, H. Leipprand, T. Lutz, and E. Krämer, "Extracting the Angle of Attack on Rotor Blades From CFD Simulations," *Wind Energy* 21, no. 10 (2018): 807–822.

31. J. Schluntz and R. H. J. Willden, "The Effect of Blockage on Tidal Turbine Rotor Design and Performance," *Renewable Energy* 81 (2015): 432–441.

32. A. Wimshurst and R. H. J. Willden, "Tidal Power Extraction on a Streamwise Bed Slope," *Ocean Engineering* 125 (2016): 70–81.

33. F. Zilic de Arcos, C. Vogel, and R. H. J. Willden, "Hydroelastic Modelling of Composite Tidal Turbine Blades," Advances in Renewable Energies Offshore: Proceedings of the 3rd International Conference on Renewable Energies Offshore (Renew 2018), October 8-10, 2018, Lisbon, Portugal, (2018): pp. 877.

34. F. Zilic de Arcos. (2021), "Hydrodynamics of Highly-Loaded Axial Flow Tidal Rotors," Ph.D. Thesis, University of Oxford.

35. F. Zilic de Arcos, C. R. Vogel, and R. H. J. Willden, "A Parametric Study on the Hydrodynamics of Tidal Turbine Blade Deformation," *Journal of Fluids and Structures* 113 (2022): 103626, https://www.sciencедirect.com/science/article/pii/S0889974622000755.

36. F. Zilic de Arcos, C. R. Vogel, and R. H. J. Willden, "Hydrodynamic Independence and Passive Control Application of Twist and Flapwise Deformations of Tidal Turbine Blades," *Journal of Fluids and Structures* 118 (2023): 103827.

37. F. Zilic de Arcos, G. Tampier, and C. R. Vogel, "Numerical Analysis of Blockage Correction Methods for Tidal Turbines," *Journal of Ocean Engineering and Marine Energy* 6, no. 2 (2020): 183–197.

38. S. W. Tucker Harvey, X. Chen, D. Rowe, et al., "Tidal Turbine Benchmarking Exercise: Geometry Specification and Environmental Characterisation," Proceedings of the Fourteenth European Wave and Tidal Energy Conference. EWTEC Plymouth, (2021).

39. S. W. Tucker Harvey, X. Chen, K. Bhavsar, et al., *Tidal Benchmarking Project Dataset: R001* (University of Oxford, 2022), https://ora.ox.ac.uk/objects/uuid:29babc14-f289-48bd-bc24-26b629380774.

40. J. Y. Luo and A. D. Gosman, "Prediction of Impeller-Induced Flow in Mixing Vessels Using Multiple Frames of Reference," Institute Of Chemical Engineers Symposium Series, (1994).

41. F. R. Menter, M. Kuntz, and R. Langtry, "Ten Years of Industrial Experience With the SST Turbulence Model," *Turbulence Heat and Mass Transfer* 44 (2003): 625–632.
42. S. A. Ning, “A Simple Solution Method for the Blade Element Momentum Equations With Guaranteed Convergence,” *Wind Energy* 17, no. 9 (2014): 1327–1345.

43. H. Rahimi, J. G. Schepers, W. Z. Shen, et al., “Evaluation of Different Methods for Determining the Angle of Attack on Wind Turbine Blades With CFD Results Under Axial Inflow Conditions,” *Renewable Energy* 125 (2018): 866–876.

44. A. Wimhurst and R. H. J. Willden, “Computational Observations of the Tip Loss Mechanism Experienced by Horizontal Axis Rotors,” *Wind Energy* 21, no. 7 (2018): 544–557.

45. A. S. Bahaj, A. F. Molland, J. R. Chaplin, and W. M. J. Batten, “Power and Thrust Measurements of Marine Current Turbines Under Various Hydrodynamic Flow Conditions in a Cavitation Tunnel and a Towing Tank,” *Renewable Energy* 32, no. 3 (2007): 407–426.

# Appendix A

## 2D Polar Coefficients

The polars used for the BEM computations are presented in Figure A1. The polar coefficients were computed with a RANS CFD solver, a target $y^+ \approx 1$ and a $k-\omega$ SST turbulence model. Turbulence intensity and length scale were adjusted to match the corresponding conditions of the testing facilities based on the spanwise-average chord lengths. Chord-based Reynolds number interpolation was used for the Sch15B ($Re \in [10.0 \times 10^6 - 40.0 \times 10^6]$) and OxB rotors ($Re \in [0.1 \times 10^6 - 0.4 \times 10^6]$), owing to the large $Re$ variations across the range of tip-speed ratios due to changes in rotational speed. A single Reynolds number was used for each foil on the MEXICO rotor given that $Re$ variations were limited due to a constant rotational speed. $Re$ for the MEXICO rotor are $0.6 \times 10^6$, $0.7 \times 10^6$, and $0.8 \times 10^6$ for the DU91-W2-250, RISØ A1-21 and NACA 64 418 aerofoils, respectively.

# Appendix B

## Impact of Blockage on Spanwise Lift and Drag Coefficient

An interesting insight derived from the numerical simulations of the Sch15B rotor at different blockage ratios is the impact of blockage on the sectional blade forces. Figure B1 shows the polar coefficients calculated with the line-average technique at three different span locations for the 5 different blockage ratios. Each simulated tip-speed ratio provides, at each spanwise location, a normal and tangential force, as well as an angle of attack. These three values are used to reconstruct the lift and drag curves, as shown in the corresponding figure. Each bullet point corresponds to a specific simulation and spanwise location, and cubic splines are used to link the extracted values.

The CFD results show a negligible effect of blockage on the lift and drag curves for all the analysed cases and throughout the entire range of available angles of attack for a given radial location, as the curves collapse on top of each other with relatively low spread. The spanwise location, however, shows a clear impact on lift and drag, a phenomenon that has been previously described in the literature [27].

The convergence between the $C_L$ and $C_D$ curves for different blockage ratios shows that the changes in thrust and torque can be explained primarily through the changes in angles of attack that arise due to the effect of blockage. This further explains the good performance of blockage corrections that are introduced as changes to the inflow speed to modify thrust, power, and tip-speed ratio, as described by Zilic de Arcos et al. [37]. These results also suggest that no further blockage-related corrections should be required by the blade-element component of a constrained BEM model such as the one presented by Vogel et al. [8].

<figure><img src="image_13.png" /><figcaption>**FIGURE A1** | Two-dimensional polar coefficients used for the BEM computations of the MEXICO, Sch15B and OxB rotors.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**FIGURE B1** | Impact of blockage on the lift and drag coefficients at different spanwise locations extracted from the Sch15B blade-resolved CFD simulations and for different blockage ratios. Each bullet point corresponds to a simulation and cubic splines are provided for reference.</figcaption></figure>