Connections between resonance and nonlinearity in swimming performance of a flexible heaving plate

TL;DR

This study clarifies how resonance influences the swimming performance of a flexible plate in viscous fluid, revealing the interplay of linear resonant effects and nonlinear behaviors at different amplitudes.

Contribution

It provides an unambiguous definition of resonance in viscous fluid-structure systems and distinguishes linear resonant effects from nonlinear effects in finite-amplitude swimming.

Findings

Resonance peaks occur near the linear stability frequency.

Non-resonant mechanisms become significant at large amplitudes.

Leading-edge separation affects performance at high amplitudes.

Abstract

We investigate the role of resonance in finite-amplitude swimming of a flexible flat plate in a viscous fluid. The role of resonance in performance remains unclear for two reasons: i) a lack of definition of resonance for the fully-coupled fluid-structure interaction system in a viscous flow, and ii) the presence of nonlinear effects, which makes it difficult to disentangle resonant and non-resonant mechanisms in finite-amplitude swimming. We address point i) and provide an unambiguous definition for system resonance by computing global linear stability modes of the fully-coupled fluid-structure interaction system that account for the viscous fluid, the plate, and the coupling between them. We then resolve point ii) by considering high-fidelity nonlinear simulations of systematically increased amplitude. By comparing the results for different amplitudes with one another and with the linear stability modes, we separate linear and/or resonant effects from nonlinear and/or non-resonant effects. Resonant behavior is observed over a wide range of plate stiffnesses, with peaks in trailing-edge motion and thrust occurring near the resonant frequency defined by the global linear analysis. The peaks broaden and weaken with increasing heave amplitude, consistent with an increased damping effect from the fluid. At the same time, non-resonant mechanisms are present at large heave amplitudes. The input power exhibits qualitatively different dynamics at large heave amplitudes compared to smaller heave amplitudes, where resonance dominates. Moreover, leading-edge separation is present for stiff plates at large heave amplitudes, which can drastically alter the performance characteristics from what one would expect through linear predictions.

Tables

Table 1. Table 1 : Range of parameter values considered for the nonlinear simulations. The frequency is varied in increments of 0.1.

$Re$	$M$	$S$	$h_0$	$f$
240	0.01	{0.02, 0.2, 2, 20}	{0.001, 0.01, 0.1}	1–3.2

Table 2. Table 2 : Definition of the oscillation amplitude, a 𝑎 a , for the quantities h T ​ E subscript ℎ 𝑇 𝐸 h_{TE} , C T subscript 𝐶 𝑇 C_{T} , C P subscript 𝐶 𝑃 C_{P} .

		$h_{TE}$	$C_T$	$C_P$
$a$		$\frac{\max \left(h_{TE}\right)}{h_0}$	$\frac{\max \left(C_T-\overline{C_T}\right)}{h_0^2}$	$\frac{\max \left(C_P-\overline{C_P}\right)}{h_0^2}$

Table 3. Table 3 : Simulation parameters for the grid convergence study. Grid 4 was used to obtain the data presented in the article.

Grid	$\Delta x$	$\Delta t$	Smallest sub-domain size	Total domain size	$\frac{{\Vert 𝝌-{𝝌}_f\Vert }_{\infty }}{{\Vert {𝝌}_f\Vert }_{\infty }}$
1	0.0083	0.0004	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	0.00337
2	0.0083	0.0002	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	0.00474
3	0.0056	0.0001	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	0.00182
4*	0.0042	0.00015	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	0.00070
5	0.0042	0.00015	$[-0.2,1.8]\times [0.5,0.5]$	$[-30.98,32.75]\times [-15.78,15.78]$	0.00019
6	0.0042	0.000075	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	0.00058
7	0.0033	0.0001	$[-0.2,1.8]\times [0.5,0.5]$	$[-15.07,16.04]\times [-7.87,7.87]$	—

Equations

$$Re=\frac{{\rho }_{f}UL}{\mu },M=\frac{{\rho }_{s}d}{{\rho }_{f}L},S=\frac{EI}{{\rho }_f{U}^2{L}^3},h_0=\frac{A}{L},f=\frac{f^{\ast }L}{U}$$

$$h_{LE}(t)=h_0\sin (\omega t)$$

$$C_T(t)=-\frac{F_x}{\frac{1}{2}{\rho }_f{U}^{2}L},C_L(t)=\frac{F_y}{\frac{1}{2}{\rho }_f{U}^{2}L},C_P(t)=-C_L(t)v_{LE}(t)$$

$$\mathbf{M}\dot{\mathbf{y}}=\mathbf{r}(\mathbf{y})$$

$${\mathbf{M}\dot{\mathbf{y}}}^p={\mathbf{Ay}}^p$$

