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 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 actuator disc flows. However, the distribution of 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.

1 Introduction

The start of rotor aerodynamics dates back more than 100 years, when the concept of the actuator disc to represent the action of a propeller was formulated by Froude (1889). In this concept the disc carries only thrust, no torque. Based on this Joukowsky (1918) published the first performance prediction that still holds today, for a hovering helicopter rotor or a propeller without forward speed. About 2 years later Joukowsky (1920) and Betz (1920) published the optimal performance of discs representing wind turbines, for which reason it is called the Betz–Joukowsky maximum (Okulov and van Kuik, 2012). The names of Betz and Joukowsky are also connected with the two concepts for actuator discs with thrust and torque. The model of Betz (1919) was similar to the vortex model of Prandtl for an elliptically loaded wing. This gives an induced velocity which is constant over the wing span, resulting in minimum induced drag. In Betz's model each rotor blade is represented by a lifting line such that the vortex sheet released by the blade has a constant axial velocity. Joukowsky (1912) developed the vortex model of a propeller based on a horseshoe vortex of a wing. In his model each blade is modelled by a lifting line with constant circulation.

The constant circulation model of Joukowsky as well as the constant velocity model of Betz represented the ideal rotor. It was not yet possible to compare the models and to conclude which was best. Both models were valid only for lightly loaded rotors as wake expansion or contraction was neglected. A solution for the wake of Betz's rotor, still restricted to lightly loaded propellers, was presented by Goldstein (1929). The non-linear solution, so including wake deformation, was published by Okulov (2014) and Wood (2015). A comparison of the models of Betz and Joukowsky for rotors was presented by Okulov et al. (2015, chap. 4) showing that Joukowsky rotors perform somewhat better than Betz rotors when both operate at the same tip speed ratio. The same conclusion was drawn for actuator discs by van Kuik (2017): at low tip speed ratio the Joukowsky disc performs somewhat better than the Betz disc. For increasing tip speed ratios, both models become the same as they converge to Froude's actuator disc.

The Joukowsky and Froude discs are still subjects for research as many modern design and performance prediction codes are based on them; see e.g. Sørensen (2015). Over the last decades the disc received the most attention from the wind energy research community, but recently new propeller research on the actuator disc concept has been published; see Bontempo and Manna (2018a, b, 2019). The performance aspects are known by many studies using momentum theory, vorticity or computational fluid dynamics (CFD) methods. Experimental verification is shown by e.g. Lignarolo et al. (2016) and Ranjbar et al. (2019). Recent research aims for deriving efficient tip corrections (see e.g. Moens and Chatelain, 2018; Zhong et al., 2019) or for configurations including a hub (Bontempo and Manna, 2016) or duct (Dighe et al., 2019).

The present paper addresses the topic which received the least attention: the velocity distribution at the disc. The paper is part of a sequence of papers, starting with van Kuik and Lignarolo (2016) concerning flows through wind turbine Froude discs calculated by a potential flow method, followed by van Kuik (2017) concerning the momentum theory and potential flow calculations for wind turbine Joukowsky discs, and the conference paper van Kuik (2018a) where the extension to propeller discs was presented. The latter paper was not yet conclusive in the explanation of the difference between wind turbine and propeller discs regarding the velocity distribution at the disc: for wind turbine discs the velocity vector in the plane containing the disc axis, the meridian plane, seems to be uniform, while it seems non-uniform for propeller discs. Compared to van Kuik (2018a) all calculations have been redone at equal, highest possible accuracy, leading to slightly different quantitative conclusions and a consistent explanation for the (non-)uniformity of the velocities at the actuator discs. Some of the content of van Kuik (2018a) regarding the average velocity at the disc is repeated, in order to make the paper readable independently of the previous papers. The open-access book van Kuik (2018b) contains the content of all papers mentioned in this paragraph.

Section 2 presents the equations of motion and the coordinate system. Section 3 discusses the average velocity at the disc and some remarkable disc flows, followed by Sect. 4 treating the velocity distribution for both actuator disc modes. Section 5 analyses the differences observed between wind turbine and propeller discs, followed by the concluding Sect. 6.

2 Equations of motion

Figure 1 shows the coordinate systems. The disc is placed perpendicular to the undisturbed velocity U₀, rotating with angular velocity Ω. All vectors are in the positive direction, apart from Γ_axis, the vortex at the axis, and γ_φ, the azimuthal component of the wake boundary vortex sheet. The steady Euler equation is valid:

$$\begin{array}{}\text{(1)}& {\displaystyle }{\displaystyle }\mathit{\rho }\left(\mathit{v}\cdot \mathbf{\nabla }\right)\mathit{v}=-\mathbf{\nabla }p+\mathit{f},\end{array}$$

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: $\mathit{v}=\mathit{\{}{v}_x,v_r,v_{\mathit{\phi }}\mathit{\}}$. ρ is the flow density and p the pressure. In some of the equations dimensionless variables for the axial velocity will be used: $u_{\mathrm{d}}=v_{x,\mathrm{d}}/U_{\mathrm{0}}$ and $u_{\mathrm{1}}=v_{x,\mathrm{1}}/U_{\mathrm{0}}$, with the subscripts ${}_{\mathrm{0},\mathrm{d},\mathrm{1}}$ denoting values far upstream, at the disc and far downstream as indicated in Fig. 2. Furthermore, $\mathit{v}=\mathit{\{}{v}_{\mathrm{s}},v_n,v_{\mathit{\phi }}\mathit{\}}$ is used, where $v_{\mathrm{s}}=\sqrt{v_x^{\mathrm{2}}+v_r^{\mathrm{2}}}$ is the velocity component along a streamline at the surface with constant Ψ, with Ψ denoting the Stokes stream function.

Figure 1 The coordinate system of an actuator disc acting extracting energy. Ψ is the Stokes stream function. All vectors are in the positive direction except Γ_axis and γ_φ.

Figure 2 The stream tube of a propeller disc from cross sections A₀, infinitely far upstream, to A₁ in the fully developed wake. Only the upper half of the stream tube is shown.

The pressure and azimuthal velocity are discontinuous across the disc when ϵ→0. For such an infinitely thin disc, integration of Eq. (1) yields

$$\begin{array}{}\text{(2)}& {\displaystyle }{\displaystyle }\mathit{F}=\underset{\mathit{ϵ}\to \mathrm{0}}{lim}\underset{\mathit{ϵ}}{\int }\mathit{f}\mathrm{d}x={\mathit{e}}_x\mathrm{\Delta }p+{\mathit{e}}_{\mathit{\phi }}\mathit{\rho }{v}_x\mathrm{\Delta }{v}_{\mathit{\phi }},\end{array}$$

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

$$\begin{array}{}\text{(3)}& {\displaystyle }{\displaystyle }{v}_{\mathit{\phi }}={\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }r}}.\end{array}$$

The Bernoulli equation reads

$$\begin{array}{}\text{(4)}& {\displaystyle }{\displaystyle }p+{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }\mathit{v}\cdot \mathit{v}=H.\end{array}$$

When this is integrated across the disc and combined with Eq. (3), the axial component of Eq. (2) becomes

$$\begin{array}{}\text{(5)}& {\displaystyle }{\displaystyle }{F}_x=\mathrm{\Delta }p=\mathrm{\Delta }H-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }\mathrm{\Delta }{v}_{\mathit{\phi }}^{\mathrm{2}}=\mathrm{\Delta }H-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\mathit{\rho }{\left({\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }r}}\right)}^{\mathrm{2}}.\end{array}$$

The power converted by an annulus dr of the actuator disc equals the torque Q times rotational speed Ω, giving ΩdQ=2πΩf_φr²dr, but also the integrated value of f⋅v with the use of Eq. (1), giving 2πr(v⋅∇)Hdr. Equating both expressions shows that

$$\begin{array}{}\text{(6)}& {\displaystyle }{\displaystyle }\mathit{f}\cdot \mathit{v}=\mathrm{\Omega }r{f}_{\mathit{\phi }}=(\mathit{v}\cdot \mathbf{\nabla })H.\end{array}$$

f_φ is expressed by the φ component of the Euler Eq. (1): $f_{\mathit{\phi }}=\mathit{\rho }(\mathit{v}\cdot \mathbf{\nabla })v_{\mathit{\phi }}$. Herewith

$$\begin{array}{}\text{(7)}& {\displaystyle }{\displaystyle }\mathrm{\Omega }r{f}_{\mathit{\phi }}=\mathit{\rho }(\mathit{v}\cdot \mathbf{\nabla })\left(\mathrm{\Omega }r{v}_{\mathit{\phi }}\right),\end{array}$$

which gives with Eq. (6)

$$\begin{array}{}\text{(8)}& {\displaystyle }{\displaystyle \frac{\mathrm{1}}{\mathit{\rho }}}\mathbf{\nabla }H=\mathbf{\nabla }\left(\mathrm{\Omega }r{v}_{\mathit{\phi }}\right)=\mathbf{\nabla }\left({\displaystyle \frac{\mathrm{\Omega }\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }}}\right).\end{array}$$

Consequently, for a Joukowsky disc

$$\begin{array}{}\text{(9)}& {\displaystyle }{\displaystyle }\mathrm{\Delta }H=\mathit{\rho }{\displaystyle \frac{\mathrm{\Omega }\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }}}=\text{constant}.\end{array}$$

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. 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

