2.1. Brief review of classical momentum theory

Consider a cylindrical CV encompassing a turbine rotor, $V_1$, as depicted in figure 1, with one planar face well upstream of the rotor in the undisturbed free stream, and one face just downstream of the rotor, $S_D$, where $z=0$. The rotor's axis of rotation is colinear with the central axis of the CV, whose radius $R_{CV} \rightarrow \infty$, allowing us to assume that the axial velocity on the curved lateral boundary is equal to $U_\infty$. Assuming the flow to be incompressible and the flow inside the CV to be steady in a rotor-fixed frame of reference, conservation of mass and momentum in the CV leads to

(2.2)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_D} \left( 2a(1-a) - p \right) {\rm d}S, \end{equation}

where $p$ is the normalized gauge pressure.

Figure 1.

Control volumes to be referenced in the present work. Here, $V_1$ is a cylindrical CV of large radius, with its downstream face just behind the blades; $S_D$ is the circular downstream face of the bounding control surface of $V_1$; $V_2$ is a mean streamtube passing through the rotor, extending from the free stream to the far wake where wake expansion is presumed to have ceased.

Typically, this CV definition is not used in axial momentum theory, e.g. Sørensen (Reference Sørensen2016), as there is no obvious and direct way to infer the behaviour of the pressure field on $S_D$. Instead, the flow is assumed to be circumferentially uniform, i.e. $a = a(r)$, and the mean streamtube passing through the rotor and extending from the free stream to the far wake is considered ($V_2$ in figure 1). This streamtube is assumed to be of circular cross-section for all $z$. It is then assumed that the gauge pressure and radial velocity are zero in the far wake, and that the Bernoulli equation holds on mean streamlines between the near and far wakes. Annular streamtubes – subsets of $V_2$ – are all assumed to expand uniformly, allowing us to infer that the axial induction in the far wake is twice that at the rotor; that is, the axial velocity in the far wake (at axial location $z = z_{w}$) is taken to be $u_z(r_{w},\theta,z_{w}) = 1-2a(r)$, where the locations $(r_w,\theta,z_w)$ and $(r,\theta,0)$ are connected by a mean streamsurface. These various assumptions allow the removal of pressure from (2.2) and lead to the following classical result:

(2.3)\begin{align} C_T &= \frac{1}{\rm \pi} \int_{S_{rot}} 4 a (1-a) \,{\rm d}S \nonumber\\ &= \frac{1}{\rm \pi} \int_0^{2{\rm \pi}} \int_{0}^{1} 4a(1-a) r \,{\rm d}r \,{\rm d}\theta, \end{align}

where $S_{rot}$ is the swept plane of the rotor, located at $z = 0^{-}$ in the chosen coordinate system.

The integrand on the right-hand side of (2.3) is maximized when $a=0.5$, and when $a=0.5$ everywhere on $S_{rot}$ we obtain $C_T = 1$. As a result, momentum theory cannot account for well-established observations of $C_T > 1$, and this high-thrust regime must be modelled empirically in BEM codes, e.g. Buhl (Reference Buhl2005). Although there are second-order modifications to momentum theory to reduce its limitations (for example, one could account for azimuthal velocities in the wake to arrive at what Bontempo & Manna (Reference Bontempo and Manna2017) call ‘generalized momentum theory’, replacing the far-wake zero gauge pressure assumption), the fundamental theoretical issue of high thrust remains unresolved. We are thus motivated to venture an entirely different approach.

2.2. Exact thrust expressions from vortical impulse theory

By applying vortical impulse theory to the CV $V_1$ in figure 1, Limacher & Wood (Reference Limacher and Wood2021) derived the following expression for the thrust on a steadily rotating turbine in an unbounded, incompressible and steady free-stream flow (neglecting viscous terms on $S_D$):

(2.4)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_D} \{u_\theta^{2} + u_r^{2} - a^{2} - r (1-a) \omega_\theta + r u_\theta \omega_z \} {\rm d}S, \end{equation}

where the two components of vorticity are defined as

(2.5)$$\begin{gather} \omega_z = \frac{1}{r} \frac{\partial}{\partial r}\left(r u_\theta\right) - \frac{1}{r} \frac{\partial u_r}{\partial \theta}, \end{gather}$$

(2.6)$$\begin{gather}\omega_\theta = \frac{\partial u_r}{\partial z} - \frac{\partial u_z}{\partial r}. \end{gather}$$

