QUICK SEARCH:   [advanced] Author: Keyword(s):
Home     Help     Feedback     Subscriptions     Archive     Search     Table of Contents

First published online October 17, 2008
Journal of Experimental Biology 211, 3490-3503 (2008)
Published by The Company of Biologists 2008
doi: 10.1242/jeb.019224

This Article Summary Figures Only Full Text (PDF) Supplementary Material Alert me when this article is cited Alert me if a correction is posted Services Email this article to a friend Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Shirgaonkar, A. A. Articles by MacIver, M. A. Search for Related Content PubMed PubMed Citation Articles by Shirgaonkar, A. A. Articles by MacIver, M. A. Social Bookmarking                         What's this?

The hydrodynamics of ribbon-fin propulsion during impulsive motion

Anup A. Shirgaonkar¹ and Oscar M. Curet¹

Neelesh A. Patankar^1,*

Malcolm A. MacIver^1,2,*

¹ Department of Mechanical Engineering, Northwestern University, Evanston, IL 60208, USA
² Department of Biomedical Engineering, R. R. McCormick School of Engineering and Applied Science and Department of Neurobiology and Physiology, Northwestern University, Evanston, IL 60208, USA

^* Author for correspondence (e-mail: n-patankar{at}northwestern.edu)

^* Author for correspondence (e-mail: maciver{at}northwestern.edu)

Accepted 14 August 2008

	Summary

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 
Weakly electric fish are extraordinarily maneuverable swimmers, able to swim as easily forward as backward and rapidly switch swim direction, among other maneuvers. The primary propulsor of gymnotid electric fish is an elongated ribbon-like anal fin. To understand the mechanical basis of their maneuverability, we examine the hydrodynamics of a non-translating ribbon fin in stationary water using computational fluid dynamics and digital particle image velocimetry (DPIV) of the flow fields around a robotic ribbon fin. Computed forces are compared with drag measurements from towing a cast of the fish and with thrust estimates for measured swim-direction reversals. We idealize the movement of the fin as a traveling sinusoidal wave, and derive scaling relationships for how thrust varies with the wavelength, frequency, amplitude of the traveling wave and fin height. We compare these scaling relationships with prior theoretical work. The primary mechanism of thrust production is the generation of a streamwise central jet and the associated attached vortex rings. Under certain traveling wave regimes, the ribbon fin also generates a heave force, which pushes the body up in the body-fixed frame. In one such regime, we show that as the number of waves along the fin decreases to approximately two-thirds, the heave force surpasses the surge force. This switch from undulatory parallel thrust to oscillatory normal thrust may be important in understanding how the orientation of median fins may vary with fin length and number of waves along them. Our results will be useful for understanding the neural basis of control in the weakly electric knifefish as well as for engineering bio-inspired vehicles with undulatory thrusters.

Key words: aquatic locomotion, vortex shedding, propulsion, vortex rings

	INTRODUCTION

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 
Weakly electric knifefish have been studied for several decades to gain insights into how vertebrates process sensory information (for reviews, see Bullock, 1986; Turner et al., 1999). Knifefish continuously emit a weak electric field, which is perturbed by objects that enter the field and distort the field because their electrical properties differ from the surrounding fluid. These perturbations are detected by thousands of electroreceptors on the surface of the body. MacIver and coworkers (Snyder et al., 2007) have shown that the knifefish is able to both sense and move omnidirectionally. The primary thruster that weakly electric fish use in achieving their remarkable maneuverability is an elongated anal fin (Fig. 1), which we generically refer to as a ribbon fin.

View larger version (71K): [in this window] [in a new window]   Fig. 1. Apteronotus albifrons, the black ghost knifefish of South America. (A) Photograph courtesy of Neil Hepworth, Practical Fishkeeping magazine. Inset shows side view of the fish to illustrate the angle of the fin with respect to the body. (B) Frame of video of the fish swimming in a water tunnel, from ventral side. (C) Model of the ribbon fin, shown together with symbols for fin length (Lfin), fin height (hfin) and two of the kinematic state variables, angular deflection () and the vector from the rotational axis to a point on the fin (rm). The Eulerian (fluid) reference frame (x, y, z; positive direction indicated by inset) is shown together with names for velocities in the body frame fixed to the rigid fin rotation axis indicated by the red line (surge, heave, sway, roll, pitch, yaw; positive direction indicated by arrows). In this work, simulations were performed with the body-fixed frame oriented to the Eulerian frame as illustrated, and forces are with respect to a traveling wave moving in the head-to-tail direction indicated.

Knifefish swim by passing traveling waves along the ribbon fin. The waveform is often similar in overall shape to a sinusoid (Fig. 1B). The body is typically held straight and semi-rigid while swimming (Fig. 1B), i.e. body deformations are small compared with fin deformations. This may facilitate sensory performance as the body houses the electric field generator, and movement of the tail causes modulations of the field that are more than a factor of 10 larger than prey-related modulations (Chen et al., 2005; Nelson and MacIver, 1999). Knifefish frequently reverse the direction of movement without turning by changing the direction of the traveling wave on the fin and are as agile swimming backward as they are swimming forward (Blake, 1983; Lannoo and Lannoo, 1993; MacIver et al., 2001; Nanjappa et al., 2000).

The ability to switch movement direction rapidly (in 100 ms) (MacIver et al., 2001) is integral to several behaviors. Previous work by MacIver et al. (MacIver et al., 2001) has shown that prey are usually detected while swimming forward, after the prey has passed the head region, and the fish then rapidly reverses the body movement to bring the mouth to the prey during the prey strike. During inspection of novel objects, the fish are observed to engage in forward–backward scanning motions (Assad et al., 1999), which may be important for increasing spatial acuity (Babineau et al., 2007).

Ribbon-fin-based swimming is commonly referred to as the gymnotiform mode by Breder (Breder, 1926). In addition to being agile, prior research on ribbon-finned swimmers has also suggested that they are highly efficient for movement at low velocities (Blake, 1983; Lighthill and Blake, 1990). This claim is supported by the discovery that these fish use half the amount of oxygen per unit time and mass as non-gymnotid teleosts (Julian et al., 2003).

Two goals motivate the current study. First, in order to advance from the mature understanding we have of sensory signal processing in weakly electric knifefishes to an understanding of how these signals are processed to control movement, we need to characterize the hydrodynamics of ribbon-fin propulsion. Second, artificial ribbon fins may provide a superior actuator for use in highly maneuverable underwater vehicles for applications such as environmental monitoring (Epstein et al., 2006; MacIver et al., 2004).

We use computational fluid dynamics to examine the flow structures and forces arising from a sinusoidally actuated ribbon fin. We compare the computed flow structures with those measured from a robotic ribbon fin using digital particle image velocimetry (DPIV) and compare the computed surge force with the drag force measured from towing a cast of the fish. Whereas tow drag can provide a useful estimate of the thrust needed during steady swimming, for the impulsive motions modeled in this study, the thrust needed to undergo typical accelerations is more directly relevant. Thus, we also compare computed forces with the thrust that we estimate is needed for two different types of swimming direction reversals: reversals that occur during prey capture strikes from kinematic data collected in a previous study (MacIver et al., 2001) and reversals that occur during refuge tracking behavior, where fish placed in a sinusoidally oscillating refuge will move to maintain constant position with respect to the refuge.

For the present study, we idealize the fin kinematics as a traveling sinusoid on an otherwise stationary (i.e. non-translating, non-rotating) membrane (Fig. 1C). As a consequence, the top edge of the fin remains fixed at all times, and all points on the fin below this edge move in a sinusoidal manner. The fin deformation is specified in Eqn 1 below. As indicated in Fig. 1C, positive surge is defined as the force on the fin from the fluid in the direction from the tail to the head. If the traveling wave passes from the tail to the head, then the force on the fin from the fluid would be from the head to tail, corresponding to negative surge. Positive heave is vertically upward.

