1 Introduction
There has been an intense effort over the past 10 years to understand the effect of fluid elasticity on micro-organism swimming. Experiments, analysis and simulations of low Reynolds number swimming of micro-organisms in complex fluids, in particular viscoelastic fluids, has led to a variety of results – some complimentary, some apparently conflicting – on the effect of fluid elasticity on swimming speed. We know that gait, body stiffness and nonlinear effects matter, but we still do not have a clear understanding of how they interact during locomotion.
Early work quantifying the effect of fluid elasticity on swimming using a linear constitutive law for the fluid and asymptotic analysis of small-amplitude motions showed that elasticity had no effect on swimming speed but increased swimming efficiency (Chaudhury Reference Chaudhury1979; Sturges Reference Sturges1981). However, in Lauga (Reference Lauga2007a ) a full analysis of the classical Taylor swimming sheet for small-amplitude undulatory motion showed that the nonlinearities of the viscoelastic fluid model must be included in a computation of swimming speed, and found that swimming speed is always hindered by fluid elasticity. Similar small-amplitude asymptotic analysis was done for waving filaments and helices (Fu, Powers & Wolgemuth Reference Fu, Powers and Wolgemuth2007; Fu, Wolgemuth & Powers Reference Fu, Wolgemuth and Powers2009) also predicting slow-downs due to fluid elasticity. Riley & Lauga (Reference Riley and Lauga2015) and Elfring & Goyal (Reference Elfring and Goyal2016) have demonstrated the importance of the details of the swimming gait in understanding the effect of fluid elasticity on swimming speed, indeed showing that elastic speed-ups are possible for some gaits, which is further highlighted by the analysis of three-sphere swimmers (Curtis & Gaffney Reference Curtis and Gaffney2013).
Biological swimmers have been shown to change their gait in response to changes in rheology (Shen & Arratia Reference Shen and Arratia2011; Gagnon, Keim & Arratia Reference Gagnon, Keim and Arratia2014; Qin et al. Reference Qin, Gopinath, Yang, Gollub and Arratia2015), making it hard to interpret the mechanisms responsible for observed changes in swimming performance. In more controlled physical models of swimmers in different fluids a variety of results have shown that fluid elasticity can boost swimming speed (Liu, Powers & Breuer Reference Liu, Powers and Breuer2011; Keim, Garcia & Arratia Reference Keim, Garcia and Arratia2012; Espinosa-Garcia et al. Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013) or retard swimming speed (Dasgupta et al. Reference Dasgupta, Liu, Fu, Berhanu, Breuer, Powers and Kudrolli2013; Godínez et al. Reference Godínez, Koens, Montenegro-Johnson, Zenit and Lauga2015). Swimmers with large amplitude motions have been theoretically investigated using numerical simulations, and have added significant information about the response of swimmers to fluid elasticity with a variety of swimming gaits (Balmforth, Coombs & Pachmann Reference Balmforth, Coombs and Pachmann2010; Teran, Fauci & Shelley Reference Teran, Fauci and Shelley2010; Zhu, Lauga & Brandt Reference Zhu, Lauga and Brandt2012; Montenegro-Johnson, Smith & Loghin Reference Montenegro-Johnson, Smith and Loghin2013; Spagnolie, Liu & Powers Reference Spagnolie, Liu and Powers2013; Li, Karimi & Ardekani Reference Li, Karimi and Ardekani2014; Thomases & Guy Reference Thomases and Guy2014; Li & Ardekani Reference Li and Ardekani2015; Salazar, Roma & Ceniceros Reference Salazar, Roma and Ceniceros2016). In addition, for recent reviews of swimming in complex fluids, see Elfring & Lauga (Reference Elfring and Lauga2015) for a theoretical view and see Sznitman & Arratia (Reference Sznitman and Arratia2015) for an experimental view. These many studies have focused on different types of swimmers, in different fluid rheologies and despite the wealth of results we still lack an understanding of the underlying principles of swimming in complex fluids.
To try to isolate physical mechanisms that are significant in a variety of biologically relevant problems, but simple enough to analyse, we focus here on undulatory swimmers in an Oldroyd-B fluid. Even in this more restrictive setting, we nevertheless still find apparently contradictory results and a lack of mechanistic explanations for those differences. Asymptotic analysis of infinitely long, prescribed shape, small-amplitude swimmers has shown that swimming is hindered by the addition of elastic stresses (Lauga Reference Lauga2007a ), although allowing for flexibility can lead to enhancements (Riley & Lauga Reference Riley and Lauga2014). Biological experiments have shown a viscoelastic slow-down for C. elegans (Shen & Arratia Reference Shen and Arratia2011), while simulations of finite-length swimmers with large tail amplitudes (Teran et al. Reference Teran, Fauci and Shelley2010; Thomases & Guy Reference Thomases and Guy2014) give a non-monotonic boost as fluid elasticity is varied. In Thomases & Guy (Reference Thomases and Guy2014) we concluded that shape changes due to body flexibility and fluid elasticity are important, but those results did not explain the results from a physical experiment which showed monotonic speed-ups due to fluid elasticity in swimmers with large tail amplitudes (Espinosa-Garcia et al. Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013). Furthermore, recent numerical simulations (Salazar et al. Reference Salazar, Roma and Ceniceros2016) appear to contradict the speed-ups reported in Teran et al. (Reference Teran, Fauci and Shelley2010) and Thomases & Guy (Reference Thomases and Guy2014).
The relevance of body elasticity in viscoelastic speed enhancements was identified for small-amplitude infinite-length swimmers in Riley & Lauga (Reference Riley and Lauga2014), where the authors attribute the speed enhancements to a viscoelastic ‘suction’ which results in an amplitude boost. However, their analysis does not extend to finite-length large-amplitude swimmers where the role of elasticity-induced shape changes has not been addressed directly. The disparity of the results in Teran et al. (Reference Teran, Fauci and Shelley2010), Espinosa-Garcia et al. (Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013), Thomases & Guy (Reference Thomases and Guy2014) and Salazar et al. (Reference Salazar, Roma and Ceniceros2016), all focusing on large-amplitude, finite-length, undulatory swimmers in Oldroyd-B fluids, indicates that something is missing in our understanding of the problem.
There remains a gap between our understanding from analysis and what we see in computational, biological and physical experiments. Here we combine analysis with numerical simulations of finite-length large-amplitude swimmers to show how fluid elasticity induces shape changes in finite-length flexible swimmers and how those shape changes can lead to speed boosts. We show how shape changes depend on both body stiffness as well as fluid elasticity and analyse the effect that shape changes alone have on swimming speed.
2 Effect of passive body dynamics
2.1 Methodology
We follow the computational framework in Thomases & Guy (Reference Thomases and Guy2014) and Guy & Thomases (Reference Guy and Thomases2015), where the swimmer is modelled as an inextensible flexible sheet of finite-length
$L$
immersed in a two-dimensional (2-D) fluid. We describe the undulatory motion of the swimmer by a curvature of the form
$$\begin{eqnarray}\unicode[STIX]{x1D705}_{0}(s,t)=\left(A_{t}(L-s)/L+A_{h}s/L\right)\sin (2\unicode[STIX]{x03C0}t/T+\unicode[STIX]{x03C0}s),\end{eqnarray}$$
where
$s\in [0,L]$
is the body coordinate. Here
$A_{t}$
is the curvature amplitude at the ‘tail’ (
$s=0$
) and
$A_{h}$
is the curvature amplitude at the ‘head’ (
$s=L$
) of the swimmer.
We use the immersed boundary method to solve for the coupled motion of the fluid and the swimmer (Fauci & Peskin Reference Fauci and Peskin1988). Both inextensibility and shape are imposed (approximately) by forces that are designed to penalize extension and deviation from the prescribed target curvature. These forces are derived from the variation of bending and extension (stretching) energy functionals. For example, the bending energy is
$$\begin{eqnarray}E_{\text{b}}=B/2\int _{\unicode[STIX]{x1D6E4}}\left(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0}\right)^{2}\,\text{d}s,\end{eqnarray}$$
where
$B$
is the bending stiffness,
$\unicode[STIX]{x1D705}$
is the curvature of the sheet,
$\unicode[STIX]{x1D6E4}$
is the swimmer body and
$\unicode[STIX]{x1D705}_{0}$
is the prescribed target curvature. One can interpret the model as an active sheet with bending stiffness
$B$
driven by an active body moment density
$B\unicode[STIX]{x1D705}_{0}$
. We scale forces relative to viscous forces so that for
$B\gg 1$
, the realized shape of the swimmer is very close to the prescribed shape. For
$B\sim 1$
, the realized shape is the result of fluid–structure interaction; i.e. passive body dynamics influences the resulting stroke.
The viscoelastic fluid is described by the Oldroyd-B model at zero Reynolds number (Bird et al. Reference Bird, Hassager, Armstrong and Curtiss1980), regularized by stress diffusion (Sureshkumar & Beris Reference Sureshkumar and Beris1995; Thomases Reference Thomases2011). The system of equations describing the fluid are
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6E5}\boldsymbol{u}-\unicode[STIX]{x1D735}p+\unicode[STIX]{x1D709}\unicode[STIX]{x1D735}\boldsymbol{\cdot }\unicode[STIX]{x1D749}_{\boldsymbol{p}}+\boldsymbol{f}=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D735}\boldsymbol{\cdot }\boldsymbol{u}=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle De\left(\unicode[STIX]{x2202}\unicode[STIX]{x1D749}_{\boldsymbol{p}}/\unicode[STIX]{x2202}t+\boldsymbol{u}\boldsymbol{\cdot }\unicode[STIX]{x1D735}\unicode[STIX]{x1D749}_{\boldsymbol{p}}-\unicode[STIX]{x1D735}\boldsymbol{u}~\unicode[STIX]{x1D749}_{\boldsymbol{p}}-\unicode[STIX]{x1D749}_{\boldsymbol{p}}~\unicode[STIX]{x1D735}\boldsymbol{u}^{\text{T}}\right)+\unicode[STIX]{x1D749}_{\boldsymbol{p}}=\dot{\unicode[STIX]{x1D738}}+De\;\unicode[STIX]{x1D700}\unicode[STIX]{x1D6E5}\unicode[STIX]{x1D749}_{\boldsymbol{p}}, & \displaystyle\end{eqnarray}$$
where
$\boldsymbol{u}$
is the fluid velocity,
$p$
is the pressure,
$\unicode[STIX]{x1D749}_{\boldsymbol{p}}$
is the viscoelastic stress,
$\dot{\unicode[STIX]{x1D738}}$
is the rate of strain tensor and
$\boldsymbol{f}$
is the elastic force density generated by the swimmer. Here
$\unicode[STIX]{x1D709}$
is the polymer to solvent viscosity ratio,
$De=\unicode[STIX]{x1D706}/T,$
the Deborah number is the ratio of elastic relaxation time to stroke period and
$\unicode[STIX]{x1D700}\ll 1$
is the stress diffusion coefficient.
The system is solved in a 2-D periodic domain of size
$[0,2]\times [0,2],$
with
$L=1,\text{d}t=10^{-3}$
and
$\text{d}x=2^{-8}.$
We fix
$\unicode[STIX]{x1D709}=0.5,$
consistent with Teran et al. (Reference Teran, Fauci and Shelley2010), and
$\unicode[STIX]{x1D700}=0.0015$
which provides a regularization to control large stress gradient growth (Thomases Reference Thomases2011). We enforce inextensibility with a dimensionless stiffness constant of
$2500$
.
2.2 Varying body stiffness
To understand the role of body elasticity, we use our simulations to calculate the Stokes-normalized swimming speed while varying
$B$
and
$De$
for a fixed period (
$T=1$
). We use a stroke defined by (2.1) with
$A_{t}=5,$
and
$A_{h}=2.$
This gives a high-amplitude stroke as in Teran et al. (Reference Teran, Fauci and Shelley2010) and Thomases & Guy (Reference Thomases and Guy2014). In figure 1(a) we plot normalized swimming speed as a function of
$De$
for three characteristic stiffness values of
$B=0.1,1.0,10.0$
, which we refer to as very soft, moderately soft and stiff, respectively. For very soft swimmers we see a monotonic boost in swimming speed, with a greater than 50 % boost for high
$De.$
Figure 1. (a) Swimming speed (normalized by Newtonian swimming speed) as a function of
$De$
for different bending stiffness
$B$
. Here we fix the period,
$T=1$
. (b–d) Normalized swimming speed as a function of
$De.$
Adaptations from the literature: (b) from Espinosa-Garcia et al. (Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013), a physical model of a swimmer with a flexible tail. (c) From Teran et al. (Reference Teran, Fauci and Shelley2010) and Thomases & Guy (Reference Thomases and Guy2014), two different numerical simulations for a soft stroke with a large-amplitude tail. (d) From Thomases & Guy (Reference Thomases and Guy2014), numerical simulations for a stiff stroke with a large-amplitude tail.
This response is similar to what was reported in Espinosa-Garcia et al. (Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013) using a physical model of a swimmer with a flexible tail (figure 1
b). For moderately soft swimmers, we see a non-monotonic speed-up, including a smaller speed boost over the Newtonian speed, followed by a slow-down at larger
$De.$
This type of non-monotonic speed-up was first reported in Teran et al. (Reference Teran, Fauci and Shelley2010) and again in Thomases & Guy (Reference Thomases and Guy2014) for a soft stroke with high amplitude (figure 1
c). Finally, for stiff swimmers we see non-monotonic behaviour but no boost over the Newtonian speed, again followed by a slow-down at larger
$De.$
This type of slow-down was reported in Thomases & Guy (Reference Thomases and Guy2014) for a stiff kicker (figure 1
d).
Figure 2. Dimensional swimming speed in a Newtonian fluid for stiff,
$B=10.0,$
and moderately soft,
$B=1.0,$
swimmers over a range of frequencies
$T^{-1}.$
Inset figures show shapes of swimmer over a period at the highest computed frequency for both soft and stiff swimmers.
In contrast to stiff, or rigid, swimmers, the dynamics of flexible swimmers involves an additional time scale. In a viscous fluid, rigid swimmers move with a velocity proportional to the beat frequency (the only time scale in the problem). The problem of a rigid swimmer in a viscoelastic fluid has two time scales, the beat frequency and the relaxation time, whose ratio is the dimensionless relaxation time
$De.$
The swimming speed of soft swimmers depends nonlinearly on the frequency because the shape changes with the frequency. Figure 2 shows the swimming speed in a Newtonian fluid for both a stiff swimmer,
$B=10.0$
, and a moderately soft swimmer,
$B=1.0$
, over a range of beating frequencies. We see a linear response in the case of a stiff swimmer and a nonlinear response for the moderately soft swimmer. Inset swimmer shapes show how the shape changes as a function of the stiffness in the high frequency case.
To illustrate the significance of multiple time scales for flexible swimmers in viscoelastic fluids we compute the Stokes-normalized swimming speed as a function of
$De$
varied two ways: by varying the relaxation time for a fixed period and by varying the period for a fixed relaxation time. Both simulations are performed with the same bending stiffness,
$B=1.0$
, where passive body dynamics is significant. Results are shown in figure 3(a) for a swimmer with the same stroke from figure 1, and the two curves show remarkable qualitative differences. For a rigid swimmer these would give equivalent results. Thus this third time scale, arising from body flexibility, needs to be explicitly included in any discussion of swimming in elastic fluids. A more complete picture of how the swimming speed depends on both the relaxation time and period when the body is soft is shown in figure 3(b). Contours of constant
$De=1-5$
are overlaid in black and the effect of body stiffness is clearly evident as the swimming speed varies significantly along any of the contours. The dashed lines denote the locations of the data in figure 3(a).
Figure 3. (a) Normalized swimming speed as a function of
$De.$
Curves generated by varying only relaxation time (for fixed period
$T=1$
) or stroke period (for fixed relaxation time
$\unicode[STIX]{x1D706}=0.5$
). Body stiffness is fixed:
$B=1.0$
. (b) Normalized swimming speed as a function of both
$\unicode[STIX]{x1D706}$
and
$T,$
with contours overlaid for constant
$De$
values. Dashed lines correspond to the locations of the data in (a).
3 Analysis of shape changes
The effect of body stiffness on swimming kinematics has been previously studied for viscous fluids (Wiggins & Goldstein Reference Wiggins and Goldstein1998; Lauga Reference Lauga2007b ). Shape changes in viscoelastic fluids have been examined (Fu, Wolgemuth & Powers Reference Fu, Wolgemuth and Powers2008), but the relationship between shape changes and swimming speed has not been examined for finite-length swimmers. Here, we review the theory and compare it with numerical simulations.
3.1 Linear theory: Newtonian fluids
We begin by considering small-amplitude displacements of a finite-length thin elastic rod in a Newtonian fluid driven by prescribed curvature,
$\unicode[STIX]{x1D705}_{0}(s,t)$
, (equivalently, prescribed moments) along the body with free ends. The theoretical analysis is similar in two and three dimensions, however, we will make note of the differences when we compare with the numerical simulations in two dimensions. We proceed with the analysis in three dimensions for simplicity. The shape of the rod is determined by the balance between elastic forces and viscous drag. The vertical displacement,
$y(s,t)$
satisfies
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D701}_{\bot }y_{t}=-{\mathcal{B}}(y_{ss}-\unicode[STIX]{x1D705}_{0})_{ss}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle y_{ss}-\unicode[STIX]{x1D705}_{0}=0,(y_{ss}-\unicode[STIX]{x1D705}_{0})_{s}=0,\quad \text{at }s=0,L. & \displaystyle\end{eqnarray}$$
Here
$\unicode[STIX]{x1D701}_{\bot }$
is the perpendicular drag coefficient and
${\mathcal{B}}$
is the bending stiffness of the rod. Note that the use of new notation
${\mathcal{B}},$
is intended to distinguish this (dimensional) bending stiffness we use in the linear theory from our previously defined (non-dimensional) bending stiffness,
$B$
which we use in our numerical simulations.
Non-dimensionalizing (3.1) using the body length
$L$
and the period of the driving force
$T$
results in the dimensionless parameter we call the body response time:
$$\begin{eqnarray}G=\frac{T}{{\mathcal{B}}^{-1}\unicode[STIX]{x1D701}_{\bot }L^{4}}=\frac{\text{period of motion}}{\text{elasto-hydrodynamic beam relaxation time}}.\end{eqnarray}$$
We note that
$G^{-1}$
could be called a body relaxation time. The same non-dimensional group has appeared previously, but has been interpreted differently. In Shelley & Ueda (Reference Shelley and Ueda2000) a quantity similar to
$G$
was considered an ‘effective viscosity’ of growing elastic filaments. The Sperm number (
$Sp=G^{-1/4}$
) is the ratio of the body length to the viscous decay length (Wiggins & Goldstein Reference Wiggins and Goldstein1998; Fu et al.
Reference Fu, Powers and Wolgemuth2007). This interpretation is natural when considering filaments driven at one end rather than along the body as we do here.
We change variables from displacement to curvature deviation
$$\begin{eqnarray}c(s,t)=\unicode[STIX]{x1D705}(s,t)-\unicode[STIX]{x1D705}_{0}(s,t)\end{eqnarray}$$
to facilitate comparing with large-amplitude simulations. For small displacements
$y_{ss}(s,t)\approx \unicode[STIX]{x1D705}(s,t)$
and (3.1) (non-dimensionalized) becomes
$$\begin{eqnarray}\displaystyle & \displaystyle c_{t}=-Gc_{ssss}-\frac{\unicode[STIX]{x2202}\unicode[STIX]{x1D705}_{0}}{\unicode[STIX]{x2202}t}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle c=0,\quad c_{s}=0\quad \text{at }s=0,1. & \displaystyle\end{eqnarray}$$
For a given
$\unicode[STIX]{x1D705}_{0},$
we use an orthogonal function expansion to solve the non-dimensional equations for
$c(s,t)$
. We let the driving curvature be given as
$$\begin{eqnarray}\unicode[STIX]{x1D705}_{0}(s,t)=\mathop{\sum }_{k=1}^{\infty }\unicode[STIX]{x1D6FC}_{k}^{\infty }\text{e}^{2\unicode[STIX]{x03C0}\text{i}\unicode[STIX]{x1D707}_{k}t}\unicode[STIX]{x1D6F9}_{k}(s),\end{eqnarray}$$
and solve the eigenvalue problem,
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D707}\unicode[STIX]{x1D6F9}(s)=-\unicode[STIX]{x1D6F9}_{ssss}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6F9}=0,\quad \unicode[STIX]{x1D6F9}_{s}=0\quad \text{at }s=0,1, & \displaystyle\end{eqnarray}$$
for eigenvalues
$\unicode[STIX]{x1D707}_{k}$
and eigenfunctions
$\unicode[STIX]{x1D6F9}_{k}(s).$
The expansion coefficients of the realized curvature,
$\unicode[STIX]{x1D705}$
, are then
$$\begin{eqnarray}\unicode[STIX]{x1D6FC}_{k}=\unicode[STIX]{x1D6FC}_{k}^{\infty }\left(1-\left(1-\frac{G\unicode[STIX]{x1D707}_{k}}{2\unicode[STIX]{x03C0}\text{i}}\right)^{-1}\right).\end{eqnarray}$$
From this solution we can see that as the rod is stiffened (
$G\rightarrow \infty$
), the resultant curvature tends to the prescribed curvature,
$\unicode[STIX]{x1D6FC}_{k}\rightarrow \unicode[STIX]{x1D6FC}_{k}^{\infty }.$
We also see that for softer rods, i.e. smaller values of the body response time
$G$
, the amplitude of the curvature decreases and there is a phase lag relative to the prescribed shape.
As mentioned above, we use intrinsic coordinates and curvature deviations, to allow us to consider large prescribed curvatures. However we note that (3.5) lacks terms coming from geometric nonlinearities and inextensibility that may not be small when the prescribed curvature is large (Goldstein & Langer Reference Goldstein and Langer1995; Camalet & Jülicher Reference Camalet and Jülicher2000). In §§ 3.3 and 4.2 we compare our simulations to theoretical analysis using (3.5), and in appendix A we show that the influence of the additional terms is in fact small for the amplitudes we consider.
3.2 Linear theory: viscoelastic fluids
We can modify the linear theory for elastic rods to include fluid elasticity. This is similar to what was done in Fulford, Katz & Powell (Reference Fulford, Katz and Powell1998) and Fu et al. (Reference Fu, Wolgemuth and Powers2008). In Fulford et al. (Reference Fulford, Katz and Powell1998) modifications to linear rod theory to include linear viscoelastic fluid effects were presented, and the authors concluded that while fluid elasticity does not change swimming speed, it reduces total work and thus can boost efficiency. However, it was pointed out in Lauga (Reference Lauga2007a ) that it is essential to use a nonlinear elasticity model in these types of calculations because the swimming speed itself is second order in amplitude, where the nonlinear effects are relevant. We note that with these higher-order terms Lauga (Reference Lauga2007a ) shows that swimming speed is always hindered by fluid elasticity for the case considered – infinite-length low-amplitude swimmer with sinusoidal undulations. In Fu et al. (Reference Fu, Wolgemuth and Powers2008) the authors analysed shape changes induced by fluid elasticity in a linearly elastic fluid. Unlike swimming speed, shape changes due to fluid elasticity come in to the asymptotic expansion at first order in amplitude, and hence it is reasonable to use a linearly elastic fluid to look at shape changes. Fu et al. (Reference Fu, Wolgemuth and Powers2008) did not make conclusions about how these shape changes affect swimming speed. Here we perform a similar analysis as in Fu et al. (Reference Fu, Wolgemuth and Powers2008), but by applying the analysis to deviations in curvature we are able to study shape changes in low- and high-amplitude finite-length flexible rods. In § 4 we discuss how these shape changes affect swimming speed.
As in (3.1) we can write a force balance relation between the force on a fluid and from the beam as
$$\begin{eqnarray}f_{fluid}-{\mathcal{B}}(y_{ss}-\unicode[STIX]{x1D705}_{0})_{ss}=0,\end{eqnarray}$$
where the
$f_{fluid}$
represents the normal force on the rod from the viscoelastic fluid. If we define the fluid force to be based on the total deviatoric stress
$\unicode[STIX]{x1D749}=\dot{\unicode[STIX]{x1D738}}+\unicode[STIX]{x1D749}_{\boldsymbol{p}}$
then (upon linearization) using (2.5):
$$\begin{eqnarray}De\;{\dot{f}}_{fluid}+f_{fluid}=(1+\unicode[STIX]{x1D709})f_{vis}+De\;{\dot{f}}_{vis},\end{eqnarray}$$
where
$f_{vis}$
is the viscous drag force. Note that given the form of the system in (2.3)–(2.5), we have assumed a total viscosity of
$1+\unicode[STIX]{x1D709}$
. The swimmer motion is time periodic so we take the Fourier transform in time of (3.12) to solve for viscoelastic modifications to the fluid drag. This yields,
$\hat{f}_{fluid}=\left((1+\unicode[STIX]{x1D709}+2\unicode[STIX]{x03C0}\text{i}De)/(1+2\unicode[STIX]{x03C0}\text{i}De)\right)\hat{f}_{vis}.$
As in the viscous theory, we can solve for modifications to the curvature from body stiffness and use the modifications to the fluid drag to account for the fluid elasticity:
$$\begin{eqnarray}\unicode[STIX]{x1D6FC}_{k}=\unicode[STIX]{x1D6FC}_{k}^{\infty }\left(1-\left(1-\frac{G\unicode[STIX]{x1D707}_{k}}{\unicode[STIX]{x1D701}_{ve}2\unicode[STIX]{x03C0}\text{i}}\right)^{-1}\right),\quad \text{with }\unicode[STIX]{x1D701}_{ve}=\frac{1+\unicode[STIX]{x1D709}+2\unicode[STIX]{x03C0}\text{i}De}{1+2\unicode[STIX]{x03C0}\text{i}De}.\end{eqnarray}$$
The coefficients in (3.13) give an analytical expression for the modifications to the rod shapes relative to the prescribed shapes as fluid and body elasticity are varied.
3.3 Elastic shape changes: theory and numerical comparison
The analysis in the previous sections made use of resistive force theory which relates the drag force and velocity on a long thin cylindrical object. More generally, for small amplitudes the vertical displacement satisfies
$$\begin{eqnarray}y_{t}=MFy,\end{eqnarray}$$
where
$M$
is the mobility operator and
$F$
is the linearized bending force operator. Resistive force theory makes the approximation
$M\approx 1/\unicode[STIX]{x1D701}_{\bot }.$
Our analysis of shapes (3.13) contains the quantity
$G\unicode[STIX]{x1D707}_{k},$
where
$G$
depends on the drag coefficient,
$\unicode[STIX]{x1D701}_{\bot },$
(3.3). To use the more general linear theory in our analysis one can identify
$\unicode[STIX]{x1D707}_{k}/\unicode[STIX]{x1D701}_{\bot }=\unicode[STIX]{x1D707}_{k}^{MF}$
, where
$\unicode[STIX]{x1D707}_{k}^{MF}$
denotes the
$k\text{th}$
eigenvalue of the operator
$MF.$
We relate
$G\unicode[STIX]{x1D707}_{k}$
to the dimensionless bending stiffness,
$B$
used in our simulations, through
$$\begin{eqnarray}G\unicode[STIX]{x1D707}_{k}=\frac{T{\mathcal{B}}}{L^{4}\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D707}_{k}=\left(\frac{T{\mathcal{B}}}{L^{4}}\right)\left(\frac{\unicode[STIX]{x1D707}_{k}}{\unicode[STIX]{x1D701}_{\bot }}\right)=B\unicode[STIX]{x1D707}_{k}^{MF}.\end{eqnarray}$$
To compare the linear analysis with our two-dimensional simulations we numerically approximate (3.14). For small deviations to the vertical displacement,
$M$
is the integral operator which is the convolution of the vertical force with the fundamental solution to Stokes equations. We approximate
$M$
using the method of regularized Stokeslets (Cortez Reference Cortez2001), which is a numerical method based on a regularized Greens function for the Stokes equations. We can also numerically approximate
$\unicode[STIX]{x1D707}_{k}^{F},$
the
$k\text{th}$
eigenvalue of the bending force operator
$F,$
using a second-order finite difference, and we find that with point spacing
$\unicode[STIX]{x0394}s=0.002,$
the eigenvalues of
$F$
are within 1 % of the eigenvalues of the continuous operator. We give the eigenvalues for the first four non-trivial modes in table 1. Note that to compute
$\unicode[STIX]{x1D707}_{k}^{MF}$
we assume a viscosity of one. Also in table 1, we give the first four (mode dependent) drag coefficients computed as
$\unicode[STIX]{x1D701}_{k}=\unicode[STIX]{x1D707}_{k}^{F}/\unicode[STIX]{x1D707}_{k}^{MF}$
.
Table 1. Eigenvalues of the discretized operators
$F$
and
$MF$
using
$\unicode[STIX]{x0394}s=0.002,$
and the effective drag coefficient
$\unicode[STIX]{x1D701}_{k}=\unicode[STIX]{x1D707}_{k}^{F}/\unicode[STIX]{x1D707}_{k}^{MF}$
for the first four non-trivial modes.
In order to compare the predicted shape changes given by (3.13) with our numerical simulations we prescribe a curvature of the form
$$\begin{eqnarray}\unicode[STIX]{x1D705}_{0}(s,t)=A\sin (2\unicode[STIX]{x03C0}t),\end{eqnarray}$$
in our model equations (2.3)–(2.1). The prescribed standing wave of constant curvature corresponds to a motion through circular arcs with peak curvature
$A$
. By symmetry, this motion does not result in any horizontal translation of the body. We refer to these non-translating ‘swimmers’ as flexors. We consider both low- and high-amplitude curvatures,
$A=0.5$
and
$A=4.0.$
The shapes are shown inset in figure 4(a).
Figure 4. (a) Normalized amplitude of the first mode,
$|\unicode[STIX]{x1D6FC}_{1}|/|\unicode[STIX]{x1D6FC}_{1}^{\infty }|,$
from (3.13), for the flexor at low
$A=0.5$
and high
$A=4.0$
amplitude over a range of bending stiffness
${\mathcal{B}}.$
Linear theory is shown in solid lines and simulation data are indicated by markers. (b,c) Data from (a) renormalized by Newtonian (
$De=0$
) data, as a function of both
${\mathcal{B}}$
and
$De.$
In figure 4 we plot the theoretical predictions from (3.13) (solid lines) along with values computed from numerical simulations; low amplitude (
$A=0.5$
) are indicated by hollow markers, and high amplitude (
$A=4.0$
) are indicated with filled markers. In figure 4(a) we plot the normalized amplitude of the first mode (
$|\unicode[STIX]{x1D6FC}_{1}|/|\unicode[STIX]{x1D6FC}_{1}^{\infty }|$
) to see how the amplitude deviates from the prescribed amplitude as a function of bending stiffness
$B$
. We see that generically the amplitude of the flexor decreases as the flexor is softened for fixed
$De$
. For sufficiently soft flexors (
$B\lesssim 1$
) viscoelasticity increases the amplitude monotonically with
$De$
, but for stiffer swimmers the amplitude changes non-monotonically with fluid elasticity.
In figure 4(b,c) we renormalize the data by the amplitude in a viscous fluid to see the effects of viscoelasticity more clearly. Again we see that fluid elasticity can increase the amplitude significantly for a soft flexor, but that effect is lost as the flexor is stiffened. When we plot the amplitude as a function of
$De$
for the very soft, moderately soft and stiff cases we see again that three qualitatively different regimes emerge. For very soft flexors the amplitude is monotonically increased by elasticity, for moderately soft flexors the response is non-monotonic, and can decrease or increase the amplitude, and for stiff flexors there is little change in the amplitude due to fluid elasticity. It is notable that the linear theory does such a good job predicting shape changes for low and high amplitude and for low and high Deborah number. In appendix A we derive the theory for both the limit of small amplitude and the limit of high stiffness. We see in figure 4(a,b) that the largest differences are for moderate stiffnesses at high amplitude. We note that we are showing results only for the first mode. For higher modes the trends are similar but the transition from stiff to soft behaviour occurs at lower values of
$B$
because the eigenvalues
$\unicode[STIX]{x1D707}_{k}$
increase with
$k$
.
4 Analysis of swimming speed
In a viscous fluid, increasing the stroke amplitude will increase the swimming speed, and we can infer from § 3 that soft swimmers in a viscoelastic fluid sometimes obtain an amplitude boost over the corresponding swimmer in a Newtonian fluid. However, when comparing swimmers in a viscoelastic fluid to those in a viscous fluid, even with an amplitude boost the viscoelastic swimmer may not swim faster than the viscous swimmer due to additional fluid elastic forces that the swimmer will encounter. Thus the effect of elasticity-induced shape changes is difficult to decouple from the overall effect of fluid elasticity. Analytical expressions for swimming speed can be obtained in certain limits, or for specialized swimmers, but even in these cases we see that the effect of fluid elasticity depends on many factors. For example infinite-length small-amplitude undulatory swimmers show that a slow-down is generically expected for stiff swimmers in a viscoelastic fluid (Lauga Reference Lauga2007a ), but allowing for body flexibility, shape changes can lead to speed boosts (Riley & Lauga Reference Riley and Lauga2014).
In regimes that are more challenging for analysis such as the large-amplitude, finite-length swimmers considered here, it is more difficult to attribute speed boosts or slow-downs to specific swimmer attributes. For large-amplitude finite-length undulatory swimmers, it was conjectured (Teran et al. Reference Teran, Fauci and Shelley2010) that speed boosts were related to large tail stresses, and in Thomases & Guy (Reference Thomases and Guy2014) stroke asymmetries were correlated with both slow-downs and speed-ups. Here we will compute a stroke-induced swimming speed that isolates the effect of fluid elasticity on shape changes, and how those shape changes affect swimming speed in a Newtonian fluid. We then compare that analysis with the full nonlinear numerical simulations where the effect of shape changes is coupled with the fluid elasticity.
4.1 Swimming speed: two-mode swimmer
To keep the analysis simple we define a gait whose swimming speed in a viscous fluid we can compute analytically. We define a ‘two-mode swimmer’ given by the curvature:
$$\begin{eqnarray}\unicode[STIX]{x1D705}(s,t)=A_{1}\cos (2\unicode[STIX]{x03C0}t/T+\unicode[STIX]{x1D719}_{1})\unicode[STIX]{x1D6F9}_{1}(s)+A_{2}\cos (2\unicode[STIX]{x03C0}t/T+\unicode[STIX]{x1D719}_{2})\unicode[STIX]{x1D6F9}_{2}(s),\end{eqnarray}$$
where the
$\unicode[STIX]{x1D6F9}_{i}(s)$
for
$i=1,2,$
are the first and second bending modes. The modulation of a single mode results in a standing wave and will not translate in a Newtonian fluid. We use a sum of the first two modes with a phase difference to generate a non-reciprocal motion. Shapes of the first, second and sum of the first and second modes are plotted in figure 5 for both low and high amplitudes.
Figure 5. Shapes for first, second and sum of first two modes for ‘two-mode’ swimmer defined by (4.1). Snapshots of one period for low- and high-amplitude strokes.
Using resistive force theory one can derive the (time-averaged) swimming speed for a given small-amplitude motion:
$$\begin{eqnarray}\left\langle U\right\rangle =\left(\frac{\unicode[STIX]{x1D701}_{\bot }}{\unicode[STIX]{x1D701}_{\Vert }}-1\right)\frac{1}{LT}\int _{0}^{\text{T}}\int _{0}^{L}y_{s}y_{t}\,\text{d}s\,\text{d}t,\end{eqnarray}$$
where
$U$
is the swimming speed,
$y(s,t)$
is the vertical displacement of the swimmer and
$\unicode[STIX]{x1D701}_{\bot }$
and
$\unicode[STIX]{x1D701}_{\Vert }$
are the perpendicular and parallel drag coefficients, respectively, (Wiggins & Goldstein Reference Wiggins and Goldstein1998; Lauga Reference Lauga2007b
).
For small amplitudes, the shape of the swimmer (up to translation and rotation) is given by integrating equation (4.1) twice in space to compute the swimming speed via (4.2). The swimming speed (in a viscous fluid) for the two-mode swimmer is proportional to the product of the amplitudes and the sine of the phase difference:
$$\begin{eqnarray}\left\langle U\right\rangle \propto A_{1}A_{2}\sin (\unicode[STIX]{x1D719}_{2}-\unicode[STIX]{x1D719}_{1}).\end{eqnarray}$$
With this expression and the theoretical prediction for shape changes, we define a stroke-induced swimming speed which is the swimming speed in a Newtonian fluid that depends on the shape changes due to fluid elasticity and body flexibility. Specifically, for a given
$De$
and
$G,$
we compute
$\unicode[STIX]{x1D6FC}_{k}$
from (3.13) (
$A_{j}=|\unicode[STIX]{x1D6FC}_{j}|,$
and
$\unicode[STIX]{x1D719}_{j}=\arg (\unicode[STIX]{x1D6FC}_{j})$
) and the stroke-induced swimming speed from (4.3). We parametrize the shape changes due to changes in
$De$
using a parameter we call the stroke Deborah number (Thomases & Guy Reference Thomases and Guy2014),
$\text{Stroke}_{De}.$
In other words,
$\text{Stroke}_{De}$
represents the value of
$De$
used in (3.13) to compute the stroke-induced swimming speed via (4.3). Our analytical expression for the shape changes is based on a linearly elastic fluid, but because the nonlinear elastic effects and swimming speed are both second order in amplitude, we do not expect the stroke-induced swimming speed to capture the true viscoelastic swimming speed. An analytical expression for the swimming speed in a nonlinear viscoelastic fluid, as was computed in Lauga (Reference Lauga2007a
), is not tractable in the finite-length case, because translational invariance, which facilitates the calculation for infinite-length swimmers, is lost.
We plot the stroke-induced swimming speed for the two-mode swimmer over a range of
$\text{Stroke}_{De},$
and as with the flexor, we see the emergence of three regimes dependent on the body stiffness, see figure 6(a). Shape changes boost the stroke-induced swimming speed if the swimmer is very soft, a smaller boost is obtained for the moderately soft swimmer and additionally there is a non-monotonic response to increasing elasticity including a regime where shape changes slow-down the swimmer, and finally if the swimmer is stiff there is a negligible effect.
4.2 Swimming speed: theory and numerical comparison
Figure 6. (a) Theoretically predicted stroke-induced swimming speed (computed using (4.3)) normalized by the Newtonian stroke-induced swimming speed. (b) Stokes-normalized swimming speed in simulations with the two-mode swimmer. (c) Ratio of speed to stroke-induced speed for low- and high-amplitude strokes. (d) Ratio of high- to low-amplitude swimming speed. (Dashed lines are for graphical interpretation.)
We simulate a two-mode swimmer of both low and high amplitude by prescribing a curvature of the form given in (4.1) with
$A_{1}=0.8A$
,
$A_{2}=0.6A$
,
$\unicode[STIX]{x1D719}_{2}-\unicode[STIX]{x1D719}_{1}=\unicode[STIX]{x03C0}/2$
, for
$A=0.5$
(low) and
$4.0$
(high). These values come from projections of the stroke used to generate figure 1. The Stokes-normalized swimming speeds for a very soft, moderately soft and stiff swimmer at both low (hollow markers) and high (filled markers) amplitude are plotted in figure 6(b). We note that the three regimes seen in figure 6(a) still emerge from these simulations, but for this two-mode swimmer the simulation swimming speeds are always slower than the stroke-induced swimming speeds. The ratio of swimming speed to stroke-induced swimming speed is shown in figure 6(c). This quantity can be interpreted as the effect of fluid elasticity that is not related to shape changes. It is notable that these curves collapse onto a single curve for the low-amplitude swimmers at all stiffnesses as well as the high-amplitude swimmer in the very soft regime. This additional elastic fluid effect on swimming speed is likely to be highly stroke dependent.
Figure 7. Amplitudes
$|\unicode[STIX]{x1D6FC}_{1}|,|\unicode[STIX]{x1D6FC}_{2}|$
and
$\sin (\unicode[STIX]{x1D719}_{2}-\unicode[STIX]{x1D719}_{1})$
, for the two-mode swimmer, equation (4.1). Linear theory is shown in solid lines (dependence on
${\mathcal{B}},$
and
$De$
coming from (3.13)) and simulation data are indicated by markers for high amplitude
$A=4.0,$
and
$De=0.5$
(a),
$De=4.0$
(b).
The additional effects of fluid elasticity are fundamentally different for the large-amplitude, large
$De$
regime. In figure 6(d) we plot the ratio of swimming speeds for the high-to-low amplitude strokes, and see that for
$De>1$
(for sufficiently stiff swimmers) a significant difference in swimming speed arises. This difference is not related to shape changes because, like the flexor, the elasticity-induced shape changes predicted by the theory for the two-mode swimmer agree very well with the simulation results, for all
$De,$
at low and high amplitudes; see figure 7. At low amplitude (not shown) the relative error between theoretical and numerical predictions is less than 1 %, and for high amplitude the error at low
$De$
is at most 5
$\%$
and at high
$De$
the error is at most 9 %. These results indicate that the theoretically predicted shape changes and their isolated effects on swimming speed can be well approximated by the analytical results for the amplitudes simulated and the range of
$De$
considered. A mechanistic understanding is lacking to explain what causes the dramatic slow-downs of swimmers in the high-amplitude, high
$De$
regime.
Figure 8. (a) Average elastic–viscous stress ratio in a two-mode swimmer over a range of bending stiffness
${\mathcal{B}}$
for both low- and high-amplitude strokes and low and high
$De.$
(b) Colour field of the magnitude of the elastic–viscous stress ratio around stiff swimmers (
${\mathcal{B}}=10$
) at low and high
$De$
and low and high amplitudes. (c) Tail amplitudes over a range of bending stiffness
${\mathcal{B}}$
, measured as the maximum displacement at the tail over a period, for both low- and high-amplitude strokes and low and high
$De.$
We conjecture that the slow-downs in the high-amplitude, high
$De$
regime must be attributed in part to the large localized stresses that accumulate near the body (Teran et al.
Reference Teran, Fauci and Shelley2010; Thomases & Guy Reference Thomases and Guy2014). To explore this conjecture, in figure 8(a) we plot the average elastic to viscous stress ratio, computed as the time average over one period of
$\Vert \unicode[STIX]{x1D749}_{\boldsymbol{p}}\Vert /\Vert \unicode[STIX]{x1D735}\boldsymbol{u}\Vert$
where
$\Vert \cdot \Vert$
is the Frobenius norm, over a range of body stiffness
$B$
for both low- (hollow markers) and high-amplitude (filled markers) strokes for
$De=0.5,4.0.$
There is a notable transition in the stress ratio in the high
$De$
, high-amplitude swimmer as the body is stiffened, while this stress ratio is flat for both low-amplitude and low
$De$
swimmers. Stiff swimmer shapes along with the elastic–viscous stress ratio are plotted on a log scale in figure 8(b). The low-amplitude strokes are surrounded by elastic stresses that are at least two orders of magnitude smaller than the high-amplitude strokes, but even at large amplitude the low
$De$
swimmer still has relatively low stress near the body. Lastly, we plot the tail amplitude, as one measure of the swimmer stroke, in figure 8(c). We see that for sufficiently soft swimmers the ‘high-amplitude’ stroke has a lower amplitude, which explains why in figure 6 the very soft high-amplitude swimmer behaves like the low-amplitude swimmers. For the high-amplitude strokes it is only in the large-amplitude and large
$De$
regime where significant stress accumulates near the swimmer.
5 Conclusions
In Thomases & Guy (Reference Thomases and Guy2014) we showed that stroke related speed-ups depend on body stiffness, and the analysis from this paper shows explicitly how the stroke changes depend on body stiffness and fluid elasticity through two dimensionless ‘relaxation times’: the fluid relaxation time,
$De$
and the body relaxation time,
$G^{-1}.$
When we look at apparently contradictory results from the literature, we see that calculating
$G$
will determine which regime the swimmer falls into.
In Espinosa-Garcia et al. (Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013) the Sperm number is reported to be between 0.5–2.5, but even with the awareness that these are soft swimmers, the authors ‘conjecture that the effect [due to shape changes] is not significant’. We use their reported parameters (
$L=25$
mm, cross-sectional radius
$a=62.5~\unicode[STIX]{x03BC}\text{m}$
, Young’s modulus
$E=80$
GPa, viscosity
$\unicode[STIX]{x1D707}=2.7$
Pa s. For moment of inertia
$I=\unicode[STIX]{x03C0}a^{4}/4,$
we get
${\mathcal{B}}=EI=9.6\times 10^{-7},$
$\unicode[STIX]{x1D701}=4\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}/\ln (L/a),$
and thus
$G=({\mathcal{B}}T/\unicode[STIX]{x1D701}L^{4})=((9.6\times 10^{-7}\ln (400))/(4\unicode[STIX]{x03C0}\cdot 2.7(25\times 10^{-3})^{4}))\approx 0.43.$
) and a characteristic frequency of
$1~\text{s}^{-1},$
and find
$G\approx 0.43.$
We cannot directly conclude that this value lies in the very soft regime (
$B\approx 0.1$
) due to differences between two and three dimensions as well as the way that the micro-swimmer is driven (it is a flexible tail with a magnetically driven head). In appendix B we repeat the calculation from § 4.1 to compute an equivalent stroke-induced swimming speed for a flexible filament driven at one end. This calculation shows that speed boosts still arise for sufficiently soft swimmers, despite the difference in driving mechanism. The most significant boost in speed from viscoelastic shape changes occurs for
$G\approx 0.1,$
but for
$G\approx 0.43$
it would be reasonable to conclude that the significant speed-ups observed in the experiment are related to shape changes.
In Salazar et al. (Reference Salazar, Roma and Ceniceros2016) the parameter reported for what they consider to be a soft swimmer is
$B=2$
. However their swimmer length is
$L=0.6$
mm (with characteristic length 1 mm) hence an equivalent dimensionless body response time
$G$
must be multiplied by
$L^{-4}\approx 7.7$
. This pushes their ‘soft’ simulations into the stiff regime where there are no speed-ups from shape changes, also agreeing with their results. Furthermore, in Salazar et al. (Reference Salazar, Roma and Ceniceros2016) it is conjectured that stress diffusion, used to regularize the simulations in Teran et al. (Reference Teran, Fauci and Shelley2010) and Thomases & Guy (Reference Thomases and Guy2014), is the source of the speed-ups, but the speed-ups we see are theoretically predicted, and realized in our simulations, even in the low-amplitude regime where no regularization is necessary.
In our analysis we quantify the effect of body and fluid elasticity-induced shape changes on swimming speed. We see that the shape change analysis holds for the amplitudes simulated and the range of
$De$
considered, and in this case we see an additional elastic slow-down that is reminiscent of the type of slow-down predicted by asymptotic analysis of infinite-length small-amplitude undulatory swimmers (Lauga Reference Lauga2007a
). It may be tractable to apply asymptotic analysis (Riley & Lauga Reference Riley and Lauga2014; Elfring & Goyal Reference Elfring and Goyal2016) to determine the form of the elastic slow-down for low-amplitude finite-length swimmers. A fundamentally different regime arises for large-amplitude swimmers in highly elastic fluids. A different approach is needed to understand the mechanisms that cause large localized stresses and their effect on swimming.
Acknowledgements
The authors would like to thank H. Fu and R. Zenit for interesting discussions and suggestions on this work, and M. Shelley for suggesting the term ‘flexors’. The authors would also like to thank the anonymous referees for suggesting useful modifications to our original manuscript. The work of R.D.G. was partially supported by NSF grants DMS-1160438 and DMS-1226386.
Appendix A. Derivation of partial differential equations (PDE) for rod dynamics
In this appendix we give the derivation for the equation of motion for a thin filament in a viscous fluid which includes terms arising from inextensibility and geometric nonlinearities. For more details on similar calculations see for example Goldstein & Langer (Reference Goldstein and Langer1995) and Camalet & Jülicher (Reference Camalet and Jülicher2000).
A.1 Geometric relations
Consider an inextensible thin rod whose centreline position is given by
$\boldsymbol{X}(s,t)$
, where
$s$
is arclength coordinate. We suppose that the deformation of the rod is planar. Let
$\hat{\unicode[STIX]{x1D749}}$
and
$\hat{\boldsymbol{n}}$
be the tangent and normal vectors,
$\unicode[STIX]{x1D713}$
be the tangent angle and
$\unicode[STIX]{x1D705}$
be the curvature. We have the following relationships between these quantities:
$$\begin{eqnarray}\unicode[STIX]{x1D713}_{s}=\unicode[STIX]{x1D705},\quad \hat{\unicode[STIX]{x1D749}}_{s}=\unicode[STIX]{x1D705}\hat{\boldsymbol{n}},\quad \hat{\boldsymbol{n}}_{s}=-\unicode[STIX]{x1D705}\hat{\unicode[STIX]{x1D749}},\quad \hat{\unicode[STIX]{x1D749}}_{t}=\hat{\boldsymbol{n}}\unicode[STIX]{x1D713}_{t},\quad \hat{\boldsymbol{n}}_{t}=-\hat{\unicode[STIX]{x1D749}}\unicode[STIX]{x1D713}_{t}.\end{eqnarray}$$
A.2 Equation of motion
Let
$\boldsymbol{F}(s,t)=F_{\bot }\hat{\boldsymbol{n}}+F_{\Vert }\hat{\unicode[STIX]{x1D749}}$
be the force density along the rod. Using resistive force theory, the motion of the rod is given by
$$\begin{eqnarray}\boldsymbol{X}_{t}=\frac{1}{\unicode[STIX]{x1D701}_{\bot }}F_{\bot }\hat{\boldsymbol{n}}+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}F_{\Vert }\hat{\unicode[STIX]{x1D749}},\end{eqnarray}$$
where
$\unicode[STIX]{x1D701}_{\bot }$
and
$\unicode[STIX]{x1D701}_{\Vert }$
are the normal and tangential drag coefficients, respectively. Taking the derivative of this equation with respect to arclength gives
$$\begin{eqnarray}\unicode[STIX]{x2202}_{t}\boldsymbol{X}_{s}=\left(\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x2202}_{s}F_{\bot }+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}F_{\Vert }\unicode[STIX]{x1D705}\right)\hat{\boldsymbol{n}}+\left(\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x2202}_{s}F_{\Vert }-\frac{1}{\unicode[STIX]{x1D701}_{\bot }}F_{\bot }\unicode[STIX]{x1D705}\right)\hat{\unicode[STIX]{x1D749}}.\end{eqnarray}$$
Because
$\boldsymbol{X}_{s}=\hat{\unicode[STIX]{x1D749}}$
, the left-hand side of the above equation can be expressed as
$\unicode[STIX]{x2202}_{t}\boldsymbol{X}_{s}=\hat{\unicode[STIX]{x1D749}}_{t}=\hat{\boldsymbol{n}}\unicode[STIX]{x1D713}_{t}$
, and thus the normal terms give
$\unicode[STIX]{x1D713}_{t}$
and the tangential terms must be zero:
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D713}_{t}=\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x2202}_{s}F_{\bot }+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}F_{\Vert }\unicode[STIX]{x1D705}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x2202}_{s}F_{\Vert }-\frac{1}{\unicode[STIX]{x1D701}_{\bot }}F_{\bot }\unicode[STIX]{x1D705}=0. & \displaystyle\end{eqnarray}$$
Equation (A 5) represents a constraint on the forces that must be satisfied to maintain inextensibility. The evolution equation for the curvature is obtained by differentiating (A 4) with respect to arclength to obtain
$$\begin{eqnarray}\unicode[STIX]{x1D705}_{t}=\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x2202}_{ss}F_{\bot }+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x2202}_{s}\left(F_{\Vert }\unicode[STIX]{x1D705}\right).\end{eqnarray}$$
A.3 Expression for elastic forces
The elastic forces are obtained from the variation of an elastic energy functional. The total elastic energy is the sum of a bending term from (2.2) with an energy associated with inextensibility:
$$\begin{eqnarray}{\mathcal{E}}=\int _{0}^{L}\frac{{\mathcal{B}}}{2}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})^{2}+\frac{\unicode[STIX]{x1D6EC}}{2}\boldsymbol{X}_{s}\boldsymbol{\cdot }\boldsymbol{X}_{s}\,\text{d}s,\end{eqnarray}$$
where
$\unicode[STIX]{x1D6EC}$
is a tension used to enforce inextensibility. The force comes from the variation of the energy
$$\begin{eqnarray}\boldsymbol{F}=-\frac{\unicode[STIX]{x1D6FF}{\mathcal{E}}}{\unicode[STIX]{x1D6FF}\boldsymbol{X}}.\end{eqnarray}$$
Using the natural free boundary conditions
$$\begin{eqnarray}\unicode[STIX]{x1D705}=\unicode[STIX]{x1D705}_{0},\quad \unicode[STIX]{x1D705}_{s}=\unicode[STIX]{x2202}_{s}\unicode[STIX]{x1D705}_{0},\quad \unicode[STIX]{x1D6EC}=0,\end{eqnarray}$$
the force is
$$\begin{eqnarray}\boldsymbol{F}=\frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}s}\left(-{\mathcal{B}}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s}\hat{\boldsymbol{n}}+\left(\unicode[STIX]{x1D6EC}+{\mathcal{B}}\unicode[STIX]{x1D705}\left(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0}\right)\right)\hat{\unicode[STIX]{x1D749}}\right).\end{eqnarray}$$
We define the total tension in the rod as
$$\begin{eqnarray}T=\unicode[STIX]{x1D6EC}+{\mathcal{B}}\unicode[STIX]{x1D705}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0}),\end{eqnarray}$$
and the expression for the force is
$$\begin{eqnarray}\displaystyle \boldsymbol{F} & = & \displaystyle \frac{\unicode[STIX]{x2202}}{\unicode[STIX]{x2202}s}\left(-{\mathcal{B}}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s}\hat{\boldsymbol{n}}+T\hat{\unicode[STIX]{x1D749}}\right)\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle \phantom{\boldsymbol{F}} & = & \displaystyle \left(-{\mathcal{B}}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{ss}+\unicode[STIX]{x1D705}T\right)\hat{\boldsymbol{n}}+\left({\mathcal{B}}\unicode[STIX]{x1D705}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s}+T_{s}\right)\hat{\unicode[STIX]{x1D749}}.\end{eqnarray}$$
The normal and tangential force densities on the rod are thus
$$\begin{eqnarray}\displaystyle F_{\bot } & = & \displaystyle -{\mathcal{B}}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{ss}+\unicode[STIX]{x1D705}T,\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle F_{\Vert } & = & \displaystyle {\mathcal{B}}\unicode[STIX]{x1D705}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s}+T_{s}.\end{eqnarray}$$
With these forces, the evolution equation for the curvature (A 6) and the inextensibility constraint which determines the tension (A 5) are
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D705}_{t}=\frac{1}{\unicode[STIX]{x1D701}_{\bot }}(-{\mathcal{B}}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{ssss}+(\unicode[STIX]{x1D705}T)_{ss})+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}({\mathcal{B}}(\unicode[STIX]{x1D705}^{2}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s})_{s}+(\unicode[STIX]{x1D705}T_{s})_{s}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{\unicode[STIX]{x1D701}_{\Vert }}T_{ss}-\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D705}^{2}T+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}({\mathcal{B}}\unicode[STIX]{x1D705}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{s})_{s}+\frac{1}{\unicode[STIX]{x1D701}_{\bot }}{\mathcal{B}}\unicode[STIX]{x1D705}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{ss}=0. & \displaystyle\end{eqnarray}$$
These equations together with the boundary conditions at the ends of the rod
$$\begin{eqnarray}\unicode[STIX]{x1D705}=\unicode[STIX]{x1D705}_{0},\quad \unicode[STIX]{x1D705}_{s}=\unicode[STIX]{x2202}_{s}\unicode[STIX]{x1D705}_{0},\quad T=0,\end{eqnarray}$$
determine the motion.
A.4 Asymptotic expansions
We examine the leading-order behaviour in two different limits: (i) small-amplitude motion in which
$\unicode[STIX]{x1D705}\rightarrow 0$
and (ii) high stiffness,
${\mathcal{B}}\rightarrow \infty$
, in which
$\unicode[STIX]{x1D705}\rightarrow \unicode[STIX]{x1D705}_{0}$
. Note that for soft bodies with large
$\unicode[STIX]{x1D705}_{0}$
the realized amplitude,
$\unicode[STIX]{x1D705}$
, is in fact small. Thus we expect the largest discrepancy between the asymptotic solutions at large
$\unicode[STIX]{x1D705}_{0}$
at intermediate stiffness. In fact this is what we see in figure 4 for the flexor and figure 7 for the two-mode swimmer.
A.4.1 Small amplitude
In the limit of small curvature, we see from (A 17) that the size of the tension is like the square of the curvature. At leading order the tension is zero, and the equation for the curvature is
$$\begin{eqnarray}\unicode[STIX]{x1D705}_{t}=-\frac{{\mathcal{B}}}{\unicode[STIX]{x1D701}_{\bot }}(\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0})_{ssss}.\end{eqnarray}$$
Changing to curvature deviation,
$c=\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0}$
, and non-dimensionalizing gives (3.5) analysed in the main text.
A.4.2 Large stiffness
In the limit
${\mathcal{B}}\rightarrow \infty$
,
$\unicode[STIX]{x1D705}\rightarrow \unicode[STIX]{x1D705}_{0}$
. We introduce the variable
$$\begin{eqnarray}c=\unicode[STIX]{x1D705}-\unicode[STIX]{x1D705}_{0}\end{eqnarray}$$
to denote the deviation from the prescribed curvature. For increasing stiffness,
$c\rightarrow 0$
. Changing variables from
$\unicode[STIX]{x1D705}$
to
$c$
and linearizing about small
$c$
gives
$$\begin{eqnarray}\displaystyle & \displaystyle c_{t}+\unicode[STIX]{x2202}_{t}\unicode[STIX]{x1D705}_{0}=\frac{1}{\unicode[STIX]{x1D701}_{\bot }}(-{\mathcal{B}}c_{ssss}+(\unicode[STIX]{x1D705}_{0}T)_{ss})+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}({\mathcal{B}}(\unicode[STIX]{x1D705}_{0}^{2}c_{s})_{s}+(\unicode[STIX]{x1D705}_{0}T_{s})_{s}), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{\unicode[STIX]{x1D701}_{\Vert }}T_{ss}-\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D705}_{0}^{2}T+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}{\mathcal{B}}\left(\unicode[STIX]{x1D705}_{0}c_{s}\right)_{s}+\frac{1}{\unicode[STIX]{x1D701}_{\bot }}{\mathcal{B}}\unicode[STIX]{x1D705}_{0}c_{ss}=0, & \displaystyle\end{eqnarray}$$
and the boundary conditions are
$$\begin{eqnarray}c=0,\quad c_{s}=0,\quad T=0.\end{eqnarray}$$
These equations contain many terms that are absent for small curvatures. However, as demonstrated in the text by comparing with numerical results, the low curvature equations appear to give a reasonable approximation at the amplitudes tested. Below we show why the two approximations are similar.
For simplicity we consider the flexor in which
$\unicode[STIX]{x1D705}_{0}$
is only a function of time. With this simplification, equations (A 21)–(A 22) become
$$\begin{eqnarray}\displaystyle & \displaystyle c_{t}+\unicode[STIX]{x2202}_{t}\unicode[STIX]{x1D705}_{0}=-\frac{{\mathcal{B}}}{\unicode[STIX]{x1D701}_{\bot }}c_{ssss}+\frac{{\mathcal{B}}\unicode[STIX]{x1D705}_{0}^{2}}{\unicode[STIX]{x1D701}_{\Vert }}c_{ss}+\unicode[STIX]{x1D705}_{0}\left(\frac{1}{\unicode[STIX]{x1D701}_{\bot }}+\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}\right)T_{ss}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{1}{\unicode[STIX]{x1D701}_{\Vert }}T_{ss}-\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D705}_{0}^{2}T+{\mathcal{B}}\unicode[STIX]{x1D705}_{0}\left(\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}+\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\right)c_{ss}=0. & \displaystyle\end{eqnarray}$$
These equations can be solved by orthogonal function expansion. First we eliminate the tension by solving the second equation. We express
$T$
as the series
$$\begin{eqnarray}T=\mathop{\sum }_{m=1}^{\infty }\unicode[STIX]{x1D6FD}_{m}(t)\sin \left(m\unicode[STIX]{x03C0}s\right).\end{eqnarray}$$
We then write the function
$T$
as
$$\begin{eqnarray}T(s,t)=\unicode[STIX]{x1D64C}^{-1}\unicode[STIX]{x1D737},\end{eqnarray}$$
where
$\unicode[STIX]{x1D737}$
represents the sequence of coefficients and
$\unicode[STIX]{x1D64C}^{-1}$
is the orthogonal operator which maps these coefficients to
$T$
, i.e.
$\unicode[STIX]{x1D64C}$
is like the Fourier transform operator. Using this expansion, the operator applied to
$T$
in (A 25) diagonalizes, and the solution is
$$\begin{eqnarray}T={\mathcal{B}}\unicode[STIX]{x1D705}_{0}\left(\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}+\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\right)\unicode[STIX]{x1D64C}^{-1}\left(\frac{1}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x1D648}^{2}+\frac{1}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D705}_{0}^{2}\right)^{-1}\unicode[STIX]{x1D64C}c_{ss},\end{eqnarray}$$
where
$\unicode[STIX]{x1D648}$
is a diagonal matrix with elements
$m\unicode[STIX]{x03C0}$
on the diagonal. After using this expression to eliminate
$T$
in (A 24), after some simplification, we get the equation for the curvature deviation as
$$\begin{eqnarray}c_{t}+\unicode[STIX]{x2202}_{t}\unicode[STIX]{x1D705}_{0}=-\frac{{\mathcal{B}}}{\unicode[STIX]{x1D701}_{\bot }}c_{ssss}+\frac{{\mathcal{B}}\unicode[STIX]{x1D705}_{0}^{2}}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x1D64C}^{-1}\unicode[STIX]{x1D63F}\unicode[STIX]{x1D64C}c_{ss},\end{eqnarray}$$
where
$\unicode[STIX]{x1D63F}$
is a diagonal matrix with elements on the diagonal
$$\begin{eqnarray}\unicode[STIX]{x1D60B}_{mm}=1-\frac{\left(1+\displaystyle \frac{\unicode[STIX]{x1D701}_{\Vert }}{\unicode[STIX]{x1D701}_{\bot }}\right)^{2}}{\left(1+\displaystyle \frac{\unicode[STIX]{x1D701}_{\Vert }}{\unicode[STIX]{x1D701}_{\bot }}\frac{\unicode[STIX]{x1D705}_{0}^{2}}{m^{2}\unicode[STIX]{x03C0}^{2}}\right)}.\end{eqnarray}$$
Because
$\unicode[STIX]{x1D701}_{\Vert }<\unicode[STIX]{x1D701}_{\bot }$
, the values of
$\unicode[STIX]{x1D60B}_{mm}$
can be bounded as
$-3\leqslant \unicode[STIX]{x1D60B}_{mm}\leqslant 1$
.
Equation (A 29) contains one additional term involving the second derivative that is not present in the corresponding low-amplitude equation (3.5). Below we argue that the additional term is small even when
$\unicode[STIX]{x1D705}_{0}$
itself is not small. Our analysis in the main text relied on performing an eigenfunction expansion using the eigenfunctions of the beam equation. We use the same expansion here
$$\begin{eqnarray}c=\mathop{\sum }_{k=1}^{\infty }\unicode[STIX]{x1D6FC}_{k}(t)\unicode[STIX]{x1D6F9}_{k}(s),\end{eqnarray}$$
and we relate
$c$
to its expansion coefficients by
$$\begin{eqnarray}c(s,t)=\unicode[STIX]{x1D64B}^{-1}\unicode[STIX]{x1D736}(t),\end{eqnarray}$$
where
$\unicode[STIX]{x1D64B}^{-1}$
is the orthogonal operator that maps the expansion coefficients to
$c$
. We can transform (A 29) into a system of differential equation for the expansion coefficients as
$$\begin{eqnarray}\frac{\text{d}\unicode[STIX]{x1D736}}{\text{d}t}+\unicode[STIX]{x1D64B}\unicode[STIX]{x2202}_{t}\unicode[STIX]{x1D705}_{0}=-\frac{{\mathcal{B}}}{\unicode[STIX]{x1D701}_{\bot }}\unicode[STIX]{x1D649}^{4}\unicode[STIX]{x1D736}+\frac{{\mathcal{B}}\unicode[STIX]{x1D705}_{0}^{2}}{\unicode[STIX]{x1D701}_{\Vert }}\unicode[STIX]{x1D64B}\unicode[STIX]{x1D64C}^{-1}\unicode[STIX]{x1D63F}\unicode[STIX]{x1D64C}\unicode[STIX]{x1D64B}^{\prime \prime }\unicode[STIX]{x1D736},\end{eqnarray}$$
where we define
$\unicode[STIX]{x1D64B}^{\prime \prime }$
so that
$\unicode[STIX]{x1D64B}^{\prime \prime }\unicode[STIX]{x1D736}=c_{ss}$
and
$-\unicode[STIX]{x1D649}^{4}$
is a diagonal matrix containing the eigenvalues of the beam equation. That is, the
$k\text{th}$
diagonal entry of
$\unicode[STIX]{x1D649}$
,
$\unicode[STIX]{x1D708}_{k}$
, is related to the
$k\text{th}$
eigenvalue,
$\unicode[STIX]{x1D707}_{k}$
, by
$\unicode[STIX]{x1D707}_{k}=-\unicode[STIX]{x1D708}_{k}^{4}$
. One expects that the contribution of
$\unicode[STIX]{x1D64B}^{\prime \prime }$
to the
$k\text{th}$
equation to scale like
$\unicode[STIX]{x1D708}_{k}^{2}$
. As argued above the norm of
$\unicode[STIX]{x1D64C}^{-1}\unicode[STIX]{x1D63F}\unicode[STIX]{x1D64C}$
is about 1, and so we expect the additional terms relative to the bending terms to contribute
$\unicode[STIX]{x1D708}_{k}^{-2}$
. The smallest eigenvalue is about
$\unicode[STIX]{x1D708}_{1}^{4}\approx 500$
, and thus we expect these additional terms to be small.
In figure 9 we compare the expansion coefficients of the first mode of the small-amplitude expansion and the high stiffness expansion for the high-amplitude flexor (
$\unicode[STIX]{x1D705}_{0}=4$
) for a range of stiffnesses. In agreement with our numerical results, the qualitative behaviour of the two expansions is the same, and at high and low stiffnesses the quantitative behaviour is the same. The largest difference occurs for moderately soft bodies where the difference is less than 25 %.
Figure 9. (a) Normalized amplitude of first mode for a high-amplitude flexor with amplitude
$\unicode[STIX]{x1D705}_{0}=4$
as a function of stiffness
${\mathcal{B}}$
for the leading-order low-amplitude expansion and the leading-order high stiffness expansion. The ratio of drag was set to
$\unicode[STIX]{x1D701}_{\Vert }/\unicode[STIX]{x1D701}_{\bot }=0.75$
for this computation. The results are relatively insensitive to this value. (b) Magnitude of the difference of the two expansion coefficients.
Appendix B. Rod driven at one end
In this appendix we give the derivation for the shape-induced swimming speed of a thin filament in a viscous or viscoelastic fluid which is driven by oscillations at one end. This type of motion is akin to the experiments of Espinosa-Garcia et al. (Reference Espinosa-Garcia, Lauga, Zenit, Diego and Jolla2013), and we perform the calculation to compute the swimming speed as a function of the dimensionless body response time
$G$
and the fluid relaxation time
$De.$
We also show the range of
$G$
where a viscoelastic speed-up is theoretically predicted for this type of motion.
B.1 Shape of the swimmer
The problem we consider is a flexible filament with one free end and one clamped end. The motion is driven by prescribing sinusoidal oscillations in the angle at the clamped end. The shape of the filament satisfies
$$\begin{eqnarray}\displaystyle & \displaystyle y_{t}=-Gy_{ssss}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle y(0,t)=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle y_{s}(0,t)=\cos (2\unicode[STIX]{x03C0}t), & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle y_{ss}(1,t)=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle y_{sss}(1,t)=0. & \displaystyle\end{eqnarray}$$
We now change variables using
$$\begin{eqnarray}y=w+s\cos (2\unicode[STIX]{x03C0}t)=w+\text{Re}(s\text{e}^{2\unicode[STIX]{x03C0}\text{i}t}),\end{eqnarray}$$
so that
$w$
represents deviations from the infinitely stiff case of a straight rod. Letting
$w$
be complex valued, the equation for
$w$
is then
$$\begin{eqnarray}w_{t}=-Gw_{ssss}-2\unicode[STIX]{x03C0}\text{i}s\text{e}^{2\unicode[STIX]{x03C0}\text{i}t},\end{eqnarray}$$
with homogeneous boundary conditions. This equation can be solved using an expansion of eigenfunctions,
$\unicode[STIX]{x1D6F9}_{k}(s)$
, which satisfy
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D707}\unicode[STIX]{x1D6F9}(s)=-\unicode[STIX]{x1D6F9}_{ssss}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6F9}(0)=\unicode[STIX]{x1D6F9}^{\prime }(0)=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D6F9}^{\prime \prime }(1)=\unicode[STIX]{x1D6F9}^{\prime \prime \prime }(1)=0. & \displaystyle\end{eqnarray}$$
We express the function
$s$
using an eigenfunction expansion
$$\begin{eqnarray}s=\sum _{k}\unicode[STIX]{x1D6FC}_{k}^{\infty }\unicode[STIX]{x1D6F9}_{k}(s),\end{eqnarray}$$
and look for a solution to the PDE of the form
$$\begin{eqnarray}w(s,t)=\text{e}^{2\unicode[STIX]{x03C0}\text{i}t}\sum _{k}\unicode[STIX]{x1D6FD}_{k}\unicode[STIX]{x1D6F9}_{k}(s).\end{eqnarray}$$
The transformation of the PDE yields
$$\begin{eqnarray}2\unicode[STIX]{x03C0}\text{i}\unicode[STIX]{x1D6FD}_{k}=G\unicode[STIX]{x1D707}_{k}\unicode[STIX]{x1D6FD}_{k}-2\unicode[STIX]{x03C0}\text{i}\unicode[STIX]{x1D6FC}_{k}^{\infty },\end{eqnarray}$$
which gives
$$\begin{eqnarray}\unicode[STIX]{x1D6FD}_{k}=\displaystyle \frac{-\unicode[STIX]{x1D6FC}_{k}^{\infty }}{1-\displaystyle \frac{G\unicode[STIX]{x1D707}_{k}}{2\unicode[STIX]{x03C0}\text{i}}}.\end{eqnarray}$$
We can then write the shape as
$$\begin{eqnarray}\displaystyle y(s,t) & = & \displaystyle \text{Re}\left\{\mathop{\sum }_{k=1}^{\infty }(\unicode[STIX]{x1D6FC}_{k}^{\infty }+\unicode[STIX]{x1D6FD}_{k})\unicode[STIX]{x1D6F9}_{k}(s)\text{e}^{2\unicode[STIX]{x03C0}\text{i}t}\right\},\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle \phantom{y(s,t)} & = & \displaystyle \text{Re}\left\{\mathop{\sum }_{k=1}^{\infty }\unicode[STIX]{x1D6FC}_{k}^{\infty }\left(1-\left(1-\frac{G\unicode[STIX]{x1D707}_{k}}{2\unicode[STIX]{x03C0}\text{i}}\right)^{-1}\right)\unicode[STIX]{x1D6F9}_{k}(s)\text{e}^{2\unicode[STIX]{x03C0}\text{i}t}\right\}.\end{eqnarray}$$
Notice that the factor multiplying
$\unicode[STIX]{x1D6FC}_{k}^{\infty }$
above is exactly the same as the one that appears in the shape analysis of swimmers driven by active moments in (3.10). We can express the shape as
$$\begin{eqnarray}y(s,t)=\text{Re}\left\{\mathop{\sum }_{k=1}^{\infty }\unicode[STIX]{x1D6FC}_{k}\unicode[STIX]{x1D6F9}_{k}(s)\text{e}^{2\unicode[STIX]{x03C0}\text{i}t}\right\},\end{eqnarray}$$
where
$\unicode[STIX]{x1D6FC}_{k}$
is defined by (3.10). As in the main body of the paper, to add viscoelastic effects, we simply use (3.13) in place of (3.10) to define the expansion coefficients. Although these expressions are the same, we note that the eigenfunctions and eigenvalues are different for this problem.
B.2 Expression for swimming speed
The expansion for the shape of the swimmer can be written as
$$\begin{eqnarray}y(s,t)=\mathop{\sum }_{k}^{\infty }A_{k}\cos (2\unicode[STIX]{x03C0}t+\unicode[STIX]{x1D719}_{k})\unicode[STIX]{x1D6F9}_{k}(s),\end{eqnarray}$$
where
$$\begin{eqnarray}\displaystyle & \displaystyle A_{k}=|\unicode[STIX]{x1D6FC}_{k}|, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \unicode[STIX]{x1D719}_{k}=\text{arg}(\unicode[STIX]{x1D6FC}_{k}). & \displaystyle\end{eqnarray}$$
To compute the swimming speed, we average in space and time the product
$$\begin{eqnarray}y_{s}y_{t}=-2\unicode[STIX]{x03C0}\sum _{n}\sum _{m}A_{n}A_{m}\cos (2\unicode[STIX]{x03C0}t+\unicode[STIX]{x1D719}_{n})\sin (2\unicode[STIX]{x03C0}t+\unicode[STIX]{x1D719}_{m})\unicode[STIX]{x1D6F9}_{n}(s)\unicode[STIX]{x1D6F9}_{m}^{\prime }(s).\end{eqnarray}$$
Averaging the above expression by integration gives the swimming speed
$$\begin{eqnarray}\left<U\right>\propto \sum _{n}\sum _{m}A_{n}A_{m}\sin (\unicode[STIX]{x1D719}_{n}-\unicode[STIX]{x1D719}_{m})\int _{0}^{1}\unicode[STIX]{x1D6F9}_{n}(s)\unicode[STIX]{x1D6F9}_{m}^{\prime }(s)\,\text{d}s.\end{eqnarray}$$
Figure 10. (a) Shape-induced Stokes-normalized swimming speed as a function of
$\text{Stroke}_{De}$
for a range of
$G$
. (b) Shape-induced Stokes-normalized swimming speed as a function of
$G$
for low and high
$\text{Stroke}_{De}$
. (c) Shape-induced swimming speed as a function of
$G$
for low and high
$\text{Stroke}_{De}$
.
We use this expression with the first six modes to compute the shape-induced Stokes-normalized swimming speed; see figure 10. As with the problem from the paper, for sufficiently soft bodies, we see an almost monotonic speed-up from the shape changes. For sufficiently stiff swimmers (
$G>1$
) we see a monotonic slow-down. There is a transition range around
$0.1<G<1$
. The location of this transition is evident in figure 10(b) where we show the swimming speed as a function of
$G$
for low and high
$\text{Stroke}_{De}$
. The qualitative results from the paper do not change in the sufficiently soft regime, but this problem has a different driving mechanism and hence there is a different effect in the stiff regime. As we see in figure 10(c) as the body is stiffened the swimming speed goes to zero and viscoelasticity always slows the swimmer.