The only additional assumption inherent to the derivation of (2.4) is that the flow field within $V_1$ appears steady in a frame of reference rotating with the blades. It is worth emphasizing that the underlying impulse theory is an exact expression of the conservation of linear momentum, and thus the stated assumptions are rather minor restrictions on the generality of (2.4). For example, no assumption of circumferential uniformity has been made, nor has any circumferential averaging of parameters been implied, in direct contrast to conventional momentum-based analyses. Rather, all parameters in (2.4) through (2.6) are free to vary on $S_D$ as functions of $r$ and $\theta$, i.e. on $S_D$ we have $u_\theta = u_\theta (r,\theta,0)$, $u_r = u_r(r, \theta,0)$, etc.

Although the generality of (2.4) is a virtue of the present approach, we are interested in making simplifications to reveal the conditions under which more well-known expressions from wind turbine theory can be recovered. A key outcome of Limacher & Wood (Reference Limacher and Wood2021) was to present a sufficient set of conditions under which (2.4) collapses to the thrust expression inherent to blade-element theory, wherein the total blade force is taken to be the spanwise integral of the sectional Kutta–Joukowsky force. One of those conditions is that the trailing vorticity be aligned with streamlines on $S_D$ in the frame of reference rotating with the blades, which may be expected when the flow over the blades remains attached at high Reynolds numbers. This will be the first physical proposition to be considered.

This proposition also happens to solve a practical measurement problem. Using stereo-PIV, all three velocity components on $S_D$ will be available, but axial velocity gradients will not, and thus the $\omega _\theta$-term in (2.4) cannot be resolved. However, when vortex lines and streamlines (in the rotating frame) are aligned on $S_D$, we have

(2.7)\begin{equation} \frac{u_\theta - \lambda r}{1-a} = \frac{\omega_\theta}{\omega_z}, \end{equation}

where $\lambda = \varOmega R/U_\infty$ is the tip speed ratio, and $\varOmega$ is the dimensional rotation rate of the rotor. This relation allows $\omega _\theta$ to be removed from (2.4), yielding

(2.8)\begin{equation} C_{T} = \frac{1}{\rm \pi} \int_{S_D} \{ u_r^{2} + u_\theta^{2} - a^{2} + \lambda r^{2} \omega_z \} {\rm d}S, \end{equation}

which can be evaluated using stereo-PIV data on $S_D$. The validity of this thrust expression, specifically in the high-thrust regime of turbine operation, will be considered in the experimental portion of this work.

As a second proposition, we may simply neglect the $\partial u_r/\partial z$-contribution to the $\omega _\theta$-term in (2.4). Setting $\partial u_r/\partial z = 0$ on $S_D$ yields

(2.9)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_D} \left\{u_\theta^{2} + u_r^{2} - a^{2} - r (1-a) \frac{\partial a}{\partial r} + r u_\theta \omega_z \right\} {\rm d}S, \end{equation}

where $\partial a/\partial r = - \partial u_z/\partial r$ has been substituted. This assumption is plausible at first glance, because one might expect the expansion rate of the mean streamtube passing through the rotor to be maximal at the rotor plane. That is, as Sørensen (Reference Sørensen2016) explains, the rate of streamtube expansion increases as one approaches the rotor from upstream, and decreases as one moves downstream from the rotor towards the non-expanding far wake. However, even for the idealized case of an actuator disk, Sørensen (Reference Sørensen2016) shows that $\partial u_r/\partial z$ takes on significant values in the rotor tip region. Whether this has a noteworthy effect on the associated force contribution is a matter for experimental validation.

Since classical thrust expressions do not involve velocity derivatives, we are motivated to manipulate the remaining derivative in (2.9) to remove explicit dependence on $\partial a/\partial r$, which is achieved by integration by parts. We express the area integral over $S_D$ as a double integral over $r$ and $\theta$, substituting ${\rm d}S = r \,{\rm d}r \,{\rm d}\theta$. For the unbounded case, we take the limits on the $r$-integral to be $0$ to $\infty$ and the $\partial a/\partial r$-term becomes

(2.10)\begin{align} -\int_0^{2{\rm \pi}} \int_0^{\infty} r^{2} (1-a) \frac{\partial a}{\partial r} {\rm d}r\, {\rm d}\theta &= \int_0^{2{\rm \pi}} \left\{ r^{2} a\left. \left(\frac{a}{2} - 1\right) \right|_0^{\infty}\right.\nonumber\\ &\quad \left. + \int_0^{\infty} \left(2a(1-a) + a^{2} \right) r \,{\rm d}r \right\} {\rm d}\theta. \end{align}

The first term inside the $\theta$-integral on the right-hand side vanishes, since $a$ decays rapidly with radius from the edge of the rotor in an unbounded flow (faster than $a \propto r^{-2}$). Substituting the remaining terms into (2.9), we obtain

