1. Introduction
Undulatory fish swimming can be classified by bending pattern, from travelling-wave-dominated anguilliform to standing-wave-dominated ostraciiform (Gray Reference Gray1933; Lindsey Reference Lindsey, Hoar and Randall1978). The bending pattern determines the flow characteristics around the fish and the resultant propulsion (Lighthill Reference Lighthill1960, Reference Lighthill1970). Standing wave propulsion produces greater thrust at high Reynolds numbers, while travelling wave propulsion is more efficient at intermediate Reynolds numbers (Borazjani & Sotiropoulos Reference Borazjani and Sotiropoulos2010; Cui et al. Reference Cui, Yang, Shen and Jiang2018).
Anguilliform swimmers use muscles along their bodies to generate travelling waves propagating towards the tails (Tytell & Lauder Reference Tytell and Lauder2004). Experiments with a robotic eel have demonstrated that travelling wave propulsion results in higher efficiency, defined as the ratio of forward thrust to power consumption (Anastasiadis et al. Reference Anastasiadis, Paez, Melo, Tytell, Ijspeert and Mulleners2023). The enhanced swimming efficiency is attributed to the ability of travelling waves to transport vortices along the eel body. As the vortices travel, they grow continuously until they are shed at the tail, thereby enhancing propulsion. In contrast, swimmers using standing waves generate thrust by moving their tails as a lever. This results in lift force perpendicular to the caudal fin generated by the circulation. While the thrust generated by this method is high, a significant amount of the thrust is directed sideways, reducing efficiency (Thandiackal et al. Reference Thandiackal, White, Bart-Smith and Lauder2021).
While travelling wave propulsion is beneficial for swimming efficiency, generating these waves in passive elastic foils that could serve as propulsors for man-made swimming devices presents a challenge (Ramananarivo, Godoy-Diana & Thiria Reference Ramananarivo, Godoy-Diana and Thiria2014). In finite-sized plates, flexural waves propagating along the plate are reflected at the free end, resulting in standing wave oscillations. Wave reflection can be suppressed by tapering the thickness of the plate near the free end. Because wave speed is proportional to stiffness, gradually decreasing stiffness causes the waves to slow down, trapping them near the free end and reducing reflection. This phenomenon is known as the acoustic black hole effect (Mironov Reference Mironov1988). Devices using this effect have been utilised for vibration absorption (Feurtado, Conlon & Semperlotti Reference Feurtado, Conlon and Semperlotti2014; Conlon, Fahnline & Semperlotti Reference Conlon, Fahnline and Semperlotti2015), sound attenuation (Conlon et al. Reference Conlon, Fahnline and Semperlotti2015) and energy harvesting (Zhao, Conlon & Semperlotti Reference Zhao, Conlon and Semperlotti2015). Furthermore, fins and rays with tapered thickness and stiffness are common in fish (Alben, Madden & Lauder Reference Alben, Madden and Lauder2007; Lauder et al. Reference Lauder, Anderson, Tangorra and Madden2007; Aiello et al. Reference Aiello, Westneat and Hale2017, Reference Aiello, Hardy, Cherian, Olsen, Ahn, Hale and Westneat2018), suggesting that thickness tapering may be beneficial for undulatory swimming.
Experiments and simulations confirm a correlation between stiffness variation and swimming performance. Oscillating fins with varying stiffness tend to exhibit larger curvatures and trailing edge deflections, redirecting fluid flow and increasing thrust (Kancharala & Philen Reference Kancharala and Philen2016). Propulsors with tapered thickness demonstrate a broader range of actuation frequencies, resulting in higher swimming efficiency (Yeh, Li & Alexeev Reference Yeh, Li and Alexeev2017). Both simulations and experiments indicate that tapered propulsors yield a high standing wave ratio
$\mathcal{S}$
(Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022; Leroy-Calatayud et al. Reference Leroy-Calatayud, Pezzulla, Keiser, Mulleners and Reis2022), a parameter characterising the ratio of travelling to standing waves that form the propulsor’s bending pattern and is correlated with swimming efficiency (Cui et al. Reference Cui, Yang, Shen and Jiang2018). Simulations show that travelling waves facilitate the transport of side edge vortices (SEVs) towards the trailing edge, enhancing the trailing edge vortices beneficial for thrust production (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). This in turn leads to a more efficient locomotion. Furthermore, the increased trailing edge vorticity yields a bifurcating reverse von Kármán wake.
In this work, we perform three-dimensional numerical simulations using a fully coupled fluid–structure interaction (FSI) computational model to probe the effects of standing and travelling waves on undulatory locomotion of an elastic propulsor. To this end, we investigate the hydrodynamic performance of an elastic plate that consists of a uniform-thickness root section and a tapered attachment with an exponentially decreasing thickness. The plate oscillates in viscous fluids at different Reynolds numbers
$\textit{Re}$
. Viscous dissipation at lower
$\textit{Re}$
facilitates the formation of travelling waves (Ramananarivo et al. Reference Ramananarivo, Godoy-Diana and Thiria2013, Reference Ramananarivo, Godoy-Diana and Thiria2014), thereby enhancing the acoustic black hole effect due to plate tapering (Mironov Reference Mironov1988). On the other hand, higher viscosity diminishes propulsion efficiency by increasing viscous loss. We compare the thrust, power and efficiency of the tapered plate with a plate of uniform thickness. We analyse the flow patterns generated by the tapered and uniform plates to gain insights into their fluid dynamics at different
$\textit{Re}$
.
The paper is organised as follows. Section 2 describes the set-up of the problem, while the computational model is described in § 3. Section 4 presents the results of this study and their discussion: § 4.1 analyses the kinematics and bending patterns of the plates; § 4.2 presents and discusses the hydrodynamic performance of the plates; and § 4.3 analyses the flow patterns generated by the plates and correlates them with their hydrodynamic performance. Our conclusions are summarised in § 5.
2. Problem set-up
We consider a thin elastic plate of length
$L$
and width
$W$
oscillating in an incompressible Newtonian fluid of density
$\rho$
and dynamic viscosity
$\mu$
. The plate is actuated at the leading edge with sinusoidal oscillations given by
$A(t) = A_0 \sin (2\pi t / \tau )$
, where
$A_0$
is the leading edge displacement amplitude,
$\tau$
is the oscillating period, and
$t$
is time. The plate material is isotropic and homogeneous with density
$\rho _s$
, Poisson’s ratio
$\nu$
, Young’s modulus
$E$
. Figure 1(a) illustrates the plate and its actuation method.
(a) Elastic plate actuated at the leading edge. (b) Thickness profile of the exponentially tapered plate with
$b=5$
and
$L_t/L = 0.5$
. (c) Computational domain with a refined inner mesh centred around the plate shown by the blue lines.
Two different plates are analysed in this study, one with uniform thickness
$h_0$
that serves as a reference case, and one with exponential tapering. Exponential tapering is chosen for the tapering shape as it enhances travelling waves in the plate bending pattern (Mironov Reference Mironov1988; Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). The thickness
$h(x)$
of the exponentially tapered plate is a piecewise function (figure 1
b)
\begin{align} h(x) = \left \{ \begin{array}{ll} h_0, & \text{for } \quad 0 \leqslant x \leqslant L_c, \\ h_0 b^{-\zeta }, & \text{for } \quad L_c \lt x \leqslant L, \end{array} \right . \end{align}
where
$b = h(L)/h_0$
is the tapering ratio, and
$\zeta = (x-L_c)/(L-L_c)$
. The local stiffness is given by the bending stiffness
$D(x) = E\,h(x)^3/12(1-\nu ^2)$
. The plate surface density is
$\rho _h = \rho _s\, h(x)$
.
We characterise the plate elasticity using
$r = \sqrt {D_r/D}$
, defined as the ratio of the resonance plate elasticity
$D_r$
, leading to the resonance oscillations for a given actuation frequency
$ {2\pi }/{\tau }$
, and the plate elasticity
$D$
. Thus
$r$
represents the vicinity of the plate oscillation to the fundamental resonance frequency where
$r =1$
. We therefore refer to
$r$
as the frequency ratio. We determine
$D_r$
computationally by varying the Young’s modulus of the plate to find the frequency that generates the largest tip displacement, signifying resonance. The value of
$D_r$
depends on the Reynolds number
$\textit{Re}$
of the plate. In our simulations, we vary
$r$
by varying plate elasticity, and keep the oscillation frequency constant.
The fluid flow generated by the oscillating plate is characterised by the Reynolds number
$\textit{Re} = {\rho U_c L}/{\mu }$
, where
$U_c = L/\tau$
is the characteristic plate velocity. The plate mass ratio is defined as
$\chi = \rho W/\rho _s h_0$
. The hydrodynamic performance of the plate is characterised using its dimensionless period-averaged thrust
$\mathcal{F} = F_x/F_0$
, power input
$\mathcal{P} = P/P_0$
, and efficiency
$\eta = \mathcal{F}/\mathcal{P}$
. Here,
$F_x=\int _{\tau } \int _{L W} g_x\, {\rm d}s\, {\rm d}t/\tau$
is the period-averaged force acting on the plate in the
$x$
-direction, and
$P=\int _{\tau } \int _{L W} \boldsymbol {g} \boldsymbol{\cdot }\boldsymbol {v}\, {\rm d}s\, {\rm d}t/\tau$
is the period-averaged power. Furthermore,
$g_x$
denotes the
$x$
-component of the local distributed hydrodynamic force
$\boldsymbol{g}$
, and
$\boldsymbol {v}$
is the local plate velocity. The characteristic force
$F_0$
is defined as
$({1}/{2})\rho U_c^2 W L$
, and the characteristic power input
$P_0$
is defined as
$( {1}/{2})\rho U_c^3 W L$
.
The flow vorticity generated by the oscillating plate is characterised using dimensionless total enstrophy averaged over one oscillation period
$\mathcal{E} = \tau L^{-3}\int _{\tau } \int _{V} |\omega (x,y,z,t)|^2\, {\rm d}V\, {\rm d}t$
, where
$\omega (x,y,z,t)$
is the instantaneous flow vorticity, and
$V$
is integration volume.
To characterise the plate bending pattern, we define the dimensionless maximum tip displacement of the plate
$\delta =d_{\textit{max}}/L$
, where
$d_{\textit{max}}$
is the maximum tip displacement. Furthermore, we define the normalised maximum bending area of the plate
$\mathcal{A}=\max(\int _0^L z(x,y=0,t)\, {\rm d}x)/A_0 L$
. Finally, we introduce the standing wave ratio
$\mathcal{S}$
defined as the ratio of travelling and standing waves in the plate bending pattern. It is computed by taking the two-dimensional fast Fourier transform of the time-dependent bending pattern (Bani-Hani et al. Reference Bani-Hani, Borazjani, Esfahani, Krovi and Karami2015). Two peaks arise at the same spatial frequency, with the taller peak amplitude
$f_1$
representing the propagating forward waves, and the shorter peak amplitude
$f_2$
representing the reflected waves. The magnitude of
$\mathcal{S}$
is defined as
$\mathcal{S}=(f_1-f_2)/f_1$
. Thus
$\mathcal{S}=0$
represents purely standing waves, whereas
$\mathcal{S}=1$
represents purely travelling waves. We use 50 plate profiles over one period of steady-state plate oscillation to evaluate
$\mathcal{S}$
.
3. Computational model
To model the oscillating plate, we use our in-house FSI computational framework. The framework integrates a lattice Boltzmann model (LBM) for modelling a viscous incompressible Newtonian fluid, and a finite difference model (FDM) for modelling the thin elastic plate. The LBM and FDM solvers are bidirectionally coupled to yield the fully coupled FSI solver.
The LBM is a mesoscale modelling method that uses fictitious fluid ‘particles’ moving along a cubic lattice with spacing
$\Delta x$
to simulate fluid flows governed by the Navier–Stokes equations. The process proceeds by these particles flowing to adjacent lattice sites, colliding with other arriving particles, and relaxing the distribution towards equilibrium. We use a D3Q19 cubic lattice with a 19-velocity continuous distribution function
$f_i(\boldsymbol {r},\boldsymbol {c}_i,t)$
that characterises the fluid flow in three spatial dimensions. The distribution function represents the mass density of the particles at lattice site
$\boldsymbol{r}$
at time
$t$
, moving in the lattice direction
$i$
with a constant lattice velocity
$\boldsymbol {c}_i$
. The time evolution of
$f_i$
is evaluated by integrating the discrete Boltzmann equation
$f_i(\boldsymbol {r}+\boldsymbol {c}_i\, \Delta t, t+\Delta t)= f_i(\boldsymbol {r} ,t) + \varOmega _i$
(Timm et al. Reference Timm, Kusumaatmaja, Kuzmin, Shardt, Silva and Viggen2016). Here,
$\varOmega _i$
is the collision operator and represents the change in
$f_i$
due to instantaneous fluid particle collisions at each lattice node, and
$\Delta t$
is the integration time step. We use a dual relaxation time collision operator introduced by Ladd & Verberg (Reference Ladd and Verberg2001) that enhances model stability in simulations of flows with larger
$\textit{Re}$
compared to a more traditional single relaxation time collision operator. The macroscopic quantities of density
$\rho$
, momentum
$\rho \boldsymbol {u}$
, and stress
$\boldsymbol{\varPi }$
are recovered from the distribution function as follows:
\begin{align} \rho = \sum _i{f_i}, \quad \rho \boldsymbol {u} = \sum _i{f_i \boldsymbol{c_i}}, \quad {\varPi } = \sum _i{f_i \boldsymbol{c_i} \boldsymbol{c_i}}. \end{align}
We consider a plate whose thickness is much smaller than its other dimensions, such that it does not significantly affect the flow. As such, the plate is modelled as a thin elastic two-dimensional plate that satisfies the Kirchhoff–Love assumptions (Timoshenko & Woinowsky-Krieger Reference Timoshenko and Woinowsky-Krieger1959). The plate governing equation is solved using a central finite difference scheme with ghost nodes at the boundary to impose the boundary conditions, and with time advancing using the explicit Euler method (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022).
The fluid and solid solvers are coupled at the fluid–solid interface using a two-way coupling method (Alexeev, Verberg & Balazs Reference Alexeev, Verberg and Balazs2005; Alexeev & Balazs Reference Alexeev and Balazs2007; Mao & Alexeev Reference Mao and Alexeev2014). A no-slip, no-penetration condition is imposed on the moving plate surface using a linearly interpolated bounce-back rule that enables subgrid resolution with second-order accuracy (Bouzidi, Firdaouss & Lallemand Reference Bouzidi, Firdaouss and Lallemand2001; Chun & Ladd Reference Chun and Ladd2007). The rule is applied to distributions that cross the plate boundary during the propagation step accounting for the instantaneous position and velocity of the plate. Hydrodynamic forces on the plate are calculated using the momentum exchange method that accounts for the momentum transfer between the plate and the fluid. To this end, the momentum of the distribution function is evaluated before and after the interaction with the plate boundary. The resultant hydrodynamic forces are distributed to the surrounding FDM nodes using a procedure that conserves force and moment (Alexeev et al. Reference Alexeev, Verberg and Balazs2005; Demirer, Oshinowo & Alexeev Reference Demirer, Oshinowo, Erturk and Alexeev2021). A more detailed description of our FSI solver as well as its validation can be found elsewhere (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022).
Figure 1(c) shows the computational domain with dimensions
$8L \times 6L \times 8L$
. The plate is surrounded by a fine grid of size
$4L \times 3L \times 3L$
, in order to better resolve the flow near the oscillating plate. The refined grid has spacing
$\Delta x_{\textit{fine}} = 1$
, and the coarse grid has spacing
$\Delta x_{\textit{coarse}} = 2$
. A periodic boundary condition is imposed at the outer boundaries of the computational domain that is sufficiently large for the wake from the oscillating plate to dissipate. We set fluid density
$\rho = 1$
, and vary the dynamic viscosity
$\mu$
to change
$\textit{Re}$
. The plate dimensions are
$L = 75$
and
$W = 30$
, leading to aspect ratio
$L/W=2.5$
. A grid sensitivity analysis was performed, indicating that further increase of grid resolution does not noticeably affect the simulation results. The thickness at the plate root corresponds to
$h_0=0.225$
. As noted above, the plate is modelled as a two-dimensional plate, in which case the plate properties are introduced in the FDM through the bending stiffness
$D$
and surface density
$\rho _h$
, which in turn depend on
$h_0$
. The tapered length is
$L_t= 0.5L$
, and the thickness ratio is
$b = 5$
. The plate is discretised with 21 and 11 FDM nodes along the plate length and width, respectively. The plate is driven at the root to oscillate with amplitude
$A_0 = 0.05L$
and period
$\tau = 2000$
. The plate mass ratio is
$\chi = 5$
, and Poisson’s ratio is
$\nu = 0.31$
. We perform the simulations for 40 periods of plate oscillations for the flow to reach a steady state.
All dimensional values used in our computational model are specified using LBM units that are defined as the lattice spacing
$\Delta x=1$
, fluid density
$\rho = 1$
, and time step
$\Delta t=1$
. In what follows, the simulation results are presented in dimensionless form to facilitate the comparison with experimental data and other computational models.
4. Results and discussion
4.1. Kinematics and bending pattern
Bending patterns for exponentially tapered and uniform plates oscillating at
$\textit{Re} = 100$
with different frequency ratios. The first column shows the exponentially tapered plates, and the second column shows the uniform plates. See also supplementary movie 1.
Figure 2 shows the bending patterns for the tapered and uniform plates oscillating at
$\textit{Re}=100$
. These plots show the bending pattern of the plate at
$10$
equally spaced time intervals over one oscillation period. The red lines show the bending profiles in the first half of the period, and the blue lines show the bending profiles during the second half of the period.
The bending patterns for both plates are similar at resonance
$r\approx 1$
, but diverge drastically for post-resonance oscillations. For
$r\approx 2$
and
$r\approx 3$
, the exponentially tapered plate develops bending patterns where tip displacement significantly exceeds root displacement. The uniform plate, on the other hand, undergoes a seesaw motion, with tip and root displacement having nearly the same amplitudes. We find that the bending patterns remain nearly unaffected by the Reynolds number in the range between 100 and 2000 (see supplementary movie 1).
(a) Maximum tip displacement
$\delta$
and (b) standing wave ratio
$\mathcal{S}$
as functions of frequency ratio
$r$
. The solid lines with solid markers represent uniform plates, while the dashed lines and empty markers represent exponentially tapered plates. Black squares, green diamonds, red triangles and blue circles represent
$\textit{Re}$
values 100, 500, 1000 and 2000, respectively.
In figure 3(a), we show the normalised maximum tip displacement
$\delta$
as a function of the frequency ratio
$r$
for the tapered and uniform plates of different
$\textit{Re}$
. At resonance, the tapered and uniform plates at
$\textit{Re}=100$
show slightly greater of tip displacement than respective plates with higher
$\textit{Re}$
that show similar magnitudes of
$\delta$
. Off-resonance, the tapered plates yield drastically higher
$\delta$
compared to the uniform plates. For both plate types,
$\delta$
weakly depends on
$\textit{Re}$
. For the tapered plates,
$\delta$
decreases to a local minimum at
$r \approx 2$
, then recovers for higher
$r$
to approximately match its magnitude at resonance. The plates with uniform thickness exhibit more than a
$50\,\%$
reduction in the tip displacement at
$r \approx 3$
compared to the resonance, then
$\delta$
gradually increases with
$r$
as the plate oscillations approach the second resonance frequency.
Figure 3(b) shows the relationship between the standing wave ratio
$\mathcal{S}$
and the frequency ratio
$r$
. We find that at resonance,
$\mathcal{S}$
for both tapered and uniform plates is close to zero, indicating that resonance oscillations are dominated by standing waves, regardless of variations in plate thickness. For tapered plates,
$\mathcal{S}$
increases with increasing
$r$
, reaching approximately 0.8 at
$r \approx 4$
. Thus travelling waves dominate the oscillations of tapered plates with higher
$r$
. For plates with uniform thickness, the standing wave ratio remains rather small for the entire range of tested
$r$
, indicating that these plates exhibit oscillations dominated by standing waves. The figure also shows that
$\mathcal{S}$
is insensitive to
$\textit{Re}$
. Decreasing
$\textit{Re}$
from 2000 to 100 results in only a minor change in
$\mathcal{S}$
. This indicates that the increase of
$\mathcal{S}$
with
$r$
for tapered plates is not a result of viscous dissipation, but can rather be attributed to the plate geometry suppressing wave reflection at the trailing edge (Mironov Reference Mironov1988).
4.2. Hydrodynamic performance
(a) Normalised thrust
$\mathcal{F}$
, (b) normalised power
$\mathcal{P}$
, and (c) plate efficiency
$\eta$
as functions of frequency ratio
$r$
. Inset in (a) shows the dependence of normalised thrust
$\mathcal{F}$
on tip displacement
$\delta$
. (d) Plate efficiency
$\eta$
as a function of plate standing wave ratio
$\mathcal{S}$
. The solid lines with solid markers represent uniform plates, while the dashed lines and empty markers represent tapered plates. Black squares, green diamonds, red triangles and blue circles represent
$\textit{Re}$
values 100, 500, 1000 and 2000, respectively.
In figure 4, we plot dimensionless thrust
$\mathcal{F}$
, power
$\mathcal{P}$
, and efficiency
$\eta$
of tapered and uniform plates as functions of the frequency ratio
$r$
. As shown in figure 4(a), for both plate types,
$\mathcal{F}$
peaks with similar magnitudes at resonance, where the tip displacement amplitude is maximised. For the uniform plate,
$\mathcal{F}$
decreases sharply to nearly zero for post-resonance frequencies. For the tapered plate,
$\mathcal{F}$
exhibits a local minimum at
$r \approx 2$
, and gradually increases for higher
$r$
, saturating for
$r\gt 3.5$
. The dependence of
$\mathcal{F}$
on
$r$
closely correlates with the dependence of
$\delta$
on
$r$
shown in figure 3(a). In the inset in figure 4(a), we show the relationship between
$\mathcal{F}$
and
$\delta$
, indicating a positive correlation between these two parameters consistent with previous studies (Ramananarivo, Godoy-Diana & Thiria Reference Ramananarivo, Godoy-Diana and Thiria2011; Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). We note, however, that there is a noticeable spread in
$\mathcal{F}$
for larger
$\delta$
, indicating that
$\delta$
is not the only parameter defining
$\mathcal{F}$
.
When comparing
$\mathcal{F}$
for plates with different Reynolds numbers, we find that higher
$\textit{Re}$
results in greater
$\mathcal{F}$
. The most significant reduction in
$\mathcal{F}$
occurs when
$\textit{Re}$
is decreased to
$\textit{Re} = 100$
. Since the bending pattern is insensitive to
$\textit{Re}$
, our result suggests that the dependence of
$\mathcal{F}$
on
$\textit{Re}$
can be related to the differences in vortex shedding and dissipation at flows with different
$\textit{Re}$
. We note also that the inset in figure 4(a) indicates that
$\mathcal{F}$
grows faster with increasing
$\delta$
for plates with greater
$\textit{Re}$
.
Figure 4(b) shows the dimensionless power consumed by the oscillating plates
$\mathcal{P}$
as a function of
$r$
. The power is largest around the resonance
$r \lt 1$
, and drastically decreases for post-resonance frequencies. Furthermore, for post-resonance frequencies, the tapered plate requires more power than the uniform plate, whereas at resonance, the power of the plate with uniform thickness exceeds that of the tapered plate. We find that power is insensitive to
$\textit{Re}$
for
$\textit{Re}\geqslant 500$
. At
$\textit{Re}=100$
the power is noticeably higher compared to the plates with larger
$\textit{Re}$
due to a greater contribution of viscous dissipation.
Figure 4(c) shows the hydrodynamic efficiency
$\eta$
defined as the ratio of the dimensionless thrust
$\mathcal{F}$
and power
$\mathcal{P}$
. For the plates with uniform thickness, the efficiency has a maximum for
$r$
slightly above the resonance. For greater
$r$
, the efficiency rapidly decreases then slightly recovers for
$r\gt 2.5$
. For the tapered plates,
$\eta$
nearly monotonically increases with
$r$
, and significantly exceeds the efficiency of uniform plates. The tapered plates operating at higher values of
$\textit{Re}$
yield greater efficiencies, with the efficiency at
$\textit{Re}=2000$
almost triple the efficiency at
$\textit{Re}=100$
for
$r\geqslant 3$
. Overall, for all frequencies tested, the efficiency of tapered plates exceeds that of plates with uniform thickness.
Figure 4(d) shows
$\eta$
as a function of
$\mathcal{S}$
. We find that for the tapered plates,
$\eta$
increases monotonically with
$\mathcal{S}$
. This implies that travelling waves facilitate more efficient undulatory swimming, which is consistent with our previous results (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). However, the increase in
$\eta$
with
$\mathcal{S}$
is less significant for plates with lower
$\textit{Re}$
. This suggests that the gains in efficiency due to travelling wave propulsion are lower in an environment with higher viscosity. We find that the highest efficiency is achieved by tapered plates with travelling wave oscillations corresponding to
$r \approx 4$
, where they produce only slightly weaker thrust than at resonance. We notice that for
$\textit{Re}=100$
, the efficiency exhibits only a minor increase with
$\mathcal{S}$
that saturates for the larger
$\mathcal{S}$
. For the uniform plate dominated by standing waves, the efficiency is not well correlated with
$\mathcal{S}$
for the tested values of
$r$
. This suggests that higher efficiency can also be achieved without relying on travelling waves.
4.3. Flow patterns and vorticity
Snapshots of
$\mathcal{Q}$
-criterion contours
$(\mathcal{Q}\tau ^2 = 5)$
coloured by the
$y$
-component of vorticity: (a–c) tapered plates at
$\textit{Re} = 100$
, (d–f) tapered plates at
$\textit{Re} = 2000$
, and (g–i) uniform plates at
$\textit{Re} = 2000$
. Plates are actuated at (a,d,g)
$r \approx 1$
, (b,e,h)
$r \approx 2$
and (c, f,i)
$r \approx 3$
. Snapshots are taken at
$t/\tau = 0$
. See also supplementary movie 2.
Figure 5 shows the contours of the
$\mathcal{Q}$
-criterion for uniform and exponential plates, coloured by the
$y$
-component of vorticity (see supplementary movie 2). The vorticity produced by the oscillating elastic plates is the most significant at the trailing edge and side edges of the plate, whereas leading edge vortices are relatively weak. This is to be expected, as the plate elasticity causes the trailing edge to undergo greater displacement than the leading edge. Furthermore, the vorticity is significantly reduced at lower
$\textit{Re}$
, where it is rapidly dissipated by the fluid.
Uniform plates exhibit high vorticity at resonance, but as the frequency ratio increases beyond resonance, the vorticity diminishes due to reduced tip displacement. Tapered plates generate vorticity at resonance comparable to that of uniform plates at
$r \approx 1$
. While the vorticity of tapered plates also decreases off-resonance, it remains substantially higher than that of uniform plates at
$r \approx 2$
and
$r \approx 3$
, due to their greater tip displacements. The magnitude of SEVs is comparable between uniform and tapered plates at higher
$r$
, although tapered plates have a significantly greater tip displacement indicating the effect of bending pattern on the strength of SEVs.
Snapshots of the
$y$
-component of the vorticity field at the plate midsection
$y=0$
: (a,c,e) tapered plates, (b,d, f) uniform plates. Plates are actuated at (a,b)
$r \approx 1$
, (c,d)
$r \approx 2$
and (e, f)
$r \approx 3$
. Snapshots are taken at
$t/\tau = 0.5$
and
$\textit{Re} = 2000$
. See also supplementary movie 3.
Figure 6 shows snapshots of the
$y$
-component of the vorticity field at a plane crossing the midsection of the plate at
$y=0$
. The snapshots compare the vorticity for tapered and uniform plates with different
$r$
. At resonance, both plates generate highly complex vortex structures with multiple differently sized vortices (figures 6
a,b). It is also noted that a slight asymmetry appears in the vortex structures of both plates.
At
$r \approx 2$
and
$3$
, the vorticity level is much lower than at resonance for both geometries. The vortex wake of the tapered plate bifurcates into two separate streams diverging with an acute angle (figures 6
c,e). The bifurcating wake was shown to be related to the longitudinal transport of SEVs by travelling waves propagating along the plate (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). Furthermore, pairs of counter-rotating vortices formed by the bifurcating wake contribute to thrust generation by imparting momentum to the plate. For uniform plates oscillating at post-resonance frequencies, reduced trailing edge displacement leads to a decrease in vortex size, resulting in the formation of a well-organised reverse von Kármán vortex street in the wake (figures 6
d, f).
Snapshots of the
$y$
-component of the vorticity field at the plate midsection
$y=0$
for tapered plates with (a)
$\textit{Re} = 100$
, (b)
$\textit{Re} = 500$
, (c)
$\textit{Re} = 1000$
, and (d)
$\textit{Re} = 2000$
. Snapshots are taken at
$t/\tau = 0.5$
and
$r \approx 3$
. See also supplementary movie 4.
The changes of vortex structures with
$\textit{Re}$
are shown in figure 7 for tapered plates oscillating at post-resonance frequency
$r \approx 3$
. We note that the bending pattern only weakly depends on
$\textit{Re}$
, therefore the changes in vortex structures are due to the differences in viscous dissipation. The figure shows that
$\textit{Re}$
has a dramatic effect on vorticity. At
$\textit{Re} = 100$
, the vortices rapidly dissipate shortly after they detach from the plate (figure 7
a). At
$\textit{Re} = 500$
, vortices persist for longer, forming a reverse von Kármán street downstream of the plate. At
$\textit{Re} = 1000$
, the trailing edge vortices intensify further, and break up into two distinct streams of vortices shed above and below the plate midplane (figure 7
c). The bifurcating wake fully develops at
$\textit{Re} = 2000$
, forming two well-separated vortex streams that diverge downstream of the plate (figure 7
d).
The changes in vortex structures with
$\textit{Re}$
correlate with the varying thrust generated by tapered plates at post-resonance
$r$
, as shown in figure 4(a). Despite exhibiting similar bending patterns, the thrust performance of plates diminishes in a low-
$\textit{Re}$
environment. This reduction can be attributed to faster vortex dissipation and the breakdown of the bifurcating wake, resulting in weaker trailing edge vortex interactions and consequently lower thrust.
4.4. Effects of standing wave ratio on propulsion
(a) Power
$\mathcal{P}$
(scaled with
$\textit{Re}$
) as a function of the total period-averaged enstrophy
$\mathcal{E}$
. (b) Total period-averaged enstrophy as a function of
$\mathcal{S}$
. Solid markers represent uniform plates, while empty markers represent tapered plates. Black squares, green diamonds, red triangles and blue circles represent
$\textit{Re}$
values 100, 500, 1000 and 2000, respectively.
(a) Normalised maximum bending area
$\mathcal{A}$
as a function of
$\mathcal{S}$
, with symbols coloured by mean power
$\mathcal{P}$
. (b) Normalised maximum tip displacement
$\delta$
as a function of
$\mathcal{S}$
, with symbols coloured by mean thrust
$\mathcal{F}$
. The solid markers represent uniform plates, while the empty markers represent tapered plates. Squares, diamonds, triangles and circles represent
$\textit{Re}$
values 100, 500, 1000 and 2000, respectively.
To examine the dependence of the power consumed by the plate on the generated vorticity, we plot in figure 8(a) power
$\mathcal{P}$
scaled with
$\textit{Re}$
as a function of dimensionless average enstrophy
$\mathcal{E}$
. The data for scaled power for uniform and tapered plates collapse into a single curve that linearly increases with
$\mathcal{E}$
. The result indicates that the power is defined by the amount of vorticity generated by the plate, regardless of its geometry or bending pattern. Furthermore, the consumed power scales linearly with
$\textit{Re}$
due to the enhanced viscous dissipation.
Figure 8(b) shows enstrophy
$\mathcal{E}$
as a function of standing wave ratio
$\mathcal{S}$
. We find that enstrophy is maximised at smaller
$\mathcal{S}$
that take place near resonance, and decreases to an approximately constant value when the standing wave ratio increases beyond
$\mathcal{S}\gt 0.3$
. This shows that higher
$\mathcal{S}$
results in bending patterns with travelling waves that generate less vorticity, thereby requiring less power for oscillations. Furthermore, at resonance, tapered plates produce significantly lower
$\mathcal{E}$
compared to uniform plates. It has been shown by Ramananarivo et al. (Reference Ramananarivo, Godoy-Diana and Thiria2011) that the deflection angle at the plate trailing edge can substantially influence hydrodynamic performance. We speculate that the high flexibility of the trailing edge in tapered plates results in a more favourable deflection angle, reducing vorticity relative to uniform plates, even though both plate types exhibit similar tip displacements.
Previous studies have shown that power consumption is proportional to the normalised maximum bending area of the plate
$\mathcal{A}$
(Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). Figure 9(a) shows
$\mathcal{A}$
as a function of
$\mathcal{S}$
, with markers coloured by
$\mathcal{P}$
. We find that indeed,
$\mathcal{P}$
tends to decrease with decreasing
$\mathcal{A}$
. Furthermore, we find that oscillations characterised by greater
$\mathcal{S}$
result in lower
$\mathcal{A}$
. Thus plate tapering leading to greater
$\mathcal{S}$
is beneficial for reducing power. Furthermore, we find that for tapered and uniform plates,
$\mathcal{A}$
is nearly independent of
$\textit{Re}$
, and therefore is defined by the thickness profile of the plates.
The reduction in power with decreasing
$\mathcal{A}$
can be attributed to the effect of bending patterns on the formation of SEVs, which arise as the oscillating plate displaces fluid during its plunging motion, causing the fluid to flow around the plate side edges and form vortices. When the leading and trailing edges of the plate move out of phase, as illustrated in figure 2 for
$r \approx 2$
and
$r \approx 3$
, the effective area
$\mathcal{A}$
is reduced. In this scenario, the fluid tends to flow along the plate length rather than around its edges, leading to weaker SEVs compared to resonance oscillations, as shown in figure 5. Since SEVs do not contribute to thrust generation, their reduced strength lowers the power required to oscillate the plate, thereby improving propulsion efficiency.
We note that the uniform plate can generate bending patterns with low
$\mathcal{A}$
at relatively small
$\mathcal{S}$
resulting in weaker SEVs compared to tapered plates with the same
$r$
(cf. figures 5(e, f) and 5(h,i)). Although such oscillations of uniform plates require lower power input
$\mathcal{P}$
, they are characterised by small
$\delta$
(figure 3
a), and thus generate weak thrust
$\mathcal{F}$
(figure 4
a). This, in turn, yields suboptimal hydrodynamic efficiency of uniform at post-resonance frequencies (figure 4
c).
Figure 9(b) shows the maximum tip displacement
$\delta$
as a function of
$\mathcal{S}$
, with markers coloured by
$\mathcal{F}$
. We find that
$\delta$
values for tapered plates are maximised for lower and higher values of
$\mathcal{S}$
. In contrast, intermediate values of
$\mathcal{S}$
that represent plate oscillation patterns with a combination of travelling and standing waves correspond to relatively low tip displacement and thrust. At lower
$\mathcal{S}$
, the bending is maximised at resonance, leading to larger
$\mathcal{A}$
, resulting in higher power consumption and inefficient hydrodynamics. Conversely, tapered plates yielding high
$\mathcal{S}$
corresponding to plate oscillations dominated by travelling waves exhibit bending patterns characterised by a combination of great
$\delta$
and low
$\mathcal{A}$
.
Thus we conclude that standing wave oscillations at resonance and travelling wave oscillations emerging at post-resonance frequencies result in high thrust. However, travelling wave oscillations suppress the production of vorticity that does not contribute to thrust generation, leading to a more efficient undulatory swimming.
5. Conclusion
We use three-dimensional computer simulations to probe the hydrodynamic performance of exponentially tapered and non-tapered elastic plates plunging in an incompressible Newtonian fluid at different Reynolds numbers. The tapering is designed to promote travelling waves in the plate bending pattern through the acoustic black hole effect.
Our results show that the ratio between travelling and standing waves, quantified by the standing wave ratio
$\mathcal{S}$
, ranges from nearly zero at resonance, indicating standing wave dominance, to approximately 0.8 at post-resonance frequencies, indicating travelling wave oscillations. Thus by altering actuation frequency, the bending pattern of tapered plates can be varied from standing wave to travelling wave oscillations. In contrast, for the tested frequencies, plates with uniform thickness exhibit oscillations dominated by standing waves. We additionally find that the bending pattern of oscillating elastic plates is largely insensitive to
$\textit{Re}$
. We note that in our simulations, the frequency ratios do not exceed 4.5, which is below the second resonance. It is therefore unclear whether these trends will persist at higher frequency ratios.
Our simulations reveal that thrust is proportional to the plate trailing edge oscillation amplitude. The maximum thrust is generated at resonance by standing wave oscillations, and its magnitude increases with
$\textit{Re}$
. While tapered and uniform plates produce similar thrust and tip displacement at resonance, plate tapering yields noticeably higher hydrodynamic efficiency than plates with uniform thickness due to the lower enstrophy produced.
Travelling wave oscillations post-resonance, while maintaining approximately the same trailing edge amplitude as standing waves at resonance, produce somewhat lower thrust. Although travelling waves generate less thrust, they offer significantly higher hydrodynamic efficiency, especially at higher
$\textit{Re}$
. Compared to standing wave oscillations, travelling waves at
$\textit{Re}=2000$
increase the efficiency by approximately
$50\,\%$
, while at
$\textit{Re}=100$
, the increase is approximately
$20\,\%$
. The enhanced hydrodynamic efficiency of travelling wave propulsion is largely attributed to the lower magnitude of side edge vortices generated by travelling wave oscillations, characterised by a reduced magnitude of the plate maximum bending area, and the increased transport of displaced fluid to the trailing edge by travelling waves. Oscillations with intermediate
$\mathcal{S}$
that mix standing and travelling waves generate less thrust than either standing or travelling waves, although they somewhat improve efficiency compared to standing wave oscillations.
Flow field analysis reveals that standing and travelling wave oscillations produce distinctly different wake structures. Travelling wave oscillations generate a characteristic bifurcating wake, with two vortex jets extending downstream at an acute angle from the plate plane. This wake bifurcation is attributed to the longitudinal transport of side edge vortices by the travelling waves, which enhances vorticity near the trailing edge (Demirer et al. Reference Demirer, Oshinowo, Erturk and Alexeev2022). The bifurcating wake becomes more pronounced at higher Reynolds numbers. Specifically, it forms at
$\textit{Re} = 1000$
and fully develops at
$\textit{Re} = 2000$
. Pairs of oppositely rotating vortices in bifurcating wake assist the thrust generation. At lower
$\textit{Re}$
, viscous dissipation suppresses flow vorticity, preventing the development of bifurcating wakes, which can be related to the decrease in thrust for plates at lower
$\textit{Re}$
.
When applied to fish locomotion, our findings support the notion that oscillatory swimming, such as thunniform and ostraciiform modes that rely on standing wave oscillations, is advantageous for generating greater thrust and faster swimming (Sfakiotakis, Lane & Davies Reference Sfakiotakis, Lane and Davies1999; Tytell et al. Reference Tytell, Borazjani, Sotiropoulos, Baker, Anderson and Lauder2010). However, this comes at the cost of reduced hydrodynamic efficiency. In contrast, travelling wave locomotion, characteristic of anguilliform swimming, produces somewhat weaker thrust but offers significantly greater efficiency. Regardless of the undulatory swimming mode, tapering is shown to improve overall hydrodynamic performance. We note, however, that our results are constrained to rectangular plates with a specific tapering geometry and a wavelength approximately twice the plate length. Additionally, the frequency ratios tested do not include the second resonance of the plates and beyond. It therefore remains an open question regarding the optimum plate shape and tapering geometry that lead to the best hydrodynamic performance, as well as how the hydrodynamics of the plates behaves at second-resonance oscillations and beyond.
Our findings extend beyond the context of fish locomotion, as oscillating elastic plates can be harnessed in engineering applications for bio-inspired undulatory propulsion in fluids. This study demonstrates that thickness tapering offers a simple, cost-effective and efficient means for drastically improving the hydrodynamic performance of biomimetic swimming robots. Furthermore, the results may be useful to guide the design of efficient energy harvesting devices by enhancing their ability to capture energy from ambient fluid flows.
Supplementary movies
Supplementary movies are available at https://doi.org/10.1017/jfm.2026.11157.
Funding
We thank the National Science Foundation for financial support (CBET 2217647).
Declaration of interests
The authors report no conflict of interest.
Open access
