On the velocity at wind turbine and propeller actuator discs

The first version of the actuator disc momentum theory is more than 100 years old. The extension towards very low rotational speeds with high torque for discs with a constant circulation, became available only recently. This theory gives the performance data like the power coefficient and average velocity at the disc. Potential flow calculations have added flow properties like the distribution of this velocity. The present paper addresses the comparison of actuator discs representing propellers and wind turbines, with emphasis on the velocity at the disc. At a low rotational speed, propeller discs have an 5 expanding wake while still energy is put into the wake. The high angular momentum of the wake, due to the high torque, creates a pressure deficit which is supplemented by the pressure added by the disc thrust. This results in a positive energy balance while the wake axial velocity has lowered. In the propeller and wind turbine flow regime the velocity at the disc is 0 for a certain minimum but non-zero rotational speed . At the disc, the distribution of the axial velocity component is non-uniform in all flow states. However, the distribution of 10 the velocity in the plane containing the axis, the meridian plane, is practically uniform (deviation < 0.2 %) for wind turbine disc flows with tip speed ratio λ > 5 , almost uniform (deviation ≈ 2 %) for wind turbine disc flows with λ= 1 and propeller flows with advance ratio J = π, and non-uniform (deviation 5 %) for the propeller disc flow with wake expansion at J = 2π. These differences in uniformity are caused by the different strengths of the singularity in the wake boundary vorticity strength at its leading edge. 15

2 Equations of motion Figure 1 shows the coordinate systems. The disc is placed perpendicular to the undisturbed velocity U 0 , rotating with angular velocity Ω. All vectors are in positive direction, apart from Γ axis , the vortex at the axis, and γ ϕ , the azimuthal component of 60 the wake boundary vortex sheet. The steady Euler equation is valid: with f the force density, in this case distributed at the disc with thickness . The velocity is presented in the cylindrical coordinate system with x pointing downstream: v = {v x , v r , v ϕ }. ρ is the flow density, p the pressure. In some of the equations dimensionless variables for the axial velocity will be used: u d = v x,d /U 0 and u 1 = v x,1 /U 0 , with the subscripts 0,d,1 denoting 65 values far upstream, at the disc and far downstream as indicated in Figure 2. Furthermore, v = {v s , v n , v ϕ } is used, where v s = v 2 x + v 2 r is the velocity component along a streamline at the surface with constant Ψ, with Ψ denoting the Stokes stream function.
The pressure and azimuthal velocity are discontinuous across the disc when → 0. For such an infinitely thin disc, integration of Eq. (1) yields: where ∆ denotes the jump across the disc, and F the applied surface load. A Joukowsky disc has a wake with swirl, induced by a vortex Γ at the axis. The vortex core radius δ is assumed to be infinitely thin. The azimuthal velocity is: The Bernoulli equation reads: s 5.1 spherecontrolvo df Figure 2. The stream tube of a propeller disc from cross sections A0, infinitely far upstream, to A1 in the fully developed wake. Only the upper half of the stream tube is shown.
When this is integrated across the disc and combined with Eq. (3), the axial component of Eq.
(2) becomes: The power converted by an annulus dr of the actuator disc equals the torque Q times rotational speed Ω giving ΩdQ = 2πΩf ϕ r 2 dr, but also the integrated value of f · v with the use of Eq. (1), giving 2πr(v.∇)Hdr. Equating both expressions 80 shows that: f ϕ is expressed by the ϕ-component of the Euler equation (1): Herewith: which gives with Eq. (6): Consequently, for a Joukowsky disc: In the wind turbine mode ∆H < 0, as energy is taken from the flow. With Ω always taken positive, Γ and v ϕ are negative in the wind turbine mode, and positive in the propeller mode. This explains why Γ axis is shown with a negative sign in Fig. 1.

90
Furthermore Eq. (9) shows that for Ω → ∞ meanwhile keeping ∆H constant, Γ vanishes and, by Eq. (7), also f ϕ . The result is the Froude disc without torque and wake swirl.
The power P converted by the disc follows by integration of Eq. (6) on the actuator disc. In dimensionless notation this becomes:

95
With λ = ΩR/U 0 and q = Γ/(2πRU 0 ), Eq. (9) becomes: and similarly Eq. (10): The thrust T is derived in the same way, based on Eq. (5). Dimensionless, the thrust coefficient is C T = T /( 1 2 ρU 2 0 A d ) = 100 C T,∆H + C T,∆vϕ according to the two terms on the right-hand side of (5): C T,∆vϕ does not contribute directly to the conversion of power, as it does not appear in Eq. (12). It is a conservative contribution to C T , delivering the radial pressure gradient balancing the swirl immediately behind the disc. For finite q and δ → 0, For a non-zero δ combined with high λ, low q, C T,∆vϕ becomes small. For typical wind turbine parameters 105 λ = 8, C T,∆H = −8/9 and δ = 0.05R with δ representing the blade root cut-out area, C T,∆vϕ ≈ −0.02.
The power and thrust have the same sign as ∆H or q: positive for propeller discs, negative for wind turbine discs. Consequently, the thrust and power (coefficients) are negative for discs extracting energy from the wake, and positive for discs adding energy to the wake.
The velocity in the far wake is characterised by v r = 0. Herewith the Bernoulli equation (4) becomes in the far wake: The radial derivative is ∂p 1 /∂r 1 = ρ(v 2 ϕ,1 /r 1 − v x,1 ∂v x,1 /∂r). When this is compared with the condition for radial pressure equilibrium in the fully developed wake, given by substitution of v r = 0 in the radial component of Eq. (1): The momentum theory presented in van Kuik (2017) is valid when a different sign convention for q is used, as in van Kuik (2017) it was defined q = −Γ/(2πRU 0 ) instead of Γ/(2πRU 0 ). This theory lacks an analytical solution. However, a numerical solution of equation 29 of this paper is possible. Expressed in λ and q, this is an implicit expression for u 1 : where q has changed sign. After solving Eq. (16) for u 1 , the wake expansion or contraction is given by van Kuik (2017, eq. 28).
The average velocity at the disc u d is given by van Kuik (2017, eq. 27), again with a change of sign of q in both equations. Figure 3 shows u d for −1 < C T,∆H ≤ +2, 0 < λ ≤ 5 as well as λ = ∞. The advance ratio J = π/λ is also given. The part of the figure with C T,∆H < 0 shows u d for wind turbine discs, with C T,∆H > 0 for propeller discs. For λ = 5 the difference 125 with λ = ∞ is smaller than 0.7 % so the Froude momentum theory results are practically recovered. Apparently, swirl has little effect when λ > 5. The flow cases a to e are defined in Table 1, together with two flow parameters: the dimensionless average velocity at the disc, u d and the dimensionless absolute velocity in the meridian plane |v| m = v 2 x,d + v 2 r,d /U 0 . |v| m is the same as the velocity along a streamline v s at the position of the disc, so |v| m = v s,d /U 0 . u d and |v| m and will be examined in the next sections. Table 1. Definition of actuator disc flow cases a to e, the average velocity at the disc u d and the absolute velocity in the meridian plane |v|m. Several particularities can be observed in Fig. 3: -For values of λ < 1.4 the minimum attainable C ∆H > −1, giving u d = 0, so the flow is blocked. Such a minimum λ exists in the wind turbine as well as propeller flow regime.
-For wind turbine discs having λ > 1.4 the minimum C ∆H is −1.0, with u d shrinking from 0.5 at λ = 5 to 0 at λ ≈ 1.4.
-For propeller discs having a very high J, u d < 1 so the wake expands. This upper boundary of the expanding wake 135 region is the line u d = 1, giving an undeformed wake. The lower boundary is defined by u d = 0, giving blocked flow.
Both boundaries put a limit to the maximum attainable C T,∆H . For low J there is no upper limit for C T,∆H : the wake can be accelerated to any value.
These particularities will be discussed in the next subsections, to start with the propeller disc.
3.2 Propeller discs having an expanding wake 140 For low rotational speed (low λ, high J), the average axial velocity at the disc u d deviates from the famous Froude result: u d < 1 2 (u 1 + 1). This happens in both flow regimes. Responsible for this is the radial pressure distribution necessary to maintain the swirl. This gives a contribution to the momentum balance, as is explained in van Kuik (2018b, chapter 6). The first term in the disc load Eq. (5) gives the contribution of ∆H to the disc load, the second term the swirl related pressure contribution. This contribution − ρ 2 (Γ/(2πr)) 2 is always < 0, while the sign of the first term depends on the actuator disc mode: for wind turbine 145 discs < 0, for propeller discs > 0. Consequently, both terms may cancel for propeller flows, resulting in a zero pressure jump at r = R. With Eqs. (5) and (9) the condition for this particular flow is derived: ΩR = − 1 2 v ϕ or: In Fig. 3 this specific flow regime is indicated by the line separating the propeller disc regime with a contracting wake, from the propeller disc regime with an expanding wake, with u d = 1 at the separation line. The resulting flow has a wake with constant 150 radius, so v x = U 0 , v r = 0 throughout the flow. In the wake v ϕ = Γ/(2πr). The vortex sheet separating the wake from the outer flow consists of axial vorticity across which ∆H = 1 2 (ΩR) 2 . For lower rotational speeds the pressure jump at the edge has become < 0, as the swirl-related pressure term in Eq. (5) overrules the ∆H term, thereby generating wake boundary vorticity as for wind turbine disc flows. Although kinetic energy in the wake is lower than outside the wake, the disc load adds potential energy (pressure) to the flow such that the total energy 155 in the wake is higher than upstream. More explanation of this remarkable flow regime is provided in van Kuik (2018b, section 6.3).