We chose to characterize the hydrodynamics of a fixed fin in a stationary flow because this is most relevant for understanding forces arising from maneuvering movements that occur when the body is at near-zero velocity with respect to the fluid far away from the body. In future work, using this approach will also allow us to compare our simulated force estimates with those obtained empirically from a robotic ribbon fin placed on a linear track pushing against a force sensor (Epstein et al., 2006). In subsequent studies, we will be examining the hydrodynamics of a stationary fin under imposed flow conditions and when the fin and an attached body are allowed to self-propel through the fluid.

Flow visualizations from computational simulations and DPIV indicate that the mechanism of thrust generation is a streamwise central jet and associated attached vortex rings. We show that, despite the lack of cylindrical symmetry in the morphology of the fin, its peculiar deformation pattern – the traveling wave – produces a jet flow often found in highly symmetric animal forms, such as jellyfish and squid. Whereas previous research focused exclusively on the surge force (Blake, 1983; Lighthill and Blake, 1990), we find that the ribbon fin is also able to generate a heave force, which pushes the body up. This arises from the generation of counter-rotating axial vortex pairs that are shed downward and laterally from the bottom edge of the fin. We hypothesize that the slanted angle of the fin base with respect to the spine observed in many gymnotids (e.g. Fig. 1A) leverages this heave force for forward translation. We also find that, in certain cases, as the number of waves on the fin decreases to below approximately two-thirds, the heave force surpasses the surge force. This switchover from an undulatory parallel thrust mode to an oscillatory normal thrust mode may provide insight into how the position and orientation of median fins varies with the length of the fin and the number of wavelengths that can be placed on it.

We show how the surge force from the fin scales as a function of a few key parameters. We found that for a stationary fin without imposed flow: (1) the surge force is proportional to (frequency)²x (angular amplitude)^3.5x(fin height/fin length)^3.9x(wavelength/fin length), where is a function that approximates the variation of surge force as a function of wavelength normalized with fin length; (2) for angular deflections above =10 deg., where is defined in Fig. 1C, previous analytical work (Lighthill and Blake, 1990) underestimated the magnitude of surge force; (3) surge force shows a peak when the wavelength is approximately half the fin length, similar to what is observed biologically (Blake, 1983) and contrary to the monotonic increase in force with wavelength predicted by the analytical results of Lighthill and Blake (Lighthill and Blake, 1990); and (4) the computed surge force compares well with empirical estimates based on fin kinematics and accelerations during swimming direction reversals, as well as with body drag measurements.

	MATERIALS AND METHODS

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 
Fin model for fluid simulations
As illustrated in Fig. 1B,C, the fin can be approximated by a membrane of infinitesimal thickness with an imposed sinusoidal traveling wave. The angular position (x,t) of any point on the fin, as illustrated in the figure, is given by:

(1)

where _max is the maximum angular deflection of the fin from the midsagittal plane, is the wavelength, f is frequency, x is the coordinate in the axial direction and t is time. Differentiation of Eqn 1 with respect to time leads to the velocity u_fin of any point on the fin. It is given by:

(2)

where r_m=|r_m| (see Fig. 1C) is the perpendicular distance of a point on the fin from the axis of rotation, and j and k are unit vectors in the y and z directions, respectively.

Computational fluid dynamics algorithm
Our objective is to solve the fluid flow around a fin that is moving with prescribed kinematics and simultaneously solve the forces exerted on the fin by the surrounding fluid. The fluid flow due to the fin is computed by using the immersed boundary technique (e.g. Fadlun et al., 2000; Mittal and Iaccarino, 2005; Peskin, 2002). To compute the velocity field at the n^th time step starting from the solution at the end of the (n–1)^th time step, first we compute an intermediate (i.e. prior to the application of the immersed boundary condition) velocity field, û, from the Navier–Stokes equations:

(3)

where p is pressure, µ is dynamic viscosity, is fluid density and is the gradient operator. The incompressibility constraint û=0 gives rise to the Poisson equation for pressure:

(4)

The entire computational domain is discretized by using a structured Cartesian uniform grid. Periodic boundary conditions are used along with a domain sufficiently large so as not to allow hydrodynamic interactions between the periodic images to contaminate the solution significantly. This is explained in the Validation and verification section in the Results. The spatial derivatives are discretized using sixth order compact finite difference schemes with high spectral resolution (Lele, 1992) with the optimized coefficients proposed by Lui and Lele (Lui and Lele, 2001). Time advancement is done with a low storage, low dissipation and dispersion Runge-Kutta (RK) scheme of fourth order (LDDRK4) (Stanescu and Habashi, 1998). The Poisson equation for pressure, Eqn 4, is solved using fast Fourier transforms.

In the second step of the algorithm, the intermediate velocity u is corrected at the location of the fin. To that end, we set:

(5)

only at the location of the fin, and:

(6)

at all other locations. The correction equation above is equivalent to imposing an external force on the fluid [the direct forcing approach discussed in Fadlun et al. (Fadlun et al., 2000)]. It imposes an internal boundary condition that ensures that the points in the flow field that are coincident with the fin location move with the imposed velocity of the fin. As u is not necessarily divergence free, it is projected onto a divergence-free velocity space by using a potential as:

(7)

where is obtained by solving the Poisson equation:

(8)

This completes the two-step algorithm at a given time step. The same procedure is repeated at each time step for the desired total simulated time, typically the time course of several traveling waves passing down the length of the fin.

Force on the fin is computed as the rate of change of momentum of the fluid in the correction step (Eqn 5). The first step of obtaining û (Eqns 3 and 4), as well as Eqns 7 and 8, conserve momentum. Hence, all of the change in the total fluid momentum comes from the immersed boundary condition of the correction step in Eqn 5. This change is computed as the difference between the total momentum before and after the correction step.

In the immersed boundary approach, a single regular grid can be used to solve the fluid and pressure equations at all times. This grid need not be body conforming. An alternative approach is to use a body-conforming method in which the grid conforms to the solid body surface. The immersed boundary approach has certain advantages over the body-conforming grids approach. In the latter case, as the fin shape is changing with time, remeshing is often necessary at every time step. This adds to the computational cost. Also, the use of fast elliptic solvers (e.g. fast Fourier transform-based) for pressure is not straightforward for complex meshes. Lastly, in the immersed boundary approach, the velocity correction step does not require significant computer time.

The correction step (Eqn 5) is implemented by using Lagrangian marker particles representing the fin. The prescribed fin deformation velocity is known at these marker particle locations. The velocity correction Eqn 5 would be exact if the grid points on which the fluid equations are solved coincided with the fin marker particle locations. In general, this is not true because the fin geometry is usually complex (Fig. 2). Therefore, when implementing Eqn 5, we interpolate the fin velocity from the marker particle location to its neighboring fluid grid points. The interpolation is carried out using a top-hat interpolation function as follows. For each Lagrangian marker particle, we determine which eight Eulerian grid nodes surround it. The velocity u_fin at the marker particle is then assigned to those Eulerian nodes. For nodes that have contributions from multiple marker particles, the arithmetic mean of the values from all contributing marker particles is assigned to that node. For Eulerian nodes in the fluid domain, no particle contribution results, hence the values at these nodes are left unchanged during the correction step.

View larger version (33K): [in this window] [in a new window]   Fig. 2. Close-up of a top view of the fluid-fin domain used for the numerical simulations. Shown here are the material points of the fin (black) and the Eulerian fluid grid (gray), with lengthwise grid spacing dx=0.375 mm and widthwise spacing of dy=0.4 mm. The inset shows the extent of the close-up with respect to the entire fluid-fin domain. Note that the material points do not conform to the fluid grid points.

The above numerical scheme gives the velocity and pressure fields in the fluid domain at each time step. This scheme is used to set up the ribbon-fin problem where the fin is fixed along the top edge. We consider a fin size that matches that of an adult black ghost knifefish (Fig. 1A).

Parameter identification
The physical parameters involved are the geometric and kinematic parameters of the fin, and the fluid properties. They are: L_fin, h_fin, f, _max, , , µ, which are, respectively, the length of the fin, the height of the fin (see Fig. 1), the frequency of the traveling wave, the maximum angular deflection of the fin from the midsagittal plane, the wavelength of the traveling wave, the fluid density and the dynamic viscosity. Typical values of the geometric parameters for an adult (15 cm in length) black ghost knifefish (Apteronotus albifrons) are L_fin10 cm and h_fin1 cm. These values are taken from photographs of live specimens in our laboratory. These geometric values were the ones used for all simulations performed in this study, except for one series of numerical experiments where we varied h_fin. Typical kinematic parameters are f3 Hz, _max30 deg. and 5 cm (Blake, 1983). The density and viscosity of water are taken as =1000 kg m^–3 and µ=8.9x10^–4 Pas.

