2. Theoretical background

The present study considers the superposition of a steady horizontal flow of velocity $\boldsymbol {u}_0=u_0 \boldsymbol {e}_x$ with an unsteady flow $\boldsymbol {u}'=u_x'(t,x,z)\boldsymbol {e}_x+u_z'(t,x,z)\boldsymbol {e}_z$ imposed by the orbital velocity of gravity surface waves. The wave dispersion relation of the gravity surface waves in the deep-water approximation, accounting for a Doppler shift, reads

(2.1) \begin{align} (\omega _r-u_0k)^2=gk, \end{align}

with $\omega _r=2\pi f$ the angular velocity and $f$ the wave frequency perceived by the hydrofoil, $k=2\pi /\lambda$ the wavenumber and $\lambda$ the wavelength. Additionally, $g=9.81$ m s⁻² is the gravitational acceleration. The wave amplitude, $A$ , is quantified by the wave steepness $\epsilon =2A/\lambda$ . The orbital motion of radius $\mathcal{R}(z)=Ae^{{2\pi z}/{\lambda }}$ and its tangential velocity $u'(z)=\mathcal{R}(z)\omega$ , with $\omega =\sqrt {gk}$ given by (2.1) for $u_0=0$ m s⁻¹, decrease exponentially with depth $z$ (here, negative).

The total velocity field $\boldsymbol {u}_0+\boldsymbol {u}'$ in the flow is estimated using linear potential flow theory and is represented in figure 1 for a selected test condition of $u_0=0.75$ m s⁻¹, $\lambda =4$ m and $\epsilon =0.04$ . The presence of a large constant velocity component $u_0$ dominates the wave-induced orbital motions, with an immersed object experiencing successive variations of relative velocity magnitude and direction about the mean flow $\boldsymbol {u}_0$ . The relative velocity magnitude is maximum at the wave peaks and minimum at their trough locations, implying a maximal and minimal acceleration on the rising and descending wave slopes, respectively. The direction of the local velocity vector that will induce an angle of attack between the flow and the hydrofoil is estimated as

(2.2) \begin{align} \phantom {.}\tan {(\alpha )}=\frac {u_z}{u_x}=\frac {\mathcal{R}\omega \sin (kx-\omega _r t)}{\mathcal{R}\omega \cos (kx-\omega _r t)+u_0}, \end{align}

where $u_x=u_0+u_x'$ and $u_z=u_z'$ are respectively the horizontal and vertical velocity components of the flow seen by the hydrofoil. The amplitude of the flow velocity fluctuation in the vertical direction $z$ is hence given by $\mathcal{R}\tan (\alpha )$ . Note that the flow field shown here is calculated through linear approximation, so waves with large amplitude (typically $\epsilon \geq 0.04$ ) will exhibit nonlinearities not taken into account in the present estimation.

Figure 1. Theoretical velocity field expected at a fixed abscissa $x$ , under surface waves calculated using potential flow theory with $u_0=0.75$ m s⁻¹, $\lambda =4$ m and $\epsilon =0.04$ . The locations of the extrema of velocity magnitude, acceleration and angle of attack are shown. Inset shows a schema of the hydrofoil used with the direction of the different forces indicated.

The expected horizontal and vertical forces $F_x$ and $F_z$ experienced by the immersed hydrofoil follow from the Morison equation (Morison et al. Reference Morison, O’BRIEN, Johnson and Schaaf1950; Chaplin & Subbiah Reference Chaplin and Subbiah1997) accounting for steady-state drag and lift, and inertia effects. They are constructed as

(2.3) \begin{align} \begin{split} F_x &=\frac {1}{2}C_d(\beta ) \rho csu^2\cos (\alpha )+\frac {1}{2}C_l(\beta ) \rho csu^2\sin (\alpha )+C_{mx}(\beta ) \rho V \frac {\partial u}{\partial t}+F_{\kern-1pt f},\\ F_z &=-\frac {1}{2}C_l(\beta ) \rho csu^2\cos (\alpha )+\frac {1}{2}C_d(\beta ) \rho csu^2\sin (\alpha )-C_{\textit{mz}}(\beta ) \rho V \frac {\partial u}{\partial t},\\ \end{split} \end{align}

where $\rho$ is the liquid density, $V$ the hydrofoil volume, and $\beta =\alpha _0-\alpha$ the angle of attack between the flow and the hydrofoil. Here, $\alpha _0$ is the static pitch angle of the hydrofoil. Additionally, $C_d$ , $C_l$ , $C_{mx}$ and $C_{\textit{mz}}$ are respectively the steady-state drag and lift, and horizontal and vertical inertia coefficients. Here, $C_d$ and $C_l$ are defined respectively parallel and perpendicular to the instantaneous local flow velocity vector. Henceforth, we consider the forces $F_x$ and $F_z$ , parallel and normal, respectively, to the mean steady velocity, $\boldsymbol {u}_0$ . A constant force $F_{\kern-1pt f}$ is added to the horizontal component to account for the friction of the end plates that are used to mount the hydrofoil in the water and enforce quasi-two dimensionality in the flow. The directions of the different forces acting on the hydrofoil are represented in the inset of figure 1. Note that because of experimental constraints, detailed later, the hydrofoil is mounted upside down, implying a negative lift and the use of $\beta =\alpha _0-\alpha$ as angle of attack. The drag and the lift are defined respectively to be aligned and perpendicular with the local flow vector, whereas the inertia forces are defined along the global $x$ and $z$ coordinates. The force coefficients depend on the angle of attack between the hydrofoil and the local flow velocity, which also varies with time. As indicated by (2.3), the flow velocity, acceleration and angle of attack will all affect the total force generation. Competition between these external parameters, which are not in phase with one another, drives variation in $F_x$ and $F_z$ .

3. Experimental set-up

Experiments were performed in a 142 m long towing tank of depth $H=2.3$ m at Delft University of Technology. The experimental set-up is shown in figure 2(a), where the hydrofoil shape is represented as a yellow outline. The hydrofoil is mounted with a baseline angle of attack of $\alpha _0=10^\circ$ to a translating platform that moves along the tank at a constant speed $u_0=0.75$ m s⁻¹. This constant angle is chosen to be less than the theoretical steady stall angle of the foil section of approximately $16^\circ$ , but sufficiently close to it to enable the observation of both attached flow and wave-triggered stall conditions across the tested parameter space and within each wave period. This constant speed imposes a steady constant flow on the hydrofoil. The hydrofoil (NACA4415), also represented in figure 2(b) and shown in figure 2(c), of chord length $c=0.16$ m $\ll \lambda$ (scaling not respected in figure 2 a) and span $s=0.48$ m, is immersed in water at a constant distance from the mean free-surface height of $h=0.4$ m. The immersion depth being 2.5 times larger than the chord length, free-surface effects are negligible and do not significantly affect hydrodynamic forces (Pernod et al. Reference Pernod, Sacher, Wackers, Augier and Bot2023). A larger depth may further reduce free surface interaction, but the wave-induced orbital flow, which decreases exponentially with depth, would be less pronounced. Elliptical end plates of principal axes 360 and 200 mm are used to attenuate any transverse flow component, mimicking a two-dimensional (2-D) flow in the $(x,z)$ plane. Although a minor transverse flow appears in steady cases, likely due to slight misalignment, the flow will be considered 2-D in the following. Two different models are tested. First, one with a uniform leading edge, as in figure 2(b–c), and then a model fitted with tubercles, shown in figure 2(d). The cross-sectional profile of the two hydrofoils can be seen in figures 5 and 8 (yellow outlines). The tubercle geometry used is spanwise sinusoidal, with amplitude $a=\, 0.1 c=16$ mm and wavelength $l=\, 0.25 c=40$ mm, following the geometry described by Wei et al. (Reference Wei, New and Cui2015) and used by Hrynuk & Bohl (Reference Hrynuk and Bohl2020) and Badia et al. (Reference Badia, Hickel, Scarano and Li2025). Note that the wavelength of the tubercles $l$ should not be confused with the wavelength of the incident waves $\lambda$ .