Minimum λ operation with blocked flow
In van Kuik (2018b, section 6.3) the operation at minimum possible λ is analysed. In this flow case u d = 0 as well as u 1 = 0, so the disc acts as a blockage to the flow. In the wake the change of axial momentum is zero, but H wake − H 0 = 0 as the The streamlines indicate stream tube values increasing with ∆Ψ = 0.1Ψ1. Table 1 show the flow cases, also indiated in Fig. 3, for which the flow field has been calculated numerically with the potential flow method used in van Kuik (2017). An assessment of the accuracy presented in van Kuik (2018b, appendix D). The highest attainable accuracy is applied: calculated values of integrated properties like wake expansion of contraction deviate less than 165 0.3 % from momentum theory values. The same holds for the local satisfaction of the boundary conditions at the wake boundary v n = 0, ∆p = 0, except within a distance 0.02R from the disc edge, where v n may deviate up to 0.02U 0 without challenging the condition Ψ = Ψ 1 and without affecting integrated flow quantities.

Flow patterns
In Fig. 4 the streamlines of flow cases a to e are shown, grouped according to their position in Fig. 3. Flow cases a looks similar to flow case e, although the latter is a propeller disc flow.

4 The velocity distribution at the disc
With v s being the velocity in the meridian plane, v s /U 0 at the upstream side of the disc equals |v| m = v 2 x,d + v 2 r,d /U 0 . Table 1 gives the numerical values of u d and |v| m for the flow cases considered. The differences between u d as calculated numerically and as resulting from the momentum theory is 0.2% or less. The value for |v| m is the value for r = 0. Figure 5 shows the distribution of the axial and radial velocity components and the meridional velocity. Most striking is the distribution 175 of this meridional velocity being practically uniform. The explanation of this is presented in section 5, but first the velocity distributions shown in Fig. 5 are analysed. -Discs with λ = ∞ but heavier disc loads: the non-uniformity in |v| d− is −0.7 % for C T,∆H = −0.97, −0.8 % for

The meridional velocity
The optimal operational regime of modern wind turbines is λ > 5 with C T ∆H > −0.9, so the non-uniformity in |v| m of flow cases representing this optimal regime is negligible.

Propeller flows
The non-uniformity in |v| m is 2 % in flow case b, J = 0. It decreases to 1.3 % in flow case d, J = π, becomes 0 for J = 1.5π when the flow case without wake deformation is reached according to Eq. (17), and becomes strongly negative for higher J as shown in flow case e: −5 % for J = 2π. Usually the advance ratio J is lower than 2.5, see for example McCormick (1994, figure 6.12). Fig. 3 shows that in this regime the impact of wake swirl is very limited, so flow case b is considered 200 representative, with a non-uniformity of ≈ 2 %.