Dynamic similarity tells us that the force on the fin from the fluid does not independently depend on these parameters but is a function of the dimensionless groups of these parameters. Using dimensional analysis, these dimensionless groups can be written as: (f²)/(2µ), _max, h_fin/L_fin, /L_fin.

In the following sections, we study the effect of these nondimensional parameters by changing the physical parameters as follows. We vary the Re = (f²)/(2µ), by changing f and keeping all the other physical variables fixed. Similarly, we vary h_fin/L_fin by varying h_fin, and vary the specific wavelength /L_fin by choosing different values of . The range of kinematic and geometric parameters numerically investigated is shown in Table 1. The Re for the baseline case is 894. Note that in all cases, the fin is simulated without a fish body (see Materials and methods) and is approximated as an infinitesimally thin membrane, within the smearing over 2–3 grid cells caused by interpolation on the computational grid (0.7–1.1 mm).

Computing environment
All of the simulations were performed using the San Diego Supercomputer Center's IA-64 Linux Cluster, which has 262 compute nodes each consisting of two 1.5 GHz Intel Itanium 2 processors running SuSE Linux. For interprocessor communication, the cluster uses the Myrinet 2000 gigabit ethernet interconnect network. The computational fluid dynamics code was written in Fortran 90 and C and uses the FFTW Library (www.fftw.org) for fast Fourier transforms.

Flow visualization
We used DPIV to visualize the flow field in the coronal (horizontal) plane of a robotic ribbon fin (Fig. 3). The ribbon fin was 23.5 cm long and 7.0 cm high. The experimental details are described in Epstein et al. (Epstein et al., 2006).

View larger version (21K): [in this window] [in a new window]   Fig. 3. Illustration of the digital particle image velocimetry (DPIV) setup used to visualize the flow field in the coronal (horizontal) plane of the robotic ribbon fin. The laser sheet was approximately 5 mm below the fin. The DPIV setup includes a charge couple device (CCD) camera to record the particle motion, a Nd:YAG laser light source, water tank and the robotic ribbon fin.

The time interval between image pairs was 750 µs, with image pairs obtained at a rate of 15 Hz. The laser sheet was approximately 5 mm below the bottom edge of the ribbon fin. The velocity field was obtained by analyzing each image pair with commercial particle velocimetry software (Insight, TSI). These velocity fields were then post-processed and analyzed using MATLAB (The Mathworks, Natick, MA, USA).

For the flow visualization, the robotic ribbon fin was operated at f=1.5 Hz, _max=25 deg., /L_fin=1.4 and h_fin/L_fin=0.3.

Drag measurements and analysis of body reversals
Drag
An accurate urethane cast of a 185-mm long Apteronotus albifrons was used from a previous study (MacIver and Nelson, 2000). It was bolted to a rigid rod suspended from a custom force balance that used three miniature beam load cells (MB-5-89, Interface, Scottsdale AZ, USA). For force balance and calibration details, see Ringuette et al. (Ringuette et al., 2007). The fish cast was then towed through a tow tank that was 450x96x78 cm in length, width and depth (GALCIT towtank, Caltech, Pasadena, CA, USA) using a gantry system driven by a speed-controlled DC servomotor above the tank. Details of the towtank and gantry system can be found in Ringuette et al. (Ringuette et al., 2007). Trials were conducted at three speeds: 10, 12 and 15 cm s^–1. Drag was measured with the cast at two different orientations: (1) head-first, with the long axis of the body (spine) parallel to the oncoming flow (pitch=0 deg.), and (2) tail-first, with the spine parallel to the oncoming flow (pitch=0 deg.). Only the data collected after the startup force transient had settled were analyzed, until just before the end of the towing distance (300 cm). The data were filtered with a digital Butterworth low-pass filter (cutoff at 5 Hz) to remove transducer transients prior to further statistical analysis.

Body reversals
We estimated the thrust needed for two types of swimming direction reversals: (1) rapid reversals made during prey strike (MacIver et al., 2001) and (2) slower reversals that occur during refuge tracking behavior (Cowan and Fortune, 2007). These thrust estimates were based on motion capture data from a prior study (MacIver et al., 2001), measurements of refuge tracking movements and fin kinematics, and the effective mass of the fish (mass and added mass).

During prey-capture reversals, the body is initially moving forward at approximately 10 cm s^–1 while searching for prey. Following prey detection, the traveling wave on the fin reverses and the fish rapidly decelerates to reverse the direction of body movement and bring the mouth to the prey. We analyzed reversal accelerations across 116 prey-capture trials using methods described previously (MacIver et al., 2001). Trials within one-half of one standard deviation (number of trials N=45) were selected for quantification of the ribbon-fin traveling wave frequency immediately after the moment of body direction reversal.

In refuge tracking behavior, a fish within a refuge (such as a transparent plastic tube) will attempt to maintain a fixed relationship with respect to the refuge. Thus, if the tube is moved sinusoidally, the fish will too. The kinematics of this behavior have been analyzed previously (Cowan and Fortune, 2007) but without measurement of traveling wave frequency and for a different species of knifefish. Thus, to collect preliminary data for comparison with our computed forces, we placed an adult Apteronotus albifrons underwater in a tube suspended from a custom XY robot described elsewhere (Solberg et al., 2008). Digital video recordings were made of the refuge and fish while the refuge was moved sinusoidally with an amplitude of 5 cm at frequencies of 0.25–0.4 Hz. The video was manually inspected to estimate the fin frequency immediately following body reversals.

Force estimates
The geometry of the fin simulated in the present study approximates that of the fish examined in the prey capture and refuge tracking behaviors. We estimate that the ribbon fin was moved with amplitudes similar to 30 deg., with a /L_fin of 0.5 and a L_fin of 10 cm. To approximate the thrust from the fin, we use Eqn 9 (to be introduced later) using these parameters and the traveling wave frequency measured from the video.

