Mitigation of Rotor Thrust Fluctuations through Passive Pitch

A rotor blade operating in an unsteady or a sheared stream experiences force fluctuations, which increase the structural requirements for both the blades and the rotor system as a whole. In this paper, we investigate whether force fluctuations can be passively mitigated without compromising the mean load. We consider a tidal turbine rotor in a shear current at a diameter-based Reynolds number of 2 × 10 7 . The blades are rigid and can passively pitch. A mass–spring system acts on the spanwise axis of the blade governing the pitch kinematics. The effectiveness of this system is demonstrated with three methodologies: an analytical model based on blade element momentum theory and Theodorsen’s theory, and two sets of Reynolds-averaged Navier–Stokes simulations performed with two independent codes. The analytical and the numerical models are validated against experiments and simulations in the literature. All methodologies predict a reduction of more than 80% of the thrust fluctuations. Furthermore, because of the more uniform thrust force exerted on the current, the wake behind a passive pitch blade does not diffuse the onset shear flow. This results in a more sheared wake and faster wake recovery. The present results demonstrate the potential benefits of passive pitch and may underpin future applications of this concept for different types of turbines and compressors, including wind turbines, propellers, helicopter rotors, etc. © 2022TheAuthor(s)


Introduction
When a blade rotates in a non-uniform or time-dependent flow field, it experiences load variations.This applies to both rotors that extract energy from the flow such as wind (Hansen and Butterfield, 1993;Sørensen, 2011) and tidal turbines (Adcock et al., 2021), and rotors that provide energy to the fluid such as ship propellers (Kerwin, 1986), helicopter rotors (Conlisk, 1997) and axial flow compressors (Hawthorne and Novak, 1969).Load fluctuations are undesirable because of their effect on the fatigue life of the blades and of the turbine structures.Furthermore, fluctuations in the power (harvested or provided), require more expensive generators and motors than if they were rated for only the mean power.
Unsteady loads are typically controlled with active systems, as they provide more flexibility than passive systems (Cattafesta and Sheplak, 2011).However, passive control might enable higher system reliability and thus it is desirable when frequent maintenance is not possible (Barlas and van Kuik, 2010).This is the case, for example, of large offshore wind turbines and, even more for underwater tidal turbines.
The most adopted unsteady load control system for wind and tidal turbines is pitch control (Bossanyi, 2003).The blades are rigid and the pitch angle (around the spanwise axis) is varied by an actuator.The pitch, however, can also be C (k) Coefficient of the force normal to the rotor plane (thrust coefficient) [-] Coefficient of the force tangential to the rotor plane [-] C M 0 Preload moment coefficient [-] controlled passively.For example, the pioneering work of Cheney and Spierings (1976) inspired technologies that exploited the centrifugal force to pitch the blade and limit the angular velocity (e.g.Chen and Shiah, 2016;Yanagihara et al., 2011).Buttrell (1981) considered a two-blade rotor where the aerodynamic pitch moment on one blade pitches the other blade.Hertel et al. (2004), instead, studied numerically a linear torsional spring and a damper acting on the span axis of the blade.Most of the passive pitch systems considered by previous authors were conceived to vary quasi-steadily the pitch angle of a wind turbine to respond to mean wind speed variations.Thus, pitch control was used to limit the maximum angular velocity and the power output.In contrast, it remains unclear whether a passive pitch system can mitigate the amplitude of fast load fluctuations for a given mean flow velocity.Tully and Viola were the first to investigate a morphing foil concept for tidal turbines.They demonstrated that a flexible foil that elastically bends when the hydrodynamic load increases, can mitigate the amplitude of the load fluctuation.Their study underpinned the work of Arredondo-Galeana et al. (2021) on trailing edge flaps equipped with a torsional spring at the hinge.Arredondo-Galeana et al. ( 2021) demonstrated that the fraction of load alleviation is given by the relative   size of the flap's and the foil's chord.This suggested that pitching the entire blade around the leading edge should enable maximum unsteady load mitigation.Such a system was studied analytically by Dai et al. (2019) and Pisetta et al. (2021).They considered the blade of a commercial-scale tidal turbine that pitched passively around an axis near the leading edge.Pisetta et al. (2021) found that the effect of the inertia of the blade and the fluid were detrimental to efficiency.However, they postulated that reasonable values of blade inertia for a tidal turbine would have not compromised the overall effectiveness of the passive pitch system.This paper aims to demonstrate with three-dimensional numerical simulations the results of Dai et al. (2019) and Pisetta et al. (2021) for a full-scale tidal turbine in a shear flow.We achieve this by using an enhanced version of the analytical model developed by Pisetta et al. (2021) to design a passive pitch system, which we then test with computational fluid dynamics (CFD) simulations.The rest of the paper is organised as follow.The flow and the operative conditions considered are described in Section 2. The analytical model by Pisetta et al. (2021) and a novel modification of the original formulation is illustrated in Section 3. The two CFD models tested in this paper are described in Section 4. Verification and validation of the numerical simulations are detailed in Section 5 and the results are discussed in Section 6.The key conclusions of the paper are summarised in Section 7.

