Spatiotemporal Asymmetry in Metachronal Rowing at Intermediate Reynolds Numbers

Abstract

In drag-based swimming, individual propulsors operating at low Reynolds numbers (where viscous forces dominate over inertial forces) must execute a spatially asymmetric stroke to produce net fluid displacement. Temporal asymmetry (that is, differing duration between the power vs. recovery stroke) does not affect the overall generated thrust in this time-reversible regime. Metachronal rowing, in which multiple appendages beat sequentially, is used by a wide variety of organisms from low to intermediate Reynolds numbers. At the upper end of this range, inertia becomes important, and increasing temporal asymmetry can be an effective way to increase thrust. However, the combined effects of spatial and temporal asymmetry are not fully understood in the context of metachronal rowing. To explore the role of spatiotemporal asymmetry in metachronal rowing, we combine laboratory experiments and reduced-order analytical modeling. We measure beat kinematics and generated flows in two species of lobate ctenophores across a range of body sizes, from 7 to 40 mm in length. We observe characteristically different flows in ctenophores of differing body size and Reynolds number, and a general decrease in spatial asymmetry and increase in temporal asymmetry with increasing Reynolds number. We also construct a one-dimensional mathematical model consisting of a row of oscillating flat plates whose flow-normal areas change with time, and use it to explore the propulsive forces generated across a range of Reynolds numbers and kinematic parameters. The model results show that while both types of asymmetry increase force production, they have different effects in different regions of the parameter space. These results may have strong biological implications, as temporal asymmetry can be actively controlled while spatial asymmetry is likely to be partially or entirely driven by passive fluid-structure interaction.

Introduction

For swimming and flying animals, the physics of how appendages interact with the surrounding fluid is a strong function of the Reynolds number |$( {Re} )$|⁠, which represents the ratio of inertial to viscous effects within a flow. At low |$Re$| (<1), viscous forces dominate; at higher |$Re$|⁠, inertial forces take on more importance. Efforts to understand the physics of swimming have mainly focused on the extremes of fully viscous-dominated or inertia-dominated flow. However, many animals swim at intermediate |$Re$|⁠, where inertia and viscosity are both important. This includes many metachronal swimmers—that is, those that move fluid via the sequential, coordinated motion of a row of closely spaced appendages (Barlow et al. 1993; Lim and DeMont 2009; Murphy et al. 2011; Caldwell et al. 2012). In this drag-based locomotion strategy, each appendage repeats a power stroke (in which it is extended and moving quickly) followed by a recovery stroke (in which it is shortened or bent and moving slowly), thus creating net fluid displacement in the power direction. A phase lag between adjacent appendages creates the appearance of a “metachronal wave.” This behavior, classically exemplified by cilia (Brennen and Winet 1977), occurs at |$Re$| up to several thousand (Lim and DeMont 2009). Here, we use ctenophores, the largest animals in the world who swim using cilia (Tamm 2014) with appendage-scale |$Re$| of 10–200, to bridge the knowledge gap between the low and intermediate Re regimes for metachronal swimming.

Ctenophores (comb jellies) swim using eight rows of paddle-like ctenes circumscribing a spheroidal body (Fig. 1). Each ctene comprises thousands of long cilia which move as a unit (Goebel et al. 2020). These cilia are similar to those found at the micron scale (low|$\,\,Re$|⁠), but are much longer, approximately 1 mm in length (Afzelius 1961). Ctenophores therefore provide a unique opportunity to study how cilia move as inertia becomes important.

Fig. 1

Morphology of Bolinopsis vitrea. (A) top view showing the eight ctene rows and their radial symmetry, (B) lateral view of ctene rows circumscribing the body, and (C) close top view showing component cilia of each ctene.

While the whole-body kinematics and fluid dynamics of ctenophore swimming are well-described (Matsumoto 1991; Colin et al. 2010; Gemmell et al. 2019), few studies directly examine how ctene kinematics affect the flow. Barlow and Sleigh (1993) showed that metachronal coordination increases efficiency, as has also been shown for low |$Re$| cilia (Elgeti and Gompper 2013). Numerical simulations have proven an effective tool to investigate the hydrodynamic interactions between ctenes responsible for this highly effective flow production (Dauptain et al. 2008). However, these studies primarily observed the effect of varying beat frequency rather than the kinematics of the individual ctenes.

The beat cycle of a metachronal appendage can be both spatially asymmetric (higher flow-normal area on the power stroke vs. the recovery stroke) and temporally asymmetric (power stroke duration shorter than recovery stroke duration). At low |$Re$|⁠, flow is time-reversible and spatial asymmetry is necessary for a single appendage to produce net fluid displacement (Purcell 1977). However, metachronal coordination can remove this constraint (Takagi 2015)—if appendages beat with a nonzero phase-lag, spatially symmetric stroking can produce net flow even at low |$Re$|⁠. However, at intermediate |$Re$| the role of spatiotemporal asymmetry is still unclear. Very recently, (Saffaraval and Goudarzi 2021) experimentally explored the roles of spatial and temporal asymmetries in a single oscillating plate at intermediate |$Re$|⁠, finding that increasing temporal asymmetry increased pumping performance and the presence of spatial asymmetry increased the pumping efficiency. However, experimental constraints in this study did not allow for continuous variation of spatial asymmetry, nor varying levels of both asymmetries simultaneously.

In this study, we explore how appendages’ spatiotemporal asymmetry affects generated flows and forces. We use high-speed video and velocimetry to measure the kinematics and fluid dynamics of the beating ctene rows of two species of lobate ctenophores, across a range of body sizes and |$Re$|⁠. We also introduce a reduced-order analytical model to help explore the parameter space (see Table 1) of spatiotemporal asymmetry in metachronal rowing. In both the experiments and the model, we quantify temporal asymmetry (Ta), as defined by Gauger et al. (2009), and propose a new non-dimensional parameter (Sa) to quantify spatial asymmetry. Using the model, we show how variations in |$Sa$| and |$Ta\,\,$| have different effects in different regions of the parameter space. We show that |$Ta$| affects force production and swimming efficiency more strongly than |$Sa$|⁠, but that these effects are not uniform across the range of observed |$Sa$| and |$Ta$|⁠. Ultimately, we find that spatiotemporal asymmetry variation plays an important role in metachronal swimming and is worthy of further investigation.

Table 1

List of kinematic, morphometric, fluid dynamic, and model parameters