To estimate the thrust needed to accomplish the reversal, we multiply the effective mass previously estimated for the knifefish (10.4 g) (MacIver et al., 2004) by the accelerations extracted from the motion capture data of MacIver et al. (MacIver et al., 2001).

	RESULTS

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 
Validation and verification
For our mixed Eulerian–Lagrangian immersed boundary approach, an array of validation tests was performed to ensure that both the Eulerian (for the fluid) and Lagrangian (for the fin) grid resolutions are sufficient to obtain accurate results. The flow solver was validated using test cases that allow comparison with an analytical solution or experimental data [for details, see Shirgaonkar and Lele (Shirgaonkar and Lele, 2007)]. Sensitivity tests with respect to the number of Lagrangian particles per grid cell (N_p) representing the fin showed that although the velocity field away from the body is relatively insensitive to N_p, N_p=8 is needed to obtain accurate fluid-solid surface forces. This is consistent with the previous observation that two particles are needed in each direction per grid cell to prevent `leakage' of fluid through the membrane (Peskin, 2002).

To validate the immersed boundary implementation, we simulated an impulsively started flat plate normal to the flow. For such a plate, flow separation occurs past the plate at both of its edges. This leads to the formation of twin vortices behind the flat plate. The development of such vortices was studied experimentally by Taneda and Honji (Taneda and Honji, 1971). For validation, we choose their case with Re=126 to compare with their experimental results. In this case, the Re is defined as Ud/v, where d is the length of the plate, U is the constant translational velocity of the plate and v is the kinematic viscosity. Fig. 4 depicts how the vortex length normalized by the plate length, s/d, grows as a function of time and compares the simulation results with those of Taneda and Honji (Taneda and Honji, 1971). The vortex length is defined as the distance of the rear stagnation point from the plate (see Fig. 4). Although this validation test was performed at Re=126, as a part of verification of the code, we performed a grid convergence study at the actual fish Re of 5600 based on the wavelength and the wave speed. This ensures that the effect of higher Re, namely the thinning of boundary layers, is captured numerically.

View larger version (5K): [in this window] [in a new window]   Fig. 4. Vortex length as a function of time for the test case of the impulsively started flat plate. s is the vortex length, d is the plate length, U is the plate velocity, t is the time and T is the non-dimensional time (Taneda and Honji, 1971).

View larger version (20K): [in this window] [in a new window]   Fig. 5. Forces on the fin from the fluid with coarse (blue), nominal (green) and fine (red) grid resolution. Colored solid lines show forces filtered to remove numerical noise at the grid-scale, and light gray lines are the raw data.

The domain size used can potentially affect the force values due to the periodic boundary conditions used. To examine whether the above domain is sufficiently large, two different simulations with domain sizes L_x, L_y, L_z=15 cm, 4 cm, 4 cm and 17 cm, 6 cm, 6 cm were carried out. The forces on the fin with the two domains differed by less than 8%. Hence, we chose the shorter of the two domains above to minimize cost while maintaining adequate accuracy.

Flow visualization
Computational results
The vortex structure around the ribbon fin is shown in Fig. 6A and B where the isosurface of the total vorticity magnitude 30 s^–1 is shown and it is colored by the x-velocity. Movies 1 and 2 in supplementary material show the time evolution of the central jet underlying these vortex ring structures. Fig. 6C depicts the same vorticity isosurface together with the x-velocity isosurface u=2 cm s^–1. This shows the correlation of the high x-velocity regions and the vortex tubes. Fig. 7 shows x-velocity contours on a horizontal slice slightly below the lower edge of the fin and also the velocity vectors in the plane of the slice. The central axial jet along the surge direction is evident from the velocity vectors. The horizontal rolls of surge vorticity are seen in Fig. 8, which in later sections will be shown to be related to the successive upward and downward y-velocity near the fin surface (Fig. 9).

View larger version (46K): [in this window] [in a new window]   Fig. 6. Vortex ring structure and x-velocity around the ribbon fin. (A) Top view, and (B) front end view of isosurface of vorticity magnitude 30 s–1, colored by the x-velocity normalized by the wave speed V (in this case V=10 cm s–1). The primary and secondary vortex rings are indicated in A (see Fig. 14A for a schematic of the primary vortex rings). (C) Bottom view, spatial correlation between an isosurface of the same vorticity magnitude as in A (yellow) and an isosurface of x-velocity at 2 cm s–1 (blue). Primary vortex rings are seen to connect across the fin surface. These results are for t=1.0 s and the baseline case: max=30 deg., f=2 Hz, /Lfin=0.5, hfin/Lfin=0.1. Movies 1 and 2 (supplementary material) show the evolution of the streamwise jet from perspective and bottom views, respectively, for the baseline case.

View larger version (64K): [in this window] [in a new window]   Fig. 7. Velocity field on a horizontal slice slightly below the lower edge of the fin (top view). The fin is shown in translucent gray color. Fin wave travels to the right. The arrows represent the velocity field and the color map indicates the x-velocity. Snapshot of flow at same time and wave parameters as for Fig. 6.

View larger version (28K): [in this window] [in a new window]   Fig. 8. Isosurfaces of positive (yellow) and negative (blue) axial vorticity (x), side view. Note the counter-rotating axial vortex pairs attached to the tip of the ribbon fin (see Fig. 14B for a simplified schematic). Snapshot of the flow at the same time for the same wave parameters as for Fig. 6.

View larger version (51K): [in this window] [in a new window]   Fig. 9. Horizontal slice of y-velocity (top view). Blue indicates velocity out of the plane of the paper, producing a net upward heave force. The fin is shown in translucent gray color. The fin wave travels to the right. Snapshot of the flow at the same time for the same wave parameters as for Fig. 6.

View larger version (117K): [in this window] [in a new window]   Fig. 10. Digital particle image velocimetry (DPIV) results for the robotic ribbon fin at three different time instances. The ribbon fin is shown in translucent gray. The colormap indicates the x-velocity normalized by the wave speed V. In this case V=0.5 m s–1. (A) DPIV image and computed flow field at t=2.67 s (after four complete wave cycles). (B1), (B2) and (B3) flow field for the robotic ribbon fin at t=2.67, 2.80 and 3.00 s, respectively, equivalent to 0, 20 and 50% of one wave cycle. Velocity values at and above 50% of the wave speed are mapped to red to highlight the streamwise jet, which includes velocities as high as 90% of the wave speed.

Forces on the ribbon fin
Fig. 11 shows the temporal behavior of the surge, heave and sway forces (see Fig. 1 for definitions) from the simulation of the baseline case where f=2 Hz, _max=30 deg., /L_fin=0.5 and h_fin/L_fin=0.1. Surge and heave forces oscillate at twice the fin frequency, and sway force varies with the fin frequency. This is because at any given axial cross-section of the fin, the fluid is being pushed by the fin in the downward and streamwise directions at both extremities of the fin oscillations. The two peaks of surge forces in one fin period are unequal. This is a result of the asymmetric wake structure observed in Fig. 6C. Asymmetric wakes have also been shown to exist in an even simpler system consisting of an oscillating cylinder in cross-flow (Williamson and Roshko, 1988). Williamson and Roshko show that the wake structure (and hence forces) can significantly change with changes in system parameters (frequency and amplitude) or the flow history (Williamson and Roshko, 1988). Hence, to examine the sensitivity of our results, we varied the initial phase of the fin displacement from 0 deg. to 180 deg. in steps of 90 deg. It was found that the time-averaged mean forces are robust against the changes in the initial phase. In addition, the forces we report later do not show any sudden change in behavior with system parameters, leading us to believe that the nature of the vortex shedding in the wake is not sensitive to the system parameters for the range of parameters studied here.

View larger version (7K): [in this window] [in a new window]   Fig. 11. Computed fin forces as a function of time for the baseline case, f=2 Hz, max=30 deg., /Lfin=0.5, hfin/Lfin=0.1.

View larger version (14K): [in this window] [in a new window]   Fig. 14. Schematic of the primary flow features constituting the mechanisms of surge and heave force generation. (A) Surge: streamwise jet with the associated vortex rings. Position of vortex rings do not necessarily follow a regular pattern. See also Fig. 6C and Fig. 7. (B) Vortex tubes created along the edge of the ribbon fin shown in Fig. 8. Yellow and blue indicate the sign of vorticity, which alternates according to the direction of the fin tip movement.

View larger version (18K): [in this window] [in a new window]   Fig. 12. Trends in surge and heave forces with respect to frequency, maximum angular deflection, fin aspect ratio and specific wavelength. In each panel, the parameter shown on the x-axis is varied while keeping other parameters of the baseline case constant (max=30 deg., f=2 Hz, /Lfin=0.5, hfin/Lfin=0.1). The gray area shows the root mean square deviation of the force for one cycle. (A) Simulation Set 1, (B) simulation Set 2, (C) simulation Set 3 and (D) simulation Set 4; surge to heave crossover is at /Lfin=1.5. The estimated error based on grid sensitivity and domain length sensitivity is 10% of the mean value.

View larger version (4K): [in this window] [in a new window]   Fig. 13. Tow drag of a 18.5 cm cast of a black ghost knifefish at 10, 12 and 15 cm s–1. Error bars indicate standard deviation. The broken line shows drag when towing the fish tail first. Values are slightly horizontally offset to show the standard deviation bars. Solid gray line indicates thrust values estimated from Blake (Blake, 1983) for a 15 cm Xenomystus nigri, a knifefish with a body shape similar to that of the black ghost knifefish.

The calculation of available thrust for this case, using the scaling laws described below, gave a value of 5.5±0.1 mN. The calculation of needed thrust using effective mass and measured acceleration gave 12.6±4 mN. We address the gap between these values in the Discussion.

For the preliminary refuge tracking data, we found reversal accelerations of 106±72 cm s^–2 (N=7). Across these accelerations, the fin traveling wave frequency was 4±0.4 Hz. The calculation of available thrust gave 2±0.4 mN. The calculation of needed thrust using effective mass and measured acceleration gave 1±0.7 mN.

	DISCUSSION

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 
General flow features
First, we study the baseline case to examine the fluid dynamics near the fin. The flow features observed in the vicinity of the ribbon fin indicate certain mechanisms for surge and heave force generation.

Fig. 6A shows an instantaneous isosurface of vorticity magnitude colored by the axial velocity. Two distinct flow features are manifested here. (1) The fin generates a series of vortex rings on both of its surfaces. These rings are attached to the fin on both sides, and their vortex axes are located near the fin's midsagittal plane, a small distance below the fin, as depicted in Fig. 6B. (2) Associated with these crab-shaped vortex rings, there is a central jet along the fin axis, which represents the momentum imparted to the fluid by the fin along the direction of wave motion. The central jet is more prominently seen in the slice of x-velocity (Fig. 7). The central jet is also observed in the flow visualization from the DPIV results for the robotic ribbon fin (see Fig. 10).

It is generally known that vortex shedding behind the posterior of a swimming animal is the core signature of thrust generation, as found in previous studies of anguilliform swimming (Hultmark et al., 2007; Kern and Koumoutsakos, 2006; Tytell and Lauder, 2004). These shed vortices can be of different types, e.g. the `2-S' vortex tubes [i.e. two single, counter-rotating vortices for each periodic cycle (see Koochesfahani, 1989; Williamson and Roshko, 1988)] shed at the trailing edge of anguilliform swimmers such as eel (Tytell and Lauder, 2004) or vortex rings shed by jellyfish (Dabiri et al., 2005a; Dabiri et al., 2005b). Traditionally, the generation of vortex rings as a means of self-propulsion has been considered prominent mainly in animals that swim using a backward jet directly emanated from their interior surfaces (e.g. jellyfish and squid). Despite lacking the high degree of cylindrical symmetry that is found in the jet-generating surfaces of such organisms, the ribbon fin is able to produce a significantly strong central propelling jet. This is observed in Fig. 6C, where the x-velocity has a very strong spatial correlation with the vortex rings. Movies 1 and 2 in supplementary material show the temporal behavior of the central jet. Such vortex ring structure was observed in experiments by Drucker and Lauder for pectoral fins of black surfperch and bluegill sunfish, where the rings are generated by cyclic flapping motion of the pectoral fins (Drucker and Lauder, 1999). Although the rings observed in our fin simulations seem similar, they are generated by a different motion pattern – a traveling wave that is perpendicular to the flapping velocity of the fin. This suggests that vortex rings can be observed along the body surface in different modes of swimming, when a predominantly axial jet is present along the body. This ring structure is distinct from the `2-S' vortex structures observed in eels (Tytell and Lauder, 2004) or the `2-P' structures (i.e. two counter-rotating vortex pairs, one on each side of the body – left and right with reference to the swimming direction) computed in eel simulations (Kern and Koumoutsakos, 2006). The difference between the current ribbon-fin flow field and the eel wakes is that the ribbon fin generates a prominent backward central jet along the length of the body.

In general, for a swimming animal in the Eulerian regime (i.e. Re

Fig. 6A also indicates that there are secondary, smaller vortex rings that are not attached to the fin but are shed into the surrounding fluid. Associated with them are the smaller jets located at the centers of these secondary rings (henceforth referred to as the secondary jets, as opposed to the primary central jet). Both the primary and secondary jets exhibit some angular deflections from the surge direction. This means that some momentum is lost in the lateral direction. Some comments can be made regarding the biological implications of this observation. For animals needing rapid maneuverability in order to capture prey (e.g. black ghost knifefish), the ability to exchange lateral momentum with the fluid is essential. The morphological characteristics required for this purpose can entail an undesirable energy loss in the lateral direction even in the cruising mode. This is also observed in numerical simulations of eel swimming with a 2-P vortex structure (Kern and Koumoutsakos, 2006), where the wake has an angular shape, reflecting energy loss in the lateral direction. By contrast, for low-maneuverability jellyfish, the jet is almost unidirectional in the surge direction (Dabiri et al., 2005a; Dabiri et al., 2005b). This would minimize the lateral energy loss but, as a consequence, the ability to generate rapid lateral forces would be diminished.

Heave generation can be understood from the streamwise vorticity (_x) shown in Fig. 8. The opposite surfaces of the fin produce opposite-signed vorticity. It then separates from the fin surface, rolls into streamwise vortex tubes at the bottom edge of the fin and is advected in the lateral direction (z direction in Fig. 1) by the fin motion. The induced velocity of such axial vortex pairs is vertically downward, imparting a downward momentum to the fluid. The resultant upward reaction on the fin causes the heave force. The downward jet from these rolls is depicted in the vertical velocity contours (Fig. 9).

The mechanisms of surge and heave generation that are inferred from Figs 6 and 8 are schematically depicted in terms of the streamwise jet and attached vortex tubes in Fig. 14A and B, respectively. The downward motion of the vortex tubes was observed only when the heave generation is significant, i.e. above _max=20 deg. Below this angle, the fin was found not to generate any significant heave (Fig. 12B). Flow at _max=20 deg. showed that these tubes are attached to the lower end of the fin, with little or no downward shedding, indicating they are weak compared with those in the higher heave case.

The total forces on the fin are shown in Fig. 11. The trends shown by these forces will be discussed in the next section.

Trends
Magnitude of surge force
Our computational results show that the surge force from the simulated ribbon fin varies between 0 and 1.85 mN over the parameter range examined. According to the measurements by Blake (Blake, 1983), a frequency of 3Hz with _max of 30 deg. and a specific wavelength of approximately 0.5 is typical of gymnotids. In our results, this regime results in a force of 1 mN (Fig. 12A). We compare this result with three empirical findings: estimated thrust needed during two kinds of body reversals and the measured drag force on the body. The most relevant comparison for our simulations, occurring with a non-translating fin in stationary fluid, is with body reversals where the body is momentarily non-translating.

Our estimate of the magnitude of surge available during the rapid body reversals typical of prey capture behavior, based on the 7 Hz traveling wave frequency that was measured (see Results), is 6 mN. Our estimate of the necessary surge force to achieve the measured accelerations, using an estimate of effective mass, is 13 mN. We hypothesize that the gap between these estimates is due to the presence, during the prey capture rapid reversals, of 1–3 rapid bilateral pectoral fin strokes. It should also be noted that our simulations do not account for the acceleration of the surrounding fluid during traveling wave reversals because we are assuming a stationary fluid.

During the refuge tracking behavior body reversals, no such bilateral pectoral fin strokes were observed. In this case, the estimated available surge force is 2 mN whereas the estimated necessary surge force based on measured acceleration and effective mass is 1 mN. In this case, the measured fin frequency was 4 Hz. This compares favorably with our computed thrust for similar kinematics.

Finally, the body tow drag measure enables us to estimate the necessary thrust capacity of the fin under steady swimming conditions. Our drag measurement, at typical steady-state swimming velocities, was 1–2 mN (Fig. 13), similar to our computed surge force. Prior theoretical estimates of the thrust of a similar ribbon-finned fish, Xenomystus nigri, are also in this range (gray line of Fig. 13) (Blake, 1983).

Effect of frequency
Fig. 12A shows the mean surge and heave forces on the ribbon fin for different frequencies. It shows that both surge and heave forces have a quadratic dependence with frequency. As the velocity scale is V=f, the force varies as V²A, where A is the characteristic area of the fin, which can be interpreted as the total wetted area. This indicates that the force generation is essentially inertial in nature at these moderate Re. Nevertheless, viscous effects may not be negligible as they change the flow in the vicinity of the fin substantially so as to affect the overall flow dynamics.

Effect of _max
The fin mean surge and heave forces as a function of _max are shown in Fig. 12B. In this case, the surge force follows a power law curve, whereas heave force becomes significant beyond around _max>20 deg. Similar trends have been shown in an oscillating airfoil study (Koochesfahani, 1989). In this work, Koochesfahani shows that an airfoil could produce drag (negative heave) or thrust (positive heave) depending on the frequency and amplitude of oscillation. In addition, he shows that the frequency at which the airfoil produces thrust depends on the amplitude of oscillation. Moreover, in previous work (see Blake, 1983), it has been reported that gymnotiform swimmers do not substantially change wave amplitude with swimming speed. Based on the results of the present study, we speculate that this may be to ensure an upward heave force. The empirically observed wave amplitude from a black ghost knifefish during swimming is approximately 30 deg. (Blake, 1983) above the 20 deg. threshold found above.

Effect of h_fin/L_fin
In the numerical simulations for Set 3, the height of the ribbon fin was varied while keeping the rest of the parameters fixed. As in all the previous cases, surge force is greater than the heave force and both forces follow a power law relationship with fin height (Fig. 12C). For the shortest h_fin, the heave force obtained was insignificant compared with the surge force, and for all other fin heights the heave force was positive.

Effect of specific wavelength
In the numerical simulations used to study the effect of specific wavelength, f and _max were kept constant. Increasing the wavelength produces two competing effects in force production: (1) an increase in the wave speed (V=f) and (2) a decrease in the total wetted area of the fin caused by a reduced number of waves on the fin (L/). Whereas increasing wave speed produces higher force, a decrease in the wetted area of the fin produces smaller forces. The interaction between these competing effects may be the basis of the peak in the surge force with specific wavelength shown in Fig. 12D.

For the range of traveling wave parameters that we considered, the optimal specific wavelength lies between 0.5 and 0.8. This is close to the empirically observed value of 0.4 in gymnotids (Blake, 1983). It should be noted that this optimal range is only with respect to maximization of surge force. It does not consider hydrodynamic energy loss due to dissipation or turbulence or physiological energy consumption in the muscles/appendages that generate the fin deformation. Our current goal is to consider the hydrodynamic effects only and relate the force generation mechanism both qualitatively and quantitatively to various kinematic parameters of the fin motion.

At /L_fin=1.5, and other variables set to the baseline case, the surge force component is equaled (and for larger values of /L_fin, surpassed) by the heave force component. This crossover point marks a distinct shift in the mode of propulsion. Below this critical value, the motion of the fin is predominantly undulatory and surge force dominates whereas above this value the motion of the fin is predominantly flapping (or oscillatory) and heave force dominates. The transition point between the undulatory and flapping propulsion modes will vary with the values of the other kinematic parameters.

Although we have not observed /L_fin near the crossover point in normal swimming of knifefish, we have observed nearly vertical translations during rare startle events (M.A.M. unpublished observations). Our results indicate that this motion may be accomplished with large values of /L_fin.

Heave force and fish morphology
During forward cruising, the sum of the force components from a fish's propulsors should ideally be parallel to the body axis. Given the relationship between surge and heave as a function of fin deformation kinematics found in the present study, we can ask what the resultant angle is and compare this with the insertion angle of the fin. For the case of f=4 Hz, /L_fin=0.5 and _max=30 deg., the angle between the base (i.e. the upper edge) of the ribbon fin and the resultant force is approximately 11 deg. This is indeed similar to the angle of insertion of the fin with respect to the spine that we observe in Apteronotus albifrons (Fig. 1A).

The surge–heave crossover point may also provide insight into the position and orientation of short median paired fins in species, such as the triggerfish (family Balistidae). As the fin becomes shorter, for a given fin ray density, the number of fin rays is reduced. This results in large values of /L_fin. Surge and heave forces will then be comparable in magnitude. In this regime, to maximize translational thrust, the fin should have an angle of insertion similar to 45 deg. This is the case for some species of triggerfish with short median paired fins [see fig. 1 in Lighthill and Blake (Lighthill and Blake, 1990)].

Force scaling
Scaling analysis was conducted to understand the relationship between surge force generation and the traveling wave variables. For the force scaling with f, _max and h_fin/L_fin, we assumed that the surge force and the independent variable have a power–law relationship of the form g()=a^b, where g() is the computed force, is one of the independent variables listed above and a and b are constants. a and b were found, respectively, from the intercept and the slope of a linear regression fit of the log-transformed data. For specific wavelength variation, the force did not follow a power law (see Fig. 12D), most likely due to the two competing effects contributing to the force behavior as explained above in `Effect of specific wavelength'. Hence, for /L_fin we used a curve fit of another form (described below).