Problem statement
We consider an 18 m diameter, 1.5 MW tidal turbine that was tested at the European Marine Energy Centre's site within the ReDAPT project (Harrison, 2015).Table 1 and Fig. 1 summarise the turbine specifications and the operative conditions.The velocity profile along the vertical axis z of the free stream velocity U ∞ is given by a 1/7 power law, where U hub = 2 m s −1 is the free stream velocity at the hub height z hub = 20 m.The rotor experiences a one-per-revolution fluctuation because of the onset shear flow (Eq.( 1)).To mitigate the load fluctuation, we consider a passive system that periodically pitches the blades over a period of revolution.

Analytical model
In this section, we present the analytical model that was originally developed by Pisetta et al. (2021) and that it is here modified to account for the blockage effect.This model is used to identify the optimal spring stiffness and preload that are used in the CFD simulations.In the Results section (Section 6), we present the force predicted with this analytical model as well as with the CFD simulations.

Blade element momentum theory
For each blade section, the axial (a) and tangential (a ′ ) induction factors are evaluated with the iterative procedure suggested by Burton et al. (2001): where r is the radial coordinate; F is the Prandtl's tip loss correction (Prandtl, 1921) to account for the finite span of the blade; γ is the inflow angle between the rotor plane and the mean velocity Ū; N b = 3 is the number of blades; c is the chord length; and C L and C D are the lift and drag coefficients, respectively.The tip loss factor is computed as where (5) For every angle of attack α, the lift and drag coefficients (C L , C D ) are taken from a database built by two-dimensional Reynolds-averaged Navier-Stokes simulations using the same numerical setup as that for the three-dimensional simulations reported in Section 4 (RKE-SM approach).The procedure is repeated until a, a ′ , C L and C D converge to steady values.

Mean loads and torsional spring preload
To ensure that the mean load is conserved, the preload of the torsional spring is such that the spring moment balances, on average, the hydrodynamic and centrifugal moments.For every ith section, the angle of attack is, on average, where β 0 is the pitch angle between the rotor plane and the time-averaged direction of the chord (Fig. 2).All angles and moments are defined positive anticlockwise.The time-averaged hydrodynamic moment is computed as the moment experienced by the blade section at a constant incidence ᾱ, that it where ρ is the fluid density, and C M,i (α) is the moment coefficient (positive anticlockwise), which is tabulated for different values of α akin to the lift and drag coefficients.For the whole blade, the hydrodynamic moment is Each blade section also experiences a centrifugal moment M c,i , which is also integrated over the blade as The torsional spring provides, on average, a moment Mk = −κθ 0 , where κ is the spring stiffness and θ 0 is the preload.
Hence, the latter is computed from the equilibrium of the moments as ) . (10)

Passive pitch formulation
Using Theodorsen's theory (Theodorsen, 1934), the instantaneous hydrodynamic pitch moment is computed as where, δα i = α i − ᾱi is the angle of attack fluctuation experienced by the ith blade section; δ αi and δ αi are their first and second time derivatives, respectively; ) ; ) C (k); b i = c i /2 is the semi-chord of the ith blade section; e = 0.8, which is constant for all blade section, is the distance in semi-chords from the pitch axis to mid chord; C (k) is Theodorsen's function for a given reduced frequency and ω T is the angular velocity of the blade.The hydrodynamic moment contributions associated with I a i and with the first term of C a i are due to the added mass.These contributions are given by, for instance, Leishman (2006), and are derived from the time derivative of the potential function of a thin flat plate with semi-chord b i pitching around a point at eb i upstream of the mid-chord.
The centrifugal moment M c depends on θ, and thus it is linearised around θ 0 as where The dynamic equation of the full blade is Each term of Eq. ( 16) can be written as a function of δθ , and its first and second time derivatives, which are the only unknowns: The inflow angle experienced by the blade varies periodically due to the onset shear flow, and its fluctuation is approximated as where ∆γ is the amplitude of the fluctuation δγ , j is the imaginary unit, and t is time.The blade is expected to oscillate with the same period, namely where ∆θ is the amplitude of the fluctuation δθ , and c θ is a phase shift.By substituting δ θ = jω T δθ, δ θ = −jω T 2 δθ, δ γ = −ω T 2 δγ , and δ γ = −jω T δγ into Eq.( 17), we find The dynamics of δθ are determined by the algebraic solution of Eq. ( 20) for each time step and the equilibrium position of the blade is computed as θ = θ 0 + δθ .
The lift (L i ) for each ith section is computed using Theodorsen's theory (Theodorsen, 1934) as where the first of the three terms of L i is the added mass computed for a flat plate with semi-chord b i (Leishman, 2006).
The Drag (D i ) is computed as The normal (N i ) and the tangential (T i ) loads to the rotor plane are where N is in the direction of the thrust and acts in the direction of U N , while T is associated with the power generation and acts in the opposite direction of U T (Fig. 2).Their coefficients are The load coefficients for the entire blade are The results are presented in the form of the non-dimensional stiffness coefficient where A s is the rotor swept area.

Blockage correction
The analytical code presented in Sections 3.1-3.3does not take into account the effect that the seabed and the free surface have on the turbine performance, i.e. the blockage effect.Whelan et al. (2009), for example, employs a correction based on Maskell (1963) and showed that this correction is significant on tidal turbines.
The streamtube through the turbine rotor expands because of the thrust that the turbine exerts on the fluid, but the blockage partially mitigates this expansion resulting in a higher flow velocity through the rotor plane compared to the unblocked condition.An estimate of the velocity increase through the rotor plane is given by the velocity increase of the bypass flow, which is the flow external to the streamtube.The bypass flow speed can be written as U 2 = τ f U ∞ , with τ f > 1.The normal and tangential force coefficients (Eq.( 25)) are thus multiplied by a factor τ 2 f to account for blockage.
To estimate τ f , we use the numerical simulations RKE-SM described in Section 4 and we consider the bypass velocity on a vertical plane through the axis of the turbine and parallel to the stream.From the average velocity in the two highlighted bypass regions in Fig. 3, we find τ f = 1.067.