Variable	Description
|$\rho $|	Seawater density (⁠|$1.027x{10^{ - 6}}\, \rm kg/m{m^3}$|⁠)
|$\nu $|	Seawater kinematic viscosity (⁠|$1.05 \, \rm m{m^2}/s$|⁠)
m	Ctenophore mass
V	Ctenophore volume
|$\alpha $|	Added mass coefficient
N	Number of appendage rows
n	Number of appendages in a row
|${x_B}$|	Animal position
|${x_A}$|	Appendage x-position
|${y_A}$|	Appendage height
r	Average radius of body's flow normal area
l	Appendage length
w	Appendage width
|${L_B}$|	Body length
f	Beat frequency
|${P_L}$|	Phase lag
|$\Phi $|	Stroke amplitude
|${C_B}$|	Drag coefficient of a prolate spheroid
|${F_D}$|	Drag force
|${F_p}$|	Propulsive force of one appendage
|${F_{net}}$|	Propulsive force of an appendage row
|$G( \cdot )$|	Relative influence operator
|$R{e_\omega }$|	Oscillatory Reynolds number
|$Sa$|	Spatial asymmetry
|${A_e}$|	Area enclosed by ctene tip trajectory
|${A_o}$|	Maximum area enclosed by tip trajectory
|$Ta$|	Temporal asymmetry
|${t_r}$|	Recovery stroke duration
|${t_p}$|	Power stroke duration
|$\eta $|	Swimming efficiency
|${P_i}$|	Power used by appendages
|${P_o}$|	Power required to move animal

Materials and methods

Animal collection and husbandry

A total of two sampling sites and species were used for this study (summarized in Table 2). Bolinopsis vitrea were collected at Flatt's Inlet, Bermuda in May 2018 and transported to the Bermuda Institute of Ocean Sciences (BIOS). Larger animals were maintained in a flow-through planktonkreisel tank with filtered seawater, while smaller specimens were kept in glass jars partially submerged in an open sea table with a continuous flow of filtered seawater. Animals were fed ad libitum on freshly hatched nauplii of Artemia sp. All experiments were carried out at ambient temperature (21–23°C) and completed within 12 h of animal collection.

Table 2

Morphometric and kinematic measurements of observed animals (mean |$ \pm $| one standard deviation)

Data set	Sample size	Body length (mm)	Ctene length (mm)	Beat frequency (Hz)
BIOS	8	18.43 ± 14.01	0.42 ± 0.11	15.11 ± 2.21
MBA	4	24.01 ± 8.12	0.87 ± 0.1	27.74 ± 1.97

Bolinopsis infundibulum were sampled from existing cultures at the Monterey Bay Aquarium (MBA), Monterey, CA, in August 2019. Cultured animals were maintained in planktonkreisels with flow-through filtered seawater at 12–13°C and fed Artemia nauplii twice daily. Animals were removed from one kreisel, immediately used in experiments, and replaced in a second kreisel to avoid re-use.

Imaging facilities

To obtain both flow velocities and ctene kinematics, we used a micro-scale Particle Shadow Velocimetry (µPSV) system following (Gemmell et al. 2014). µPSV is also called “brightfield Particle Image Velocimetry” (Khodaparast et al. 2013). Collimated white light illuminates the entire filming vessel and produces sharply-defined shadows of tracer particles and the target organism, which are then imaged via an extra-long working distance (ELWD) microscope objective (Fig. 2B). After inverting and preprocessing raw images, they can be input to standard PIV algorithms to produce two-dimensional velocity fields (Goss et al. 2007). The imaging volume is controlled by the optical depth of correlation (Koutsiaris 2012), such that strongly out-of-focus tracers do not contribute to the calculated flow field. The use of ELWD objectives permits the focal plane to be centered in the filming vessel, thus avoiding wall effects (Gemmell et al. 2014).

Fig. 2

(A) Schematic of experimental μPSV set up showing in-line configuration of camera, ELWD objective, filming vessel, and light source. (B) Example image collected via μPSV.