Figure 2. (a) Experimental set-up, implemented in a 142 m long towing tank of $H=2.3$ m water depth. A hydrofoil (yellow outline) is immersed in water between two end plates and attached to a translating platform that moves along the tank at a speed $u_0$ . Waves are generated by a wavemaker, implying an orbital flow $u'$ . (b) Model of the set-up in which the PIV system (cameras and laser sheet) is fully represented. The location of the force transducers is indicated. The hexapod to which the model is mounted does not impose any motion in the currently reported experiment. (c) Photograph of the immersed model, with the PIV system in operation. (d) Photograph of hydrofoil with tubercles at its leading edge. The three planes in which PIV measurements are taken are represented by green vertical lines. The baseline straight leading edge location is depicted, for reference, as a red dashed line. (e) Photograph of the towing tank with waves of wavelength $\lambda =4$ m and steepness $\epsilon =0.04$ .

At one extremity of the tank, a wave generator creates sinusoidal gravity waves of chosen wave length $\lambda \in [2,6]$ m and amplitude, i.e. of chosen steepness $\epsilon \in [0.02, 0.06]$ . At the other extremity of the tank, a parabolic beach is used to dissipate the waves, thus avoiding most reflections. An example of typical waves generated, with $\lambda =4$ m and $\epsilon =0.04$ , is shown in figure 2(e). The wave surface elevation $\eta (t)$ is measured using an acoustic sensor, which is traversed with the hydrofoil, at an acquisition rate of $1$ kHz. The waves, as observed in a reference frame co-translating with the hydrofoil, impose an orbital velocity field $u'(z,t)$ that makes the total flow field $u=u_0+u'$ experienced by the hydrofoil unsteady (see figure 1). To quantify the flow variation, a stereoscopic PIV system is suspended from the translating platform. The PIV set-up, shown in figure 2(b) and in operation in figure 2(c), includes two cameras mounted in a torpedo structure to the side of the hydrofoil and a laser sheet projected from the downstream direction in the $(x,z)$ plane. Because the end plates obscure the flow around the hydrofoil, measurements are made from 31 $^\circ$ below the horizontal, with 58 $^\circ$ subtended between camera optical axes about the $z$ -coordinate; the hydrofoil was mounted upside down for this purpose. The focal plane of the cameras was aligned with the measurement plane using Scheimpflug adapters and the measurement plane was accurately mapped to the camera sensors via a pinhole calibration, a typical approach in PIV experiments (Adrian & Westerweel Reference Adrian and Westerweel2011). Frame pairs, separated by 1 ms, are recorded by the two cameras at an acquisition rate of 45 Hz. Velocity fields are computed using DaVis 10.2.1.90519, via cross-correlation using an initial interrogation window of $96\times 96$ pixels and a final of $64\times 64$ pixels, with 50 % overlap, yielding a final vector spacing of $3.4$ mm, or approximately 50 vectors per chord length of the hydrofoil. Although both in-plane and spanwise velocity components are measured, we focus in the present paper on in-plane phenomena, while the $y$ -direction component is left for future treatment. The horizontal, $F_x$ , and vertical, $F_z$ , forces acting on the hydrofoil are obtained from a six-component force balance at an acquisition rate of $1$ kHz. A low-pass filter, with cut-off frequency 5 Hz, is used to remove noise that mainly originates from high-frequency vibrations of the towing platform. Measurements are performed on both hydrofoil models for a range of different wave conditions.

The parameters chosen here result in a Reynolds number $Re=u_0c/\nu \approx 1.2\times 10^5$ , with $\nu =10^{-6}$ m $^2$ s⁻¹ the water kinematic viscosity. The Reynolds number is a consequence of the value, $u_0$ , of the steady translational velocity, which was chosen to ensure the desired wave-induced unsteadiness relative to the steady translation, inducing stall and reattachment each wave period. To quantify the unsteady effects of waves on the hydrofoil, the Keulegan–Carpenter number, $\textit{KC}=\varDelta u_x/(\textit{fc})\in [0.15,3.97]$ , is used, with $\varDelta u_x$ the amplitude of velocity variation and $f$ the wave encounter frequency of the hydrofoil. The values of $\textit{KC}$ in this experiment indicate that inertia forces (i.e. added mass) dominate or are of similar order to steady drag forces (Keulegan & Carpenter Reference Keulegan and Carpenter1958). The rate of change of angle of attack can be expressed as a pitch rate (Hrynuk & Bohl Reference Hrynuk and Bohl2020; Badia et al. Reference Badia, Hickel, Scarano and Li2025), non-dimensionalised as $\varOmega ^*=\varDelta \alpha fc/U_0\in [0.01,0.06]$ , with $\varDelta \alpha$ the amplitude of variation of the angle of attack. Unsteady effects associated with vortex formation are expected to emerge from $\varOmega ^*\gt 0.005$ , whereas from $\varOmega ^*\geq 0.05$ , a quasi-static description is no longer valid and a dynamic stall regime dominated by unsteady aerodynamic effects is achieved (Choudhry et al. Reference Choudhry, Leknys, Arjomandi and Kelso2014). Considering the experimental parameters chosen here, both regimes should be observable. The various experimental conditions used in the present work are summarised in table 1. The maximum and minimum phase-averaged values of the imposed velocity components, $u_x$ and $u_z$ , and of the flow orientation, $\alpha$ , are estimated via PIV in an area upstream of the leading edge at one chord length below the hydrofoil. This location is considered to represent the external flow field, as it is far enough from the hydrofoil so as to barely be influenced by the hydrofoil’s presence. The case without waves (first row of table 1) shows this to not be precisely true for the horizontal velocity component $u_x=0.8$ m s⁻¹, which is slightly larger than the constant value $u_0=0.75$ m s⁻¹ expected for this case. The value of $u_x$ is thus slightly overestimated in our measurements. As the hydrofoil is upside down with a baseline pitch angle $\alpha _0=10^\circ$ , the angle of attack is given by $\beta =\alpha _0-\alpha$ . Note that the amplitude of variation of the angle of attack is here at most $25^\circ$ , which is smaller than the values of $30^\circ$ and $50^\circ$ used by Hrynuk & Bohl (Reference Hrynuk and Bohl2020) and Badia et al. (Reference Badia, Hickel, Scarano and Li2025), respectively. Note also that the case $\lambda =6$ m, $\epsilon =0.06$ is not tested because the large wave amplitude led to indications of wave-breaking phenomena in our experimental facility.

Table 1. Experimental conditions.

Figure 3. Time history of (a) the surface elevation $\eta (t)$ , (b) the horizontal force component $F_x$ and (c) the vertical force component $F_z$ for $\lambda =4$ m and various wave amplitudes: $\epsilon =0$ (black dashed lines), 0.02 (blue solid lines), 0.04 (green solid lines) and 0.06 (red solid lines). (d,e,f) Phase-averaged data over 40 wave periods of the same. Coloured dashed lines correspond to the standard deviation of the average. Vertical black dashed lines correspond to time instants shown in figure 5.

4. Effects of waves on a regular hydrofoil

We first focus on how waves interact with a hydrofoil with a straight-leading-edge geometry. Simultaneous time histories of the surface elevation $\eta (t)$ (panel a), and of the horizontal and vertical forces $F_x(t)$ and $F_z(t)$ (panels b and c) are plotted in figure 3, for $\lambda =4$ m and various wave amplitudes, $\epsilon =0$ , i.e. no waves, (black dashed lines), 0.02 (blue solid lines), 0.04 (green solid lines) and 0.06 (red solid lines). Phase-averaged data over 40 wave periods of the same quantities are plotted in figure 3(d–f) (solid lines) with their standard deviation (dashed lines), where the time $t^*$ is non-dimensionalised by the wave period $T=2\pi /\omega _r$ . Time histories of other incident wave lengths are plotted in figures S1 and S2 of the supplementary material available at https://doi.org/10.1017/jfm.2025.1112, and evidence similar results. As previously indicated, precise definition of the external flow velocity is difficult given the limited upstream field of view. For this reason, we present the measured forces in dimensional form, thereby avoiding additional uncertainty in non-dimensionalisation.