Selected values for stiffness and inertia
The value of the spring stiffness that minimises the amplitude of the thrust fluctuations is found with an exhaustive search algorithm.The performances are almost constant for C κ < 10 −2 and thus we select C κ = 10 −3 .For this value of the spring stiffness coefficient, the preload moment coefficient is C M 0 = 8 × 10 −3 .It is noted that, as the stiffness tends to zero, the preload tends to infinity and the spring exerts a constant moment (Viola et al., 2021).
The effectiveness of the passive pitch decreases with the blade inertia.This is due to the delayed pitch response of the blade as the onset flow varies.Fig. 4 shows the amplitude of the normal force fluctuations ∆C N versus the moment of inertia of the blade.We normalise ∆C N with the amplitude of the normal force fluctuations of a rigid, fixed-pitch blade.
Therefore, the ordinate takes the value of zero when the passive pitch allows perfect unsteady load cancellation, and one when the fluctuations are unchanged with respect to a fixed pitch blade.The moment of inertia I is representative of a hollow blade with a shell made of composite layers then filled with water, and it is normalised with the moment of inertia I 0 of an equal volume of water.The composite density is set as 600 kg m −3 and, for a blade shell thickness varying from 14.5 mm at the tip to 56.6 mm at the root (Bir et al., 2011), we find I = 2400 kg m 2 and I/I 0 = 0.99.Fig. 4 shows that the performance of the passive pitch does not vary significantly for I/I 0 between 0.99 and 2.5.To improve the convergence of the weakly-couple fluidstructure-interaction approach, we set I/I 0 = 2.5, which corresponds to a conservative estimate of the blade inertia, about one order of magnitude higher than estimated.

Hydroelastic stability
The stability of the passive pitch system depends on the position of the pitch axis and the force and moment coefficients of the blade.For each section, we fix the pitching axis near the leading edge, at x = 0.1c.The forward position of the pitch axis with respect to the aerodynamic centre, which is around 0.25c, ensures that the hydrodynamic pitch moment M h has the same sign as the lift.Note that the pitch moment is here defined positive anticlockwise, thus opposite to the typical aeronautic convention.Hence, M h > 0 for a wide range of angles of attack including from 0 • to 15 • .Furthermore, for the chosen geometry and pitch axis location, ∂M h /∂α > 0 for 0 • ≤ α ≤ 15 • .These two conditions ensure that the passive pitch system is stable.In fact, a small increase in the angle of attack due to, for example, a gust, For completeness, it is also noted that the pitch axis could be chosen such that M h < 0 as, for example, at the trailing edge.In this case, the system would be stable as long as ∂M h /∂α > −κ.Finally, it is noted that the modelled system has only one degree of freedom in pitch, and thus flutter is inhibited.However, if flexural vibrations occurred, then flutter instabilities should also be investigated.

Numerical method
The unsteady incompressible Reynolds-Averaged Navier-Stokes (RANS) equations are solved with two independent codes: (1) the commercial tool Star-CCM+ with the Realisable k − ϵ (RKE) turbulence model and a Sliding Mesh (SM) strategy, hereafter denoted as RKE-SM; and (2) the in-house code χ navis with the Spalart-Allmaras (SA) turbulence model and an Overlapping Grid (OG) approach, hereafter denoted as SA-OG.An overview of the key features of the two numerical models is provided in Table 2, while the details of each model are provided in Sections 4.2 and 4.3, respectively.In both cases, the RANS equations are solved in the frame of reference fixed with the channel, while the turbine rotates at a constant angular velocity.

Fluid-structure interaction
A weakly coupled approach is adopted as shown in Fig. 5.At every time step, the fluid solver computes the hydrodynamic moment M n around the pitch axis.The structural analysis consists in solving the equation of motion for the pitch.A second-order implicit finite difference Newmark-Wilson scheme is used to compute the blade angular velocity β and its angular position β at time-step n, based on the values computed at the previous time step n − 1: where ∆t is the time step and βn is computed solving the mass-spring system equation: where 2).Once β n is computed, the mesh is updated and the new hydrodynamic moment M n+1 is computed.

Computational setup of RKE-SM
The volume containing the blades and the hub constitutes an inner domain that rotates inside a stationary cuboid domain.The free surface is modelled as rigid.The two domains are coupled by a sliding internal interface, recomputed at every time step.The stationary domain is made of 12.7 million cells, discretising the 17D T × 11D T × 2D T fluid region.Fig. 6a presents the boundary conditions.The velocity profile modelled with Eq. ( 1) is enforced at the inlet with zero turbulence level.A reference zero static pressure is set at the outlet.The slip condition is selected for the rest of the outer surfaces of the stationary domain.A no-slip condition is set on the blade and hub surfaces.
The inner rotating domain is discretised with a hexahedral structured mesh of approximately 23 million cells with y + = 70, ∆x/c = 0.015 and ∆z/∆x = 1, where ∆x and ∆z are the linear dimensions of the first layer cell adjacent to the blade in the azimuthal and radial direction, respectively.Because the first cell centre is in the logarithmic layer of the blade boundary layer (y + = 70), wall functions are employed to compute the friction velocity.The tower is not included in the model.To generate the volume mesh for the three-bladed tidal turbine, the 120 • periodicity of the rotor is exploited by only meshing the grids around one blade and the remaining two blades are included in the computational domain using the same grid and connected with internal interface boundary conditions.It should be noted that two relative motions are considered, that of the external stationary domain with the internal domain of the rotor, as well as that of the pitching blades with respect to the rotor.
Fig. 6b shows that the grid is divided into eight blocks to control the mesh along the blade.A block with a mesh of 12 layers and a total thickness of 0.1% of the blade chord at the root is generated in the gap between the blade and the hub to enable the blade to pitch (Fig. 6c).A C-grid is adopted surrounding the blade and an O-grid is employed to build the mesh on the collapsed surface of the tip (Fig. 6d-f).