Our simulations give a surge force power–law relationship with exponent 2, 3.5 and 3.9 for f, _max and h_fin/L_fin, respectively (see Fig. 15A–C). Now we can express the propulsive fin force as:

(9)

where C₁ is a constant and (/L_fin) is a function of the specific wavelength which can be approximated by:

(10)

The constant 0.6 in Eqn 10 is found by curve fitting (see Fig. 15D). This function was motivated by the velocity profile of the Lamb–Oseen vortex used to model aircraft trailing vortices (Shirgaonkar and Lele, 2006), which has a similar shape to the force behavior in this case and the property of blending two physical effects from two different regimes. Based on the simulations, the constant C₁ in Eqn 9 is approximately 86.03. Note that this equation is only verified for the parameter space considered in this study and _max is in radians, is in kg m^–3, f is in s^–1, L_fin is in m and F_surge is in N.

View larger version (14K): [in this window] [in a new window]   Fig. 15. Scaling of the surge force. Computational results are shown with a dot, and the curve fit with a black line. The error is defined as the root mean square error in mN for each curve fit. Estimation of surge force for a stationary fin using eqn 11 in Lighthill and Blake (Lighthill and Blake, 1990) are show in gray. In each panel, the factor shown on the x-axis is varied while keeping other parameters of the baseline case constant (max=30 deg., f=2 Hz, /Lfin=0.5, hfin/Lfin=0.1). (A) linear fit of the natural log-transformed data for frequency (Hz), (B) linear fit of the natural log-transformed data for max (deg.) (C) linear fit of the log-transformed data for hfin/Lfin, (D) Curve fit of surge force as a function of specific wavelength.