(2.11)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_D} \{u_\theta^{2} + u_r^{2} + 2a(1-a) + r u_\theta \omega_z \} {\rm d}S. \end{equation}

This result is of clear theoretical interest, as the appearance of the $2a(1-a)$ term prompts us to ask: Under what conditions will (2.11) reduce to conventional results of momentum theory?

To start, we may note that the $\omega _z$-term vanishes in the special case of circumferential uniformity (i.e. $\partial /\partial \theta \equiv 0$), which is typically assumed in momentum theory. This is shown as follows:

(2.12)\begin{equation} \int_{S_D} r u_\theta \omega_z \,{\rm d}S = 2{\rm \pi} \int_0^{\infty} r u_\theta \frac{\partial }{\partial r}(r u_\theta)\,{\rm d}r = \left. {\rm \pi}r^{2} u_\theta^{2} \right|_0^{\infty} = 0. \end{equation}

One obtains the third expression from the second by integration by parts, and the whole expression vanishes because the circulation around any contour outside the wake must be zero. When the flow is not circumferentially uniform, this term cannot be so easily dismissed. In the experiment to be presented, however, it will be evaluated and will be shown to be negligible for the cases considered.

Although the validity of neglecting $\partial u_r/\partial z$ on $S_D$ is yet to be tested, it is less dubious further downstream. Note that (2.9) and (2.11) are both instantiations of our second physical proposition that $\partial u_r/\partial z = 0$ on the downstream CS, but nothing in their derivation requires the CS to be directly behind the rotor. If we move it to the far wake and call this new surface $S_W$, we can invoke the typical assumption in momentum theory that $u_r \equiv 0$ on $S_W$, from which $\partial u_r/\partial z \equiv 0$ follows. We then obtain

(2.13)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_W} \{u_\theta^{2} + 2u_z(1-u_z) \} {\rm d}S, \end{equation}

which is identical to what one finds through a conventional momentum balance. The second term is the momentum flux from the large cylindrical CV with $S_W$ as its downstream CS, and the first term is equivalent to the pressure contribution on $S_W$ when wake swirl is considered (McCutchen Reference McCutchen1985; Wood Reference Wood2007).

Our success in recovering a result of a classical momentum balance, expressed in terms of far-wake velocities, corroborates the approach taken thus far. However, we must return our attention to the near wake, with the goal of expressing the thrust in terms of velocities just downstream of the rotor. If we are to recover the familiar $4a(1-a)$ integrand in (2.3) under some circumstances, we might presume that the $\partial u_r/\partial z$-term contribution is, in fact, significant, and that it is of similar magnitude to (and of the same sign as) the $\partial a/\partial r$-term in (2.9), from which the $2a(1-a)$ term in (2.11) was extracted. This leads us to our third proposition, to be presented in the next subsection.

2.3. A third physical proposition, its thrust expression and a comparison with classical momentum theory

The thrust expressions thus far derived are exact, meaning that the expressions follow directly from the conservation of linear momentum when the corresponding physical propositions strictly hold (alignment of streamlines and vortex lines, or $\partial u_r/\partial z \rightarrow 0$ on $S_D$). The following derivation, based on a third physical proposition, does not attain to quite the same rigour as the first two, and some speculation must be introduced. This third proposition nonetheless warrants consideration, because the resulting thrust expression offers a possible explanation for the high-thrust regime, and it also results in an experimentally based thrust estimate which more closely matches the measured thrust than any other expression thus far presented.

Outside of the wake, we expect the flow to be irrotational, meaning $\partial u_r/\partial z = -\partial a/\partial r$. We now hypothesize that within the wake the opposite is true, i.e. $\partial u_r/\partial z = \partial a/\partial r$. This corresponds to local strain rates of zero, suggesting that the trailing vortex cores are purely ‘rolling’ structures that exhibit no shear. This zero-strain-rate hypothesis is our third physical proposition about the near wake.

Since we are making opposite assertions over different regions of the flow, we must express the transition between them. Let us define $k(r,\theta )$ as a smooth function such that

(2.14)\begin{equation} \omega_\theta = \frac{\partial u_r}{\partial z} + \frac{\partial a}{\partial r} = k(r,\theta) \frac{\partial a}{\partial r}, \end{equation}

everywhere on $S_D$, and $k \in [0,2]$. Where $k=0$, the flow is irrotational; where $k=2$, the flow is strain-rate free.

The $\omega _\theta$-contribution in the general thrust equation (2.4) then becomes