Computational setup of SA-OG
Computations for both fixed and passive pitch blades are also performed with the χnavis solver, which is a generalpurpose simulation code developed at the Institute for Marine Engineering, National Research Council (CNR-INM).It is based on a finite volume formulation with variables co-located at the cell centres.The spatial discretisation of the convective terms is achieved by means of a third-order upwind biased scheme, whereas, a second-order scheme is used for the approximation of the diffusive terms.Time derivatives are approximated by a second-order implicit formula (three points backwards).The divergence-free solution at each time step is computed iteratively by a dual time stepping integration.Convergence to steady-state in pseudo-time is accelerated by an Euler implicit scheme with approximate factorisation, local pseudo-time step and multi-grid iteration (Favini et al. 1996).The turbulent viscosity is calculated by means of the one-equation model by Spalart and Allmaras (1992).Further details on the numerics can be found in (Di Mascio et al. 2001, 2009, Broglia and Durante 2018, Posa and Broglia 2019).The computational domain and the boundary conditions are the same as for RKE-SM, and the tower is also not considered in these simulations.A dynamic overlapping grid approach is used as in Di Mascio et al. ( 2006) and Zaghi    2015).An overview of the computational mesh on the surface of the tidal turbine is shown in Fig. 7a.The inner domain is discretised by body-fitted patched and overlapped blocks (Fig. 7b-c).The overset approach allows a finer grid near the blade and an overall lower number of cells.In fact, the thickness of the first cell on the wall is always below one wall unit, thus the near-wall flow is resolved without the need of wall functions, while ∆x/c and ∆z/∆x have the same values as the RKE-SM simulations.The total number of cells is about 5M cells.Fig. 8a compares the grids for the two numerical methods on a meridian plane (i.e. through the vertical axis of the rotor disk) parallel to the stream.It can be seen that the grids are well refined in the near wake region, although one diameter downstream of the turbine, the grid resolution is lower for RKE-SM than for SA-OG.Figs.8b-c show the grids on the cross-section of the wake at x/D T = 0.3 and 3.

Verification and validation
The numerical simulations are verified and validated to ensure that the numerical and modelling errors on the force computation are sufficiently low compared to the difference between the forces on a fixed-pitch and passive-pitch blade.Verification is performed by testing different spatial and time resolutions, and then computing the numerical uncertainty on the forces.Validation, instead, is based on the comparison between computed and measured forces.As no experimental data is available for a full-scale tidal turbine, RKE-SM results are validated modelling relevant wind tunnel experiments using a mesh and a numerical setting as similar as possible to that of the full-scale tidal turbine.We use high spatial and temporal resolution to validate the forces to a low uncertainty level.Hence, RKE-SM simulations will be used to investigate in detail the effect of the passive pitch.Conversely, the SA-OG results will be used to demonstrate that the predicted unsteady load mitigation is independent of the numerical method.For this purpose, we use a lower spatial resolution to reduce the computational cost.The numerical uncertainty of the SA-OG is computed by testing different spatial and time resolutions for the full-scale fixed-pitch rotor.

Verification and validation test case
We model the experiments of Piziali (1994), who tested a periodically pitching cantilever blade with NACA 0015 sections, a span of 1.5 m, a chord of 0.3 m, zero twist and a square tip cap (Figs.9a-b).The wing was pitched around the quarter chord point with a mean angle of attack ᾱ = 4 • , a mean pitch angle β = 7.4 • , a pitching amplitude δβ = 4 • and a reduced frequency k = 0.038 at Re = 2 × 10 6 .The normal force N and the tangential force T at 50% (mid-span) and 99.5% (tip) of the span were measured by pressure integration.Hence, in this section, N and T are sectional forces, i.e. forces per unit span.The angle between T and the flow velocity U ∞ is the same as the blade section at 75% span of the tidal turbine blade modelled in Section 6 (11.4 • ).
We define the sectional force coefficients C N and C T as the ratio between the measured force per unit span and 1/2ρU 2 ∞ c.We focus on the values of these force coefficients at two-phase angles: φ = 0.4π and φ = 1.6π , where φ = 0 corresponds to the mean pitch angle.We compare the measured coefficients with the pressure force coefficients computed with CFD, consistently with the method adopted to measure the forces.
The domain and the mesh are as similar as possible to those used with the RKE-SM for the full turbine, including the grid resolution of the first layer of cells near the blade surface, and the grid growth ratio as the distance from the blade increases (Fig. 9c-d).The grid is made of 2 million cells with y + = 50.The velocity is uniform at the inlet and a reference static pressure is set at the outlet.A no-slip condition is applied at both the plane at the root of the blade and the blade surface, while a slip condition is set on the remaining boundaries.The Realisable k − ϵ turbulence model is used with the same numerical setting described in Section 4 to generate the results presented in Section 6.