In the case $\epsilon =0$ , no wave is generated, so the angle of attack is only given by the baseline pitch angle $\beta =\alpha _0=10^\circ$ . The quantities $\eta =0$ mm, $F_x\approx 2.1$ N and $F_z\approx -23.8$ N remain constant over time as the flow field is only driven by the constant velocity $u_0$ . Note that $F_z$ is negative because the hydrofoil is mounted upside down. Here, $F_x$ , respectively $F_z$ , are coincident with the steady drag, respectively lift, i.e. by (2.3) with $\alpha =0^\circ$ , $\beta =\alpha _0=10^\circ$ and $u=u_0$ . The coefficients $C_d=0.04$ and $C_l=1.17$ are obtained using Xfoil (Drela Reference Drela1989) with a steady angle of attack $\alpha _0$ . To calculate these coefficients, a transition to turbulence in the boundary layer is assumed to occur near the leading edge (specifically, at 5 $\,\%$ of the chord length) due to the roughness of the hydrofoil surface. Knowing the value of $C_d$ for the hydrofoil enables, by subtraction in this steady test case, an estimate of the viscous friction drag $F_{\kern-1pt f}=1.3$ N added by the end plates. This value agrees with the order of magnitude of the friction drag generated by rough plates of similar dimension parallel to a flow (White Reference White1991). For simplicity, we assume that the viscous drag on the end plates varies minimally under wave conditions, and take its value to be constant. The lift coefficient $C_l$ obtained using XFoil agrees well with the experimental value of 1.11 found here at $\alpha _0=10^\circ$ .

Adding waves to the steady flow generates, as expected, variations over time of the force components $F_x$ and $F_z$ that follow the wave periodicity. The larger the wave amplitude is, the greater the magnitude of the force fluctuation. Only small linear force variations appear for small incident waves (blue curves), whereas large nonlinear variations with the emergence of local peaks are observed for waves of large amplitude (red curves). During these variations, a competition occurs between drag, lift and inertia forces, which are driven by the velocity magnitude, the angle of attack and the flow acceleration. As shown in figure 1, the flow velocity is in phase with the wave amplitude, whereas both the flow acceleration and the angle of attack have an offset of $\pi /2$ . In our case, as the hydrofoil is upside down, the maximum angle of attack is achieved on the descending slope of the waves where the deceleration of the flow is maximum.

We observe that $F_x$ has a $\pi /2$ offset with $\eta$ ; its minimal values are achieved when the deceleration and the angle of attack are maximum, and its maximal values are when the acceleration is maximum and the angle of attack minimum. This implies that the fluctuations are mostly driven by inertial effects, i.e. added-mass, consistent with the small Keulegan–Carpenter number $\textit{KC}\in [0.15,3.97]$ . In particular, the importance of inertia makes the horizontal force component highly negative during the flow deceleration, despite a flow moving in the positive direction. We can also note that $F_x$ does not vary symmetrically around its stationary value of $F_x\approx 2.1$ N and deviates more at its minima compared with its maxima. This is due to the angle of attack that, when it is maximal, increases the value of $C_{mx}$ and hence inertia effects, which here occur during the flow deceleration, and decreases $F_x$ . Temporal variation of the vertical force component $F_z$ is close to being in phase opposition with $\eta$ , meaning that a larger flow velocity will result in a smaller $F_z$ . The shift observed between $F_z$ and $\eta$ is due to the angle of attack that increases lift effects on the descendant wave slope. This indicates that $F_z$ is mostly driven by the lift (which is negative) and influenced by the angle of attack. Note also that for the highest wave steepness presented, $\epsilon =0.06$ , $F_z$ does not experience a single maximum each wave period, but remains elevated (i.e. at its minimal lift) for a longer fraction of the wave period. The hydrofoil responses to different wave forcing are also quantified via probability density functions (p.d.f.) in figure S3.

Figure 4. Theoretical model of the force variations (a) $F_x(t)$ and (b) $F_z(t)$ . Experimentally measured values are plotted as solid lines and the theoretical prediction as dashed lines, for selected cases $\lambda =4$ m and $\epsilon =0$ (black), $\epsilon =0.02$ (blue) and $\epsilon =0.06$ (red). (c,d) Variation of the steady-state drag and lift coefficients obtained using Xfoil, as well as the added mass coefficients obtained via potential theory, with the angle of attack $\beta$ . The vertical dashed lines are indicative of the minimum and maximum angles of attack encountered in the different cases $\epsilon =0.02$ (blue), $\epsilon =0.04$ (green) and $\epsilon =0.06$ (red).

A theoretical estimate, using (2.3), of the time-varying force components $F_x(t)$ and $F_z(t)$ is shown in figure 4. For the purposes of this estimate, the values of the external velocity $u$ and of its orientation $\alpha$ are measured by PIV as previously described in § 3. Then, the steady time-averaged coefficients of drag and lift, $C_d$ and $C_l$ , that vary with $\alpha$ are obtained using Xfoil (Drela Reference Drela1989) and plotted in figure 4(c). The coefficients of inertia, $C_{mx}$ and $C_{\textit{mz}}$ , which include Froude–Krylov and added mass forces (Faltinsen Reference Faltinsen1993), are estimated using potential flow theory around a cylinder and the Kármán–Trefftz transform (Milne-Thomson Reference Milne-Thomson1973) to approximate the NACA4415 closely. They are plotted as a function of the angle of attack $\beta$ in figure 4(d). The theoretically predicted value of $C_{mx}$ at the baseline angle of attack is close to the value found experimentally during an acceleration of the hydrofoil along the water tank. The theoretical force predictions obtained from (2.3) are then plotted in figure 4(a,b) (dashed lines) and compared with the experimental data (solid lines) for weak and strong forcing cases, $\epsilon =0.02$ and 0.06. The order of magnitude of the force variation is well captured by the model, especially for the weak forcing, but is not fully in agreement with the measurements. The discrepancy stems mainly from the estimation of the coefficients. The numerical estimations of $C_d$ and $C_l$ are indeed based on steady cases, and do not take into account the dynamics present in our system. Furthermore, the estimation of the inertia coefficients is based on linear potential flow theory, which is known to not be valid at long time and to not accurately represent the physics in stalled conditions (Chaplin & Subbiah Reference Chaplin and Subbiah1997; Grift et al. Reference Grift, Vijayaragavan, Tummers and Westerweel2019; Reijtenbagh, Tummers & Westerweel Reference Reijtenbagh, Tummers and Westerweel2023) where nonlinear effects occur, which is expected under large amplitude wave forcing. The different components (drag + lift and inertia) from (2.3) are plotted in figure S4. These evidence how inertia effects dominate the minimum value of $F_x$ and how lift (and a small drag component) drive the variation of $F_z$ . The mean influence of the flow-aligned end-plates has been accounted for here, and their role in unsteady flow assumed to be equivalent to that in steady flow, i.e. the presence of a roughly constant friction force in the horizontal direction. This assumption, while not precise, is supported by the larger velocity component in the horizontal compared with the vertical direction, and does not directly influence comparison between test cases. Further numerical, experimental and theoretical works, which are beyond the scope of the present study, would be needed to properly define the coefficients and accurately model the force variations in this way.

Figure 5. Phase-averaged spanwise vorticity $\varOmega _y$ at selected non-dimensional times $t^*=0.3$ , 0.45, 0.6 and 0.75 (from left to right), corresponding to the vertical dashed black lines in figure 3(d–f), for $\lambda =4$ m, and wave amplitudes $\epsilon =0.02$ , 0.04 and 0.06 (from top to bottom). The hydrofoil here has a straight leading-edge and is not equipped with tubercles. Black solid lines represent isocontours of $Q$ -criterion at an arbitrarily chosen value $Q=100$ s $^{-2}$ , which thus enclose regions of high vorticity, which we consider as vortices. Black cross, centre of vortex rotation (maximum of $\varGamma _1$ , following Michard et al. (Reference Michard, Graftieaux, Lollini and Grosjean1997)); black dot, geometric centroid of the vortex. Every second velocity vector is skipped, for clarity.