The axial velocity
In all flow cases the axial velocity is far from uniform, as was already shown by e.g. Sørensen et al. (1998), Madsen et al. (2010). For Froude wind turbine discs, the cause of this has been addressed in van Kuik and Lignarolo (2016), for Joukowsky disc flows in van Kuik (2017). In terms of the momentum balance, the source of this non-uniformity is the pressure acting 205 on the sides of a stream annulus used as control volume. When the stream tube boundary is used as boundary of the control volume, the pressure at this boundary does not give a contribution in axial direction, but for stream annuli this is not the case.
When this pressure is calculated and included in the momentum balance, the prediction of u d per annulus by the momentum theory matches the calculated, non-uniform distribution of the u d . This may serve as the explanation of the non-uniformity of u d , but cannot be used as a prediction model, as the pressure is not known a priori. For Froude discs the u d distribution has 210 been calculated for −1 < C T,∆H < 0 enabling a surface-fit engineering approximation for u d ( r R , C T,∆H ), see van Kuik and Lignarolo (2016, section 5.2).

The radial velocity
The radial velocity receives little attention in actuator disc and rotor publications compared to the axial velocity. Some exceptions are Madsen et al. (2010), presenting an engineering model for the decreased axial velocity close to the disc or rotor edge 215 based on the radial velocity, Micallef et al. (2013), comparing calculated and measured radial velocity near rotor blade tips to asses blade bound chordwise vorticity in order to explain the initially inward motion of the tip vortex, and van Kuik et al.
(2014, section 4), quantifying this chordwise vorticity and the associated tip load responsible for this inward tip vortex motion, and Sørensen (2015, section 3.2), analysing ∂v r /∂x at the plane of the disc.
Recently Limacher and Wood (2019) found a relation between the axial and radial velocity component in the rotor or disc 220 plane: where S d is the plane of the disc or rotor from r = 0 to r = ∞. The quantity 1 − v x,d /U 0 is known as the induction a. Based on (18) Limacher and Wood (2019) conclude that v r,d /U 0 and a have to be equal close to the disc edge or rotor tip, so: Eqs. (18) and (19) have been evaluated using the velocity distributions of Fig. 5. For flow case a the left-hand side of Eq. (18) indeed approaches 0 for increasing radius of S d . Table 2 gives the radial coordinate where (19) is satisfied: almost at the disc edge for the flow states with an expanding wake a,c,e, while flow states b,d with a contracting wake show this property at a smaller radius. The expanding flows exhibit steep changes in v x and v r close to r = R,. An exact assessment the radial position where Eq. (19) is satisfied is difficult for which reason a range of r is given.
230 Eq. (19) provides a second relation between v x,d and v r,d , besides the conclusion of section 4.1 that |v| m is practically constant for r < R. This allows an engineering estimate of the wake expansion at the disc for wind turbine flows, when it is assumed that |v| m = v 2 x,d + v 2 r,d /U 0 = constant and v x,d +v r,d = U 0 at r/R = 1. As an example the flow with v x,d = v r,d is evaluated, giving |v| m = 0.707 and v x,d = v r,d = 0.5U 0 at r = R. This gives a slope of the vortex sheet shape of 45 0 at r = R. This is close to flow state a, where the numerically calculated slope is 46 0 , and |v| m = 0.684 which is 3, 3% lower than the 235 estimate. Further exploration of such an engineering estimate is left for future work.

Explanation of the (non-)uniformity of |v| m
The Euler equation of motion (1) offers a first-order explanation for the observation that |v| m is practically uniform for λ ≥ 1.
The radial component of Eq.