Verification and validation methodology
We follow the verification and validation procedure proposed by Viola et al. (2013).The error δ Φ in the estimate of a computed quantity is defined as the difference between the computed value Φ and the true value Φ TRUE , The error is broken down into a modelling error δ Φ mod and a numerical error δ Φnum , δ Φ = δ Φ mod + δ Φnum . (31) The modelling error is due to the difference between the reality and the model.For example, it is due to the underlying assumptions of the turbulence model, to the values set at the boundaries (which are unlikely to be the true values because the latter are unknown), etc. Conversely, the numerical error is due only to the discrete nature of the numerical model.
Verification is the process that estimates the numerical uncertainty U Φ associated with the computed solution Φ due to the numerical error.The uncertainties due to the different sources of numerical errors are typically which are associated with the grid size, the time step, the round-off, the input parameters, and the convergence, respectively.For a well-converged simulation in double precision, the most significant uncertainties are U Φg and U Φt .These uncertainties for a reference, base simulation, are estimated by undertaking a set of simulations where the spatial and the time resolutions are systematically varied.
Let h be the step size of the source of error.For example, for the grid size, h is the ratio between the linear grid size (e.g.∆x) of the current and the reference grid.Similarly, for the time step, h is the ratio between the current and the reference time step.Consider where Φ(base) is the value of the computed solution with the base setting (grid size, time step, etc.), for which the uncertainty is computed.Fit ϕ(h) with where the coefficients η, p and ϕ 0 are computed with the least-squares method.
The uncertainty in the computation of Φ due to a step size h is, for the base step size, where σ is the standard deviation of the fit of Eq. ( 33), ϕ max and ϕ min are the maximum and the minimum of ϕ(h), whilst h max and h min and the maximum and the minimum of h.
The uncertainties due to the different sources of numerical errors are combined to give one numerical uncertainty, which is It is noted that the convergence uncertainty U Φc is not under the square root because it is not considered independent from the other sources of error.The absolute uncertainty is Validation examines the error due to the difference between the model and the truth, which can be estimated from experimental results.The absolute numerical and experimental uncertainties are combined in the absolute validation uncertainty, where E Φ TRUE is the absolute experimental uncertainty.Then, the following cases are considered.
It must be emphasised that the objective is to validate the results to the lowest possible level of the absolute validation uncertainty, and not just to validate the results for any level.In fact, the criteria for validation can always be satisfied by increasing the numerical uncertainty.

Verification and validation analysis
We use a double-precision solver and set the number of inner iterations such that all residuals decrease by at least three orders of magnitude.We monitor the forces to ensure that the differences between each period are smaller than 1/1000.Hence the uncertainties due to the round-off and convergence are negligible.The uncertainties due to the grid and the time step are assessed as follows.The uncertainty due to the time step is considered in Fig. 10.We run simulations with time steps ∆t = 0.0006U ∞ c −1 , 0.003U ∞ c −1 , and 0.015U ∞ c −1 .The step size h is the ratio between ∆t and the base time step 0.003U ∞ c −1 , which is the value used for the turbine simulations presented in Section 6.For each time step ∆t, ϕ is the ratio between the mid-span C N computed with ∆t and that computed with the base time step.The solution converges asymptotically towards ϕ 0 when h → 0. The parameters of Eq. ( 33) are η = 3.35 × 10 −4 , p = 1.54, and ϕ 0 = 1.00.Using Eq. ( 34), we find U Φt = 0.04%.Fig. 10b shows that the uncertainty on the tangential force coefficient C T is one order of magnitude higher, U Φt = 0.3%.The grid study is shown in Fig. 11 for three grids with 0.25M, 2M(base) and 16k cells.We find U Φg = 0.83% and 15% for C N and C T , respectively.Finally, the same verification process is performed for mid-span C T , tip C N and tip C T at φ = 0.4π and φ = 1.6π .
The combined relative numerical uncertainties due to the grid and the time step together are shown in Table 3 for all tested cases.The overall uncertainty in the sectional force coefficients C N and tip C T are within 1% and 8%, respectively, with higher peaks near the tip where the absolute value of the forces drops significantly.As the low sectional force in the region very close to the tip has a small contribution to the overall blade force, we estimate that the uncertainty in the full-blade force coefficients computed with RKE-SM are within 10%.
To verify that the modelling error is small, we validate the predicted forces with those measured by Piziali (1994).
The measurement error is E Φ TRUE = 0.01.Fig. 12 shows the computed and measured C T and C N at both the mid-span and the tip.Table 3 summarises the results of the validation for each force coefficient at the two different phase angles.All computed forces are validated at a level of the validation uncertainty.We conclude that the numerical error of the forces computed with RKE-SM is relatively small compared to the validation uncertainty The uncertainty of SA-OG is computed following the same procedure as in Section 5.2 for simulations of the tidal turbine with fixed-pitch blades at full-scale conditions.We test three grids (by halving and doubling the distance between the grid nodes of the reference grid) and three time steps (by doubling and quadruplicating the reference time step).The power (C P ) and thrust (C N ) coefficients converge monotonically for both the three grids and the three time steps with order p ranging between 1.48 and 1.75 (results not presented for brevity).The numerical uncertainties due to the grid are U Φg = 23% and 9.8% for C P and C N , respectively, while those due to the time step are U Φt = 1.6% and 0.25% for C P and C N , respectively.