What is clear from the discrepancies observed in figure 4 is that linear flow theory fails to capture the intricacies of the unsteady hydrodynamics of a hydrofoil under wave forcing. This motivates an in-depth experimental characterisation of the unsteady flow which produces these forces. Contours of the phase-averaged spanwise vorticity, $\varOmega _y$ , are plotted for several times $t^*\in [0.3,0.6]$ for the case of $\lambda =4$ m, and wave amplitudes $\epsilon =0.02$ , 0.04 and 0.06 in figure 5. The corresponding phase-averaged velocity fields are represented by black arrows. Areas of high positive (anticlockwise) vorticity are shown in red (negative in blue) and the $Q$ -criterion is used to discriminate between shear and vortices. Black solid curves represent the contours of the $Q$ -criterion, set as $Q=100$ s $^{-2}$ (threshold chosen to avoid noise influence), and enclose regions of positive $Q$ -criterion, which we refer to as vortices. Any areas of $Q\gt 100$ s $^{-2}$ smaller than 2.3 cm $^2$ are attributed to measurement noise and not taken into account. The centre of vortex rotation (black cross) as well as the geometric centroid of the vortex (black dot) are plotted for every vortex and discussed in § 6. The time instants represented in this figure correspond to key moments in the development of stall, indicated by vertical dashed lines in figure 3(d–f). The full temporal evolution of the velocity and vorticity fields for this case and other forcing can be seen in the supplementary movies 1–8, in which all time instants $t^*\in [0,1]$ are shown. In the interval $t^*\in [0,0.3]\cup [0.6,1]$ , no relevant variation of the flow occurs as the angle of attack is close to 0 $^\circ$ . When $\epsilon =0.02$ (top row), no signs of boundary layer separation appear, as the angle of attack variations are small. Increasing the wave steepness to $\epsilon =0.04$ (middle row), generates angles of attack of greater amplitude and thus flow velocity fluctuations that induce a thickening of the vortical region near the trailing edge corresponding to the early stage of stall. Such trailing-edge boundary layer thickening and incipient separation is typical of quasi-steady stall phenomena. Finally, when the forcing amplitude is large, $\epsilon =0.06$ (last row), a compact accumulation of vorticity is observed near the leading edge. Significant boundary layer separation and roll-up is simultaneously observed near the trailing edge. The two vortices emerge at approximately $t^*=0.45$ , when the angle of attack is large, and are then advected along the hydrofoil while remaining close to its surface. These observations indicate an intermediate, quasi-static type of stall, but with significant roll-up of the boundary layer near the trailing edge of the suction surface. This is in contrast to the deep dynamic stall regime, characterised by a single large leading-edge vortex, as reported by Hrynuk & Bohl (Reference Hrynuk and Bohl2020) and Badia et al. (Reference Badia, Hickel, Scarano and Li2025). The deep dynamic stall regime is not observed on the straight-leading-edge hydrofoil in any of the wave forcing cases studied here, because the angle variation imposed by the waves is relatively limited.

5. Effects of leading-edge tubercles

In this section, the hydrodynamic behaviour of our second model, i.e. a hydrofoil equipped with leading edge tubercles, is tested and compared with the straight-leading-edge case. The tubercled geometry is expected to modify the flow behaviour, including the dynamic stall characteristics (Hrynuk & Bohl Reference Hrynuk and Bohl2020; Badia et al. Reference Badia, Hickel, Scarano and Li2025). The effect of tubercles on the generation of hydrodynamic force is quantified in figure 6, where the surface elevation, $\eta (t)$ , and the horizontal and vertical force components, $F_x(t)$ and $F_z(t)$ , are plotted over time for each of the two geometries. The steady case with no waves (dashed lines) and a selected wave-forced case, with $\lambda =4$ m and $\epsilon =0.06$ (solid lines), are shown. The same quantities, phase averaged over 40 wave periods (solid lines) along with the corresponding standard deviation (dashed lines), are plotted in figure 6(d–f). This large amplitude case is plotted here because it corresponds to conditions where the effects of tubercles are most visible. Other weaker wave forcing cases are plotted in figures S5 and S6. For the steady case, with no waves, $F_x\approx 2.1\pm 0.1$ N for the straight-leading-edge case and $F_x^{\textit{tub}}\approx 2.4\pm 0.1$ N for the tubercled hydrofoil, meaning that tubercles increase the steady-state drag force by approximately 14 $\,\%$ . This drag increase could preclude the use of tubercles in some contexts, but is relatively small compared with the observed magnitude of the $F_x$ fluctuations under wave forcing. Regarding the vertical force component, $F_z\approx -23.8\pm 0.2$ N and $F_z^{\textit{tub}}\approx -22.1\pm 0.2$ N, so an increase (i.e. decreased lift) of approximately $7\,\%$ is observed with tubercles. These two observations, i.e. small increase of drag and small decrease of lift, correspond to what has been observed before in steady cases for angles of attack smaller than the stall angle (Johari et al. Reference Johari, Henoch, Custodio and Levshin2007; Fan et al. Reference Fan, Dong, Li, Sun and Feng2022).

Figure 6. Time history of (a) the surface elevation $\eta (t)$ , (b) the horizontal force component $F_x$ and (c) the vertical force component $F_z$ for waves of wavelength $\lambda =4$ m and steepness $\epsilon =0.06$ . Data are given for a hydrofoil with straight leading-edge (blue solid lines) and a hydrofoil with leading-edge tubercles (red solid lines). (d,e,f) Phase-averaged data over 40 wave periods of the same. Coloured dashed lines correspond to the standard deviation of the average. Vertical black dashed lines correspond to time instants shown in figure 8.

Under unsteady wave forcing, just as in the straight-leading-edge case, variation of the force components $F_x$ and $F_z$ emerges, following the wave periodicity. A comparison with the hydrodynamic forces generated in the straight-leading-edge case reveals that the primary effect of tubercles appears in the $F_x(t)$ component (figure 6 b,e). Fluctuations of $F_x(t)$ are attenuated by approximately $35\,\%$ when leading-edge tubercles are introduced. This effect is quantified for several wave conditions in figure 7, where the difference $F_x^{\textit{str}}-F_x^{\textit{tub}}$ between the horizontal forces for the two hydrofoils is plotted over the non-dimensional time $t^*$ . The force difference becomes non-negligible (deviation above $15\,\%$ of the initial value) when $t^*\in [0.2,0.8]$ , which corresponds to the period during which the hydrofoil experiences a large angle of attack and where stall occurs. Outside of this interval, where the angle of attack is small, the flow remains attached to the hydrofoil and $F_x^{\textit{str}}-F_x^{\textit{tub}}\approx 0$ N. This implies that the hydrodynamic effect of tubercles under wave-induced forcing is primarily driven by their mediation of hydrofoil stall dynamics. It is of practical interest to note that the observed reduction in force fluctuation magnitude could be beneficial with regards to reducing structural fatigue, but at the cost of an increased mean drag. The impact of tubercles on horizontal forces generated during dynamic stall has, to our knowledge, never been previously investigated, as former studies focused primarily on vertical forces (Hrynuk & Bohl Reference Hrynuk and Bohl2020; Badia et al. Reference Badia, Hickel, Scarano and Li2025).

Figure 7. Variation over the non-dimensional time $t^*$ of the difference between the horizontal hydrodynamic force component in the cases with $F_x^{\textit{tub}}$ and without $F_x^{\textit{str}}$ tubercles, for (a) $\lambda =2$ m, (b) $\lambda =4$ m and (c) $\lambda =6$ m, for wave steepness $\epsilon =0.02$ (blue), $\epsilon =0.04$ (green) and $\epsilon =0.06$ (red). Note the differing scales of the vertical axes.

We observe only a minor impact of the presence of tubercles on the vertical force component, $F_z$ . The strong increase in lift expected at large angle of attack in steady cases (Johari et al. Reference Johari, Henoch, Custodio and Levshin2007; Fan et al. Reference Fan, Dong, Li, Sun and Feng2022) is not observed here. The effect of tubercles on lift generation under dynamic stall has not previously been measured. Hrynuk & Bohl (Reference Hrynuk and Bohl2020) did not directly measure forces and the authors only assumed a tubercle-induced increase of lift based on a greater measured circulation in the leading-edge stall vortex. Although our measurements agree with Hrynuk & Bohl (Reference Hrynuk and Bohl2020) that circulation in the leading-edge vortex is increased by the presence of tubercles, we show by direct force measurement that this does not translate to an increased lift force. Indeed, it is the total associated with the foil, both bound and shed, that contributes to lift force generation and this is not well captured under any usual definition of vortex. Indeed, the numerical simulations of Badia et al. (Reference Badia, Hickel, Scarano and Li2025) indicate that tubercles decrease the peak lift experienced by a foil undergoing dynamic stall induced by a linear pitch-up manoeuvre. This, they attribute to a decrease in circulation which was experimentally observed in the same paper. The apparent inconsistency between the two papers might be due to the differing pitch rate ranges they investigated ( $\varOmega ^*\gt 0.1$ for Hrynuk & Bohl Reference Hrynuk and Bohl2020 and $\varOmega ^*=0.05$ for Badia et al. Reference Badia, Hickel, Scarano and Li2025), although the absence of direct force measurements in either study prevents definitive conclusions. In the present study, the effect of tubercles on the vertical force is minimal, with the exception of a small deviation near the maximal values of $F_z$ . The vertical force in the straight-leading-edge case peaks twice per wave period, while with tubercles, the first peak is attenuated and the second enhanced, resulting in a single, asymmetric peak in $F_z$ . Because of the slight period-to-period variability in the wave forcing, this deviation in $F_z$ is best observed via the phase-averaged force curves in figure 6(d–f) as well as in probability density functions (p.d.f.s) of the measured force signals, which are included in figure S7. The apparent discrepancy between the present observations on peak lift and those of Badia et al. (Reference Badia, Hickel, Scarano and Li2025) can likely be attributed to the relatively small amplitude of the angle of attack variations investigated here. A divergence between the lift coefficient computed in the straight-leading-edge and tubercled cases is not observed by Badia et al. (Reference Badia, Hickel, Scarano and Li2025) until the foil is beyond the maximum angle of attack we achieve in the present study. Although the foil section used in the two studies differs, this strongly suggests that we should indeed not see a significant effect of tubercles on lift at these angles.