(2.15)\begin{align} C_{T,\omega_\theta} &={-}\frac{1}{\rm \pi} \int_0^{2{\rm \pi}} \int_0^{\infty} r^{2}(1-a)\omega_\theta \,{\rm d}r \,{\rm d}\theta, \nonumber\\ &={-}\frac{1}{\rm \pi} \int_0^{2{\rm \pi}} \int_0^{\infty} r^{2}(1-a)k \frac{\partial a}{\partial r} {\rm d}r \,{\rm d}\theta. \end{align}

Integrating by parts, the inner $r$-integral becomes

(2.16)\begin{align} -\int_0^{\infty} r^{2}(1-a)\omega_\theta \,{\rm d}r &= \left. -r^{2}a(1-a) \right|_0^{\infty} + \int_0^{\infty} 2 r a(1-a) \,{\rm d}r - \int_0^{\infty} r^{2} k a \frac{\partial a}{\partial r} {\rm d}r \nonumber\\ &\quad + \int_0^{\infty} r^{2} a(1-a)\frac{{\rm d} k}{{\rm d} r} {\rm d}r. \end{align}

The second-last integral can be integrated by parts again, and after algebraic manipulation we obtain

(2.17)\begin{align} -\int_0^{\infty} r^{2}(1-a)\omega_\theta \,{\rm d}r &= \left.-r^{2} k a \left( 1 - \frac{a}{2} \right) \right |_0^{\infty} + \int_0^{\infty} ( 2a (1-a) + a^{2} ) r k \,{\rm d}r \nonumber\\ &\quad + \int_0^{\infty} r^{2} a \left(1 - \frac{a}{2} \right) \frac{{\rm d} k}{{\rm d} r} {\rm d}r. \end{align}

The first term on the right-hand side vanishes because the flow is irrotational at infinity; that is, $k=0$ as $r \rightarrow \infty$.

Evaluation of the second and third integrals, however, would require precise knowledge of how the flow transitions from the irrotational to strain-rate-free conditions. In absence of such information, we lean on theoretical expectations to proceed. Our express intention is to recover something which approaches the classical thrust expression in (2.3) for low thrust, but which is capable of yielding thrust estimates in excess of one for highly loaded rotors.

With that aim in mind, let us make the rather simplistic assumption that $k=2$ everywhere within the wake, and $k=0$ everywhere outside it. Defining the radius at the outer edge of the tip vortices to be $r = R_e$ – with the expectation that $R_e \gtrapprox 1$ – we set $k = 2$ when $r < R_e$, and $k = 0$ when $r > R_e$, and assume that $k$ varies smoothly between $R_e^{-}$ and $R_e^{+}$. The second integral in (2.17) then reduces to

(2.18) \begin{equation} \int_0^{\infty} (2a (1-a) + a^{2}) r k \,{\rm d}r = \int_0^{R_e} (4a (1-a) + 2a^{2}) r \,{\rm d}r, \end{equation}

recovering the familiar $4a(1-a)$ term from classical analyses.

The limits on the last integral in (2.17) can changed to be from $r = R_e^{-}$ to $r = R_e^{+}$, because $dk/dr = 0$ has been assumed everywhere else. Taking $a$ to be constant on $R_e^{-} \leqslant r \leqslant R_e^{+}$ (i.e. —$\partial a/\partial r| \ll |{\rm d}k/{\rm d}r|$ at $r = R_e$), $a (1- a/2)$ can be moved out of the integral, which then reduces to

(2.19)\begin{equation} \int_{R_e^{-}}^{R_e^{+}} r^{2} a \left(1 - \frac{a}{2} \right) \frac{{\rm d} k}{{\rm d} r} {\rm d}r = R_e^{2} a_e \left( 1 - \frac{a_e}{2} \right)\int_{R_e^{-}}^{R_e^{+}} \frac{{\rm d} k}{{\rm d} r} {\rm d}r ={-} 2 R_e^{2} a_e \left( 1 - \frac{a_e}{2} \right), \end{equation}

where $a_e = a(R_e, \theta )$. Combining the preceding five equations, we obtain

(2.20) \begin{equation} C_{T,\omega_\theta} = \frac{1}{\rm \pi}\int_0^{2{\rm \pi}} \left\{ \int_0^{R_e} ( 4a(1-a) + 2a^{2} ) k r \,{\rm d}r - 2 R_e^{2} a_e \left( 1 - \frac{a_e}{2} \right) \right\} {\rm d}\theta . \end{equation}

It is usually assumed that $a=0$ for $r>1$ in classical treatments, and, to proceed, we may assume that $a_e = 0$, causing the second term on the right-hand side of (2.20) to vanish. Realistically, the obstruction of the rotor causes some acceleration of fluid beyond the free-stream velocity around the periphery of the swept disk, but this gives a zero crossing of $a$ somewhere near the wake's edge, once again making the assumption of $a_e = 0$ plausible.