In a series of papers, Lighthill and Blake (Lighthill and Blake, 1990) and Lighthill (Lighthill, 1990a; Lighthill, 1990b; Lighthill, 1990c) analyzed the fluid dynamics and force production for ribbon-fin swimmers. Lighthill and Blake did a theoretical analysis of a ribbon fin to obtain an analytical expression for the fin's mean propulsive force assuming a two-dimensional irrotational flow in the plane transverse to the fin surface. The expression for the mean propulsive force consists of two terms: (1) the rate of backward momentum transport and (2) the pressure force.

The rate of backward momentum transport is proportional to Umv calculated at the trailing edge of the fin, where the overbar denotes a mean over one period of the traveling wave on the fin, U is the swimming velocity, m is the added mass per unit length of an infinitely thin fin and v is the lateral fin velocity component that is perpendicular to the fin surface. Because in our simulations the fin is stationary (U=0), this term vanishes. This is exactly what will occur in a swimming motion starting from rest, or during the switch from forward to reverse swimming direction, a stereotypical maneuver for a black ghost knifefish during prey capture (MacIver et al., 2001). In these cases, there are instants when the force from the momentum transport term will tend to zero and the propulsive force generated by the pressure term will be dominant.

Lighthill and Blake state that the surge force due to pressure is proportional to mv²/2 (Lighthill and Blake, 1990). In their analysis, the added mass per unit length, m, of the infinitely thin fin is proportional to . At zero swimming speed (U=0), their expression for v is proportional to f_maxh_fincos_av. Here, _av is the mean angle made by the fin profile (e.g. see Fig. 1B) at the trailing end with respect to the baseline of the fin (Fig. 1C, red line). is zero at the fin base and maximum at the fin tip. Using terms for the added mass and fin tip velocity, the following equation is derived for the surge component of the mean pressure force at U=0:

(11)