Figure 8. Phase-averaged spanwise vorticity $\varOmega _y$ below a hydrofoil fitted with leading-edge tubercles, at selected non-dimensional times $t^*=0.15$ , 0.3, 0.45 and 0.6 (from left to right), corresponding to the vertical dashed black lines in figure 6(d–f), for $\lambda =4$ m and wave amplitude $\epsilon =0.06$ . Three separate measurement planes are shown: at a tubercle peak, halfway between tubercle peak and trough, and at the trough between two tubercles (from top to bottom). Black solid lines represent isocontours of $Q$ -criterion at an arbitrarily chosen value $Q=100$ s $^{-2}$ , which thus enclose regions of high vorticity, which we consider as vortices. Black cross, centre of vortex rotation (maximum of $\varGamma _1$ , following Michard et al. (Reference Michard, Graftieaux, Lollini and Grosjean1997)); black dot, geometric centroid of the vortex. Every second velocity vector is skipped, for clarity. The black double arrow illustrates the typical length scale, $R_v$ , discussed in § 7.

To better understand how the flow is affected by tubercles, the phase-averaged spanwise vorticity $\varOmega _y$ around the hydrofoil for the case of wavelength $\lambda =4$ m and large wave steepness $\epsilon =0.06$ is plotted in figure 8 for selected dimensionless times $t^*\in [0.15,0.6]$ , corresponding to the instants of interest indicated by vertical dashed lines in figure 6(d–f). The three planes in which measurements were taken (at a tubercle peak, half-way between tubercle peak and trough, and at the trough between two tubercles) are represented in the figure (from top to bottom). Comparing the vorticity contours in this figure to those in figure 5 demonstrates that the flow is strongly altered by the presence of tubercles. The corresponding videos over the complete wave period $t^*\in [0,1]$ , including for cases of weaker wave forcing amplitude and larger wavelength, can be found in supplementary movies 1–8. A concentrated region of spanwise vorticity first emerges near the leading edge in the planes of measurement located between the tubercles and half-way between tubercle peak and trough (figure 8, bottom two rows). At $t^*=0.3$ , the leading-edge vortex has not yet appeared in the tubercle peak measurement plane, indicating that the initial boundary layer separation and stall vortex roll-up are localised phenomena associated with the tubercle troughs and the stall vortex is, in its early stages, spanwise discontinuous. By $t^*=0.45$ , the leading-edge vortex appears in all measurement planes and remains at this stage attached to the leading edge via a strong shear layer. This large vortex continues to grow over time, both spatially and in strength, and then detaches and is advected along the hydrofoil from approximately $t^*=0.6$ . The leading-edge vortex generated in the tubercled case is much greater in spatial extent than that observed in the straight-leading-edge case. Because of the limited angle-of-attack variation amplitude ( $20^\circ$ for $\lambda =4$ m and $\epsilon =0.06$ in our most extreme forcing case, against $30^\circ$ of Badia et al. Reference Badia, Hickel, Scarano and Li2025 and $50^\circ$ of Hrynuk & Bohl Reference Hrynuk and Bohl2020) in the present work, deep stall is not fully achieved in the case of straight-leading-edge geometry, placing our experiments close to a quasi-static stall regime. When tubercles are applied, however, behaviour similar to the other studies occurs, i.e. vortex generation initiates between tubercles at smaller angle of attack and leads to the generation of a more diffuse large leading-edge vortex structure. This means that, although the angle-of-attack variation is here too small to trigger a deep stall regime on the straight-leading-edge hydrofoil, the presence of tubercles initiates early dynamic stall at the leading edge, leading eventually to a large leading-edge vortex present across the span. This indicates that, for relatively small angle variations and rates, tubercles’ primary effect is that of enforcing a transition to a qualitatively different flow regime, which could, depending on application and design objective, be deleterious to the operation of such a hydrofoil. It would be of interest to ascertain whether stall is also triggered early in the case of a pitching foil and to quantify comprehensively for which unsteady motion classes tubercles can be beneficial or not.

6. Analysis of the vortex motion

This section presents a more quantitative in-depth analysis of dynamic stall vortex initiation and evolution on a hydrofoil under action of waves, and examines the effect of leading-edge tubercles on the dynamic stall process. We first examine the tangential velocity of the flow along the suction surface of each hydrofoil, defined as $u_t=\boldsymbol {u}\boldsymbol{\cdot }\boldsymbol {t}$ , where $\boldsymbol {t}$ are the vectors locally tangent to the hydrofoil suction surface. The values of $u_t$ are plotted as a function of space and time in figure 9 for the straight-leading-edge hydrofoil (panel a) and the tubercled hydrofoil (panels b–d) , for the case of $\lambda =4$ m and $\epsilon =0.06$ . Black solid contours delimit the spatio-temporal locations at which reverse flow occurs, i.e. where $u_t\lt 0$ . Discrete instances of reverse, or near-reverse, flow are marked using magenta ellipses and labelled alphabetically to ease description. Corresponding figures for other wave forcing cases, $\lambda =4$ m, $\epsilon =0.04$ and $\lambda =6$ m, $\epsilon =0.04$ are given in figures S8 and S9, where similar behaviour is displayed.

During the time intervals $t^*\in [0,0.2]\cup [0.8,1]$ , which correspond to angles of attack close to $0^\circ$ , $u_t$ is similarly distributed (correlation coefficient approximately 0.8) for the two hydrofoils. A streamwise acceleration from the leading edge to a maximum value of $u_t\approx 1.25$ m s⁻¹ at approximately $x/c\approx 0.25$ and then subsequent deceleration. No reverse flow is present during these time intervals. In the straight-leading-edge case, a reverse flow first appears near the trailing edge, $x/c\approx 0.75$ , in the zone marked ‘A’. This is followed by the emergence of a nearly reversed flow (zone B) near the leading edge, which subsequently propagates towards the hydrofoil’s trailing edge. It should be noted that, because of the finite resolution of the PIV measurement, $u_t$ is obtained at approximately 4 mm offset from the actual surface of the hydrofoil, and the boundary layer is not well resolved. For comparison, this offset is of similar order to the thickness $0.16c/Re^{1/7}\approx 4.7$ mm of a fully turbulent boundary layer over a flat plate of length $c$ (Schlichting & Kestin Reference Schlichting and Kestin1961). As such, values of $u_t$ are indicative of the approximate magnitude of the outer flow velocity and its direction, but may not accurately capture very localised, small-scale, flow reversal within the boundary layer. This means that there is uncertainty as to whether flow reversal actually occurs or not in zone B. The initial onset of flow reversal near the trailing edge is typical of static, or quasi-static, stall phenomena (Buchner et al. Reference Buchner, Soria, Honnery and Smits2018), which is not unexpected due to the low dimensionless pitch rates $\varOmega ^*\in [0.01,0.06]$ studied here. The subsequent leading edge (near-)reversal and associated vortex roll-up (figure 5) indicate, however, that incipient dynamic stall behaviour is also present.

Figure 9. Spatio-temporal variations of the non-dimensional phase-averaged velocity, $u_t/u_0$ , tangent to the hydrofoil suction surface at non-dimensional times $t^*\in [0,1]$ for $\lambda =4$ m and $\epsilon =0.06$ . The two cases of (a) straight leading edge and (b–d) with tubercles are shown. Panel (b) shows the measurement plane on a tubercle peak, panel (c) half-way between tubercle peak and trough, and panel (d) at the trough between tubercles. Black solid curves delimit spatio-temporal regions in which $u_t\lt 0$ . Alphabetically labelled magenta ellipses indicate zones of interest to ease description.