System components (MBA//BIOS) include a collimated white LED light source (Thorlabs, Newton, NJ and Dolan-Jenner Industries, Lawrence, MA), a transparent filming vessel, an ELWD objective (Achrovid, Centennial, CO and Mitutoyo, Aurora, IL), and a high-speed camera (Photron AX200, Photron USA Inc., San Diego, CA and Phantom VEO640S, Vision Research Inc., Wayne, NJ; Fig. 2A). In both experiments, images were collected at 1000 frames per second with an exposure time of 900 µs. To avoid adverse animal reactions to plastic or glass microbeads, the water was seeded with phytoplankton, which served as approximately passive flow tracers since their swimming speed was much lower than the measured flow speeds (Gemmell et al. 2014). To accommodate differences in resolution between the two imaging systems, we used Tetraselmis chuii (average diameter of 10–12 µm) for the MBA dataset and Nannochloropsis oculata (2–3 µm) for the BIOS dataset. Tracer shadows were recorded in digital images with particle image diameters of 2–3 pixels, as appropriate for particle velocimetry (Raffel et al. 2018). We note that the brightfield imaging system used for the BIOS dataset is described in (Karakas et al. 2020), and that those data are contemporary with the results presented here.

Velocimetry and kinematic analysis

We recorded image sequences when animals rested against the bottom of the filming vessel (zero forward swimming velocity). In some cases, we constrained animals using transparent acrylic plates. In all cases, the focal plane remained far from solid surfaces to avoid wall effects. We pursued kinematics analysis only if the ctenes were normal to the focal plane, and velocimetry analysis only if the plane of focus was roughly bisecting the ctene row. These restrictions yielded a total of 12 image sequences for kinematics analysis and two sequences for velocimetry analysis. The Measurement Plane Width (MPW) or thickness of the useful measurement volume, equal to twice the depth of correlation, is typically around 10 µm for µPSV, using magnifications of 10–40x (Gemmell et al. 2014). Our system required a 5x magnification in order to visualize several ctenes simultaneously (with a field of view of approximately 4 × 2.5 mm), leading to an MPW of ∼130 µm (Koutsiaris 2012). At this magnification, standard PIV image prefilters retain significant noise from out-of-focus tracers. We therefore developed a custom prefiltering algorithm that finds each focused particle and generates a new binary image containing only the focused particles (Herrera-Amaya 2020). In these binarized images, the average seeding density though the focal plane is |$1.3x{10^{ - 3}}$| particles per pixel (Raffel et al. 2018). We performed vector computation (two-pass iteration with subwindows of 128 × 128 px and 64 × 64 px, with 50% overlap) using the MATLAB-based tool PIVLab (Thielicke and Stamhuis 2014). All image preprocessing and vector field postprocessing was conducted with MATLAB (Mathworks Inc., Natick, MA).

Kinematic analysis was based on two-dimensional tracking data performed via the MATLAB-based tool DLTdv8 (Hedrick 2008), including deep learning features. We tracked the base and tip of all fully visible ctenes for each image sequence. Figure 3 shows examples of two tracked ctene rows from two different animals. From the measured tip trajectory, we calculated the spatial asymmetry (degree of difference in flow-normal area between the power and recovery stroke) and temporal asymmetry (difference in duration of the power and recovery stroke). These parameters are defined in the next section.

Fig. 3

Snapshots showing the tracked ctene base (filled circles) and ctene tip trajectories (solid white line). White dots overlaid on tip trajectory indicate equal time intervals, showing the difference in tip speed between the power and recovery strokes. Ctene length (l) and stroke amplitude (Φ) are also marked. (A) shows B. vitrea with |$L_B = 11.56$| mm and |$l = 0.62\pm 0.02$| mm, with high spatial asymmetry; (B) shows B. vitrea with |$L_B = 40$| mm and |$l = 0.76\pm0.01$| mm, with lower spatial asymmetry.

Morphological and metachronal parameters

We describe the motion of a ctene row using seven parameters: two morphometric variables, and five kinematic/coordination variables. Morphometric parameters include the ctene length |$( l )$|⁠, defined as the longest measured distance between the ctene base and tip during a beat cycle (Fig. 3) and the ctene width |$( w )$|⁠, which was assumed to be |$\,\,l/2$| (Afzelius 1961). Coordination parameters include beat frequency |$( f )$|⁠, phase lag |$( {{P_L}} )$|⁠, stroke amplitude |$( \Phi )$|⁠, and temporal asymmetry (⁠|$Ta$|⁠). Phase lag is defined as the percentage of the cycle time that a ctene's beat cycle leads a neighboring ctene. Stroke amplitude |$( \Phi )$| is defined as the angle swept out by the ctene tip over one power stroke (Fig. 3B). Temporal asymmetry (also known as the kinematic asymmetry parameter (Gauger et al. 2009)) is defined as:

$$\begin{equation*} Ta\,\, = \frac{{{t_r} - {t_p}}}{{{t_r} + {t_p}}}\,\, \end{equation*}$$

(1)

where |${t_r}$| is the recovery stroke duration and |${t_p}$| the power stroke duration. This parameter is zero for a time-symmetric beat cycle and one for an infinitely fast power stroke (Fig. 4B).

Fig. 4

(A) Schematic of ctene tip speed over one beat cycle. The temporal asymmetry parameter |$Ta$| is a function of |$t_p$| (time between the first two minima of the tip speed curve) and |$t_r$| (remaining cycle time). (B) Graphical definition of spatial asymmetry parameter |$Sa$| overlaid on a timeseries of ctene kinematics (gray lines). Black half-circle is the reachable space of a ctene with length l; larger ellipse |$A_o$| (red dotted line) estimates maximum area that could be enclosed by ctene tip over a complete cycle; inner shape |$A_e$| (blue dotted line) is the actual ctene tip trajectory.

To quantify the level of spatial asymmetry in each beat cycle for each ctene, we define a new parameter:

$$\begin{equation*} Sa\,\, = \frac{{{A_e}}}{{{A_o}}}\,\, \end{equation*}$$

(2)

where |${A_e}$| is area enclosed by the tip trajectory and |${A_o}$| is the largest possible area of an ellipse inscribed within the reachable space of the ctene (that is, a half-circle with radius l), as shown in Fig. 4A. Because of the ctene's flexibility, this ellipse represents a rough estimate of the theoretical maximum of the area that could be traced out by its tip. The ratio |$Sa\,\, = \,\,\frac{{{A_e}}}{{{A_o}}}$| would therefore tend to one for highly asymmetric beat cycles, and zero for a rigid ctene moving symmetrically.

Experimental results

Velocimetry

To examine the flows generated by ctene rows beating across a range of sizes, Reynolds numbers, and spatiotemporal asymmetries, we used µPSV as described in the previous section. Figure 5A and D show sample velocity fields from two specimens of B. vitrea. Flow shown in Fig. 5A is similar to previous observations of ctenophores (Barlow et al. 1993; Colin et al. 2020), with no apparent vortex shedding—the peak flow speed is co-located with the ctene tip as it reaches its own maximum speed. Flow reversal is only present in the sub-ctene layer (i.e., space between the ctenes); the outer flow is primarily unidirectional and relatively steady. The animal shown in Fig. 5D is larger in both body size and ctene length, and generates a slightly more complex flow. The maximum flow speed is no longer attached to the ctene tip. Vortex shedding is present, an indication of the growing importance of fluid inertia. Vortex shedding such as that seen in Fig. 5D has not been previously seen in ctenophores, but is consistent with observations of metachronal swimmers operating at slightly larger |$Re$| (250–1000) ( Murphy et al. 2013; Colin et al. 2020; Ford and Santhanakrishnan 2021). Fig. 5B and E show the radial velocity magnitude, averaged over seven beat cycles, while   Fig. 5C and F show the time-averaged tangential velocity magnitude, with the angular division based on the center of a circular arc fitted to the body wall. The flow is tangentially dominated in both animals; however, the radial velocity close to the ctenes is relatively higher for the larger animal, indicating increased two-dimensionality.

Fig. 5

(A) Instantaneous velocity field for a ctenophore with body length |$L_B = 11.56 $| mm and average ctene length |$l = 0.62\pm0.02$| mm. (B) Time averaged radial velocity magnitude. (C) Time averaged tangential velocity magnitude. (D) Instantaneous velocity field for a ctenophore with body length |$L_B = 38.99$| mm and average ctene length |$l = 0.93\pm0.03$| mm. (E) Time averaged radial velocity magnitude. (F) Time averaged tangential velocity magnitude. Panels B, C, E, and F are averaged over seven beat cycles and share the color scale displayed on the right.

To quantify the viscous-inertial balance for a single ctene, we use the oscillatory Reynolds number, |$R{e_\omega } \equiv \frac{{2\pi f{l^2}}}{\nu }$| (where |$\nu $| is the kinematic viscosity of seawater). The animal shown in Fig. 5D has a slightly higher |$R{e_\omega }$| (see Table 3) and shows a lower degree of spatial asymmetry (⁠|$Sa\,\, = \,\,0.29$|⁠) relative to the animal shown in Fig. 5A (⁠|$Sa\,\, = \,\,0.53$|⁠), which has a slightly lower |$R{e_\omega }$| (see Table 3). It is likely that both the increased |$R{e_\omega }$| and the lower degree of spatial asymmetry in the stroke contribute to the presence of vortices in Fig. 5D. Though the Reynolds number is lower and the appendage is more deformable, the measured flows are similar to the well-characterized vortex shedding behavior observed in flexible flapping foils (Marais et al. 2012; Dewey et al. 2013; Quinn et al. 2014).

Table 3

Flow parameters for Fig. 5 (mean |$ \pm $| one standard deviation of 30 (small animal) and 5 (large animal) cycles)

|${L_b}( {mm} )$|	|$R{e_\omega }$|	|$l( {mm} )$|	|$f( {Hz} )$|	|$Sa$|	|$Ta$|	|${P_L}( {\rm{\% }} )$|	|$\Phi ( ^\circ )$|
11.56	23.93 ± 0.5	0.62 ± 0.02	16.74 ± 0.98	0.53 ± 0.01	0.04 ± 0.07	13.11 ± 1.22	101.76 ± 4.83
38.99	68.27 ± 3.39	0.93 ± 0.03	11.23 ± 0.23	0.29 ± 0.03	0.15 ± 0.04	13.04 ± 1.24	104.15 ± 2.49

We note that obtaining high quality flowfields is technically challenging, as it requires delicate alignment between the ctene row and the focal plane of the imaging system. Furthermore, the ctene row must remain aligned for a nontrivial period of time in order to draw conclusions about the beat cycle. This alignment is difficult to achieve without strong external constraints immobilizing the animal (such as pins (Tamm and Moss 1985) or a vice grip (Barlow et al. 1993)), which we wished to avoid. These quality standards resulted in only two video sequences suitable for µPSV analysis. While these sequences do show animals of different sizes and ctenes beating at slightly different |$R{e_\omega }$|⁠, the sample size is not large enough to draw definitive conclusions about the role of spatiotemporal asymmetry in flow generation across scales. We include these data here because they provide a helpful visualization of concepts we discuss in subsequent sections, and because they exemplify how we expect the flow to become more two-dimensional at higher Reynolds numbers. We also note that these animals are somewhat constrained and measured with zero swimming speed as discussed in the previous section. The velocity fields, particularly in the far field, are likely to be different for a freely swimming animal. The ctene kinematics, particularly those that are driven by passive fluid–structure interactions, may therefore be different for this case as well. However, our qualitative observations of the general character of the velocity fields—a more two-dimensional flowfield for larger animals with higher beat frequencies—are likely to hold for the free-swimming case as well.

Appendage kinematics

Of the previously defined kinematic/coordination parameters, several are commonly reported (stroke amplitude |$\Phi $|⁠, phase lag |${P_L}$|⁠, and beat frequency f). Table 4 shows these parameters for several ctenophore species; our results are consistent with previous observations. However, the degree of spatiotemporal asymmetry is not typically reported. In the current study, we measure spatial asymmetry |$Sa$| and temporal asymmetry |$Ta$| for ctenes engaged in idle, steady beating as previously defined (Fig. 4).

Table 4

Measured values for stroke amplitude,|$\,\,( \Phi )$|⁠, phase-lag |$( {{P_L}} )$|⁠, and the beat frequency |$( f )$| of ctene rows. Where margins are indicated, they are the mean plus or minus one standard deviation

Species	|$\Phi \,\,( ^\circ )$|	|${P_L}\,\,( \% )$|	|$f( {Hz} )$|	Beating mode	Source
B. vitrea	112.1 ± 8.9	13.2 ± 1.7	15.1 ± 2.2	Idle/Tethered	Current study
B. infundibulum	102.1 ± 5.2	13.9 ± 4.9	27.7 ± 1.9	Idle/Tethered	Current study
M. ovum  *	−	−	7	Active swimming	Matsumoto (1991)
B. artic  *	−	−	15	Active swimming	Matsumoto (1991)
L. pulchra  *	−	−	13	Active swimming	Matsumoto (1991)
P. pileus+	135 ± 15.5	−	10, 27	Tethered	Barlow et al. (1993)
P. pileus +	−	7.07 ± 3.5	5 to 25	Tethered	Barlow and Sleigh (1993)
P. bachei	−	−	10.1 ± 5.1	Active swimming	Goebel et al. (2020)

*

Average of > 25 measurements (σ not reported), + Beat frequency controlled by a mechanical actuator.

Our own results show that for B. vitrea, as |$R{e_\omega }$| increases, |$Sa$| declines and |$Ta$| is nearly constant; for B. infundibulum, |$Sa$| is nearly constant, and |$Ta$| increases (Fig. 6). All specimens of B. infundibulum were observed to have higher |$R{e_\omega }$| than that of B. vitrea. While these data are somewhat sparse, they indicate that the roles of spatial and temporal asymmetry may interact with one another, and that this interaction may change depending on the relative importance of viscosity and inertia. The drop in |$Sa$| is non-intuitive: since higher spatial asymmetry should correspond to greater force production, it is unclear why it would decrease with increasing |$R{e_\omega }$|⁠. To further investigate the relationship between |$Sa$| and |$Ta$|⁠, we turn to a simple reduced-order analytical model. Using this simple model, we can qualitatively investigate how increasing or decreasing spatiotemporal asymmetry can affect swimming efficiency. The model is simple by design, meant to predict trends in force and power generation (and how they vary with spatiotemporal asymmetry) rather than the true values of these quantities. Due to its necessary simplifying assumptions, this model will not capture the full complexity of an actual ctene row; however, it serves as a useful tool for investigating the role of spatiotemporal asymmetry across a large parameter space, which would be computationally prohibitive for a higher-fidelity model.

Fig. 6

(A) Spatial asymmetry |$Sa$| vs. oscillatory Reynolds number |$Re_\omega$|⁠. Lower |$Re_\omega$| data (shown in blue) indicates B. vitrea (BIOS); higher |$Re_\omega$| data (shown in red) indicates B. infundibulum (MBA). Linear regressions for each dataset (solid lines) have regression coefficients |$R_{BIOS} = 0.4$| and |$R_{MBA} = 0.01$|⁠. (B) Temporal asymmetry |$Ta$| vs. oscillatory Reynolds number |$Re_\omega$|⁠. Linear regressions for each dataset (solid lines) have regression coefficients |$R_{BIOS} = 0.01$| and |$R_{MBA} = 0.79$|⁠. Open circles represent the mean value for a single ctene; errorbars show 95% confidence interval based on measured cycles (>20). Dotted lines represent 95% confidence bands on the linear regressions.

Ctene row propulsion model

A limitation of behavioral experiments is that parameters cannot be varied independently; without a prohibitively large dataset, it is difficult to fully explore the parameter space and identify key functional dependencies. Computational Fluid Dynamics (CFD) modeling can partially relieve these constraints, but is computationally expensive to apply over a large number of parametric permutations; CFD models are also typically more difficult to generalize and apply to similar problems in other taxa. Therefore, we create a simple one-dimensional analytical model: a row of flat, oscillating plates whose area is a function of time (Fig. 7) and whose oscillatory Reynolds number |$R{e_\omega }$| can be varied from low to intermediate values. The plates do not rotate; they simply move back and forth along a fixed one-dimensional trajectory. Hydrodynamic interactions between plates are not considered. This model is similar to that considered by (Alben et al. 2010), but in our case the drag coefficient of each plate is a function of |$R{e_\omega }$|⁠, allowing intermediate Reynolds number effects to be considered for the first time.

Fig. 7

(A) Sample trajectory of a real ctene and (B) simplified elliptical trajectory for the analytical model. Thin solid lines (dark blue) depict tip trajectories; thick solid line (light blue) represents oscillating plate, whose tip position |$(x_A,y_A)$| is a function of θ and ellipse geometry. Dotted lines denote the stroke amplitude Φ. The spatiotemporal asymmetry parameters, |$Sa$| and |$Ta$|⁠, are prescribed for each model run; these determine the time-varying tip position |$(x_A,y_A)$|⁠, which subsequently determine the oscillating plate's speed and flow-normal area on the power and recovery strokes.

This modeled ctene row, which can also be considered as a general model for metachronal appendage motion, can be examined in isolation (to investigate propulsive force generation) or coupled to a spheroidal body (to estimate swimming efficiency). We stress that this reduced-order model is not meant to predict the actual propulsive forces and efficiencies of a swimming ctenophore, but to provide a tractable way to explore the interplay between spatial and temporal asymmetry throughout the kinematic parameter space.

Model parameters

Drag-based models have been successfully used to model metachronal rowing, even while neglecting the hydrodynamic interactions between paddles (Alben et al. 2010; Takagi 2015). These interactions do affect the overall system performance (Dauptain et al. 2008; Ford et al. 2019); however, by neglecting them, we can quickly and easily sweep the parameter space to gain insight into the role of spatiotemporal asymmetry without intensive CFD modeling or hard-to-control laboratory experiments. Our model has three notable features: (i) we directly prescribe the level of spatiotemporal asymmetry (Sa and Ta), (ii) we model the phase lag between appendages as a continuous (rather than discrete) variable, and (iii) we explicitly include both viscous and inertial effects by using an intermediate-|$Re$| drag coefficient formulation (Shih and Buchanan 1971).

Each ctene is treated as an oscillating flat plate, whose distal end traces an ellipse (Fig. 7). The ellipse's major axis is set by the prescribed stroke amplitude (Φ), and its minor axis is set by the prescribed spatial asymmetry |$Sa$|⁠. The plate reciprocates horizontally over a total cycle period |$T\,\,$|equal to the inverse of the prescribed frequency f. Time spent in the power (⁠|$0 \le \theta \le \pi $|⁠, where |$\theta $| is the angle swept out by the tip as it moves around the ellipse) and recovery (⁠|$\pi < \theta < 2\pi $|⁠) strokes is set by the prescribed temporal asymmetry |$Ta$|⁠. This yields the parametric tip position (⁠|${x_A}( t ),\,\,{y_A}( t )$|⁠) as a function of the prescribed variables (for further details of the modeled kinematics, see supplement).

Propulsive force generated by a single plate as a function of time is given by:

$$\begin{eqnarray*}{F_p} = - \frac{1}{2}\rho {C_A}\left( {\left| {{{\dot x}_A} + {{\dot x}_B}} \right|,f} \right)w{y_A} \left[ {\left( {{{\dot x}_A} + {{\dot x}_B}} \right)\left| {{{\dot x}_A} + {{\dot x}_B}} \right|} \right] \end{eqnarray*}$$

(3)

where |$\rho $| is the fluid density,|$\,\,w$| is the width of the appendage, |${\dot x_B}$| is the body swimming speed, and |${C_A}$| is the plate drag coefficient. |${C_A}$| is a function of both the plate absolute speed |${\dot x_A} + {\dot x_B}$| and frequency f, following (Shih and Buchanan 1971). Further details are given in the supplement.

To consider multiple appendages operating at a phase lag, we sum equation (3) over n plates:

$$\begin{eqnarray*} {F_{net}} = - \frac{1}{2}\mathop \sum \limits_{k = 1}^n \rho{C_A}\left( {\left| {{{\dot x}_A} + {{\dot x}_B}}\right|,f}\right)w{y_A}\left( {t + \left( {k - 1} \right)\tau } \right)\nonumber \\ \times \big[\left( {{{\dot x}_A}\left({t + \left( {k - 1} \right)\tau } \right) + {{\dot x}_B}} \right)\left| {{{\dot x}_A}\left( {t + \left( {k - 1} \right)\tau } \right) + {{\dot x}_B}} \right|\big]\ \end{eqnarray*}$$

(4)

By dephasing the position variables |${x_A}$| and |${y_A}$| by an amount |$( {k - 1} )\tau $|⁠, we introduce metachronal coordination (where |$\tau \,\, = {P_L}\,\, \cdot T$|⁠). The model does not consider inter-appendage interactions, so it is not affected by spacing between appendages.

Spatiotemporal asymmetry in a single appendage

We compute the total propulsive force (averaged over one cycle) for a single ctene on a fixed substrate |$( {\,\,{{\dot x}_B} = \,\,0} )$|⁠, modeled by Equation (3). Because of its dependence on |${x_A}$| and |${y_A}$|⁠, |${F_p}$| is implicitly a function of Sa, Ta, Φ, and f. Increasing Φ increases the magnitude of |${F_p}$| but does not change |$R{e_\omega }$|⁠; we therefore set Φ to a fixed value (⁠|$\Phi \,\, = \,\,100^\circ $|⁠). We compute |${F_P}$| for a ctene with |$l\,\, = \,\,1\,\, \rm mm$|⁠, |$w\,\, = \,\,0.5\,\, \rm mm$|⁠, and beating at frequencies corresponding to |$R{e_\omega }$| between 5 and 200 (Afzelius 1961; Matsumoto 1991). For each frequency, we continuously vary spatial and temporal asymmetry over a range of spatiotemporal asymmetry (⁠|$0.1 < Sa < 0.6$| and |$0.1 < Ta < 0.6$|⁠). To examine the effect of varying spatiotemporal asymmetry throughout the parameter space, we calculate the gradient of |${F_p}$| with respect to |$Sa$| and with respect to |$Ta$| throughout the tested |$Sa - Ta$| space. That is, we calculate two partial derivatives: |$\frac{{\partial {F_p}}}{{\partial Sa\,\,}}$|⁠, the rate of change of |${F_p}$| with respect to |$Sa$|⁠, and |$\frac{{\partial {F_p}}}{{\partial Ta\,\,}}$| the rate of change of |${F_p}$| with respect to |$Ta$|⁠. We compare the value of these two quantities over the full range of asymmetry observed in the experiments. These derivatives can measure how strongly |${F_P}$| depends on either |$Sa$| or |$Ta$|⁠; for example, if |$\frac{{\partial {F_p}}}{{\partial Sa\,\,}}$|is small but |$\frac{{\partial {F_p}}}{{\partial Ta\,\,}}$| is large, this means that force production is more sensitive to changes in temporal asymmetry than changes in spatial asymmetry.

We find that |$\frac{{\partial {F_p}}}{{\partial Sa\,\,}} \gt 0$| and |$\frac{{\partial {F_p}}}{{\partial Ta\,\,}} > 0$| always, indicating that an increase in either |$Sa$| or |$Ta$| increases propulsive force. However, |${F_p}$| is influenced differently by changes in |$Ta$| or |$Sa$| in different regions of the parameter space, depending on |$R{e_\omega }$|⁠. To measure their relative influence, we define the operator |$G( \cdot ) \equiv \frac{{{{\partial ( \cdot )}}/{{\partial Ta}}}}{{{{\partial ( \cdot )}}/{{\partial Sa}}}}$|⁠, so that |$G\,\,( {{F_p}} ) = \frac{{{{\partial {F_p}}}/{{\partial Ta}}}}{{{{\partial {F_p}}}/{{\partial Sa}}}}\,\,$|⁠. A high value of |$G( {{F_p}} )$| indicates that propulsive force is more sensitive to changes in |$Ta$|⁠; a low value of |$G( {{F_p}} )$| indicates that propulsive force is more sensitive to changes in |$Sa$|⁠. |$G( {{F_p}} )$| is therefore a summary variable that quantifies the behavior of |${F_p}$| with respect to |$Sa$| and |$Ta$|⁠.

Figure 8A shows |$G( {{F_p}} )$| averaged over |$Sa - Ta$| space as a function of |$R{e_\omega }$|for a single appendage. At lower |$R{e_\omega }$|⁠, the influence of |$Sa$| is higher than that of |$Ta$|⁠, as expected from stronger viscous effects. However, as we increase |$R{e_\omega }$|⁠, |$G( {{F_p}} )$| increases as well, making|$\,\,Ta$| the most influential asymmetry as inertia starts to dominate. Recall from Fig. 6 that in our experimental observations, |$Sa$| decreases and |$Ta$| increases as |$R{e_\omega }$| increases. This trend is corroborated by the value of |$G( {{F_p}} )$| in Fig. 8A. An increasing value of |$G( {{F_p}} )$| corresponds to a higher sensitivity of |${F_p}$| to changes in |$Ta$|⁠; in other words, it is easier to boost force production by increasing |$Ta$| than by increasing |$Sa$| as |$R{e_\omega }$| increases.

Fig. 8

(A) The ratio |$G(F_p)$| of the gradients of |$F_p$| in |$Sa-Ta$| space as a function of |$Re_\omega$|⁠, averaged over the tested parameter space. On the right, contour plots show the behavior of |$F_p$| for several values of |$Re_\omega$|⁠: (B) |$Re_\omega=5$|⁠, (C) |$Re_\omega = 60$|⁠, (D) |$Re_\omega = 145$|⁠, and (E) |$Re_\omega = 200$|⁠. High |$Re_\omega$| leads to |$Ta$| being more influential in force production than |$Sa$|⁠, an effect that increases as |$Re_\omega$| increases; for lower |$Re_\omega$|⁠, |$Sa$| more strongly influences force production.

Fig. 8B–E shows the value of |${F_P}$| for several discrete values of |$R{e_\omega }$|⁠, over the full range of |$Sa$| and |$Ta$| considered. The role of |$Sa$| vs. |$Ta$| can be inferred from the slope of the contours: a very shallow slope indicates that |${F_P}$| is very sensitive to |$Sa$|⁠, and a very steep slope indicates that |${F_P}\,\,$|is more sensitive to |$Ta$|⁠. At low |$R{e_\omega }$|⁠, spatial asymmetry is at its highest level of influence; the slope of the contours is fairly uniform across |$Sa - Ta$| space, indicating that |${F_P}\,\,$| changes slightly favor |$Sa$| over|$\,\,Ta$| in the same manner for all points in this space. As |$R{e_\omega }$| increases, the contours steepen: a change in |$\,\,Ta$| produces a higher |${F_P}$| than a similar change in |$Sa$|⁠. As is expected, the influence of |$Ta$| is weaker at lower |$R{e_\omega }$|⁠, where the flow approaches time-reversible conditions (however, the drag coefficient in our model is not appropriate for these conditions so we do not attempt to extend the model to |$R{e_\omega } \ll 1$| where the flow would be time-reversible). In general, the contours are shallowest (indicating maximal sensitivity to |$Sa)$| for low |$Ta$|⁠; the contours are steepest (indicating maximal sensitivity to |$Ta$|⁠) at high |$Ta$|⁠.

The analysis above only considers a single appendage. We may gain further insight by considering a row of appendages operating at a phase lag, and by allowing several rows to propel a simple spheroidal body. This approach also allows for the consideration of swimming efficiency. The following section describes how we incorporate Equation (4) into a simple model of a moving ctenophore.

Spatiotemporal asymmetry in metachronal rowing and one-dimensional swimming

To explore how spatiotemporal asymmetry can affect swimming efficiency, and to investigate the role of metachronal coordination, we consider several rows of oscillating plates as described by Equation (4), propelling a spheroidal body. Drag on the spheroid is modeled as

$$\begin{equation*} {F_D} = \,\, - \frac{1}{2}\rho {C_B}\left( {{{\dot x}_B}} \right)\pi {r^2}{\dot x_B}\left| {{{\dot x}_B}} \right| \end{equation*}$$

(5)

where |${C_B}$| is the intermediate-|$Re$| drag coefficient for a prolate spheroid aligned with the flow, and is a function of its speed and geometry (Loth 2008; see supplement); r is the average radius of the body's flow-normal area; and |${\dot x_B}$| is the overall swimming speed. Summing forces (propulsion, body drag, and acceleration reaction force) in the x-direction yields :

$$\begin{eqnarray*} {\ddot x_B} &+& \frac{{\rho \pi {r^2}}}{{2\left( {m + \alpha \rho V} \right)}}{C_B}\left( {{{\dot x}_B}} \right){\dot x_B}\left| {{{\dot x}_B}} \right| \nonumber \\ &+& \frac{N}{{\left( {m + \alpha \rho V} \right)}}{F_{net}}\left( {Sa,Ta,\Phi ,f,{P_L},{{\dot x}_B}} \right) = 0 \end{eqnarray*}$$

(6)

where N is the number of appendage rows (in ctenophores, |$N\,\, = \,\,8$|⁠), m is the mass of the body,|$\,\,\alpha $| is the added mass coefficient (Daniel 1983), and V is the spheroid volume. The “added mass” term (⁠|$\alpha \rho V$|⁠) accounts for displaced water that must be accelerated along with the body (Vogel 1994). This approach allows us to compute the swimming efficiency |$\eta $|⁠, which we calculate as the ratio of the output power (⁠|${P_o}$|⁠) to the input power (⁠|${P_i}$|⁠) (⁠|$\eta \equiv {P_o}/{P_i}$|⁠, see supplement). Appendix C in the supplementary material includes a verification of the predicted swimming velocities against observations of freely-swimming ctenophores; we find that despite its simplicity, the model is able to reproduce the time-varying swimming speed of a ctenophore accelerating unidirectionally from near rest to a steady speed.

We consider a modeled ctenophore whose morphometric parameters are based on the experimental measurements used to verify the predicted swimming speeds (see supplement). The geometry for the model is |${L_b} = \,\,8.18\,\, \rm mm$|⁠, |$l\,\, = \,\,0.41\,\, \rm mm,\,\,n\,\, = \,\,9,$| and |$r\,\, = \,\,3.12\,\, \rm mm$|⁠. We solve Equation (6) for a total of 40,000 combinations of |$Sa$|⁠, |$Ta$|⁠, |$R{e_\omega }$|⁠, and |${P_L}$| (see Table 5), separating each solution into two sequential events: unsteady swimming (acceleration from rest) and steady swimming. Once the model reaches a steady speed, we continue to run the calculations until the model “swims” for three times its body length (Fig. 9A). We note that oscillations in |${\dot x_b}$| are still present in the steady swimming portion for some values of |${P_L}$| due to the time-dependence of the metachronal motion.

Fig. 9

(A) Modeled swimming speed of a mid-size ctenophore |$(L_B = 8.18 mm, l = 0.41 mm, n = 9, and \ r = 3.12 mm)$| accelerating from rest at a beat frequency of |$f = 20 Hz \ (Re_\omega = 20.1)$| with |$Sa = 0.4$| and |$Ta = 0.5$|⁠. The black line shows the highest oscillations in speed |$(P_L = 0 \%)$|⁠, the blue line shows the phase lag that gives optimal efficiency for unsteady swimming (⁠|$P_L = 6.12 \%$| for this specific case) and the magenta line shows the phase lag that reduces body speed oscillations to a minimum (⁠|$P_L = 19.59 \%$| for this specific case). (B) Swimming efficiency |$\eta \equiv P_o/P_i$| (ratio of power output to power input) as a function of phase lag for both the unsteady and steady swimming period, as demarcated in panel (A). Efficiency is calculated over the entire unsteady/steady periods to obtain the results shown. Maximum efficiency (black dot) occurs at a phase lag of 6.12% during the unsteady period for this case; there is no clear maximum efficiency for the steady period. (C) Magnitude of swimming speed oscillations |$\dot{x_b}'$| as a function of phase lag. Oscillations in swimming speed are minimized for this case at a phase lag of 19.59% (black dot).

Table 5

Range and resolution of the parameter sweep study

	|$Sa$|	|$Ta$|	|$f( {Hz} )$|	|$R{e_\omega }$|	|$PL\,\,( \% )$|
Range	0.1–0.6	0.1–0.6	10–30	10–30.1	0–30
Resolution	0.1	0.025	1	1	1.5

We find that swimming efficiency |$\eta \,\,$|depends only weakly on |$R{e_\omega }$| during both unsteady and steady swimming throughout the tested parameter space; it varies less than 1% for each |$Sa - Ta$| combination across the tested range of |$R{e_\omega }$|⁠. However, efficiency does depend on the other parametrically swept variables (⁠|${P_L}$|⁠, |$Sa$|⁠, and |$Ta$|⁠). Fig. 9B shows a representative example of |$\eta \,\,( {{P_L}} )$| (for |$Sa\,\, = \,\,0.4$| and |$Ta\,\, = \,\,0.3$|⁠), where |$\eta $| is calculated over the full duration of steady swimming or unsteady swimming.

For unsteady swimming, there is a clear maximum efficiency at a slight phase lag; this peak appears for every |$Sa - Ta$| combination at a phase lag of approximately |$7.04 \pm 1.12\% $| (mean |$ \pm $| standard deviation). For steady swimming, the efficiency is largely independent of |${P_L}$| (Fig. 9B). We attribute this difference to the stronger influence of the added mass term (see Equation (6)) during the acceleration, since this term is proportional to |${\ddot x_b}$||$.$| Any energetic advantage of metachronal coordination during steady swimming is likely to be due to hydrodynamic interactions between ctenes (de Brouwer et al. 2013; Ford et al. 2019; Ford and Santhanakrishnan 2020). While our model does not consider these interactions, we can consider how phase lag |${P_L}$| affects oscillations in body speed |$\dot x_b^{\prime}$|⁠. Figure 9C shows that zero phase lag produces strong oscillations in body speed; however, |$\dot x_b^{\prime}$| is minimized for |${P_L} = \,\,19.59\% $| (for the specific combination of |$Sa\,\, = \,\,0.4$| and |$Ta\,\, = \,\,0.3$|⁠). Averaged across the tested |$Sa - Ta$| space, this minimum appears at |${P_L} = \,\,20.11 \pm 0.46\% $| (mean |$ \pm $| standard deviation). This phase lag is close to what we observe experimentally (see Table 4) and also within the optimal |${P_L}$| range reported by Ford et al. (2019) for a krill-inspired robot (⁠|$16.7\% < {P_L} < 25\% $|⁠). In our model, this corresponds to the range that produces the smallest values of |$\dot x_b^{\prime}$| (i.e., the smoothest and least “jerky” swimming), though these are not necessarily the values that maximizes swimming efficiency (⁠|$\eta $| has no clear maximum for the steady-swimming portion of the modeled trajectory). These two cases (unsteady swimming at maximum |$\eta $| and steady swimming with minimum |$\dot x_b^{\prime}$|⁠) represent cases of particular interest, which we will further explore in the next section.

We note that the swimming efficiencies shown here are artificially low, since the model does not consider hydrodynamic interactions between ctenes nor any two-dimensional aspects of the generated flows. However, this analysis highlights the importance of considering the added mass force as a significant component of propulsion, especially during bouts of whole-body acceleration; it also highlights a correlation between the range of optimal phase lags found in other studies and the range of phase lags which produced the smoothest modeled swimming speeds.

Variation of swimming efficiency over S_a - T_a space

Previously, we used a gradient-based analysis to examine the sensitivity of propulsive force |${F_P}$| to both |$Sa$| and |$Ta$| throughout |$Sa - Ta$| space. However, this analysis considered a ctene on a fixed substrate, and thus could not include swimming efficiency. We extend this gradient-based analysis here to consider how |$Sa$| and |$Ta$| affect swimming efficiency |$\eta $|⁠. The quantity |$G( \eta ) = \frac{{\frac{{\partial \eta }}{{\partial Ta}}}}{{\frac{{\partial \eta }}{{\partial Sa}}}}\,\,$|⁠, defined for all points in |$Sa - Ta$| space, represents the relative sensitivity of |$\eta $| to |$Ta$| vs. |$Sa$|⁠. Recall that |$\eta $| is weakly dependent on beat frequency f; for unsteady swimming, |$\eta $| is maximized at a moderate phase lag |${P_L}$|⁠, and for steady swimming |$\eta $| depends only weakly on phase lag |${P_L}.$| We consider a map of |$\eta $| over |$Sa - Ta$| space for two cases of interest: (A) unsteady swimming at |$f\,\, = \,\,20Hz$| and |${P_L} = \,\,6.12\% $|⁠, which maximizes |$\eta $|⁠, and (B) steady swimming at |$f\,\, = \,\,20Hz$| and |${P_L} = \,\,19.59\% $|⁠, which minimizes |${\dot x_b}^{\prime}$| (Fig. 10).

Fig. 10

Contours of swimming efficiency at |$f = 20 Hz \ (Re_\omega = 20.1)$| for (A) unsteady swimming, at |$P_L = 6.12 \%$| (which maximizes efficiency η) and (B) steady swimming, at |$P_L = 19.59 \%$| (which minimizes body speed oscillations |$\dot{x_b}'$|⁠), throughout |$Sa-Ta$| space.

In general, efficiency |$\eta $| increases if either |$Sa$| or |$Ta$| increase, with a maximum at the edge of the considered parameter space where |$Sa\,\, = \,\,0.6$| and |$Ta\,\, = \,\,0.6$| (Fig. 10A and B). However, efficiency is not dependent on |$Sa$| and |$Ta$| in the same way. This is apparent from the slope of the contours in Fig. 10, and how they vary over |$Sa - Ta$| space: very steep contours indicate higher sensitivity to |$Ta$|⁠, and very shallow contours indicate higher sensitivity to |$Sa$|⁠.

We find that |$\eta $| is more responsive to |$Ta$| than it is to |$Sa$| throughout |$Sa - Ta$| space. The gradient ratio |$G( \eta )$| is greater than 1 over nearly the entire parameter space considered, indicating that |$\eta $| is always more sensitive to changes in |$Ta$| vs. changes in |$Sa$|⁠. Interestingly, |$\eta $| is most sensitive to |$Sa$| at very high values of |$Sa$| for both unsteady and steady swimming. This means that when the stroke is very spatially asymmetric, further increases in |$Sa$| can still produce nontrivial increases in |$\eta $|⁠; there is no saturation point.

We note that the analysis in Fig. 10 is representative of all tested frequencies (⁠|$10Hz < f < 30Hz$| and |$10 < R{e_\omega } < 30$|⁠), and that |$G( \eta )$| can be inferred from the shape of the contours of |$\eta $|⁠.

Discussion and conclusions

For drag-based paddling (rowing) at intermediate Reynolds numbers, both spatial and temporal asymmetry of the rowing appendage play important roles. These asymmetries are particularly salient for metachronal rowing, which is used across a range of Reynolds numbers from the viscous-dominated to a viscous-inertial regime. Laboratory experiments on lobate ctenophores showed both qualitative and quantitative differences between the flows generated by two animals whose ctenes beat at different Reynolds numbers and different spatiotemporal asymmetries, including the presence of vortex shedding and increased two-dimensionality for the higher |$R{e_\omega }$|/lower |$Sa$| case. The higher |$R{e_\omega }$|/lower |$Sa$| case also showed increase radial momentum production, which may impact turning. Observation of ctene kinematics revealed that as |$R{e_\omega }$| increases, spatial asymmetry decreases and temporal asymmetry increases. These trends are intuitive: viscous-dominated flows are time-reversible, so a fast power stroke and a short recovery stroke would be energetically expensive without increasing the net fluid displacement. A simple mathematical model of a single appendage reinforces this conclusion, showing that force production is more sensitive to temporal asymmetry at higher Reynolds numbers. In general, the model results showed that the relative sensitivity to spatial vs. temporal asymmetry is a function of both Reynolds number and the position in the |$Sa - Ta$| space.

This simple analysis would also dictate that increasing the spatiotemporal asymmetry always increases propulsive force generation. Since the experimental data do not suggest that animals always maximize both |$Ta$| and |$Sa$|⁠, we expand the mathematical model to consider multiple rows of appendages propelling a freely-swimming model ctenophore to consider the energetic costs. Results showed that swimming efficiency |$\eta $| increases with both |$Sa$| and |$Ta$|⁠, but that |$\eta $| is always more sensitive to changes in |$Ta$| than changes in |$Sa$| (and that |$\eta $| is most sensitive to |$Sa$| when |$Sa$| is already high). It is likely that |$Sa$| is somewhat passive—that is, some ctene deformation results from the fluid-structure interaction between the beating ctene, its material properties, and the surrounding flow (Colin et al. 2020). Ctenes become longer and thicker as an animal grows (Afzelius 1961; Hernandez-Nicaise 1991; Tamm 2014), so the material properties of larger ctenes (which tend to beat at higher |$R{e_\omega }$|⁠) are different from those of smaller ctenes. It is therefore likely that some changes in |$Sa$| passively arise from variations in size and beat frequency. In other words, increasing beat frequency is correlated with increasing size, so a larger and faster animal may have lower |$Sa$| simply due to a higher bending modulus; |$Sa$| may remain low simply because it is not mechanically possible to increase it, even if doing so would increase force generation and/or efficiency. For these larger animals which are observed to beat with lower |$Sa$|⁠, increasing |$Ta$| represents a pathway to compensate for efficiency and speed losses (this is what we observe in the experimental data). The model shows that force production is generally more sensitive to |$Ta$|⁠, making it easier to compensate for |$Sa$| losses as |$R{e_\omega }$| increases. This is highly relevant for swimming animals since it is likely that |$Ta$| is more actively controlled than |$Sa$|⁠. Further study is necessary to understand if and how ctenophores may actively control both the spatial and temporal asymmetry of the ctene beat cycle, as well as how the material properties of ctenes might dictate how |$Sa$| changes passively with |$R{e_\omega }$|⁠.

Because our model does not consider hydrodynamic interactions between ctenes, it only partially captures the effects of phase lag. Nonetheless, it predicts enhanced efficiency while the animal is accelerating with a slight phase lag, which we attribute to the role of the added mass force. For steady swimming, the model results show that phase lag affects neither the overall magnitude of the propulsive force nor the swimming efficiency. We know that this is not true, and that hydrodynamic interactions lead to increased efficiency for metachronal coordination relative to synchronous stroking (Ford et al. 2019). However, even this simple model shows that there is an optimal phase lag to minimize unsteadiness in the swimming trajectory. This “smoothest trajectory” optimal phase lag was found to be similar to those which maximize efficiency in swimming animals and robots (Barlow and Sleigh 1993; Murphy et al. 2011; Ford et al. 2019).

Further exploration of the role of spatiotemporal asymmetry in metachronal rowing would greatly benefit from a computational fluid dynamics (CFD) approach and from the use of flexible robotic models for which frequency and ctene material properties can be varied. This treatment will more fully quantify the role of phase lag, and also provide a tool for investigating the impact of |$Sa$| and |$Ta$| on radial and tangential flow production, as is hinted by the results shown in Fig. 5. A better understanding of the directionality of generated flows could hold the key to understanding the remarkable maneuverability of ctenophores, as well as other metachronal swimmers—and potentially lead to new developments in bioinspired technology.

Author contribution statement

A.H. created the analytical model, analyzed the experimental data, and developed the conceptual framework for the study. M.L.B. supervised the project and conceived and designed the experiments. M.L.B. and E.S. collected the experimental data. M.M.B., W.P., and T.K. cultured and cared for animals at the Monterey Bay Aquarium and provided guidance in data collection. A.M. facilitated animal collection and maintenance at the Bermuda Institute of Ocean Sciences and provided guidance in data collection. D.M. provided equipment for experiments and guidance on initial data analysis. M.L.B. and A.H. wrote the manuscript.

Acknowledgments

The authors gratefully acknowledge Andrea Miccoli for his assistance with data collection at BIOS.

Funding

This work was supported in part by the National Academies Keck Futures Initiative #DBS5 well as the Consejo Nacional de Ciencia y Tecnologia (CONACYT) of Mexico #2019–000021–01EXTF-00633.

Data availability

All data are available upon request from the corresponding author. Associated software products can be found at https://sites.psu.edu/byronlab/software/.

Notes

From the symposium “Metachronal coordination of multiple appendages for swimming and pumping” presented at the annual meeting of the Society for Integrative and Comparative Biology virtual annual meeting, January 3–February 28, 2021.