Given the crudeness of the assumption that $k$ abruptly changes from $k=2$ to $k=0$, there is little to be gained in arguing about the precise location of $R_e$ so as to evaluate $a_e = a(R_e,\theta )$. Rather, we will proceed by neglecting the second term in (2.20), satisfied that the form of the resulting expression conforms to our theoretical expectations.

When this is done, and the remaining terms on the right-hand side of (2.20) are substituted for the $\omega _\theta$-term in (2.4), we obtain

(2.21) \begin{align} C_T &= \frac{1}{\rm \pi} \int_0^{2{\rm \pi}} \int_0^{\infty} \{u_\theta^{2} + u_r^{2} - a^{2} + r u_\theta \omega_z \} r \,{\rm d}r \,{\rm d}\theta \nonumber\\ &\quad + \int_0^{2{\rm \pi}} \int_0^{R_e} \{ 4a(1-a) + 2a^{2}\} r \,{\rm d}r \,{\rm d}\theta. \end{align}

Having recovered the familiar $4a(1-a)$ integrand, the squared velocity terms in (2.21) appear as second-order corrections to the classical results of momentum theory, allowing for a thrust coefficient that grows monotonically with increasing $a$ and allowing $C_T >1$ for $a > 0.5$. This is perhaps made clearer if we approximate $R_e \approx 1$, and define $S_d \subset S_D$ as the projection of the swept rotor plane on $S_D$, defined by $\theta \in [0, 2{\rm \pi} ]$ and $r \in [0,1]$ at $z = 0$. Under the assumption that $a \approx 0$ for $r>1$, we then obtain

(2.22)\begin{equation} C_T = \frac{1}{\rm \pi} \int_{S_D} \{u_\theta^{2} + u_r^{2} + a^{2} + r u_\theta \omega_z \} {\rm d}S + \frac{1}{\rm \pi} \int_{S_d} 4a(1-a) \,{\rm d}S. \end{equation}

The classical thrust expression in (2.3) is recovered if the first integral on the right-hand side of (2.22) is neglected. (The integral over $S_d$ in (2.22), located at $z = 0$, is expected to closely approximate the integral over $S_{rot}$ in (2.3), located at $z = 0^{-}$, as no discontinuity in $a$ is expected in the neighbourhood of $z=0$.) If valid, (2.22) suggests that momentum theory, and the BEM codes based on it, will underpredict thrust at typical operating values of $a$. In fact, underprediction of experimentally measured thrust has been observed for BEM codes across a wide range of tip speed ratios in Krogstad & Eriksen (Reference Krogstad and Eriksen2013).

The novel corrective terms in (2.22) may be interpreted as accounting for wake expansion in a manner not captured by classical momentum theory. The contribution of $u_r^{2}$ is obviously associated with expansion, but so is the $a^{2}$-contribution, since the latter was derived from a term involving $\partial u_r/\partial z$.

The various thrust equations derived here in § 2 stand as competing hypotheses. We are unaware of any experiments that have specifically tested their inherent physical assumptions, and neither will we attempt such direct validation in the experimental portion of this work. Rather, we will employ the thrust expressions resulting from each physical proposition, using stereo-PIV data to furnish thrust estimates therefrom, and we will compare those estimates with direct thrust measurements. In this way, we will attempt to corroborate, if not fully validate, at least one of the three physical propositions.

4.1. Experimental apparatus

A stereo-PIV investigation on a stream-normal plane just behind a turbine rotor was undertaken in a water channel facility at Princeton University. Figure 2 shows top-, front- and side-view sketches and a photograph of the experimental set-up. The illuminated plane was $0.09D$ behind the trailing edge of the turbine blades, where $D = 127$ mm is the turbine diameter. This was judged to be as close as practically possible without encountering severe reflections from the blades in the acquired images.

Figure 2.

Experimental set-up; (a) front-view sketch, showing the FoV for stereo-PIV data capture; (b) photograph of the set-up, showing the laser sheet intersecting the transparent hub just behind the blades; (c) top-view sketch, showing arrangement of mirrors and cameras, and the $60^{\circ }$ angle between the cameras’ optical axes; (d) side view with cutaways, showing internal shafts driving the rotor, and showing the mounting of the gearmotor and force transducer above the waterline. Sketches are not to scale.

The sampled field of view (FoV) is below the axis of the rotor where both cameras have optical access, as shown in figure 2. Individual sectors of the velocity field are measured sequentially, and a complete $360^{\circ }$ velocity field is reconstructed using the stitching method described in § 4.5.