As seen in § 5, the addition of tubercles modifies the flow. First, for a short instant, $u_t$ is observed to increase on tubercles at $t^*=0.25$ and $x/c=0.2$ (figure 9 b). At the same streamwise location, a reverse flow emerges very shortly thereafter between the tubercles (zone C in figure 9 d). The emergence of this reverse flow is presumably triggered by three-dimensional flows that correspond to streamwise vortices (Badia et al. Reference Badia, Hickel, Scarano and Li2025) and corresponds temporally with the moment at which the horizontal force in the tubercled case deviates most significantly from that in the straight-leading-edge case (figures 6 e and 7 b). The reverse flow of zone C corresponds to the accumulation of vorticity observed previously in figure 8 and remains at an approximately constant $x/c \approx 0.3$ location until $t^* \approx 0.6$ . This relatively stationary leading-edge-attached vorticity accumulation, and associated reverse flow region, is distinct from the primary leading-edge vortex structure and appears related to the shear layer connecting the leading-edge vortex to the leading edge. The signature of the leading-edge vortex on the surface-tangential velocity is observed in a distinct region (zone D) of reverse flow, which presents at all measured planes along the spanwise direction from $t^*\in [0.4,0.75]$ . In this time interval, the reverse flow in zone D is observed to advect downstream at an approximately constant rate. The observed behaviour of zone D is qualitatively similar to the slightly earlier evolution of zone B in the straight-leading-edge case, but of a clearly larger magnitude. The advection of these separated flow structures past the trailing edge at $t^* \approx 0.75$ corresponds roughly with the moment at which the horizontal force measurements re-converge in figures 6(e) and 7(b). Note also the presence of a trailing-edge separation in figure 9(d), $x/c\approx 0.8$ , which is qualitatively similar to the quasi-static stall observed in the straight-leading-edge case (zone A in figure 9 a). This separation does not appear in the on-tubercle measurement plane (figure 9 b), as the streamwise vortices created by the tubercles induce a wall-normal flux of momentum towards the surface in that plane, but not in the plane situated between tubercles (Badia et al. Reference Badia, Hickel, Scarano and Li2025).

Figure 10. Location of the (a–c) geometric centres and (d–f) centres of rotation of each identified vortex, for $t^*\in [0,1]$ . Three cases of wave forcing are shown: (a,d) $\lambda =4$ m, $\epsilon =0.04$ ; (b,e) $\lambda =4$ m, $\epsilon =0.06$ ; and (c,f) $\lambda =6$ m, $\epsilon =0.04$ . Symbol shape and colour indicate different planes of measurement, i.e. straight leading edge (black diamonds), measurement plane coincident with a tubercle peak (red dots), measurement plane coincident with a tubercle trough (blue dots), and measurement plane halfway between tubercle peak and trough (green dots). Empty symbols correspond to the structures from zones A and C, and full symbols to zones B and D. Hydrofoil shape and location are indicated via black solid outline.

Individual phase-averaged vortices are identified and tracked over time from the PIV measurement data, using a method similar to that used by Fiscaletti, Kat & Ganapathisubramani (Reference Fiscaletti, de Kat and Ganapathisubramani2018) with two separate definitions: vortex locations can be extracted as the geometric centres of regions in which $Q\gt 100$ s $^{-2}$ (black dots in figures 5 and 8), or the locations can alternatively be found using the $\varGamma _1$ method, described by Michard et al. (Reference Michard, Graftieaux, Lollini and Grosjean1997), which computes the centres of rotation of the region in which $Q\gt 100$ s $^{-2}$ (black crosses in figures 5 and 8). We focus here on tracking and characterising those vortex structures with an origin at leading edge or in the suction-surface boundary layer. The clockwise-rotating trailing edge vortex that appears in some large-amplitude test cases is not considered in this analysis. The tracked vortex locations, defined by the geometric centres, are plotted in figure 10(a–c), and as defined by the centres of rotation in figure 10(d–f), for three separate wave forcing cases: (a,d) $\lambda =4$ m, $\epsilon =0.04$ , (b,d) $\lambda =4$ m, $\epsilon =0.06$ and (c, f) $\lambda =6$ m, $\epsilon =0.04$ . These three cases correspond to wave forcing of sufficient magnitude such that stall vortices can form and be detected in our measurement. Each marker’s shape and colour indicates the leading-edge geometry and measurement plane to which it relates. Empty symbols relate to the vortex structures associated with zones A and C, and full symbols to the vortex structures associated with zones B and D that are defined in figure 9. In the case of a straight leading edge (black diamonds), almost no vortex emerges for $\lambda =4$ m, $\epsilon =0.04$ (panels a, d), as the maximum achieved angle of attack of $17^\circ$ in this case is too small to generate measurable dynamic stall over the time scale of wave forcing, only the boundary layer thickening and roll-up near the trailing-edge (zone A, empty diamonds) tends to emerge. Increasing the wave steepness to $\epsilon =0.06$ (panels b, e) produces the quasi-static stall case where flow separation, boundary layer roll-up and resulting vortex structures occur near both the trailing edge (zone A, empty diamonds in $x/c\in [0.8,1]$ ) and the leading edge (zone B, full diamonds in $x/c\in [0.2,1]$ ). Wave forcing at longer wavelength, $\lambda =6$ m (panels c, f), implies larger variation of the angle of attack leading to stronger reverse flow and therefore stall. When tubercles are present at the leading edge, vorticity accumulation near the leading edge at $x/c\approx 0.2$ (zone C) is detected, in all plotted cases, in between tubercles (blue empty markers) and in the measurement plane half-way between tubercle trough and peak (green empty markers). The emergence of the large leading-edge vortex, indicating a dynamic stall regime, is detected in all measurement planes (full blue, green and red markers) at $x/c\in [0.4,1]$ . This vortex forms farther from the hydrofoil leading edge, at $x/c\approx 0.45$ , than in the straight-leading-edge hydrofoil wherein it forms at $x/c\approx 0.2$ . The two methods of vortex detection employed here yield similar results to one another, with the only difference being that the vortices’ centres of rotation are located closer to the hydrofoil surface than their geometric centres. This indicates a vortex asymmetry, which is also clear in figure 8. We select only the geometric centroid for use in subsequent analysis.

Figure 11. Variation of the (a–c) $x$ and (d–f) $z$ location of vortex geometric centroid with the non-dimensional time $tu_x/c$ , for the three wave forcing cases (a,d) $\lambda =4$ m, $\epsilon =0.04$ ; (b,e) $\lambda =4$ m, $\epsilon =0.06$ ; and (c,f) $\lambda =6$ m, $\epsilon =0.04$ . Symbol shape and colour indicate different planes of measurement, i.e. straight leading edge (black diamonds), measurement plane coincident with a tubercle peak (red dots), measurement plane coincident with a tubercle trough (blue dots), and measurement plane half-way between tubercle peak and trough (green dots). Black dash-dotted lines, linear best fits. Alphabetically labelled ellipses indicate zones of interest to ease description, with labelling identical to that used in figure 9. Black ellipses relate to the case of the straight-leading-edge hydrofoil, while data for tubercled hydrofoils are delineated using magenta ellipses. Empty symbols correspond to the vortices from zones A and C and full symbols to zones B and D. Vertical dashed magenta lines refer to the instants of maximal and minimal horizontal velocity $u_x$ (wave peak and trough) and the blue vertical dashed lines to the instants of maximal and minimal angle of attack $\beta$ .