$$\begin{array}{}\text{(10)}& {\displaystyle }\begin{array}{rl}{C}_p& ={\displaystyle \frac{\mathrm{1}}{\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }{U}_{\mathrm{0}}^{\mathrm{3}}{A}_{\mathrm{d}}}}\underset{A}{\int }\mathit{f}\cdot \mathit{v}\mathrm{d}{A}_{\mathrm{d}}={\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}{\displaystyle \frac{\mathrm{\Delta }H}{\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }{U}_{\mathrm{0}}^{\mathrm{2}}}}\\ & =\mathrm{2}{\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}{\displaystyle \frac{\mathrm{\Omega }R}{U_{\mathrm{0}}}}{\displaystyle \frac{\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }R{U}_{\mathrm{0}}}}.\end{array}\end{array}$$

With $\mathit{\lambda }=\mathrm{\Omega }R/U_{\mathrm{0}}$ and $q=\mathrm{\Gamma }/\left(\mathrm{2}\mathit{\pi }R{U}_{\mathrm{0}}\right)$, Eq. (9) becomes

$$\begin{array}{}\text{(11)}& {\displaystyle \frac{\mathrm{\Delta }H}{\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }{U}_{\mathrm{0}}^{\mathrm{2}}}}=\mathrm{2}q\mathit{\lambda },\end{array}$$

and similarly Eq. (10) becomes

$$\begin{array}{}\text{(12)}& {\displaystyle }{\displaystyle }{C}_p=\mathrm{2}q\mathit{\lambda }{\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}.\end{array}$$

The thrust T is derived in the same way, based on Eq. (5). Dimensionless, the thrust coefficient is $C_T=T/\left(\frac{\mathrm{1}}{\mathrm{2}}\mathit{\rho }{U}_{\mathrm{0}}^{\mathrm{2}}{A}_{\mathrm{d}}\right)=C_{T,\mathrm{\Delta }H}+C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}$ according to the two terms on the right-hand side of (5):

$$\begin{array}{}\text{(13)}& {\displaystyle }{\displaystyle }\left.\begin{array}{lcl}{\displaystyle }{C}_T& =& {\displaystyle }{C}_{T,\mathrm{\Delta }H}+C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}\\ {\displaystyle }{C}_{T,\mathrm{\Delta }H}& =& {\displaystyle }\mathrm{2}\mathit{\lambda }q\\ {\displaystyle }{C}_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}& =& -q^{\mathrm{2}}\ln {\left(\frac{R}{\mathit{\delta }}\right)}^{\mathrm{2}}\end{array}\right\}.\end{array}$$

$C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}$ 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, $C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}\to \mathrm{\infty }$. For a non-zero δ combined with high λ and low q, $C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}$ becomes small. For typical wind turbine parameters λ=8, $C_{T,\mathrm{\Delta }H}=-\mathrm{8}/\mathrm{9}$ and δ=0.05R with δ representing the blade root cut-out area, $C_{T,\mathrm{\Delta }{v}_{\mathit{\phi }}}\approx -\mathrm{0.02}$.

The power and thrust have the same sign as ΔH or q: positive for propeller discs and 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 characterized by v_r=0. Herewith the Bernoulli Eq. (4) becomes in the far wake