Results
We first consider the turbine with fixed pitch and compute the power C P and thrust coefficient C N for different tipspeed ratios.The same turbine was numerically modelled by Gretton and Ingram (2011).Fig. 13 shows that the maximum C P is achieved at the tip-speed ratio λ = 4.5, which is, therefore, the condition considered hereafter.

Unsteady load mitigation
We now consider the effectiveness of the passive pitch system in mitigating the force fluctuations is assessed.Fig. 14a shows the thrust coefficients C N calculated with the analytical model and with the RKE-SM simulations for a tidal turbine blade with fixed and passive pitch.Simulations are run for about 20 cycles, of which the first 10 cycles are not considered to allow a periodic solution to develop, and data is phase-averaged over the latter 10 cycles.The angular position φ of the blade is defined from the horizontal position, such that φ = π/2 when the blade tip is near the top boundary and φ = 3/2π when it is closed to the seabed.This definition is consistent with the phase angle φ used in Section 5.
The analytical predictions and the simulations agree within 4% of the average, and both show a significant reduction of the amplitude of the force oscillation for the passive pitch.This is quantified in Table 4, which shows the amplitude and the time average of C N for a single blade.The results of the analytical model and of both the RKE-SM and the SA-OG simulations are presented.All three methods indicate that the passive pitch enables a reduction of the amplitude of the thrust fluctuation between 75% and 81%, while the mean thrust remains constant within less than 1%.Viola et al. (2021) showed that a blade that pitches around an axis at x p = 0.1c can enable a maximum unsteady load mitigation of 90%.This occurs when each blade section is equipped with a separate spring locally optimised, and the blade inertia is negligible.In contrast, here the spring is optimised in a global fashion, i.e. the whole blade shares one spring.Furthermore, as shown in Fig. 4, at I/I 0 = 10, the performance is not as good as for I/I 0 → 0. This is due to the phase lag associated with the slow response of the blade to the hydrodynamic forcing, as shown in Fig. 14b, where the

Table 4
Amplitude of normal force coefficients ∆C N and time-averaged force coefficients CN on a single tidal turbine blade computed with RKE-SM and SA-OG, and the ANL.

∆C N,Fixed
∆C N,Passive pitch angle variation ∆β is calculated analytically and with the RKE-SM simulations as shown.The trends in the pitch angles show some differences but both demonstrate that the blade passively pitches by approximately ±1 • with a small phase lag with respect to the forcing frequency, whose maxima and minima are at φ = π/2 and 3/2π .Overall, despite the global spring optimisation and the non-negligible inertia, the unsteady load mitigation enabled by the passive pitch is within 15% of the theoretical maximum (90%).

Flow separation and dynamic stall
In this section, the effect of the passive pitch on flow separation is considered.Flow separation can be detected from the sign of the azimuthal component of the wall shear stresses.Fig. 15 shows the azimuthal component of the phasedaveraged wall shear stress on the suction side of the blade.Twelve different phase positions are shown.Regions with negative values on average are marked in grey and are associated with flow reversal and thus separated flow.Solid lines with different colours show the boundary of the region of separated flow at different phase angles.
The area of the grey region for the fixed (left) and passive (right) pitch are not dissimilar, revealing that the passive pitch has a marginal effect on flow separation.For both fixed and passive pitch, all bounding lines collapse roughly on the same contour, revealing that the area of separated flow remains almost constant over a period of revolution.Hence, at this tip speed ratio, dynamic stall does not occur for the fixed nor for the passive pitch, and the latter does not alter the extent of the blade area where the flow is separated.

Kinetic energy flux
We showed in Section 6.1 that the fixed pitch blade experiences high and low thrust when the onset flow velocity is high and low, respectively.Therefore, the blade exerts on the fluid a higher resistance when the flow velocity is higher, and a lower resistance when the flow velocity is lower.Consequently, the effect of a fixed pitch blade is to smooth the onset velocity gradient and to make the downstream flow more uniform.This is akin to the effect of a mesh in a wind/water tunnel.In contrast, the passive pitch blade experiences a more constant thrust and, therefore, it results in a less uniform velocity distribution in the wake.This is shown in Fig. 16, which shows the non-dimensional instantaneous kinetic energy flux (U 3 /U 3 ∞ ) through a cross-section of the wake at 0.3 diameters downstream of the turbine.The kinetic energy flux ρU 3 /2 is non-dimensionalised with the far-field kinetic energy flux ρU 3 ∞ /2.Note that the latter varies with the vertical coordinate.The contours clearly show that the passive pitch blade results in a higher kinetic energy flux on the upper  region of the wake than on the lower region.Differences between the RKE-SM and the SA-OG simulations are believed to be due to the differences between the two meshes (see Fig. 8b).On the cross-section presented in Fig. 16, the averaged velocity computed with SA-OG is about 5% higher than that computed with RKE-SM.