Figure 10 indicates detected vortex locations in space, for various measurement times, but without explicitly including temporal information. In figure 11, the horizontal $x$ (panels a–c) and vertical $z$ (panels d–f) locations of the vortices are plotted over time, indicating similar dynamics for the three different wave forcing cases that are shown here, despite different flow conditions. Note that here, the time axis is non-dimensionalised using the convective time, $c/u_x$ , with $u_x$ the horizontal flow velocity measured by PIV, instead of the wave period, $T=1/f_r$ . This non-dimensionalisation is appropriate here as motion of the vortices post-detachment is driven by the external flow velocity. It should be noted, however, that the use of the non-dimensionalisation $t^* = t/T$ yields similar results, as the horizontal velocity $u_x$ depends on the wave period as $u_x\sim u_0\pm \mathcal{R}\omega$ . The horizontal advective velocity of the vortices is inferred by a linear fit to the experimentally observed variation of $x$ (dash-dotted lines). In the same way as in figure 9, various features of interest in figure 11 are indicated by black (straight leading edge) and magenta (tubercled leading edge) ellipses, keeping consistent alphabetic notation. As was observed previously in the case of a straight leading edge, except for in the case $\lambda =4$ m and $\epsilon =0.04$ , vortices are formed both near the trailing edge (zone A) and the leading edge (zone B), and are advected following roughly the hydrofoil surface. With tubercles present, the vorticity accumulation (zone C) near the leading edge appears early and remains closely associated with the leading edge for the duration of the measurement. Its motion is characterised only by a small horizontal velocity of magnitude between $0.04u_x$ and $0.11u_x$ , and a slight decrease over time of its $z$ coordinate, related to its spatial growth. Some time later, the large leading-edge vortex (zone D) forms and is observed to travel horizontally at a velocity of approximately $0.47\pm 0.09 u_x$ , similar to the advective velocity in straight-leading-edge case (zone B) despite a much earlier emergence. The vortex convection speed and its earlier emergence compared with the straight leading-edge case are found to be similar for the different cases. This similarity demonstrates the quantitative repeatability of the phenomenon despite the different flow unsteadiness applied. Except for in the case of $\lambda =4$ m and $\epsilon =0.04$ , this large leading-edge vortex moves only negligibly in the vertical direction, remaining located at $z/c\approx 0.1$ , and does not follow the hydrofoil shape, which is in contrast to the behaviour observed in the straight-leading-edge case. The difference is also easily observed in figure 10. This observation differs from that of Hrynuk & Bohl (Reference Hrynuk and Bohl2020), where tubercles make the vortex travel closer to the hydrofoil surface compared with the straight-leading-edge case. This difference could be due to the difference between using wave-induced flow motion rather than a pitching movement to trigger dynamic stall. Note also that from $tu_x/c\gt 3.5$ , the leading-edge vortex (zone D) has been advected beyond the tubercled hydrofoil’s trailing edge, whereas it is still present in the straight-leading-edge case (zone B). This difference in the length of time over which the leading-edge vortex resides in the vicinity of the hydrofoil suction surface may relate to the small qualitative difference in the vertical force $F_z$ time history creating the two peaks observed previously in figure 6(c) at each period.

The circulation contained within the identified vortex structures is calculated as

(6.1) \begin{align} \varGamma =\int _{Q\gt 100} \varOmega _y \,{\rm d}x\,{\rm d}z \end{align}

and plotted for the three different wave forcing conditions over time in figure 12. The $Q$ -criterion threshold, $Q\gt 100$ s $^{-2}$ , was found to qualitatively capture the main features of the vortex structures in our flow well. It should also be noted that, by construction, $Q$ -criterion only captures regions in which rotation dominates over strain, thus omitting shear and boundary layers, and serving here as a quantifiable definition delineating the primary vortex structures associated with flow separation. Care must be taken in the interpretation of these structures, as any chosen threshold value of the $Q$ -criterion is necessarily arbitrary and will never perfectly capture all relevant vorticity (see e.g. figures 5 and 8 at $t^*=0.45$ on the plane halfway between tubercle peak and trough). This source of uncertainty, linked to experimental noise, can cause scatter in computed circulation values when applying the analysis to experimental data, as is done here. Circulation contributions from the various vortex structures are plotted separately, zone A (panels a, b, c), zone C (panels d, e, f), and zone B and D (panels g, h, i), whereas the total circulation $\varGamma _{\textit{tot}}$ , with every contribution included, is reported in panels (j, k, l). Vortices that start to leave the measurement domain are no longer taken into account.

Figure 12. Evolution of the circulation $\varGamma$ as function of non-dimensional time $t^*$ for the three wave forcing cases: (a,d,g,j) $\lambda =4$ m, $\epsilon =0.04$ ; (b,e,h,k) $\lambda =4$ m, $\epsilon =0.06$ ; and (c,f,i,l) $\lambda =6$ m, $\epsilon =0.04$ . The different contributions are represented separately: (a,b,c) zone A; (d,e,f) zone C; and (g,h,i) zones B and D; and (j,k,l) the total circulation $\varGamma _{\textit{tot}}$ including all contributions. Symbol shape and colour indicate different planes of measurement, i.e. straight leading edge (black diamonds), measurement plane coincident with a tubercle peak (red dots), measurement plane coincident with a tubercle trough (blue dots), and measurement plane half-way between tubercle peak and trough (green dots). Vertical dashed magenta lines refer to the instants of maximal and minimal horizontal velocity $u_x$ (wave peak and trough) and the blue vertical dashed lines to the instants of maximal and minimal angle of attack $\beta$ . Note the varying y-axis scale.

The time-evolution of shed circulation in figure 12 is qualitatively similar across the various cases; only the magnitude changes, increasing with wave forcing intensity (wavelength and steepness). In the straight-leading-edge case (black diamonds), circulation increases first because of the near-trailing-edge boundary layer roll-up (zone A), then subsequently because of the formation of the leading-edge vortex (zone B). The two structures, once they are fully formed, have a similar circulation magnitude to one another, and so both contribute substantially to the total circulation $\varGamma _{\textit{tot}}$ . In the case of the tubercled hydrofoil, circulation appears first at the measurement planes between tubercles and half-way between tubercle peaks and troughs because of the leading-edge vorticity accumulation (zone C). The circulation contained within this structure then decreases from $t^*\approx 0.35$ as the angle of attack decreases and circulation is transferred to the large leading-edge vortex (zone D). The rate of increase of circulation in the leading-edge vortex is similar in all three planes of measurement. This leading-edge vortex contains a much greater circulation (approximately four times greater magnitude) compared with that observed in the straight-leading-edge case (zone B) (figure 12 h,i). The increase in total circulation in the plane coincident with the tubercle peaks lags behind that in the other measurement planes, as the circulation in the leading-edge vortex primarily originates from flow separation in the tubercle troughs, where a stronger shear layer is observed. Upon reduction in the strength of the feeding shear layer in between tubercles and in the half-tubercle plane from $t^*\approx 0.5$ , the total circulation measured in those planes reduces gradually towards parity with the circulation in the tubercle-coincident plane. Thereafter, $\varGamma _{\textit{tot}}$ is approximately identical for the three planes, as it consists almost entirely of circulation in the large diffuse leading-edge vortex. Direct comparison with previous studies on the effect of tubercles on dynamic stall (Hrynuk & Bohl Reference Hrynuk and Bohl2020; Badia et al. Reference Badia, Hickel, Scarano and Li2025) is difficult due to differing parameter ranges. We nevertheless observe concurrence with Hrynuk & Bohl (Reference Hrynuk and Bohl2020), in that the circulation is of greater magnitude between tubercles than on a tubercle. They observe that, as time progresses, the circulation in the tubercle peak plane becomes greater, a discrepancy with our observations which might be attributed to the non-dimensional pitching rate $\varOmega ^*$ being more than twice in that work than in the present study. As is common in the study of steady flow around a foil, it has previously been assumed that circulation in the separated vortex structures under dynamic stall is directly linked to lift force generation (Hrynuk & Bohl Reference Hrynuk and Bohl2020; Badia et al. Reference Badia, Hickel, Scarano and Li2025). In the present work, the measured circulation in the leading-edge vortex structures is greater in the presence of tubercles, but the time history of the vertical force, $F_z$ , is similar for both leading-edge geometries (figure 6 c). This seemingly contradicts earlier works’ assumptions that tubercles’ effects on leading-edge stall structures would significantly impact lift force, pointing to nuance in the role of shed circulation on force generation. Indeed, by Kelvin’s circulation theorem, identical kinematics will generate identical circulation, differing only in spatial distribution, and allocation between shed and bound vorticity. So long as this total circulation remains closely associated with the hydrofoil, the lift force should be unaffected. One must recognise however that not all generated circulation is included in computing the strength of what we call vortex structures, depending on the ‘vortex’ definition used. Such definitions, including $Q$ -criterion, typically exclude the circulation-containing boundary and shear layers, and are subject to noise, and therefore computed values’ relation to lift force generation should be interpreted with care.

To fully characterise the vortices’ formation, we must also consider their spatial dimension over time. We define vortex size by computing an effective radius $r$ , defined as $r=\sqrt {\mathcal{A}/\pi }$ , where $\mathcal{A}$ is the measured area in which $Q\gt 100$ s $^{-2}$ . This radius can be compared with $R$ , the theoretical typical length scale of the wave-induced flow, which is defined as the vertical displacement,

(6.2) \begin{align} \phantom {,}R=\mathcal{R}\tan \alpha _m, \end{align}