$$\alpha (t)=\{\begin{array}{ll}0& t\le 5\\ 6{\left(\frac{t-5}{15}\right)}^5-15{\left(\frac{t-5}{15}\right)}^4+10{\left(\frac{t-5}{15}\right)}^3& 5\le t<20\\ 1& t\ge 20\end{array}$$

$$E_r=\frac{\text{plate energy}}{\text{plate energy + flow kinetic energy}}$$

Full text

1 Introduction

The dual aims of understanding fish swimming and developing agile, efficient underwater vehicles have propelled research on the canonical problem of flow past a flapping flexible plate. Experimental (Ramananarivo et al., 2011; Alben et al., 2012; Dewey et al., 2013; Quinn et al., 2014, 2015), computational (Vanella et al., 2009; Hua et al., 2013; Zhu et al., 2014; Zhang et al., 2017), and theoretical work (Alben, 2008; Ramananarivo et al., 2011; Alben et al., 2012; Floryan & Rowley, 2018) on this problem has revealed propulsive benefits to flexibility over a range of Reynolds numbers (and for inviscid flows), plate parameters, and amplitudes and frequencies associated with the plate’s kinematics. Many studies considered actuation frequencies and material parameters for which the plate exhibits ‘first-mode’ flapping; i.e., with the plate shape similar to that of the first mode of a clamped Euler-Bernoulli beam. However, propulsive benefits were shown to persist for parameters that led to ‘higher mode’ shapes (Alben, 2008; Alben et al., 2012; Quinn et al., 2014; Floryan & Rowley, 2018). Because this canonical problem involves the passive motion of deformable structures, a prevailing question is the extent to which performance is connected to structural resonance. This question is the focus of the current work.

There is a lack of consensus in the literature about the role of resonance in propulsive performance: some investigations have found actuation at or near a resonant frequency to provide benefits to thrust and/or efficiency (Alben, 2008; Alben et al., 2012; Dewey et al., 2013; Quinn et al., 2014; Zhu et al., 2014; Floryan & Rowley, 2018), whereas others have observed optimality through off-resonant actuation (Vanella et al., 2009; Ramananarivo et al., 2011; Zhu et al., 2014). Connecting resonance to swimming performance is challenging for a number of reasons. First, the most relevant definition of resonance is not agreed upon in the literature for viscous flows. For an inviscid flow, Michelin & Llewellyn Smith (2009) computed eigenvalues of the fully coupled fluid-structure interaction (FSI) system to provide a meaningful definition of resonance. Such computations have thus far not been performed for viscous flows, and resonance has instead been defined using either the scaling for an Euler-Bernoulli beam in a vacuum (Vanella et al., 2009; Hua et al., 2013) or the plate response to finite-amplitude actuation (Quinn et al., 2014; Zhang et al., 2017). These approaches leave open questions about the role of resonance, as the former neglects fluid effects and the latter provides the possibility that nonlinearities are playing a role in the observed performance peaks.

This latter effect of nonlinearity provides the second challenge to connecting resonance to performance, as the vast majority of studies have been performed at finite amplitudes. A variety of nonlinear mechanisms that impact performance have been identified. Michelin & Llewellyn Smith (2009) used inviscid theory and computation to show that while thrust peaks aligned well with the eigenvalues of the fully-coupled FSI system, input power and efficiency exhibited off-resonant peaks (though this may be due in part to the fact that the authors did not incorporate instantaneously negative values of input power in their calculations). At Reynolds numbers $O(10^{3}-10^{4})$, Ramananarivo et al. (2011) experimentally identified superharmonic resonance for large actuation amplitudes that the authors attributed to nonlinear fluid damping. Zhu et al. (2014) demonstrated through high-fidelity computations at Reynolds numbers of $O(10^{2}-10^{3})$ that wake asymmetry can occur for certain parameters, and observed that this effect was largely deleterious for efficiency. Using experimental data at $O(10^{4})$, Moored et al. (2014) used a local linearization222This study is included in the collection of nonlinear results because the analysis was performed about a mean flow that was generated through geometrical nonlinearities in the structure as well as flow nonlinearities such as flow separation, vortex interactions, etc. of the wake behind the plate to argue that optimal efficiency was caused by selecting parameters that were closest to the wake’s natural frequency. Quinn et al. (2015) showed through experiments at Reynolds numbers of $O(10^{3}-10^{4})$ that sufficiently large swimming motions led to flow separation that adversely affected efficiency.

Given these nonlinear mechanisms, a remaining question is how to identify the role of resonance in performance of finite-amplitude swimming. Floryan & Rowley (2018) clarified the role of resonant behavior in an inviscid setting by comparing eigenvalues of the fully-coupled FSI system to linear (infinitesimal amplitude) dynamics. The authors showed that actuating at the resonant frequency corresponds to peaks in thrust and input power of the linear system. The authors additionally showed that the relative growth in magnitude of these quantities at resonance is comparable, resulting in no gains in efficiency.

The goal of the present work is to address the two aforementioned challenges in the context of a viscous fluid by considering the two-dimensional flow-plate system at a Reynolds number of $240$. First, global linear stability modes of the fully coupled FSI system are investigated to unambiguously identify the natural frequencies of the system. Next, to separate the role of resonance from nonlinearities, high-fidelity nonlinear simulations of increasing amplitude are studied and compared to the linear stability modes. Resonant and non-resonant mechanisms driving performance are characterized and, where appropriate, connections are drawn to previous studies performed at other Reynolds numbers and swimming parameters.

We note that although this study is performed in the context of swimming, nothing precludes the approach and methods employed here from being applied to scenarios relevant to flexible-wing flight, where the physical mechanisms are distinct due to a greater plate inertia with respect to the flow.

2 Problem setup and parameters considered

We consider here uniform flow past a passively deformable plate that is heaved sinusoidally at its leading edge. The fluid is modeled via the two-dimensional Navier-Stokes equations, and the plate is treated as a clamped geometrically nonlinear Euler-Bernoulli beam. The important dimensionless parameters are the Reynolds number ($Re$), dimensionless mass ratio ($M$) and stiffness ratio ($S$), and dimensionless heave amplitude ($h_{0}$) and frequency ($f$), defined as

[TABLE]

where $\rho_{f}$ ($\rho_{s}$) is the 2D fluid (structure) density, $U$ is the freestream fluid velocity, $L$ is the length of the plate, $\mu$ is the dynamic fluid viscosity, $d$ is the plate thickness, $EI$ is the flexural rigidity composed from Young’s elasticity modulus ($E$) and the second moment of area ($I$), $A$ is the dimensional heave amplitude, and $f^{*}$ is the dimensional frequency related to the dimensional angular frequency ($\omega^{*}$) via $\omega^{*}=2\pi f^{*}$. See figure 1 for a schematic of the dimensional variables defining the FSI system.

In terms of these dimensionless variables, the heave motion is written as

[TABLE]

where $\omega=2\pi f$ is the dimensionless angular frequency and $t$ is the dimensionless time defined using the convective time scale $L/U$.

In this study, the Reynolds number and mass ratio are fixed at $Re=240$ and $M=0.01$, respectively, whereas the stiffness, actuation frequency, and heave amplitude are varied. The moderate Reynolds number is representative of swimming of small fish and larvae (Webb, 1988), and simultaneously enables a comparison with higher-Reynolds-number studies and allows for a relatively large sweep of $S$ and $f$. The small mass ratio is a representative value for fish, which are neutrally buoyant and much thinner than they are long. The heave amplitudes are chosen to enable a systematic comparison of resonant and nonlinear effects in performance. Table 1 summarizes the range of parameter values considered in this work.

In this article, the response of the plate will often be summarized using the transverse amplitude of the plate’s trailing-edge, $h_{TE}(t)$. Quantities relevant to swimming performance are the dimensionless thrust, lift, and input power, defined here as

[TABLE]

where $F_{x}$ and $F_{y}$ are the dimensional drag and lift defined as the horizontal and vertical forces integrated over the plate, respectively, and $v_{LE}=dh_{LE}/dt$ is the heave velocity. Performance is measured with respect to the temporal mean of these quantities—written as $\overline{C_{T}}$, $\overline{C_{L}}$, and $\overline{C_{P}}$, respectively—and in terms of some notion of efficiency ($\eta$) defined here as the Froude efficiency: $\eta=\overline{C_{T}}/\overline{C_{P}}$.

3 Numerical methods: global linear stability modes and nonlinear simulations

The nonlinear fully-coupled FSI equations are written symbolically here in a spatially discrete, temporally continuous setting as (see Goza & Colonius (2017) for details)

[TABLE]

where $\bm{y}$ represents the full state of the system (i.e., all flow and structure variables required to characterize the full system dynamics). Equation (4) contains the Navier-Stokes equations modeling the viscous flow, the geometrically nonlinear structural equations, and the no-slip boundary condition that couples the fluid and structure. Note the distinction between the matrix $\bm{M}$ and the mass ratio $M$. The matrix $\bm{M}$ is singular, which reflects the differential-algebraic nature of the governing equations (due to the no-slip condition on the plate).

The global linear system is defined from the nonlinear system (4) by linearizing about a base state, $\bm{y}^{b}$. Writing $\bm{y}=\bm{y}^{b}+\bm{y}^{p}$ in (4), where $\bm{y}^{p}$ is a small perturbation from the base state, and retaining only linear terms in $\bm{y}^{p}$ leads to the linear system

[TABLE]

where $\bm{A=}d\bm{r}/d\bm{y}|_{\bm{y}^{b}}$. The specific form of the various terms are provided in Goza et al. (2018).

The global eigenvalues ($\lambda_{i}$) that give the natural frequency of the fully-coupled system and the corresponding modes ($\bm{\phi}_{i}$) satisfy the generalized eigenvalue problem $\lambda_{i}\bm{M\phi}_{i}=\bm{A\phi}_{i}$, $i=1,\dots,n$ (with $n$ defined as the state dimension). These eigenvalues and modes are computed by an implicitly restarted Arnoldi algorithm (Lehoucq et al., 1998), and are converged to within $||\lambda_{i}\bm{M\phi}_{i}-\bm{A\phi}_{i}||_{2}<10^{-8}$ for all results shown below.

The nonlinear simulations are performed by solving (4) using the immersed-boundary algorithm of Goza & Colonius (2017). The method treats the fluid equations using a discrete streamfunction formulation (Colonius & Taira, 2008) and the plate (modeled as a geometrically nonlinear Euler-Bernoulli beam) with a corotational finite element formulation (Criesfield, 1991). The method is strongly coupled; i.e., the nonlinear coupling between the structure and the fluid is enforced at each time instance, resulting in a stable algorithm in the presence of large displacements and rotations. This fluid-structure coupling is enforced by the stresses on the immersed surface, and immersed-boundary methods are well known to produce spurious computations of these stresses. These unphysical stress computations were remedied by Goza et al. (2016), and the techniques described therein are incorporated into the FSI algorithm of Goza & Colonius (2017) to ensure appropriate treatment of the fluid-structure coupling. The FSI solver has been validated on several flapping flag problems for flags in both the conventional configuration (pinned at the leading edge) and the inverted configuration (clamped at the trailing edge) (Goza & Colonius, 2017).

The simulation parameters used in this work are described in appendix A. To facilitate studies of steady swimming, for all nonlinear simulations the heave amplitude was ramped up over a period of twenty convective time units via $h^{ramp}_{LE}(t)=\alpha(t)h_{LE}(t)$, where $\alpha(t)$ is a scaling factor defined by

[TABLE]

The zero amplitude used for $0\leq t\leq 5$ is to allow the impulsive start from the fluid to decay. A shifted and scaled smoothstep function (Ebert et al., 2003) is used over $5\leq t\leq 20$ because of its zero first and second derivative at $t=5,20$. Analysis from the nonlinear simulations was performed after $t=20$, using data from a minimum of 10 heave periods.

4 Resonance for viscous plates: global modes of the FSI system

To provide a framework by which to define resonance in the simulations, we first present the eigenvalues of the fully coupled plate-fluid system in figure 2 for several different stiffnesses. The eigenvalues are color-coded by an energy ratio, $E_{r}$, defined here as

[TABLE]

In (7), the plate energy is defined as the sum of plate’s bending strain energy and kinetic energy. All energies are computed from the eigenvector associated with the eigenvalue being considered. Note that because the denominator of the energy ratio is the total system energy, $E_{r}\leq 1$.

Also provided for each stiffness in figure 2 is the eigenvector associated with the largest value of $E_{r}$. The eigenvector is plotted using the fact that the time response associated with a given eigenvalue-eigenvector pair $(\lambda_{j},\bm{\phi}_{j})$ is given by the real part of $\bm{\phi}_{j}e^{\lambda_{j}t}$. The inserts contain plots of this response at the time $t$ for which the trailing-edge displacement is maximal. Note that the there is no amplitude scale in the inserts, as the amplitude is irrelevant in this linear setting.

To explain the plots in figure 2, we first consider the specific stiffness $S=20$ in figure 2. The top-left plot of the figure shows that there is one complex eigenvalue pair associated with a large value of $E_{r}$333There is also a zero frequency mode of non-negligible energy that coincides with divergence of the plate in a given transverse direction, depending on the nature of the perturbation. This mode is not expected to be significant in swimming, which inherently involves flapping behavior of non-zero frequency., and the insert demonstrates that the associated eigenvector has a plate shape reflective of the first mode of an Euler-Bernoulli beam. There is minimal response in the fluid for this eigenvector (also shown in the insert). Despite the similarity of plate shape to the first mode of an Euler-Bernoulli beam, note that the frequency is smaller by two orders of magnitude than for that of an Euler-Bernoulli beam in a vacuum. This discrepancy attests to the significant role of added mass and viscosity in defining resonance for plates of low $M$ values. In addition to this structure-driven mode, there are several flow-driven modes that correspond to a small value of $E_{r}$. Indeed, these modes are similar to those associated with flow past a rigid, stationary body; cf., figure 3.

As the stiffness is decreased from $S=20$, the structure-driven mode decreases in frequency and there is an increasingly pronounced wake response as the eigenvalue approaches those of the flow-driven modes. For $S=2$ and $S=0.5$, the proximity of the structure- and flow-driven modes leads to multiple eigenvectors with nontrivial value of $E_{r}$. By $S=0.2$, the mode reminiscent of the first mode of an Euler-Bernoulli beam is associated with a much smaller value of $E_{r}$. Moreover, a new structure-driven mode appears that is similar to the second mode of an Euler-Bernoulli beam (though again, the natural frequency is different from that of an Euler-Bernoulli beam in a vacuum by two orders of magnitude). As occurred previously, the mode is initially associated with isolated wake structures, though as the stiffness is decreased the eigenvalue approaches those of the flow-driven modes and the mode’s wake becomes progressively expansive. Finally, for $S=0.02$, there is no signature of a mode reminiscent of the second mode of an Euler-Bernoulli beam with a large value of $E_{r}$. Instead, a new structure-driven mode appears that resembles the third mode of an Euler-Bernoulli beam and has a wake structure localized near the plate. One expects that the pattern would persist for lower stiffnesses, with the eigenvalue approaching those of the flow-driven modes (corresponding to a more extensive wake signature), and eventually disappearing and giving way to a higher-order structure-driven mode444Though for very low stiffnesses, one would expect flutter modes to appear and affect the trends described in the main text. The appearance of these flutter modes was, for example, seen by Floryan & Rowley (2018) for inviscid flows. .

The remainder of this article will investigate the extent to which this resonant behavior plays a role in finite-amplitude swimming. The nonlinear results will show that for a given stiffness, the ‘resonant response’ that leads to performance peaks coincides with actuating at a frequency near the natural frequency associated with the mode for which $E_{r}$ is largest. Actuating at frequencies corresponding to the other eigenvalues is demonstrated to have a negligible effect on performance. At the same time, the results will show that there are parameters for which nonlinear mechanisms obscure this resonant behavior. An example of such a nonlinear process is the vortex shedding and interaction that is a key feature of finite-amplitude swimming. This effect is notably absent in several of the modes depicted in figure 2, where vorticity is isolated near the plate.

5 Connections between nonlinearity and resonance

We now compare results from nonlinear high-fidelity simulations for $S=0.02,0.2,2,20$ to the resonant behavior provided by the linear stability modes described in section 4. Results are first presented at relatively large amplitudes relevant to fish swimming, $h_{0}=0.1$. Connections between this large-amplitude behavior and resonance are then drawn by considering nonlinear simulations of increasing amplitude.

5.1 General observations for large-amplitude motions

Figure 4 provides various quantities related to performance. Some features of these performance plots appear reflective of resonant behavior. For example, for every stiffness except the stiffest case ($S=20$), there are peaks in trailing-edge amplitude, mean thrust and efficiency at roughly the natural frequency of the global linear system (though note that in the case of mean thrust, the peak near the resonant frequency is only a local maximum for $S=2$).

However, there are also differences between the nonlinear and linear results that suggest the presence of non-resonant mechanisms. First, even when there are peaks in trailing-edge amplitude, mean thrust, and efficiency for a given stiffness, the peaks occur at slightly different frequencies for each quantity. For example, for $S=2$ the maximum trailing-edge amplitude, mean thrust, and efficiency occur at $f=1,1.3,$ and $1.1$, respectively. Second, the swimming performance for $S=20$ is entirely different from what one would expect from resonance-based arguments; e.g., there are broad peaks in maximum trailing-edge amplitude and efficiency, with maxima occurring at $f=2.4$ and $f=2.1$, respectively (well below the resonant frequency of $f\approx 3.1$). Third, for $S=2$, there is an increase in maximum trailing-edge amplitude and mean thrust for $f\geq 2.4$ despite the fact that there are no eigenvalues corresponding to this frequency range (or to a mode at higher frequency up to $f=4$) for this stiffness.

Connections to and distinctions from resonant behavior are also suggested by the snapshots of large-amplitude swimming presented in figure 5. One connection to the resonant predictions is that near the resonant frequency, the plate shape in the nonlinear simulations resembles that of the global mode. However, the vortex dynamics in the nonlinear simulations are distinct from those observed in the global linear analysis: the finite-amplitude motion of the plate generates vortices at the leading and trailing edges, and these vortices subsequently interact in a manner not seen in the linear setting.

These observations motivate several questions about the relevance of resonance in large amplitude swimming. First, for stiffnesses where there are performance peaks near the natural frequency of the global linear system, how similar are the large-amplitude dynamics to the linear stability modes? That is, to what extent is resonance responsible for optimal performance at finite amplitudes? Second, what dynamics are associated with the boosts in trailing-edge amplitude, mean thrust, and mean input power at high frequency for $S=2$? Finally, for $S=20$ where there is a noticeable departure from resonance-based performance predictions, what dynamics are associated with disappearance in performance peaks? We probe these questions in the remainder of this section by considering nonlinear simulations of increasing heave amplitude, $h_{0}=0.001,0.01,$ and $0.1$.

5.2 Probing the effect of increasing nonlinearity

In this section, we consider stiffnesses $S=0.02,0.2,2,$ and $20$ separately. For each stiffness, we investigate the effect of increasing the heave amplitude from $h_{0}=0.001$, where the linear analysis of section 4 is relevant, to $h_{0}=0.1$, which better represents fish swimming. Through this systematic variation in amplitude, we identify connections to and distinctions from resonance in swimming performance at finite heave amplitudes.

To give context for the results below, we note that for the majority of parameters considered, the plate-fluid system exhibits limit-cycle dynamics (exceptions to this will be discussed further below). Figure 6 depicts time traces of various quantities for $S=0.2$ and three different frequencies, which are qualitatively representative of the typical behavior of the system. To enable a meaningful comparison across frequencies, time is shifted so that $t=0$ coincides with the instance when the leading-edge displacement has maximal derivative (i.e., velocity), and scaled by $T$ so that the leading-edge motion has a period of unity. Note that the thrust and input power have twice the frequency of the trailing and leading-edge displacement. This fact is intuitive, since one would expect a peak in thrust or input power for both the top and bottom portions of a flapping cycle.

In the shifted and scaled representation of time, any of the quantities in figure 6 can be written as $aw((t+\varphi)/T)$, where $a$ is the peak amplitude of the function, $w(t/T)$ is a periodic function (with a period of unity for $h_{TE}$ and $h_{LE}$ and of one half for $C_{T}$ and $C_{P}$) constructed so that $w(0)$ corresponds to the signal having a maximal derivative value, and $\varphi$ is the phase shift to account for the fact that the signal does not, in general, exhibit its maximal derivative at $t/T=0$. The expression is left in terms of a generic periodic function $w(t/T)$, since the signals are not generally sinusoids for every frequency and stiffness considered (due to a nonzero offset and the possible presence of additional harmonics). The convention of $w(0)$ coinciding with the time instant for which the derivative $dw/dt$ is maximal is used to enable an intuitive interpretation of $\varphi$: $\varphi=0$ corresponds to motion that is in phase with the prescribed leading-edge motion characterized by $h_{TE}(t)$. The prevalence of the periodic behavior in this FSI system forms the basis for much of the analysis below.

5.2.1 Minimal stiffness: $S=0.02$

Here we study the minimally stiff case ($S=0.02$) in more detail. Figure 7 shows the phase shift, $\varphi$, and amplitude, $a$, associated with the periodic signals $h_{TE},C_{T},$ and $C_{P}$ at a given frequency. Note that the phase shift exhibits a discontinuous jump at $\varphi=1$ for $h_{TE}$ and $\varphi=0.5$ for $C_{T}$ and $C_{P}$, which is the period (in terms of the scaled time $t/T$) of each of these quantities.

To enable a meaningful comparison across the various heave amplitudes, the amplitude, $a$, of a quantity is scaled as described in table 2. If the dynamics were entirely linear, the scaled amplitudes in table 2 would not change with increasing heave amplitude $h_{0}$. This fact is readily apparent for the trailing-edge displacement, which scales linearly with $h_{0}$ under linear dynamics; the quadratic scaling of $C_{T}$ and $C_{P}$ under linear dynamics is explained in, e.g., Wu (1961); Floryan & Rowley (2018).

Figure 7 demonstrates that resonance plays a significant role in swimming performance for $h_{0}=0.001$, where the swimming behavior is essentially linear. In particular, $f=2.6$ (the frequency nearest to the natural frequency of the mode that maximizes the energy ratio $E_{r}$) is associated with the largest peaks in trailing-edge amplitude, thrust, and input power. It is intuitive from the linear stability analysis of section 4 that actuating near this natural frequency leads to maximal energy put into the plate, resulting in the largest displacements and forces on the plate. We note that there is a local maximum in trailing-edge amplitude at $f=1$, which is near the frequency range associated with the flow-driven modes. Yet, consistent with the fact that these modes correspond to low values of $E_{r}$, these large plate displacements do not coincide with peaks in performance.

It is also clear from figure 7 that the linear dynamics associated with $h_{0}=0.001$ largely persist at the larger heave amplitude of $h_{0}=0.01$. The largest peak in trailing-edge amplitude, thrust, and input power continues to occur at $f=2.6$. Moreover, the dynamics of a given quantity at $h_{0}=0.01$ are of extremely similar phase to those at $h_{0}=0.001$. At the same time, nonlinear effects are also apparent: the peak value of $a$ for $C_{T}$ and $C_{P}$ is smaller than that associated with $h_{0}=0.001$.

For $h_{0}=0.1$, small peaks in the amplitude $a$ persist for $h_{TE}$ and $C_{T}$ near—but no longer at—the natural frequency. The slight change in the frequency corresponding to maximal trailing-edge amplitude and thrust does not preclude the possibility of resonant-type behavior; that is, a modified natural frequency could occur due to finite-amplitude effects of added mass and viscosity. At the same time, the smaller, broader peaks in trailing-edge amplitude and thrust at $h_{0}=0.1$ compared with $h_{0}=0.001,0.01$, as well as the qualitatively different phase behavior at high frequencies for $h_{0}=0.1$, suggest non-resonant mechanisms.

Whereas the trailing-edge displacement and thrust appear to exhibit a mixture of resonant and non-resonant behavior for $h_{0}=0.1$, the input power is markedly different in both $a$ and $\varphi$ at this higher heave amplitude. The phase is roughly constant for $h_{0}=0.1$, and the maximum amplitude of the input power does not exhibit a peak.

To continue to explore the role of resonance, we show in figure 8 the maximum circulation in the wake of the plate. The positive wake circulation, $\gamma^{p}_{wake}$, is defined as the integrated positive vorticity within a rectangle defined by $1\leq x\leq 5/3$ and $-1/2\leq y\leq 1/2$ (where the leading edge of the plate is positioned at $(x,y)=(0,0)$ when $h_{LE}=0$). The maximum circulation, $\max(\gamma^{p}_{wake})$, is defined as the maximum value of $\gamma^{p}_{wake}$ over 5 flapping periods. The analogous plots associated with the negative circulation are nearly identical to the ones provided here for the positive circulation. The figure demonstrates that for all amplitudes, the maximal wake circulation occurs near the resonant frequency associated with the mode of largest $E_{r}$ value. The fact that maximal circulation occurs near the resonant frequency for all heave amplitudes suggests that, at this stiffness, resonance plays a role in the circulation formed in the plate’s wake. The mean thrust is maximal at this same frequency (cf., figure 4), indicating that this connection between resonance and circulation extends to mean swimming performance.

To summarize the results for this stiffness, the role of resonance is less straightforward than what one might conclude from figure 4, even for $S=0.02$, which exhibits peaks near the natural frequency for several mean performance quantities. The thrust dynamics are different at large heave amplitudes from their small-amplitude counterparts, particularly at high frequencies. Moreover, the input power transitions to qualitatively different dynamics at large amplitudes. Together, these facts indicate that non-resonant, nonlinear effects play a nontrivial role in performance, though the persistence of a peak in mean thrust and maximal wake circulation near the resonant frequency suggests that features of resonance may continue to be important in finite-amplitude swimming at this stiffness. Regarding this resonant phenomena, the general trend is that the peaks broaden and shorten as heave amplitude is increased. This could be reflective of a number of phenomena not accounted for by the global linear analysis and small-amplitude simulations, including added-mass effects from finite-amplitude heaving or formation and interaction of vortical structures along the plate’s length and in the plate’s wake. Interestingly, introducing damping into a spring-mass system has a similar effect of broadening and weakening the frequency response of the system, and it is possible that the aforementioned effects could be providing similar damping mechanisms. Probing this possibility in more depth is left to future study.

5.2.2 Moderate stiffness: $S=0.2$

We now consider increasing the stiffness to $S=0.2$. Figure 9 shows the phase and maximum amplitude associated with $h_{TE},C_{T},$ and $C_{P}$ at a given frequency. As was observed for $S=0.02$, the figure suggests the presence of both resonant and non-resonant mechanisms.

For small heave amplitudes ($h_{0}\leq 0.01$), the dynamics appear to be primarily linear: the behavior of $a$ and $\phi$ is similar for all three quantities (though the smaller peaks in amplitude for $h_{0}=0.01$ than for $h_{0}=0.001$ demonstrate the appearance of nonlinear mechanisms). Moreover, resonance appears to play the primary role for these values of $h_{0}$, as there are peaks in all three quantities near the resonant frequency associated with the mode of largest $E_{r}$ value.

For $h_{0}=0.1$, either global or local peaks appear in all quantities near the resonant frequency. Moreover, the phase associated with the thrust is similar to that of the lower heave amplitudes near the resonant frequency. These facts suggest a nontrivial role of resonance in dictating performance of finite-amplitude swimming at this stiffness. At the same time, the peaks are broader and smaller at this larger heave amplitude, and the phase behavior of all quantities reflect significant differences from the smaller-heave-amplitude cases (where resonance dominates). Moreover, while there is a local peak in $a$ for the input power, the overall shape of the curve is different from those associated with $h_{0}=0.001,0.01$.

To establish connections between vortex formation and performance, we provide in figure 10 the maximum wake circulation (see the text surrounding figure 8 for details on how this is defined) as a function of frequency for the various heave amplitudes. As with $S=0.02$, the maximum wake circulation occurs near the resonant frequency for all heave amplitudes. The maximum mean thrust also occurs near this frequency, which again suggests that resonance could be connected to maximal mean thrust production at this stiffness.

5.2.3 Moderate stiffness: $S=2$

Figure 11 provides the phase and maximum amplitude associated with $h_{TE},C_{T},$ and $C_{P}$ at a given frequency for a stiffness $S=2$. The figure indicates that the dynamics are essentially linear for the small heave amplitudes $h_{0}=0.001,0.01$ (both $a$ and $\varphi$ exhibit similar behavior at these smaller values of $h_{0}$), and that resonance continues to play a primary role for these values of $h_{0}$ (resonant peaks are present in all quantities). However, the behavior of $a$ for $C_{P}$ suggests that there are non-resonant mechanisms even in this linear regime: the maximum value of $a$ occurs at $f=3.2$—well away from the natural frequency for this stiffness555The global maximum may occur at a higher frequency, as the amplitude is still increasing at $f=3.2$, though in either case this is away from the natural frequency..

For the largest heave amplitude, $h_{0}=0.1$, there are both similarities and differences to the cases with smaller heave amplitudes. The main similarity is that the thrust and trailing-edge amplitude exhibit either a local or global peak near the resonant frequency corresponding to the mode with largest $E_{r}$ value. This again points to resonant features present in the dynamics exhibited by this stiffness, particularly for low frequencies.

There are also several differences between the large-heave-amplitude dynamics and those of smaller amplitude. First, the phase associated with $h_{0}=0.1$ is different for all quantities to that associated with smaller heave amplitudes. An even more conspicuous difference is the significant variations in $a$ for $h_{TE}$ and $C_{T}$ when $f\geq 2.4$. This variability is due to the fact that the signals are not periodic at this stiffness for these frequencies, rendering this amplitude-phase representation inappropriate.

Figure 12 provides time traces of $h_{TE}$, $C_{T}$, and $C_{P}$ for $f=1.4$ (for which the quantities exhibit periodic behavior), and for $f=2.6$ and $3.2$ (for which the quantities are aperiodic). These aperiodic dynamics are consistent with the pronounced leading-edge vortical structures (and associated separation at the leading edge) and asymmetric wake structures observed in figure 5 for $S=2$ at $f=2.4$ and $3$. These nonlinear, non-resonant mechanisms cause a significant departure from the essentially linear dynamics displayed in figure 11 for small heave amplitudes.

To connect performance to wake circulation, we present in figure 13 the maximum circulation in the wake of the plate (see the text surrounding figure 8 for details on how this is defined) as a function of frequency for the various heave amplitudes. For $h_{0}=0.001,0.01$ the maximum circulation is observed near the resonant frequency associated with the mode for which $E_{r}$ is maximal. For the larger heave amplitude of $h_{0}=0.1$, there is a local peak near this resonant frequency before the maximum circulation increases significantly at higher frequencies.

Together, the results contained in figures 11–13 indicate that for large-amplitude swimming at $S=2$, there is a distinction between low frequencies ($f<2.4$), for which leading-edge separation is sufficiently small, and high frequencies ($f\geq 2.4$), for which leading-edge separation and vortex dynamics are prominent. At low frequencies, the results are analogous to those found for the lower stiffnesses. Many features of linear dynamics persist at large heave amplitudes, and thrust peaks are associated with actuating at the resonant frequency to obtain maximal wake circulation. At high frequencies, the results are substantially different from those associated with the more flexible plates. There are prominent leading-edge vortex dynamics that are associated with significantly larger thrust peaks and wake circulation than what is created by resonant mechanisms at lower frequencies.

5.2.4 Maximal stiffness: $S=20$

We now turn to the stiffest case considered, $S=20$. Figure 14 provides the phase and maximum amplitude associated with $h_{TE},C_{T},$ and $C_{P}$ at a given frequency. Note that for all results shown in this subsection (corresponding to $S=20$), the frequency range is $f\in[1,3.5]$ rather than $f\in[1,3.2]$. This is done to provide enough data points beyond the natural frequency of $f\approx 3.1$ to enable meaningful conclusions about the role of resonance to be drawn.

As was observed for the more flexible cases, resonance appears to be the dominant indicator of performance at the smallest heave amplitude ($h_{0}=0.001$), with peaks in $a$ occurring for all quantities near the natural frequency associated with the linear stability mode of largest $E_{r}$ value. These resonant peaks largely persist at the larger value of $h_{0}=0.01$, and the similar phases across $h_{0}=0.001$ and $h_{0}=0.01$ suggest that the dynamics remain primarily linear up to this larger amplitude.

The trends at the larger heave amplitude are similar to those observed for the other stiffness values. The peaks broaden and weaken significantly compared with the smaller heave amplitudes, but persist near the resonant frequency. In addition, the phase behavior is distinct from that corresponding to the lower heave amplitudes.

We again probe the connection between thrust and wake circulation by showing the maximum wake circulation (see the text surrounding figure 8 for details on how this is defined) as a function of frequency for the various heave amplitudes. For $h_{0}=0.001,0.01$, and $0.1$, the maximum circulation is observed near the resonant frequency associated with the mode for which $E_{r}$ is maximal (though it is difficult to arrive at concrete conclusions about the role of resonance in wake circulation at this stiffness because the resonant frequency is near the maximal frequency considered in this study). Irrespective of the connection to resonance, the frequency associated with maximal wake circulation is also associated with maximal mean thrust (cf., figure 4), as was observed for more flexible plates.

Finally, we note that for $h_{0}=0.1$ at this larger stiffness, the plate-fluid system contains significant leading-edge vorticity that indicates separation; cf., figure 5. This leading-edge vorticity is most pronounced for $S=20$ when $1.8\leq f\leq 2.4$. From figure 4, it is clear that these frequencies coincide with small increases in efficiency. Moreover, unlike for $S=2$, the separated flow behavior is not associated with aperiodic system dynamics for $S=20$, nor does it correspond to substantial changes in performance. This demonstrates that leading-edge separation and vorticity generation, while not accounted for in the linear analysis, do not necessarily indicate significant changes in performance from what one would predict through resonance-based arguments.

6 Conclusions and connections to other work

In this article, we investigated the role of resonance in finite-amplitude swimming of a flexible flat plate in a viscous fluid ($Re=240$). First, we unambiguously defined resonance using linear stability modes of the fully coupled fluid-structure interaction system. These modes represent eigenvalues and eigenvectors that account for the viscous fluid, the plate, and the coupling between them. The eigenvalues associated with modes that maximized the energy put into the plate, quantified in terms of an energy ratio ($E_{r}$), were shown to provide an appropriate measure for system resonance.

We then probed the connections between finite-amplitude swimming and resonant behavior through high-fidelity nonlinear simulations of systematically increased heave amplitude. By comparing the results for different heave amplitudes with one another and with the linear stability modes, connections to and distinctions from resonant behavior were established for finite-amplitude swimming over a wide range of stiffnesses ($S=0.02,0.2,2,20$). The primary observations of this article (with connections to other work provided where appropriate) are:

Small heave amplitudes. For sufficiently small heave amplitudes ($h_{0}\leq 0.01$ in our studies) the dynamics were similar across heave amplitude, suggesting primarily linear behavior. Beyond this, peaks in trailing-edge displacement, thrust, and input power all occurred near the resonant frequency corresponding to the mode with largest $E_{r}$ value, indicating that resonance is a dominant contributor to performance in this small-amplitude regime.

Large heave amplitudes. At larger heave amplitudes more reflective of fish swimming ($h_{0}=0.1$ in our studies), the results indicated that both resonant and non-resonant mechanisms played a role.

Regarding resonant behavior, different plate shapes were observed depending on the plate stiffness and heave frequency, and these shapes were consistent with the shapes predicted by the global linear analysis. More flexible plates were more likely to exhibit shapes that are reminiscent of higher-modes of an Euler-Bernoulli beam. This phenomenon was also observed by Quinn et al. (2014) through experiments at higher Reynolds numbers.

In addition, peaks in trailing-edge motion and thrust persisted near the resonant frequency for all stiffnesses considered, with the general trend that the peaks weakened, broadened, and occurred at slightly shifted frequencies from the resonant frequency identified through the linear stability analysis. This distinction from the linear behavior can be due to a number of phenomena, including added mass effects from finite amplitude heaving, leading-edge separation, and nonlinear vortex generation and interaction not present for small-amplitude kinematics. Interestingly, the broader and weaker frequency response is consistent with the response of a spring-mass system subjected to damping. The aforementioned nonlinear mechanisms could therefore be providing an effective damping mechanism, which in turn could be related to the nonlinear fluid damping effect argued by Ramananarivo et al. (2011) (though the authors considered larger mass ratios and heave amplitudes than those of this study that almost certainly led to differences in the detailed system response). Clarifying the specific nonlinear mechanisms that dictate the broader, shorter peaks associated with finite-amplitude heaving is an avenue for future research efforts.

We note that connections to resonance have been argued for at higher Reynolds numbers as well. Quinn et al. (2014) identified peaks in thrust and efficiency at specific frequencies from experimental data at $Re\sim O(10^{4})$, and used this to claim that resonance was a contributor to performance. Our results clarify that performance peaks in finite-amplitude swimming are indeed connected to resonant behavior at lower Reynolds numbers. This observation indicates that performance peaks observed at higher Reynolds may also be related to resonance of the FSI system, though a global stability analysis at these higher Reynolds is necessary to make conclusive statements.

Despite these peaks near the resonant frequency, non-resonant behavior was observed as well. For all stiffnesses considered, the input power exhibited a qualitatively different amplitude and phase response from the essentially linear dynamics at lower heave amplitudes. In addition, leading-edge separation (and associated vorticity generation) was present for the larger stiffness values, $S=2,20$. These cases were associated with either marginal improvements (for $S=20$, $1.8\leq f\leq 2.4$) or substantial decreases (for $S=2,f\geq 2.4$) in efficiency. This is consistent with the observations of Quinn et al. (2015) that leading-edge separation was generally deleterious to efficiency. Quinn et al. (2015) used experiments at Reynolds numbers of $O(10^{3}-10^{4})$ to arrive at their conclusions, which suggests potential extensions of our conclusions at lower Reynolds numbers to higher Reynolds numbers.

We also note that leading-edge separation corresponded to an asymmetric wake. Zhu et al. (2014) did an extensive study of the parameters that led to asymmetric wake responses, and arrived at the conclusion that the wake circulation had to be above a threshold value for asymmetry to occur. Our studies suggest that sufficiently large leading-edge vorticity is another prerequisite.

Irrespective of the connections to resonance, the maximal mean thrust in large-amplitude swimming occurs at or very near the frequency for which the wake circulation is maximized. This was observed for all stiffnesses considered in this study, and is consistent with the observations of Moored et al. (2014), who used a local linear analysis of the wake to demonstrate that performance gains were obtained by actuating at a frequency that maximally destabilized the wake. Moored et al. (2014) considered experimental data of $O(10^{4})$, which again indicates that there may be extensions of the conclusions drawn here for the moderate Reynolds number of 240 to higher Reynolds numbers.

7 Acknowledgments

The authors gratefully acknowledge funding from the Office of Naval Research (grant # 00014-14-1-0533).

Appendix A Simulation parameters and grid convergence study

The flow equations are treated using a multidomain approach: the finest grid surrounds the body and grids of increasing coarseness are used at progressively larger distances (Colonius & Taira, 2008). In all computations performed for this article, the domain size of the finest sub-domain is $[-0.2,1.8]\times[-0.5,0.5]$ and the total domain size is $[-15.07,16.67]\times[-7.87,7.87]$. The grid spacing on the finest domain is $\Delta x=2/480\approx 0.0042$ and the grid spacing for the flag is $\Delta s=2\Delta x\approx 0.0083$. The time step is $\Delta t=0.00015$, which gives a maximum Courant-Friedrichs-Levy number of approximately $0.2$.

The suitability of these parameters is demonstrated in this appendix for the case when $S=20$, $f=3.2$, and $h_{0}=0.1$. Note that this is a strenuous test case, as this choice of parameters leads to one of the largest vorticity and plate accelerations of the physical parameters considered in this study. Table 3 contains details of the various grids used in this convergence study. The final column of the table provides the quantity $||\bm{\chi}-\bm{\chi}_{f}||_{\infty}/||\bm{\chi}_{f}||_{\infty}$, where $\bm{\chi}$ is the plate’s position computed on the grid associated with the current row at $t=30$ and $\bm{\chi}_{f}$ is the plate’s position computed on grid 7 at $t=30$. The position of the plate is well-converged for all grids considered, which demonstrates that the chosen grid is adequate for the simulations used in this article.