Velocity profiles
To further investigate the effect of the passive pitch on the wake shear, Figs. 17 and 18 show the phased-averaged streamwise velocity profile along a vertical line at x = 0.3D T and y = 0.The velocity is non-dimensionalised with the free stream velocity at the hub height U hub .Results from RKE-SM and SA-OG for four phase angles are reported: φ = π/2, 2π /3, 5π /6 and π.The velocity is generally lower for negative z/D T , and higher for positive z/D T , showing a jump near z/D T = 0 associated with the hub vortex.This trend is accentuated for the passive pitch (note the dotted line being on the right-hand side of the solid line for z/D T > 0, and vice versa for z/D T < 0).

Vortex structures
The vortex system generated by the turbine is presented in Figs.19a and 19b for the fixed and passive pitch blades, respectively.The flow field is characterised by strong tip and hub vortices, shown by the isosurfaces of the Q -criterion.The figure focuses on the near wake region, where there are no obvious differences between the vortex systems generated by the fixed and the passive pitch blades.As discussed in Section 6.5, however, the wake of the passive pitch is less uniform than that of the fixed pitch, with the flow being faster on the upper part of the wake and slower on the bottom part.Hence, the tip vortices are more tilted for the passive pitch blade.This is shown in Fig. 20, where the trace of the tip vortex is visualised through the out-of-plane vorticity on the meridian plane parallel to the stream.Data is shown for the two phase angles: φ = π/2 and 5/6π .The tip vortices near   the top boundary convect faster than the vortices near the seabed.This difference is accentuated for the passive pitch compared to the fixed pitch, as shown by the dotted vertical lines.

Effect of passive pitch on the far wake
In this section, the effect of the passive pitch on the far wake is investigated.We showed in Fig. 20 that the traces of the tip vortices are accurately resolved up to approximately one diameter downstream of the turbine.Further downstream, coherent vortical structures are dissipated because of the coarser grid.Vortex breakdown occurs around two diameters downstream of the turbine (Lignarolo et al. 2014(Lignarolo et al. , 2015)), hence this complex physics is not resolved by the present simulations.However, the results showed in Sections 6.4 and 6.5 suggest that the higher wake shear enabled by the passive pitch blades could result in earlier vortex breakdown and higher turbulent mixing compared to the fixed pitch blades.This is shown in Fig. 21a, where the non-dimensional kinetic energy flux averaged over a circular area parallel to the rotor disk is plotted versus the streamwise coordinate for up to four rotor diameters downstream of the turbine.The mean kinetic energy flux decreases in the near wake because of the wake expansion.The SA-OG simulations predict a higher kinetic energy flux in the wake, i.e. a higher wake recovery, because of the greater numerical diffusion associated with the coarser grid.Both the SA-OG and the RKE-SM simulations, however, show that the passive pitch enables a higher wake recovery than the fixed pitch.The differences between fixed and passive pitch are further shown in Fig. 21b.Overall, these results suggest that the passive pitch would not only mitigate unsteady load fluctuations but also improve wake recovery and thus the power available to downstream turbines.

Conclusions
A passive pitch mechanism to reduce the fatigue loads associated with unsteady thrust fluctuations on rotor blades is assessed.Three different methodologies are used: an analytical model based on blade element momentum theory and Theodorsen's theory, and two independent numerical codes solving the unsteady Reynolds-averaged Navier-Stokes equations.The passive pitch system consists of a torsional spring acting on the spanwise axis of the blade.The blade passively pitches to feather when the flow velocity increases, and vice versa.This is implemented in all numerical codes as a mass-spring system.
The system is demonstrated on a one-megawatt tidal turbine in a shear current, with a diameter Reynolds number of 2 × 10 7 .The unsteady loads are due to the blade rotation through the sheared current.The computed thrust coefficient is validated at a level of 1% by modelling a pitching blade tested in a wind tunnel at a blade-span Reynolds number of 2 × 10 6 .Furthermore, the power and thrust curves versus the tip speed ratio for the fixed pitch conditions are compared with the numerical simulations of previous authors showing differences within 4%.The passive pitch system is assessed at the optimal tip-speed ratio of 4.5, which maximises the power output.
All three methodologies show that the passively pitching blade allows an unsteady load reduction of about 80%.Flow separation occurs on the suction side of the blades near the root.The passive pitch has a negligible effect on the variations in the area of the separated region.The most significant difference between the wakes downstream of a rotor with fixed and passive pitch blades is in the vertical shear.In fact, by providing a more uniform thrust force over the rotor disks, the passive pitch blades are less efficient in smoothing the onset sheared velocity profile.The grid resolution in the far wake is insufficient to accurately resolve tip vortex instabilities, but the averaged kinetic energy flux in the wake is found to be consistently higher for the passive pitch than for the fixed pitch by both numerical models.
The consistency between the three different methodologies adopted gives confidence in the above conclusions.However, the following limitations of this work should be considered.Firstly, while the hydrodynamic damping is considered analytically and resolved numerically, it remains to be assessed whether a sufficiently low level of mechanical friction can be practically achieved to be negligible.Secondly, it remains to be verified the effectiveness of the passive blade system to mitigate the wide spectrum of load fluctuations experienced by a tidal turbine (Scarlett and Viola, 2019).Future work should also consider how to decrease power fluctuations.In fact, the design of the present passive pitch mechanism is aimed at mitigating thrust fluctuations, but the design can be adapted to mitigate power instead of thrust fluctuations.Notwithstanding all of the above, these results are a strong indication that a passive pitch mechanism would be effective in mitigating unsteady load fluctuations on tidal turbines.Furthermore, they may underpin future research in other types of turbines and compressors, where it is desirable to mitigate unsteady loads.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