A Litron NANO L50-50PIV dual-pulse laser was used for illumination, with the beam manipulated into a $2$ mm-thick sheet using standard sheet optics. Two LaVision sCMOS CLHS cameras were mounted about $2.76D$ downstream of the illuminated plane, oriented with their viewing axes normal to the acrylic channel walls to prevent optical distortion. The viewing axes were redirected upstream using $76$ mm square mirrors mounted inside the channel, achieving a $60^{\circ }$ angle between the two views. Image capture and laser firing were synchronized using a LaVision PTU X programmable timing unit. Particle image-pair acquisition was triggered once per rotation based on the motor encoder signal using a programmable OSOYOO Due R3 microcontroller. The flow was seeded with silver-coated hollow glass spheres of mean diameter $14\ \mathrm {\mu }{\rm m}$ (Potters Industries, S5000-S3).

The rotor was driven by a gearmotor with a 3.7 : 1 gear ratio (Phidgets 3265_0) mounted above the waterline using shafts internal to the supporting apparatus. This apparatus consisted of a cylindrical fuselage ($0.23D$ in diameter) with an elliptical aft body, and a vertical sting with an airfoil cross-section (NACA0024 near the fuselage, transitioning to a NACA0021 profile above) whose chord length was $60$ mm, or $0.47D$. Both the fuselage and the sting, as well as auxiliary mounting components above the waterline, were fabricated using the fusion deposition modelling (FDM) method of three-dimensional (3-D) printing.

The investigated turbine rotor is modelled after that of McTavish, Feszty & Nitzsche (Reference McTavish, Feszty and Nitzsche2014), which they designed specifically for the low-Reynolds-number environment of a water channel. In the present study, the rotor hub diameter has been increased from their value of $0.1D$ to $0.23D$ to accommodate the internal components needed to drive the rotor. The chord and twist distributions of the blades are shown in figure 3(a).

Figure 3.

(a) Chord and twist distributions – $c/R$ and $\theta$, respectively – of the rotor blades; (b) standard Eppler E387 airfoil cross-section (black markers) and the modified E387 airfoil with a trailing edge of finite thickness (red markers), truncated for ease of manufacture by 3-D printing.