$$\begin{array}{}\text{(14)}& {\displaystyle }{\displaystyle \frac{\mathrm{1}}{\mathit{\rho }}}\left(p_{\mathrm{0}}-p_{\mathrm{1}}\right)={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\left(v_{x,\mathrm{1}}^{\mathrm{2}}-U_{\mathrm{0}}^{\mathrm{2}}+v_{\mathit{\phi },\mathrm{1}}^{\mathrm{2}}\right)-\mathrm{\Delta }H.\end{array}$$

The radial derivative is $\partial p_{\mathrm{1}}/\partial r_{\mathrm{1}}=\mathit{\rho }(v_{\mathit{\phi },\mathrm{1}}^{\mathrm{2}}/r_{\mathrm{1}}-v_{x,\mathrm{1}}\partial v_{x,\mathrm{1}}/\partial 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),

$$\begin{array}{}\text{(15)}& {\displaystyle }{\displaystyle \frac{\partial p_{\mathrm{1}}}{\partial r_{\mathrm{1}}}}=\mathit{\rho }{\displaystyle \frac{v_{\mathit{\phi },\mathrm{1}}^{\mathrm{2}}}{r_{\mathrm{1}}}},\end{array}$$

the result is $v_{x,\mathrm{1}}=\text{constant}$ or, dimensionless, $u_{\mathrm{1}}=v_{x,\mathrm{1}}/U_{\mathrm{0}}=\text{constant}$.

3 Flow pattern and average velocity

3.1 Momentum theory results for propeller and wind turbine discs

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=-\mathrm{\Gamma }/\left(\mathrm{2}\mathit{\pi }R{U}_{\mathrm{0}}\right)$ instead of Γ∕(2πRU₀). This theory lacks an analytical solution. However, a numerical solution of Eq. (19) of van Kuik (2017) is possible. Expressed in λ and q, this is an implicit expression for u₁:

$$\begin{array}{ll}{\displaystyle }& {\displaystyle \frac{\left(\mathrm{1}-u_{\mathrm{1}}\right)u_{\mathrm{1}}^{\mathrm{2}}{q}^{\mathrm{2}}}{\mathrm{1}-\mathrm{2}\mathit{\lambda }q-u_{\mathrm{1}}^{\mathrm{2}}}}\\ \text{(16)}& {\displaystyle }& {\displaystyle }=\left(-q\mathit{\lambda }-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{q}^{\mathrm{2}}\left(\mathrm{1}-\ln \left({\displaystyle \frac{q^{\mathrm{2}}}{\mathrm{1}-\mathrm{2}\mathit{\lambda }q-u_{\mathrm{1}}^{\mathrm{2}}}}\right)\right)\right),\end{array}$$

where q has changed sign. After solving Eq. (16) for u₁, the wake expansion or contraction is given by van Kuik (2017, Eq. 28). The average velocity at the disc ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ is given by van Kuik (2017, Eq. 27), again with a change of sign of q in both equations.

Figure 3 The axial velocity ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ for wind turbine discs ($-\mathrm{1}<C_{T,\mathrm{\Delta }H}\le \mathrm{0}$) and propeller discs ($\mathrm{0}\ge C_{T,\mathrm{\Delta }H}<\mathrm{2}$) for $\mathrm{0}\le \mathit{\lambda }\le \mathrm{5}$ and for λ=∞. The white markers a to e refer to flow cases defined in Table 1 and analysed in the next sections. The figure is a modified version of van Kuik (2018b, Fig. 6.2).

Figure 3 shows ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ for $\mathrm{0}<\mathit{\lambda }\le \mathrm{5}$ as well as λ=∞ and $-\mathrm{1}<C_{T,\mathrm{\Delta }H}\le \mathrm{2}$. The advance ratio $J=\mathit{\pi }/\mathit{\lambda }$ is also given. The part of the figure with $C_{T,\mathrm{\Delta }H}<\mathrm{0}$ shows ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ for wind turbine discs and with $C_{T,\mathrm{\Delta }H}>\mathrm{0}$ for propeller discs. For λ=5 the difference 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, ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$, and the dimensionless absolute velocity in the meridian plane $|\mathit{v}{|}_{\mathrm{m}}=\sqrt{v_{x,\mathrm{d}}^{\mathrm{2}}+v_{r,\mathrm{d}}^{\mathrm{2}}}/U_{\mathrm{0}}$. $|\mathit{v}{|}_{\mathrm{m}}$ is the same as the velocity along a streamline v_s at the position of the disc, so $|\mathit{v}{|}_{\mathrm{m}}=v_{\text{s,d}}/U_{\mathrm{0}}$. ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ and $|\mathit{v}{|}_{\mathrm{m}}$ will be examined in the next sections.

Table 1 Definition of actuator disc flow cases a to e, the average velocity at the disc ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ and the absolute velocity in the meridian plane $|\mathit{v}{|}_{\mathrm{m}}$.

Several particularities can be observed in Fig. 3:

• For values of λ<1.4 the minimum attainable $C_{\mathrm{\Delta }H}>-\mathrm{1}$, giving ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}=\mathrm{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 ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ shrinking from 0.5 at λ=5 to 0 at λ≈1.4.

• For propeller discs having a very high J, ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}<\mathrm{1}$ so the wake expands. This upper boundary of the expanding wake region is the line ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}=\mathrm{1}$, giving an undeformed wake. The lower boundary is defined by ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}=\mathrm{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

For low rotational speed (low λ, high J), the average axial velocity at the disc ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ deviates from the famous Froude result: ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}<\frac{\mathrm{1}}{\mathrm{2}}(u_{\mathrm{1}}+\mathrm{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, Chapt. 6). The first term in the disc load Eq. (5) gives the contribution of ΔH to the disc load and the second term the swirl related pressure contribution. This contribution $-\frac{\mathit{\rho }}{\mathrm{2}}(\mathrm{\Gamma }/(\mathrm{2}\mathit{\pi }r){)}^{\mathrm{2}}$ is always <0, while the sign of the first term depends on the actuator disc mode: for wind turbine 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: $\mathrm{\Omega }R=-\frac{\mathrm{1}}{\mathrm{2}}{v}_{\mathit{\phi }}$ or

$$\begin{array}{}\text{(17)}& {\displaystyle }\mathit{\lambda }=q/\mathrm{2}.\end{array}$$

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 ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}=\mathrm{1}$ at the separation line. The resulting flow has a wake with a constant radius, so v_x=U₀, v_r=0 throughout the flow. In the wake $v_{\mathit{\phi }}=\mathrm{\Gamma }/\left(\mathrm{2}\mathit{\pi }r\right)$. The vortex sheet separating the wake from the outer flow consists of axial vorticity across which $\mathrm{\Delta }H=\frac{\mathrm{1}}{\mathrm{2}}(\mathrm{\Omega }R{)}^{\mathrm{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 in the wake is higher than upstream. More explanation of this remarkable flow regime is provided in van Kuik (2018b, Sect. 6.3).

3.3 Minimum λ operation with blocked flow

In van Kuik (2018b, Sect. 6.3) the operation at minimum possible λ is analysed. In this flow case ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}=\mathrm{0}$ as well as u₁=0, so the disc acts as a blockage to the flow. In the wake the change of axial momentum is zero, but $H_{\text{wake}}-H_{\mathrm{0}}\ne \mathrm{0}$ as the azimuthal velocity is non-zero. Lower values of λ are not possible.

3.4 Flow patterns

Table 1 shows the flow cases, also indicated 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 0.3 % from momentum theory values. The same holds for the local satisfaction of the boundary conditions at the wake boundary $v_n=\mathrm{0},\mathrm{\Delta }p=\mathrm{0}$, except within a distance 0.02R from the disc edge, where v_n may deviate up to 0.02U₀ without challenging the condition Ψ=Ψ₁ and without affecting integrated flow quantities.

Figure 4 The flow patterns of wind turbine discs (a) and (c) and propeller discs (b) and (d) with a contracting wake, (e) with an expanding wake. The streamlines indicate stream tube values increasing with ΔΨ=0.1Ψ₁.

In Fig. 4 the streamlines of flow case 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₀ at the upstream side of the disc equals $|\mathit{v}{|}_{\mathrm{m}}=\sqrt{v_{x,\mathrm{d}}^{\mathrm{2}}+v_{r,\mathrm{d}}^{\mathrm{2}}}/U_{\mathrm{0}}$. Table 1 gives the numerical values of ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ and $|\mathit{v}{|}_{\mathrm{m}}$ for the flow cases considered. The differences between ${\stackrel{\mathrm{‾}}{u}}_{\mathrm{d}}$ as calculated numerically and as resulting from the momentum theory are 0.2 % or less. The $|\mathit{v}{|}_{\mathrm{m}}$ value in the table 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 of this meridional velocity being practically uniform. The explanation of this is presented in Sect. 5, but first the velocity distributions shown in Fig. 5 are analysed.

Figure 5 The velocity distribution at the disc, for flow cases (a) to (e) defined in Fig. 4. Black line: $|\mathit{v}{|}_{\mathrm{m}}=\sqrt{v_x^{\mathrm{2}}+v_r^{\mathrm{2}}}/U_{\mathrm{0}}$; red line: $u=v_x/U_{\mathrm{0}}$; blue line v_r∕U₀. All vertical axes have the same scale. The percentages denoting the non-uniformity of $|\mathit{v}{|}_{\mathrm{m}}$ are explained in Sect. 4.1.

4.1 The meridional velocity

Figure 5 shows the amount of non-uniformity in $|\mathit{v}{|}_{\mathrm{m}}$. This non-uniformity is defined as $\left|\mathit{v}{|}_{\mathrm{m}}\right(\mathrm{0.97})/|\mathit{v}{|}_{\mathrm{m}}\left(\mathrm{0}\right)-\mathrm{1}$, expressed in percentages, except for flow case a. In all flow cases except a, $|\mathit{v}{|}_{\mathrm{m}}$ increases or decreases monotonically from r=0 towards r=R. In flow case a, $|\mathit{v}{|}_{\mathrm{m}}$ increases with increasing r, with the maximum, 0.2 %, reached at $r/R=\mathrm{0.8}$ after which it decreases towards the disc edge. At $r/R=\mathrm{0.97}$, $|\mathit{v}{|}_{\mathrm{m}}$ differs −0.1 % from its value at $r/R=\mathrm{0}$, so for a the non-uniformity number indicates $\left|\mathit{v}{|}_{\mathrm{m}}\right(\mathrm{0.8})/|\mathit{v}{|}_{\mathrm{m}}\left(\mathrm{0}\right)-\mathrm{1}$. These numbers for a are within the uncertainty range of the calculations, so their significance is not clear. The choice for r=0.97R in the other flow cases is somewhat arbitrary but is motivated by the argument that the sharp transition at $r/R=\mathrm{1}$ shown in Fig. 5 is not physically realistic. Viscosity will smooth this transition depending on the Reynolds number used, as shown in Sørensen et al. (1998).

4.1.1 Wind turbine flows

As shown in Fig. 5 $|\mathit{v}{|}_{\mathrm{m}}$ is practically uniform in flow case a: the non-uniformity is −0.2 %. For low λ operation the non-uniformity is stronger: −1.8 % for flow case c. The non-uniformity is checked (but not shown in a figure) for several other flow cases.

• Disc load $C_{\mathrm{\Delta }H}=-\mathrm{8}/\mathrm{9}$, λ=5 instead of ∞: the result differs less than 0.1 %.

• Discs with λ=∞ but heavier disc loads: the non-uniformity in $|\mathit{v}{|}_{\mathrm{m}}$ is −0.7 % for $C_{T,\mathrm{\Delta }H}=-\mathrm{0.97}$ and −0.8 % for $C_{T,\mathrm{\Delta }H}=-\mathrm{0.995}$.

The optimal operational regime of modern wind turbines is λ>5 with $C_{T_{\mathrm{\Delta }H}}>-\mathrm{0.9}$, so the non-uniformity in $|\mathit{v}{|}_{\mathrm{m}}$ of flow cases representing this optimal regime is negligible.

4.1.2 Propeller flows

The non-uniformity in $|\mathit{v}{|}_{\mathrm{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, Fig. 6.12). Figure 3 shows that in this regime the impact of wake swirl is very limited, so flow case b is considered representative, with a non-uniformity of ≈2 %.

4.2 The axial velocity

In all flow cases the axial velocity is far from uniform, as was already shown by Sørensen et al. (1998) and Madsen et al. (2010), for example. For Froude wind turbine discs, the cause of this has been addressed in van Kuik and Lignarolo (2016) and for Joukowsky disc flows in van Kuik (2017). In terms of the momentum balance, the source of this non-uniformity is the pressure acting on the sides of a stream annulus used as control volume. When the stream tube boundary is used as the boundary of the control volume, the pressure at this boundary does not give a contribution in the 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 been calculated for $-\mathrm{1}<C_{T,\mathrm{\Delta }H}<\mathrm{0}$, enabling a surface-fit engineering approximation for $u_{\mathrm{d}}(\frac{r}{R},C_{T,\mathrm{\Delta }H})$; see van Kuik and Lignarolo (2016, Sect. 5.2).

4.3 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 based on the radial velocity, Micallef et al. (2013) comparing calculated and measured radial velocity near rotor blade tips to assess blade bound chordwise vorticity in order to explain the initially inward motion of the tip vortex, van Kuik et al. (2014, Sect. 4) quantifying this chordwise vorticity and the associated tip load responsible for this inward tip vortex motion, and Sørensen (2015, Sect. 3.2) analysing $\partial v_r/\partial 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 plane:

$$\begin{array}{ll}\text{(18)}& {\displaystyle }& {\displaystyle }\underset{S_{\mathrm{d}}}{\int }\left({\left({\displaystyle \frac{v_{r,\mathrm{d}}}{U_{\mathrm{0}}}}\right)}^{\mathrm{2}}-a^{\mathrm{2}}\right)\mathrm{d}S=\mathrm{0},{\displaystyle }& {\displaystyle }\text{from Limacher and Wood (2019)},\end{array}$$

where a is the induction $\mathrm{1}-v_{x,\mathrm{d}}/U_{\mathrm{0}}$ and S_d is the plane of the disc or rotor from r=0 to r=∞. Based on Eq. (18), the authors conclude that $v_{r,\mathrm{d}}/U_{\mathrm{0}}$ and a have to be equal close to the disc edge or rotor tip, so

$$\begin{array}{ll}\text{(19)}& {\displaystyle }& {\displaystyle \frac{v_{x,\mathrm{d}}}{U_{\mathrm{0}}}}+{\displaystyle \frac{v_{r,\mathrm{d}}}{U_{\mathrm{0}}}}=\mathrm{1}\text{ at }r\approx R,{\displaystyle }& {\displaystyle }\text{ adapted from Limacher and Wood (2019)}.\end{array}$$

Equations (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 Eq. (19) is satisfied: almost at the disc edge for the flow cases with an expanding wake a, c and e, while flow cases b and 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 accurate assessment of the radial position where Eq. (19) is satisfied is difficult for which reason a range is indicated.

Table 2 The radial position where Eq. (19) is satisfied, for flow cases a to e.

Equation (19) provides a second relation between v_x,d and v_r,d, besides the conclusion of Sect. 4.1 that $|\mathit{v}{|}_{\mathrm{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 $|\mathit{v}{|}_{\mathrm{m}}=\sqrt{v_{x,\mathrm{d}}^{\mathrm{2}}+v_{r,\mathrm{d}}^{\mathrm{2}}}/U_{\mathrm{0}}=\text{constant}$ and $v_{x,\mathrm{d}}+v_{r,\mathrm{d}}=U_{\mathrm{0}}$ at $r/R=\mathrm{1}$. As an example the flow with $v_{x,\mathrm{d}}=v_{r,\mathrm{d}}=\mathrm{0.5}{U}_{\mathrm{0}}$ at r=R is evaluated, giving $|\mathit{v}{|}_{\mathrm{m}}=\mathrm{0.707}$. This gives a slope of the vortex sheet shape of 45^∘ at r=R. This is close to flow state a, where the numerically calculated slope is 46^∘, and $|\mathit{v}{|}_{\mathrm{m}}=\mathrm{0.684}$ which is 3.3 % lower than the estimate. Further exploration of such an engineering estimate is left for future work.

5 Explanation of the (non-)uniformity of $|\mathit{v}{|}_{\mathrm{m}}$

The Euler equation of motion (Eq. 1) offers a first-order explanation for the observation that $|\mathit{v}{|}_{\mathrm{m}}$ is practically uniform for λ≥1. The radial component of Eq. (1) reads

$$\begin{array}{}\text{(20)}& {\displaystyle }{\displaystyle \frac{\partial p}{\partial r}}=-\mathit{\rho }{v}_{\mathrm{s}}{\displaystyle \frac{\partial v_r}{\partial s}}+{\displaystyle \frac{v_{\mathit{\phi }}^{\mathrm{2}}}{r}}.\end{array}$$

Equation (3) for v_φ combined with Bernoulli's Eq. (4) gives a second equation for $\partial p/\partial r$:

$$\begin{array}{}\text{(21)}& {\displaystyle }{\displaystyle \frac{\partial p}{\partial r}}=-\mathit{\rho }{v}_{\mathrm{s}}{\displaystyle \frac{\partial v_{\mathrm{s}}}{\partial r}}+{\displaystyle \frac{v_{\mathit{\phi }}^{\mathrm{2}}}{r}},\end{array}$$

so the result is $\partial v_{\mathrm{s}}/\partial r=\partial v_r/\partial s$, or at the disc

$$\begin{array}{}\text{(22)}& {\displaystyle }{\displaystyle \frac{\partial v_{\text{s,d}}}{\partial r}}={\displaystyle \frac{\partial v_{r,\mathrm{d}}}{\partial s}}.\end{array}$$

Consequently, the distribution of v_s,d is determined by the derivative $\partial v_r/\partial s$ along the streamline. In case v_r has a maximum or minimum at the disc, $v_{\text{s,d}}/U_{\mathrm{0}}=|\mathit{v}{|}_{\mathrm{m}}$ is uniform.

Qualitative observations regarding the increase or decrease in 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 observation 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 $\partial v_r/\partial s>\mathrm{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 $\partial v_{r,\text{upstream}}/\partial s<\mathrm{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 $\partial v_{r,\text{downstream}}/\partial s\approx \mathrm{0}$. This gives for the total induction in the wake $\partial v_r/\partial s<\mathrm{0}$.

Consequently, according to (a) and (b) $\partial v_r/\partial s=\mathrm{0}$ at the disc and with Eq. (22) $\partial v_{\text{s,d}}/\partial r=\mathrm{0}$ so $|\mathit{v}{|}_{\mathrm{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 $|\mathit{v}{|}_{\mathrm{m}}$.

However, these qualitative considerations miss the effect that a vortex ring induces a non-zero $\partial v_r/\partial s$ in the plane of the ring. Figure 6 shows the calculated radial velocity induced by a vortex ring positioned at $x=\mathrm{0},R=\mathrm{1}$ along the lines $r/R=\mathrm{0.8},\mathrm{0.9},\mathrm{0.97}$. The shape of the plot resembles the induction $v_r=\frac{\mathrm{\Gamma }}{\mathrm{2}\mathit{\pi }}\frac{\cos \mathit{\alpha }}{\text{dist}}$ by a point vortex in a 2−D plane, where dist is the distance to the vortex, and α is the angle of 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 $\partial v_r/\partial x\to \mathrm{\infty }$ for $r/R\to \mathrm{1}$. Apart from the distance to the ring, the strength of the ring determines the local value of $\partial v_r/\partial x$, as its value is linear in this strength.

Figure 6 The radial velocity induced by a unit vortex ring positioned at $x=\mathrm{0},$ R=1, at the lines $r/R=\mathrm{0.8},\mathrm{0.9},\mathrm{0.97}$.

Figure 7 The distribution of the vortex sheet strength ${\mathit{\gamma }}_{\mathit{\phi }}\left(x\right)/{\mathit{\gamma }}_{\mathit{\phi },\mathrm{1}}$, for flow cases (a) to (e), defined in Fig. 4. The vertical axes have the same scale, except the axis of (e), which covers a range of γ 4 times larger.

Figure 8 Curved lines: the radial velocity along the streamline passing the disc at $r/R=\mathrm{0.97}$; straight lines: the tangent of the distribution v_s,d(r)∕U₀ at $r/R=\mathrm{0.97}$, plotted through the s=0 position at the curved line.

For a vorticity tube things are slightly different, as is easily shown by the example of a tube of constant strength with a semi-finite length. Each elementary vortex ring γdx induces a non-zero $\partial v_r/\partial 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 $\partial v_r/\partial 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 $\partial v_r/\partial r$ at x=0 is opposite to the sign of $\partial v_r/\partial r$ far upstream, as is clear from the line r=0.97R in Fig. 6.

The argument of non-zero $\partial v_r/\partial s$ at $s^{*}=s_{\mathrm{d}}$ due to the vortex sheet leading edge has to be added to the arguments (a) plus (b).

• c.

  At $s^{*}=s_{\mathrm{d}}$ the induction by the leading edge vorticity at the disc edge adds a contribution to $\partial v_r/\partial s$ depending on the local vorticity strength and the inverse of the distance to the disc edge. The sign of the contribution is opposite to the sign of $\partial v_r/\partial s$ upstream of s_d.

• d.

  According to (a) and (b), the position where $\partial v_r/\partial s=\mathrm{0}$ is at the disc. With (c) it moves upstream of the disc, for all disc flows. How far it moves upstream depends on the strength of the leading edge vorticity. For discs with an expanding wake, using Eq. (22), $\partial v_{r,\mathrm{d}}/\partial s=\partial v_{\text{s,d}}/\partial r<\mathrm{0}$, and for discs with a contracting wake $\partial v_{r,\mathrm{d}}/\partial s=\partial v_{\text{s,d}}/\partial r>\mathrm{0}$. This is in agreement with Fig. 5, showing that v_s,d diminishes towards r=R for flow cases a, c and e, while it increases for flow cases b and d.

This qualitative line of arguments (a)–(d) requires a numerical validation and quantification. The calculated wake vorticity ${\mathit{\gamma }}_{\mathit{\phi }}\left(x\right)/{\mathit{\gamma }}_{\mathit{\phi },\mathrm{1}}$ is shown in Fig. 7, with γ_φ,1 being the azimuthal vorticity in the far wake: ${\mathit{\gamma }}_{\mathit{\phi },\mathrm{1}}=v_{x,\mathrm{1}}-U_{\mathrm{0}}$. In all flow cases the distributions have a singularity at the leading edge. Flow case a has the weakest singularity and flow case e the strongest. Figure 8 shows the calculated v_r along a streamline passing the disc at $r/R=\mathrm{0.97}$ (curved lines) and the tangent at $r/R=\mathrm{0.97}$ of the distribution v_m(r) (straight line), plotted through the s=0 position at the curved line. As is clear from the graphs, these straight lines coincide with the tangents to the v_r(s) distribution, confirming Eq. (22). Furthermore, downstream of the disc v_r decreases for flow cases a, c and e and increases for b and d, thereby confirming the assumption made in (b).

The absolute value of the slope of the tangents is lowest in flow case a and highest in e. This is in agreement with the strength of the leading edge singularity of ${\mathit{\gamma }}_{\mathit{\phi }}\left(x\right)/{\mathit{\gamma }}_{\mathit{\phi },\mathrm{1}}$ and the non-uniformity of v_m. In all flow cases v_r(s) reaches a maximum or minimum just upstream of the disc: at $s/R=-\mathrm{0.00155}$ for a and −0.00252 for e, with the values for other flow cases in between these positions.

6 Conclusions

With respect to the average velocity at the actuator disc, the following applies.

• For Joukowsky discs in wind turbine and propeller mode, the average velocity has been found, from λ=0 up to λ→∞ 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 disc flows without wake expansion or contraction are possible for specific values of J, marking the transition from the contracting wake operational mode at low J to the expanding wake mode at high J.

• In the propeller as well as wind turbine 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,

• $|\mathit{v}{|}_{\mathrm{m}}$ is practically uniform for wind turbine disc flows with λ>5 (deviation on the order of a few per mille).

• $|\mathit{v}{|}_{\mathrm{m}}$ is almost uniform for wind turbine disc flows with low λ and propeller flows with J≈π (deviation on the order of a few percent).

• $|\mathit{v}{|}_{\mathrm{m}}$ is non-uniform for the propeller disc flow with wake expansion at very high J (deviation on the order of several percent).

• the differences in uniformity are caused by the different strengths of the leading edge singularity in the wake boundary vorticity strength.