D
Drag [N] d Distance between the pitch axis and midchord [m] D T Diameter of the rotor [m] E Φnum Absolute numerical uncertainty [same units as Φ] E Φ val Absolute validation uncertainty [same units as Φ] f ′ Prandtl's tip loss factor exponent [-] h Relative step size [-] I Moment of inertia [kg m 2 rad −1 ] I 0 Moment of inertia of a volume of water enclosed by a surface equal to that of the blade [kg m 2 rad −1 ]I a Added inertia [kg m rad −1 ] k Reduced frequency [-] K a Added stiffness [kg m s −2 rad −1 ] Spring moment [N m] M qs Quasi-steady hydrodynamic moment [N m] M c,0 Centrifugal moment at θ = 0 [N m] M c,θ Slope of the centrifugal moment with θ [N m] N Force normal to the rotor plane [N] Inflow velocity in the direction normal to the rotor plane [m s −1 ] U T Inflow velocity in the direction tangential to the rotor plane [m s −1 ] U ∞ Free stream velocity [m s −1] U Φ h Relative uncertainty due to the step size h [-] U Φnum Relative numerical uncertainty [-] U Φc Relative uncertainty due to the convergence [-] U Φg Relative uncertainty due to the grid [-] U Φp Relative uncertainty due to other parameters [-] U Φr Relative uncertainty due to the round-off [-] U Φt Relative uncertainty due to the time step [-] U hub Free stream velocity at the hub height [m s −1 ] x Streamwise coordinate [m] x p Chordwise coordinate of the pitch axis [m] y Crossflow horizontal coordinate [m] y + Non-dimensional wall unit [-] z Vertical coordinate [m] z hub Hub height from seabed [m]

Fig. 1 .
Fig. 1.Schematic diagram of the tidal turbine operating in a sheared turbulent flow.

Fig. 2 .
Fig. 2. Schematic diagram of a blade section with a torsional spring for passive pitch.

Fig. 3 .Fig. 4 .
Fig. 3. Schematic diagram of the bypass flow on an azimuthal plane parallel to the stream.Shaded regions show where the bypass flow velocity is averaged to compute the blockage correction coefficient τ f .

Fig. 5 .
Fig. 5. Flowchart of the weakly coupled fluid-structure interaction algorithm.The fluid solution is computed either with the SA-OG solver or with the RKE-SM solver.

Fig. 6 .
Fig. 6.Computational domain and mesh of the RKE-SM simulations.(a) Domain, sub-domain and boundary conditions.(b) Blocks and grid around the rotor.(c) Grids near the blade root and the gap between the hub and the blade.(d) C-grid around the blade section.(e) Grid on the surface of the tip.(f) Detail of the O-grid near the tip.

Fig. 7 .
Fig. 7. Computational domain and mesh for the SA-OG simulation.(a) Grid on the surface of the hub and the blade.(b) Body-fitted patched blades.(c) Inner domain.(d) Overall computational domain.

Fig. 8 .
Fig. 8.Comparison between the RKE-SM and SA-OG mesh.(a) Grid on a meridian plane parallel to the stream (y = 0), (b) on a cross-section of the wake at x = 0.3D T and (c) at x = 4D T .

Fig. 9 .
Fig. 9. Computational domain and mesh of the wind tunnel tests modelled for validation of the RKE-SM simulations: (a) domain and boundary conditions, (b) pressure taps locations, (c) grid on the surface of the wall and (d) grid around wingtip.

Fig. 10 .Fig. 11 .
Fig. 10.Fitting of ϕ(h) with Eq. (33) when h is the time step and Φ is (a) C N and (b) C T .

Fig. 12 .
Fig. 12.(a) C N and (b) C T versus the instantaneous angle of attack at mid-span and at the tip computed with the RKE-SM simulations and comparison with experimental data by Piziali (1994).

Fig. 13 .
Fig. 13.(a) Power coefficient C P and (b) thrust coefficient C N versus the tip speed ratio λ computed with the analytical model, RKE-SM, SA-OG, and the simulations of Gretton and Ingram (2011).

Fig. 14 .
Fig. 14.(a) Thrust force coefficients C N versus the angular position of a blade with fixed and passive pitch, computed with the RKE-SM simulations (phase-averaged data) and with the analytical model (ANL).(b) Pitch angle of the blade with passive pitch computed with the RKE-SM simulations (phase-averaged data) and with the analytical model (ANL).

Fig. 15 .
Fig. 15.Phase-averaged azimuthal component of the wall shear stress on the suction side of the blade with fixed (left) and passive (right) pitch: the grey region shows the area with negative values and thus where the flow is separated.The boundary of the grey region is plotted for 12 equally spaced angular positions of the blade with different colours.Data is computed with the RKE-SM simulations.

Fig. 19 .
Fig. 19.Isosurfaces of non-dimensional Q -Criterion, Q = 20, coloured by non-dimensional streamwise velocity for the rotor equipped with (a) fixed and (b) passive pitch blades, computed with the RKE-SM simulations.

Table 1
Tidal turbine specifications and operative conditions.

Table 2
Comparison between the simulation setup of the two numerical simulations, RKE-SM and SA-OG.
would result in a higher hydrodynamic moment that would pitch the foil towards the original angle of attack.It is noted, however, that the system would turn unstable if the blade stalled because ∂M h /∂α would change sign.

Table 3
Validation of force coefficients.