Lighthill and Blake give the `lack of enhancement coefficient' =0.53 for a fin without body (Lighthill and Blake, 1990). To compare the scaling in Eqn 11 with that from our simulations (Eqn 9), we note that _av is a function of _max, and h_fin. Hence, in Fig. 15 we plot log-transformed data for the mean propulsive force from Eqn 11 for the same parameters used in our simulations, with no attached body. Assuming a power–law form, this gives exponents 2, 1.9 and 3.6 for f, _max and h_fin, respectively. We now compare these scaling values with our own values.

The frequency scaling from simulations and Eqn 11 is f². This contribution comes from the v² dependence of the pressure force as discussed by Lighthill and Blake (Lighthill and Blake, 1990).

The scaling for the fin height parameter, h_fin/L_fin, is also in good agreement and is close to the fourth power. The h_fin scaling arises from two distinct factors in the pressure force: (1) the added mass of the fin, proportional to , and (2) a contribution from v², which is also proportional to . The effect of fin twisting (the cos³_av term in Eqn 11) reduces the power of h_fin seen in Fig. 15C (gray line) from 4 to 3.6.

For _max, our results show a scaling power of 3.5 vs 1.9 of Lighthill and Blake (Lighthill and Blake, 1990). Note that ²_max scaling arises from the v² contribution to the pressure force. The difference between this scaling and the 1.9 scaling shown in Fig. 15B (gray line) is due to the effect of the fin twisting term (cos³_av). Our simulations give an extra factor of ^1.5_max, which we speculate more accurately resolves the effect of the fin twist on _max scaling. This could also be because of the local two-dimensional assumption by LB, which does not capture the full three-dimensionality of the flow-field around the fin.

For the specific wavelength variation, the estimation of the propulsive force with Lighthill and Blake's theory (Eqn 11) shows an initial rapid increase in surge force with specific wavelength followed by an asymptotic value at large specific wavelengths (see Fig. 15D). However, our computational results show a decrease in the propulsive force after a certain specific wavelength, contrary to the results of Lighthill and Blake (Lighthill and Blake, 1990). Such a peak in the surge force has been observed experimentally for a robotic ribbon fin mimicking a gymnotiform swimmer (Epstein et al., 2006). We attribute this mismatch between our results and those of LB to the fact that their inviscid theory does not capture the effect of decreased thrust due to decreased wetted area for increasing wavelength.

For very small angular deflections of the fin, the analytical estimates of Lighthill and Blake (Lighthill and Blake, 1990) are larger than our computed propulsive forces. But for the majority of the parameter space, Lighthill and Blake's (Lighthill and Blake, 1990) forces are lower than our results by as much as a factor of five (Fig. 15). There are several simplifying assumptions in LB's theory that may account for this mismatch: (1) a lack of hydrodynamic interaction between different parts of the fin along the surge direction, (2) the two-dimensional flow approximation and (3) the inviscid flow assumption. Our results show a clear flow interaction and axial propagation of vorticity along the ribbon fin (see Fig. 8). This interaction will affect the force-generation capability of the ribbon fin.

Furthermore, it is possible that in start-up motions, or when the swimming speed is very small, viscous forces may contribute to thrust, something that is neglected in the analytical solution of the force estimate, which relies on the inviscid flow assumption.

According to Lighthill and Blake's theory, in cruise mode, the pressure term is small compared with the overall thrust produced (Lighthill and Blake, 1990). However, in regimes with U0 (e.g. start-up motions and during rapid forward-to-reverse maneuvers), the propulsive force is dominated by the pressure contribution. Our computations show that this pressure contribution is large enough to provide sufficient thrust for refuge tracking body reversals and can play a significant role in pectoral-fin enhanced reversals in prey capture strikes. This suggests that real swimmers could in fact use this mechanism during impulsive motion.

Conclusions
The goal of this study was to understand the mechanical basis of propulsive force generation by the anal ribbon fin of a black ghost knifefish during impulsive starting movements. Using computational fluid dynamics and DPIV, we studied the key flow mechanisms underlying surge and heave force generation by an undulating but non-translating ribbon fin. We found that the mechanism of thrust generation is a streamwise central jet associated with a crab-like vortex ring structure attached to the fin. Heave force is generated for certain wave parameters by vortex tubes attached to the tip of the ribbon fin. The simulated surge force shows good agreement with the estimates for refuge tracking body reversals, as well as our drag measurements from a cast model of the knifefish. The flow features from simulations reproduce those from DPIV. Our results also provide quantitative scaling of the surge force with amplitude, frequency and wavelength of the traveling wave, and the fin height. The present high-fidelity simulations account for many of the fluid dynamical aspects not captured in the inviscid theory of Lighthill and Blake (Lighthill and Blake, 1990) and provide force mechanisms based on flow structure, qualitative trends in forces and better quantitative accuracy desirable for subsequent modeling and engineering efforts.

LIST OF ABBREVIATIONS

a, b

constants in the power–law equation

A

fin area

C₁

constant

d

flat plate length

f

frequency of the fin

F_surge

surge force

h_fin

fin height

i, j, k

unit vectors in the x, y and z directions

L_fin

fin length

L_x, L_y, L_z

computational domain size in x, y and z

m

added mass

N

number of trials

N_p

number of Lagrangian particles per grid cell

N_x, N_y, N_z

number of grid points in x, y and z

p

pressure

r_m

vector from the rotational axis to a point on the fin

r_m

distance of a point on the fin from the axis of rotation

Re

Reynolds number

s

vortex length

t

time

u

velocity vector

u_fin

velocity at any point on the fin

û

intermediate velocity field

intermediate velocity solution after velocity correction at the fin location

u^n–1

velocity field at previous time step

U

constant translational velocity of the plate

v

fin velocity perpendicular to the fin surface

V

wave speed

x, y, z

Cartesian coordinate directions

av

mean angle of a fin section with respect to the sagittal (vertical) plane

t

time step

lack of enhancement coefficient

angular position of any point on the fin

max

maximum angular deflection of the fin from the midsagittal plane

wavelength

µ

dynamic viscosity

kinematic viscosity

independent variable

fluid density

function that approximates the variation of surge force with respect to /L_fin

x

vorticity in the axial direction

divergence operator

	Acknowledgments

 
N.A.P. led the computational and theoretical fluid analysis effort; M.A.M. led the experimental fluids, robotics and all other efforts. O.M.C. and A.A.S. contributed equally. A.A.S. wrote the computational fluid dynamics code. We thank Richard Lueptow and Nick Pohlman for the DPIV setup and help with data collection and analysis. The drag measurements were made by M.A.M. with Ann Marie Polsenberg, using Joel Burdick's and Mory Gharib's facilities at Caltech, with Matt Ringuette's generous assistance. We thank two anonymous reviewers for their very helpful suggestions. This work was supported by NSF grant IOB-0517683 to M.A.M. Partial support for N.A.P and A.A.S. was provided by NSF CAREER CTS-0134546 to N.A.P. Computational resources were provided by the San Diego Supercomputer Center through NSF's TeraGrid project grant CTS-070056T to A.A.S., N.A.P. and M.A.M., and in part by a developmental grant from the Argonne National Laboratory to N.A.P.

	Footnotes

 
Supplementary material available online at http://jeb.biologists.org/cgi/content/full/211/21/3490/DC1

	References

TOP Summary INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION References

 

Assad, C., Rasnow, B. and Stoddard, P. K. (1999). Electric organ discharges and electric images during electrolocation. J. Exp. Biol. 202,1185 -1193.[Abstract]

Axelrod, H. R. (1989). Dr Axelrod's Atlas of Freshwater Aquarium Fishes. Neptune City, NJ: TFH Publications.

Babineau, D., Lewis, J. E. and Longtin, A. (2007). Spatial acuity and prey detection in weakly electric fish. PLoS Comput. Biol. 3, 402-411.

Blake, R. W. (1983). Swimming in the electric eels and knifefishes. Can. J. Zool. 61,1432 -1441.

Breder, C. M. (1926). The locomotion of fishes. Zoologica 4,159 -297.

Bullock, T. H. (1986). Electroreception. New York: Wiley.

Chen, L., House, J. H., Krahe, R. and Nelson, M. E. (2005). Modeling signal and background components of electrosensory scenes. J. Comp. Physiol. 191,331 -345.[CrossRef]

Cowan, N. J. and Fortune, E. S. (2007). The critical role of locomotion dynamics in decoding sensory systems. J. Neurosci. 27,1123 -1128.[Abstract/Free Full Text]

Dabiri, J. O., Colin, S. P., Costello, J. H. and Gharib, M. (2005a). Flow patterns generated by oblate medusan jellyfish: field measurements and laboratory analyses. J. Exp. Biol. 208,1257 -1265.[Abstract/Free Full Text]

Dabiri, J. O., Colin, S. P., Costello, J. H. and Gharib, M. (2005b). Vortex motion in the ocean: in situ visualization of jellyfish swimming and feeding flows. Phys. Fluids 17,1257 -1265.

Drucker, E. G. and Lauder, G. V. (1999). Locomotor forces on a swimming fish: three-dimensional vortex wake dynamics quantified using digital particle image velocimetry. J. Exp. Biol. 202,2393 -2412.[Abstract]

Epstein, M., Colgate, J. E. and MacIver, M. A. (2006). Generating thrust with a biologically-inspired robotic ribbon fin. Proceedings of the 2006 IEEE/RSJ International Conference on Intelligent Robots and Systems, 2006, pp.2412 -2417. October 9-15. DOI.10.1109/IROS.2006.281681

Fadlun, E. A., Verzicco, R., Orlandi, P. and Mohd-Yusof, J. (2000). Combined immersed-boundary finite-difference methods for three-dimensional complex flow simulations. J. Comput. Phys. 161,35 -60.[CrossRef]

Hultmark, M., Leftwich, M. C. and Smis, A. J. (2007). Flowfield measurements in the wake of a robotic lamprey. Exp. Fluids 43,683 -690.[CrossRef]

Julian, D., Crampton, W. G. R., Wohlgemuth, S. E. and Albert, J. S. (2003). Oxygen consumption in weakly electric Neotropical fishes. Oecologia 137,502 -511.[CrossRef][Medline]

Kern, S. and Koumoutsakos, P. (2006). Simulations of optimized anguilliform swimming. J. Exp. Biol. 209,4841 -4857.[Abstract/Free Full Text]

Koochesfahani, M. M. (1989). Vortical patterns in the wake of an oscillating airfoil. AIAA Stud. J. 27,1200 -1205.[CrossRef]

Lannoo, M. J. and Lannoo, S. J. (1993). Why do electric fishes swim backwards? An hypothesis based on gymnotiform foraging behavior interpreted through sensory constraints. Environ. Biol. Fishes 36,157 -165.[CrossRef]

Lele, S. K. (1992). Compact finite difference schemes with spectral-like resolution. J. Comput. Phys. 103,16 -42.[CrossRef]

Lighthill, J. (1990a). Biofluiddynamics of balistiform and gymnotiform locomotion 2. the pressure distribution arising in 2-dimensional irrotational flow from a general symmetrical motion of a flexible flat-plate normal to itself. J. Fluid Mech. 213, 1-10.[CrossRef]

Lighthill, J. (1990b). Biofluiddynamics of balistiform and gymnotiform locomotion 3. momentum enhancement in the presence of a body of elliptic cross-section. J. Fluid Mech. 213, 11-20.[CrossRef]

Lighthill, J. (1990c). Biofluiddynamics of balistiform and gymnotiform locomotion 4. short-wavelength limitations on momentum enhancement. J. Fluid Mech. 213, 21-28.[CrossRef]

Lighthill, J. and Blake, R. (1990). Biofluiddynamics of balistiform and gymnotiform locomotion 1. Biological background, and analysis by elongated-body theory. J. Fluid Mech. 212,183 -207.[CrossRef]

Lui, C. and Lele, S. K. (2001). Direct numerical simulation of spatially developing, compressible, turbulent mixing layers, AIAA Paper 2001-0291, Amer. Inst. Aero. Astro. www.aiaa.org.

MacIver, M. A., Fontaine, E. and Burdick, J. W. (2004). Designing future underwater vehicles: principles and mechanisms of the weakly electric fish. IEEE J. Ocean. Eng. 29,651 -659.[CrossRef]

MacIver, M. A. and Nelson, M. E. (2000). Body modeling and model-based tracking for neuroethology. J. Neurosci. Methods 95,133 -143.[CrossRef][Medline]

MacIver, M. A., Sharabash, N. M. and Nelson, M. E. (2001). Prey-capture behavior in gymnotid electric fish: motion analysis and effects of water conductivity. J. Exp. Biol. 204,543 -557.[Abstract]

Mittal, R. and Iaccarino, G. (2005). Immersed boundary methods. Annu. Rev. Fluid Mech. 37,239 -261.[CrossRef]

Nanjappa, P., Brand, L. and Lannoo, M. J. (2000). Swimming patterns associated with foraging in phylogenetically and ecologically diverse American weakly electric teleosts (Gymnotiformes). Environ. Biol. Fishes 58, 97-104.[CrossRef]

Nelson, M. E. and MacIver, M. A. (1999). Prey capture in the weakly electric fish Apteronotus albifrons: sensory acquisition strategies and electrosensory consequences. J. Exp. Biol. 202,1195 -1203.[Abstract]

Nelson, M. E., MacIver, M. A. and Coombs, S. (2002). Modeling electrosensory and mechanosensory images during the predatory behavior of weakly electric fish. Brain Behav. Evol. 59,199 -210.[CrossRef][Medline]

Peskin, C. S. (2002). The immersed boundary method. Acta Numer. 11,479 -517.

Ringuette, M. J., Milano, M. and Gharib, M. (2007). Role of the tip vortex in the force generation of low-aspect-ratio normal flat plates. J. Fluid Mech. 581,453 -468.[CrossRef]

Shirgaonkar, A. A. and Lele, S. K. (2006). Large eddy simulation of contrails: effect of atmospheric properties, AIAA Paper 2006-1414, Amer. Inst. Aero. Astro. www.aiaa.org.

Shirgaonkar, A. A. and Lele, S. K. (2007).Large eddy simulation of early stage aircraft contrails . Technical Report TF-100, Flow Physics and Computation Division, Department of Mech. Eng., Stanford University.

Snyder, J. B., Nelson, M. E., Burdick, J. W. and MacIver, M. A. (2007). Omnidirectional sensory and motor volumes in an electric fish. PLoS Biol. 5,2671 -2683.

Solberg, J. R., Lynch, K. M. and MacIver, M. A. (2008). Active electrolocation for underwater target localization. Int. J. Robotics Res. 27,529 -548.[CrossRef]

Stanescu, D. and Habashi, W. G. (1998). 2N-storage low dissipation and dispersion Runge-Kutta schemes for computational acoustics. J. Comput. Phys. 143,674 -681.[CrossRef]

Taneda, S. and Honji, H. (1971). Unsteady flow past a flat plate normal to the direction of motion. J. Phys. Soc. Jpn. 30,262 -272.[CrossRef]

Turner, R. W., Maler, L. and Burrows, M. (1999). Special issue on electroreception and electrocommunication. J. Exp. Biol. 202,1167 -1458.[Abstract]

Tytell, E. D. and Lauder, G. V. (2004). The hydrodynamics of eel swimming: I. Wake structure. J. Exp. Biol. 207,1825 -1841.[Abstract/Free Full Text]

Uhlmann, M. (2003). First experiments with the simulation of particulate flows. Madrid: Technical Report 1020, CIEMAT, Madrid, Spain. ISSN 1135-9420.

Williamson, C. H. K. and Roshko, A. (1988). Vortex formation in the wake of an oscillating cylinder. J. Fluids Struct. 2,355 -381.[CrossRef]

CiteULike    Complore    Connotea    Del.icio.us    Digg    Reddit    Technorati    Twitter    What's this?

This article has been cited by other articles:

				C. M Postlethwaite, T. M Psemeneki, J. Selimkhanov, M. Silber, and M. A MacIver Optimal movement in the prey strikes of weakly electric fish: a case study of the interplay of body plan and movement capability J R Soc Interface, May 6, 2009; 6(34): 417 - 433. [Abstract] [Full Text] [PDF]

This Article Summary Figures Only Full Text (PDF) Supplementary Material Alert me when this article is cited Alert me if a correction is posted Services Email this article to a friend Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Shirgaonkar, A. A. Articles by MacIver, M. A. Search for Related Content PubMed PubMed Citation Articles by Shirgaonkar, A. A. Articles by MacIver, M. A. Social Bookmarking                         What's this?