made by a fluid particle over a half-period of the wave forcing, at the hydrofoil nominal immersion depth, $z=-h$ , with $\mathcal{R}=Ae^{- {2\pi h}/{\lambda }}$ and $\tan (\alpha _m)={\mathcal{R}\omega }/{u_0}$ given by the amplitude of variation of the angle $\alpha$ from (2.2). The variation of $r/R$ with non-dimensional time $t^*$ is plotted in figure 13 for the different vortex structures: zone A(panels a, b, c), zone C (panels d, e, f), and zones B and D (panels g, h, i) and for the different wave forcing conditions. Data are given for each vortex up to the moment that the vortex begins to leave the measurement domain. For the straight-leading-edge hydrofoil (black diamonds), the two vortices (boundary layer roll up at trailing edge, zone A, and leading-edge vortex, zone B) quickly grow until approximately equivalent radii of $r\approx 0.6R$ . For the tubercled hydrofoil, the leading-edge region of vorticity accumulation (zone C) first grows in spatial extent, before shrinking from $t^*\approx 0.4$ both between tubercles and in the half-tubercle measurement plane (blue and green, respectively). The large leading-edge vortex (zone D) emerges subsequently in every measurement plane at the same time. For $\lambda =4$ m and $\epsilon =0.04$ (figure 13 g), the vortex radius increases over time until $r/R\approx 1.5$ , with experimental data being relatively spread. For the two other cases (panels band c), the large diffused leading-edge vortex grows from $t^*\approx 0.35$ to $t^*\approx 0.5$ in every measurement plane. Its radius grows until $r\approx R$ and tends to maintain this value until the vortex leaves the measurement domain. As seen in figure 10(b,c), in the tubercled cases, the leading-edge vortex is advected almost horizontally, without following the hydrofoil surface. The tendency of the vortex effective radius to converge to $r=R$ supports an interpretation that the spatial extent of leading-edge vortex formation is linked to the vertical displacement length scale of the wave induced fluctuating flow. This length scale appears to play a more direct role in the case of the tubercled leading edge, where the wave-induced vertical velocity component can penetrate between tubercles.

Figure 13. Evolution of the non-dimensional effective radius, $r/R$ , of the different vortex structures as a function of the non-dimensional time $t^*$ for the three wave forcing cases: (a,d,g) $\lambda =4$ m, $\epsilon =0.04$ ; (b,e,h) $\lambda =4$ m, $\epsilon =0.06$ ; and (c,f,i) $\lambda =6$ m, $\epsilon =0.04$ . The different vortices are represented separately: (a,b,c) zone A; (d,e,f) zone C; and (g,h,i) zones B and D. Symbol shape and colour indicate different planes of measurement, i.e. straight leading edge (black diamonds), measurement plane coincident with a tubercle peak (red dots), measurement plane coincident with a tubercle trough (blue dots), and measurement plane half-way between tubercle peak and trough (green dots). Vertical dashed magenta lines refer to the instants of maximal and minimal horizontal velocity $u_x$ (wave peak and trough) and the blue vertical dashed lines to the instants of maximal and minimal angle of attack $\beta$ .

7. Quantification of the increased drag due to tubercles

We can now use our observations of the flow behaviour to mechanistically explain the effect of tubercles on the generation of horizontal force $F_x$ , which is plotted in figures 6 and 7. The observed difference in $F_x$ occurs primarily in the time interval $t^*\in [0.2,0.6]$ , which corresponds to conditions where the angle of attack is maximum and the vortex structures on the suction side of the hydrofoil are in their early stages of formation. As observed in figure 6, the effect of tubercles on the value of $F_x$ in the steady case (without waves) is small compared with the magnitude of the force fluctuations in the unsteady cases. We, therefore, for expediency in modelling the effect of tubercles, assume that the steady-state drag coefficient is approximately identical for the two leading-edge geometries. We hypothesise the same for added mass, given that the planform areas of the two hydrofoil models are identical and the tubercles represent a relatively minor perturbation to the overall geometry. Under these assumptions, the observed force difference between straight- and tubercled-leading-edge cases can primarily be attributed to modified stall behaviour in the presence of tubercles. The primary qualitative flow difference between the straight-leading-edge and tubercled cases investigated here is an accumulation of vorticity near the leading edge, within the first third of the chord length. This accumulation relates to the separated shear layer which will later form the leading-edge vortex; however, in the tubercled case, it emerges earlier within each wave period, remains distinct from the leading-edge vortex as indicated by the $Q$ -criterion and exhibits minimal downstream advection. The downstream advection rate $u_v$ observed is only between $0.04u_x$ and $0.11u_x$ in $t^*\in [0.2,0.6]$ (figure 11 a–c). The motion of this vortex structure with the hydrofoil, near its leading edge, presents an additional obstruction to the horizontal passage of the external flow. We evaluate the consequence for the drag force if we make the simplistic assumption that this structure’s effect can be considered analogous to that of a bluff body of equivalent dimension, given by its effective radius $R_v=r(x/c\lt 1/3)$ , moving at the observed advective velocity $u_v$ . As an illustrative example, the length scale $R_v$ is shown in figure 8 between tubercles at $t^*=0.45$ , and corresponds to zone C discussed in § 6. To model this assumption, we introduce an additional steady-state drag term,

(7.1) \begin{align} \phantom {,}F_x^{\textit{tub}}-F_x^{\textit{str}}=C_{dv}R_v \rho (u-u_v)^2s, \end{align}

with $C_{dv}=1$ , according to the order of magnitude expected for a bluff body, such as a circular cylinder (Heddleson, Brown & Cliffe Reference Heddleson, Brown and Cliffe1957). If these assumptions are sufficiently accurate to explain the observed horizontal force difference between straight-leading-edge and tubercled cases, the constructed length scale $F_x^{\textit{tub}}-F_x^{\textit{str}}/[\rho (u-u_v)^2s]$ should be equal to the observed effective radius $R_v$ of the accumulation of leading-edge vorticity. Figure 14 shows the relationship between the experimentally determined values of these two variables. A linear regression on the data for the three cases, $\lambda =4$ m, $\epsilon =0.04$ (black dots); $\lambda =4$ m, $\epsilon =0.06$ (red squares); and $\lambda =6$ m, $\epsilon =0.04$ (green diamonds) yields a slope of $0.98\pm 0.18$ and vertical-axis intercept of $0.35\pm 2.39$ , with R-square of 0.74 (dashed line). The solid line of slope 1 represents the prediction of (7.1), and lies well within the 95 % confidence interval (blue shading) of the linear regression. The few clear outliers (markers with magenta borders) originate from the sensitivity of the $Q=100$ s $^{-2}$ contour to experimental measurement noise, which in occasional cases can artificially displace the estimated vortex centre by a large amount (see supplementary movies 5, 7 and 8). These outliers are not taken into account in the linear fit. Including these outliers in the linear regression yields a similar slope of $0.89\pm 0.26$ , intercept of $1.73\pm 3.47$ , with R-square of 0.48. The prediction of (7.1) remains in that case well within the 95 $\,\%$ confidence interval of the fit. This surprisingly simplistic model explains much of the variation in force caused by stall in the presence of large amplitude wave forcing, using only knowledge of the dynamic stall vortex dimensions and velocity. Force variation in unstalled conditions, $\lambda =2$ m in figure 7(a) for example, cannot of course be accounted for by this method. Further study is needed to better understand the bounds of the assumptions made in this work, and extend it across a greater range of forcing conditions and stall regimes. Nevertheless, the present work quantifies for the first time with direct force measurements how tubercles modify horizontal forces on a foil undergoing dynamic stall. Furthermore, we propose a simple model, which explains to a first approximation tubercles’ effect on horizontal force generation by reference to experimentally observed flow behaviour.

Figure 14. Effective length scale of the leading-edge vorticity accumulation region, as predicted by the difference in observed horizontal force between straight leading edge and tubercles hydrofoils, following (7.1) and assuming $C_{dv}=1$ . The three stall cases are represented here as: $\lambda =4$ m, $\epsilon =0.04$ (black dots); $\lambda =4$ m, $\epsilon =0.06$ (red squares); and $\lambda =6$ m, $\epsilon =0.04$ (green diamonds). The black dashed line is a best linear fit to the experimental data and yields to a slope of $0.98\pm 0.18$ and vertical-axis intercept of $0.35\pm 2.39$ , with R-square of 0.74. The clear outliers outlined in magenta are not taken into account in the fit. Blue shading indicates the 95 $\,\%$ confidence interval of the fit. The black solid line (slope 1) is the prediction of (7.1).