240
Equation (3) for v ϕ combined with Bernoulli's equation (4) gives a second equation for ∂p/∂r: so the result is ∂v s /∂r = ∂v r /∂s, or at the disc : Consequently, the distribution of v s,d is determined by the derivative ∂v r /∂s along the streamline. In case v r has a maximum 245 or minimum at the disc, v s,d /U 0 = |v| m is uniform.
Qualitative observations regarding the in-or decrease of v r are possible when moving the position along a streamline in the meridian plane. The radial velocity depends only on the vorticity γ ϕ distributed along the wake boundary, and the position of obervation s * . For a disc with an expanding wake, the following relations hold: a) At the upwind side of the streamline: when moving towards the disc, the distance to γ ϕ decreases, so v r increases, and 250 ∂v r /∂s > 0.
b) At the downwind side of the disc the streamline is to be distinguished in two parts: upstream and downstream of s * .
The upstream vorticity induces a negative v r,upstream , becoming more negative when s * moves downstream, leading to ∂v r,upstream /∂s < 0. The part of the wake downstream of s * remains a semi-infinite wake, so v r,downstream is expected to vary only little for increasing s * (this is to be verified later), leading to ∂v r,downstream /∂s ≈ 0. This gives for the total 255 induction in the wake ∂v r /∂s < 0.
Consequently, according to a) and b) ∂v r /∂s = 0 at the disc and with Eq. (22) ∂v s,d /∂r = 0 so |v| m is uniform.
For flow cases with a contracting wake the same reasoning is valid, with an appropriate change of signs, leading to a minimum v r at the disc and a uniform |v| m .
However, these qualitative considerations miss the effect that a vortex ring induces a non-zero ∂v r /∂s in the plane of 260 the ring. Figure 6 shows the calculated radial velocity induced by a vortex ring positioned at x = 0, R = 1 along the lines r/R = 0.8, 0.9, 0.97. The shape of the plot resembles the induction v r = Γ 2π cosα dist by a point vortex in a 2−D plane, where dist is the distance to the vortex, and α the angle the angular coordinate around the vortex position. As is clear by Fig. 6, this effect is strongest close to the position of the ring, as ∂v r /∂x → ∞ for r/R → 1. Apart from the distance to the ring, the strength of the ring determines the local value of ∂v r /∂x, as its value is linear in this strength.

265
For a vorticity tube things are slightly different, as is easily shown by the example of a tube of constant strength with a semifinite length. Each elementary vortex ring γdx induces a non-zero ∂v r /∂r in its own plane, but due to symmetry considerations this is annihilated except near and at the beginning of the tube. Also for the vorticity tube surrounding the actuator disc wake, the singular behaviour of ∂v r /∂s is annihilated everywhere by the induction of upstream and downstream vorticity, except at the leading edge of the wake. There the sign of the contribution to ∂v r /∂r at x = 0 is opposite to the sign of ∂v r /∂r at x < 0, 270 as is clear from Fig. 6.
The considerations a) plus b) mentioned above include the effect of increasing or decreasing distance to the vorticity and the change of sign of v r for vorticity upstream and downstream of s * , but the non-zero ∂v r /∂s at s * = s d has to be added: The absolute value of the slope of the tangents is lowest in flow case a, highest in e. This is in agreement with the strength of the leading edge singularity of γ ϕ (x)/γ ϕ,1 and the non-uniformity of v m . In all flow cases v r (s) reaches a maximum or 290 minimum just upstream of the disc: at s/R = −0.00155 for a, and −0.00252 for e, with the values for other flow cases in between these positions.

Conclusions
With respect to the average velocity at the actuator disc: -For Joukowsky discs in wind turbines and propellers mode, the average velocity has been found, from λ = 0 up to 295 λ → ∞, or J → ∞ to J = 0.
-For a very high J, propeller disc flows have an expanding wake while still energy is put into the wake. The high angular momentum of the wake flow creates a pressure deficit in the wake, which is supplemented by the pressure added by the disc. This results in a positive energy balance while the wake axial velocity has gone down.
-Propeller discs flows without wake expansion or contraction are possible for specific values of J, marking the transition 300 from the contracting wake operational mode at low J, to the expanding wake mode at high J.
-In the propeller as well as wind tubine flow regimes the velocity at the disc becomes 0 for very low rotational speed, resulting in a flow with a blocked disc.
With respect to the distribution of the velocity in the meridian plane at the disc position: -|v| m is practically uniform for wind turbine disc flows with λ > 5 (deviation in the order of a few ‰)

305
-|v| m is almost uniform for wind turbine disc flows with low λ and propeller flows with J ≈ π (deviation in the order of a few %) -|v| m is non-uniform for the propeller disc flow with wake expansion at very high J (deviation in the order of several %).
the differences in uniformity are caused by the different strengths of the leading edge singularity in the wake boundary vorticity strength.