The rotor was fabricated using the stereolithography apparatus method of 3-D printing to achieve a smoother surface finish than is possible using FDM. The three blades and the hub were printed separately and then friction-fit together. A polymer known by the tradename WaterShed 11122XC was chosen for its translucence to minimize reflections at the surface that may saturate particle images and prevent the acquisition of adequate particle-image correlations near the hub and blades. A modified Eppler E387 airfoil cross-section was used for the blades, where the trailing edge has been truncated to have sufficient thickness for 3-D printing. The modified airfoil cross-section is shown in figure 3(b) alongside the standard Eppler E387 airfoil, as obtained from the online airfoil database maintained by the Applied Aerodynamics Group at University of Illinois at Urbana-Champaign (UIUC). (https://m-selig.ae.illinois.edu/ads/coord_database.html.)

The entire rotor, fuselage and sting assembly was mounted to a six-axis force transducer (ATI Mini40) above the waterline, which itself was rigidly mounted to the water channel frame.

4.2. Motor control

In the present work, tip speed ratios of $\lambda = \{4.37, 5.25, 6.12, 7.00, 7.87, 8.75, 9.62, 10.50\}$ were investigated. The free-steam velocity was held fixed at $U_\infty = 0.19\ {\rm m}\ {\rm s}^{-1}$ for all tip speed ratios and the rotor speed was set to $\varOmega = \{125, 150, 175, 200, 225, 250, 275, 300\}$ r.p.m.

The gearmotor speed was controlled using closed-loop feedback of an encoder signal (Anaheim Automation ENC-A4TS-0360-197-H-M). Extant motor speed variations can be characterized as the superposition of an oscillation at the rotor frequency, possibly caused by uneven friction in the drivetrain, and a non-stationary slow drift. This superposition can be expressed as

(4.1)\begin{equation} \varOmega = \varOmega_{target} + \varOmega_{osc} \sin (2{\rm \pi} f_{rot} t) \pm \varOmega_d, \end{equation}

where $\varOmega _{target}$ is the target speed, $\varOmega _{osc}$ is the oscillating component at the frequency of rotation, and $\varOmega _d$ is treated as a random error due to the slow drift.

Examples of the r.p.m. signals acquired by numerical differentiation of the encoder signal (using the method of central differences) are shown as red lines in figure 4(a) and 4(b). Applying a modified-Gaussian low-pass filter, as defined in Limacher et al. (Reference Limacher, Morton and Wood2019), with passband and stopband frequencies of $0.2f_{rot}$ and $0.3f_{rot}$, we obtain estimates of the slow drift, shown as blue lines on the same figure.

Figure 4.

Examples of the motor-speed signal, acquired by numerical differentiation of the raw encoder signal. The red lines are the raw signals, the dashed black lines are the target r.p.m., and the blue lines are the low-pass filtered signal highlighting the non-stationary long-term drift: (a) 150 r.p.m. target; (b) 300 r.p.m. target.

Calculating the standard deviation of the r.p.m. signal over ten rotor cycles, $\sigma _{10cyc}$, we estimate the oscillatory amplitude as $\varOmega _{osc} = \sqrt {2} \sigma _{10cyc}$, and the instantaneous uncertainty associated with the long-term drift is characterized as twice the standard deviation of the low-pass filtered signal over a 20 s sample. Values of $\varOmega _{osc}$ and $\varOmega _d$ are reported in table 1.

Table 1.

Characterization of the motor-speed variability. Here, $\varOmega _{osc}$ represents the amplitude of the oscillatory component at the frequency of rotation; $\varOmega _d$ represents the instantaneous uncertainty associated with the non-stationary slow drift of the mean motor speed.

4.3. Direct force measurements

The measured rotor thrusts were small, amounting to approximately $0.2N$ for a thrust coefficient of $C_T = 1$. This is at the low end of the $0N$ to $20N$ measurement range of the ATI Mini40 force transducer, and it was deemed necessary to perform an in-house calibration over the range from $0N$ to $0.4N$. In addition, to obtain repeatable results, it was found necessary to apply a constant pre-load in the downstream direction of a few Newtons using an elastic band. Unrepeatable results in absence of this preload were attributed to static friction internal to the transducer. Over this range of forces, linear regression revealed a linearity error of about $0.6\,\%$, and a zero-offset error of $0.0008N$. This value is less than the measurement resolution of $0.005N$, indicating that the zero offset is effectively zero. After normalization, the thrust coefficient resolution is $0.022$.

Another challenge to the use of the Mini40 transducer is the non-stationary drift in the signal over long times. Leaving the unloaded transducer over 30 minutes, this drift could be as high as $0.05N$, which is significant for our purposes. Only relative changes in thrust over short times were deemed trustworthy. To address this challenge, we performed two types of thrust test. In the first, which we call a ‘ramp-up test’, we started recording thrust with the water channel and the rotor off; after 45 s, we ramped up the water channel speed to $0.19\ {\rm m}\ {\rm s}^{-1}$ and turned the rotor on at $\varOmega = 250$ r.p.m.; taking the mean thrust over at least 25 s, and repeating this test four times, we established a repeatable reference thrust, $T_{250}$, for subsequent trials. To determine the thrust at all other rotor speeds, we began trials with the channel on and the rotor at $250$ r.p.m., then changed the rotor set point to the new desired value and measured the change in mean thrust, $\Delta T_\varOmega$, allowing the thrust at each speed to be calculated as $T_\varOmega = T_{250} + \Delta T_\varOmega$. Due to this summation, the zero-order uncertainty due to sensor resolution is increased by a factor of $\sqrt {2}$ (Moffat Reference Moffat1988), resulting in an uncertainty contribution of $\delta C_{T,a} = 0.031$ in the thrust coefficients. In the reported results, this uncertainty will be added to two standard deviations of individual thrust measurements to yield a total uncertainty estimate.

Ten repeated trials for each rotor speed were conducted, five of which were conducted after the rotor apparatus was disassembled and reassembled. The thrust on the fuselage and sting were also characterized using five ramp-up tests with a bladeless rotor (a hub only), and this thrust has been subtracted from all measured thrust values to give the measured rotor thrust coefficients reported in § 5.

All force data were acquired at a sampling rate of $200$ Hz.

4.4. The PIV methodology

For each sector, 200 image pairs were captured with both cameras, with an interframe time of $2$ ms. Mean velocity fields were calculated from particle-image pairs in LaVision's DaVis 10 software using the sum-of-correlations method (Meinhart, Wereley & Santiago Reference Meinhart, Wereley and Santiago2000), wherein all 200 cross-correlation maps are averaged before the correlation peak is identified and used to calculate velocity. This method reduces random error in the calculated mean velocities by effectively weighting sharper correlation peaks more strongly. By contrast, if each correlation map were used to calculate a particle displacement and these displacements were averaged, a single erroneous vector from a diffuse correlation peak could greatly skew the mean value.

By using 200 correlation maps ($N = 200$) to compute the mean velocity field, the uncertainty in local mean velocity components is made acceptably low. If instantaneous deviations from the true mean (whether physical or spurious) are normally distributed with a standard deviation of $\sigma$, then repeated experiments would yield calculated means with a standard deviation about the true mean of $\sigma /\sqrt {N}$ (Moffat Reference Moffat1988) – in this case, $0.07 \sigma$ – and we may presume the error in our measured mean velocities to be of the order of $1\,\%$. When those velocity data are used as inputs to the thrust expressions given in § 3, the integration operations reduce uncertainty in the thrust estimates even further and this source of error can be neglected.

PIV calibration was performed using a 3-D calibration plate (LaVision 106-10), with $2.2$ mm-diameter marks separated by $10$ mm on each plane, and the two planes separated by $2$ mm. The calibration was refined using the planar self-calibration function in LaVision DaVis 10 (Wieneke Reference Wieneke2005), applied on over $900$ particle-image pairs for the free-stream flow with the experimental apparatus removed. The mean free-stream velocity of $U_\infty = 0.19\ {\rm m}\ {\rm s}^{-1}$ was also obtained from this dataset.

To correct for variations in laser sheet brightness over the field of view, each image underwent local brightness normalization over a $300$ pixel region. The brightness of the second image of each image pair was normalized to the brightness of the first to correct for any difference in energy between the two laser pulses. A min/max filter over a length of four pixels was also applied.

Multi-pass correlations were employed, performing one pass with a window size of $64 \times 64$ pixels and $50\,\%$ overlap, and a second pass with $32 \times 32$ pixels and $75\,\%$ overlap. The resulting vectors are spaced in the $x$- and $y$-directions by $0.00655R$.

In vector post-processing, vectors with a correlation value of less than $0.5$ were removed. Universal outlier detection (Westerweel & Scarano Reference Westerweel and Scarano2005) was also employed on $5 \times 5$ vector regions to remove spurious vectors. Vectors rejected in these ways were replaced by interpolation. No additional smoothing was applied, as Limacher et al. (Reference Limacher, McClure, Yarusevych and Morton2020) found that even modest spatial filtering can cause impulse-based force expressions to underestimate drag.

4.5. Velocity field processing

As shown in figure 2, only the lower portion of the rotor is visible to both cameras, allowing a $36^{\circ }$-sector with an outer radius of $1.2R$ to be captured at any given moment. To reconstruct the complete $360^{\circ }$ flow field, 15 of these sectors were captured, each separated by $24^{\circ }$ from its neighbours, yielding $12^{\circ }$ of overlap on either side of each sector. In the reconstruction, the velocity fields in the overlapping regions were taken as a linearly varying weighted average of the two adjacent sectors.

For ease of plotting and thrust calculations, the PIV data, acquired on a Cartesian grid, are linearly interpolated onto a circular grid, with equal radial and angular spacing of $\Delta r = 0.0039R$ and $\Delta \theta = 0.0031$ rad between data points. Using these stitched and interpolated velocity fields, the equations presented in § 3 will be used to estimate thrust. In all of those equations, the integrand takes the form $r f(r,\theta )$, and the integrals will be approximated using a midpoint integration method

(4.2)\begin{equation} \int_0^{2{\rm \pi}} \int_{r_1}^{r_2} f(r,\theta) r \,{\rm d}r \,{\rm d}\theta \approx \Delta \theta \Delta r \sum_i \sum_j f_{ij} r_{ij}. \end{equation}

Although the turbine extracts energy from the flow for most of the investigated tip speed ratios, its rotation is driven by a motor. This is made possible by the small magnitude of power extracted, which is readily dissipated by the friction of the drivetrain, or by the electrical resistance internal to the motor. To confirm that the system extracts positive energy over at least some of the investigated tip speed ratio range, we will calculate the power coefficient, $C_P$, using the near-wake velocity data. Applying the conservation of angular momentum to the volume $V_1$, and calculating power as torque multiplied by rotational speed, we obtain

(4.3)\begin{equation} C_P ={-}\frac{2\lambda}{\rm \pi} \int_0^{2{\rm \pi}} \int_{R_i}^{R_{CV}} (1-a)u_\theta r^{2} \,{\rm d}r \,{\rm d}\theta, \end{equation}

which is numerically approximated using discrete data according to (4.2). The calculated $C_P$ curve serves as a check on the reasonability of the velocity data, and we can verify that the measured thrust curve is consistent with the power regime in which the rotor operates.

Where the laser sheet intersects the translucent rotor hub, the background brightness in the images is significant. To avoid biased correlations that will introduce error into the calculations of force and power, all data within $r \leqslant R_i = 0.28 R$ will be omitted.