We carry out fluid–structure interaction simulations of self-propelled virtual swimmers to investigate the effects of body shape (form) and kinematics on the hydrodynamics of undulatory swimming. To separate the effects of form and kinematics, we employ four different virtual swimmers: a carangiform swimmer (i.e. a mackerel swimming like mackerel do in nature); an anguilliform swimmer (i.e. a lamprey swimming like lampreys do in nature); a hybrid swimmer with anguilliform kinematics but carangiform body shape (a mackerel swimming like a lamprey); and another hybrid swimmer with carangiform kinematics but anguilliform body shape (a lamprey swimming like a mackerel). By comparing the performance of swimmers with different kinematics but similar body shapes we study the effects of kinematics whereas by comparing swimmers with similar kinematics but different body shapes we study the effects of form. We show that the anguilliform kinematics not only reaches higher velocities but is also more efficient in the viscous (Re∼102) and transitional (Re∼103) regimes. However, in the inertial regime (Re=∞) carangiform kinematics achieves higher velocities and is also more efficient than the anguilliform kinematics. The mackerel body achieves higher swimming speeds in all cases but is more efficient in the inertial regime only whereas the lamprey body is more efficient in the transitional regime. We also show that form and kinematics have little overall effect on the 3-D structure of the wake (i.e. single vs double row vortex streets), which mainly depends on the Strouhal number. Nevertheless, body shape is found to somewhat affect the small-scale features and complexity of the vortex rings shed by the various swimmers.
Fish that mainly use body/caudal fin (BCF) undulations for locomotion are classified into five different types depending on the manner they swim: anguilliform, sub-carangiform, carangiform, thunniform and ostraciform (Breder, 1926; Lindsey, 1978; Videler, 1993; Webb, 1975) [for more detail on the various modes of swimming, see Sfakiotakis et al. (Sfakiotakis et al., 1999)]. Anguilliform and carangiform swimmers, which are the focus of this paper, differ from each other in both body morphology (form) and body undulations (kinematics). Anguilliform swimmers, such as the eel (Anguilla anguilla), typically have long narrow bodies and the width of the body remains almost constant from head to tail. By contrast, carangiform swimmers, such as the mackerel (Scomber scombrus), have thicker bodies with their body width decreasing at the peduncle where the body attaches to a V-shaped and typically symmetrical caudal fin of relatively high aspect ratio. The amplitude of body undulations in anguilliform swimmers is large over the whole body length whereas for carangiform swimmers the large amplitude body undulations are restricted to one-half or even one-third of the posterior part of the body and the undulation amplitude increases sharply in the caudal area. The wavelength of the body undulation wave is usually lower for anguilliform swimmers (∼70% of body length) than for carangiform swimmers (∼ one body length) (Videler and Wardle, 1991).
Two important non-dimensional parameters that characterize steady inline undulatory swimming are the Reynolds number (Re*) of the flow and the Strouhal number (St*) of the undulatory body motion, which can be defined as follows (Lauder and Tytell, 2006; Triantafyllou et al., 2000): (1) (2) In the above equations, L (m) is the fish length, U (m s–1) is the steady inline swimming speed, ν (m2 s–1) is the kinematic viscosity of the water, A (m) is the width of the wake, which is approximated by the maximum lateral excursion of the tail over a cycle and f (1 s–1) is the tail-beat frequency.
Previous experiments have provided a wealth of data in terms of body kinematics and wake flow field for undulatory swimming [for reviews, see Fish and Lauder, 2006; Triantafyllou et al., 2000]. A single row of vortices has been observed for the carangiform swimmers (Barrett et al., 1999; Muller et al., 1997; Nauen and Lauder, 2001; Nauen and Lauder, 2002; Wolfgang et al., 1999) while a double row of vortices has been observed for anguilliform swimmers (Hultmark et al., 2007; Liu and Wassersug, 1997; Muller et al., 2001; Tytell and Lauder, 2004). Comparing the performance of anguilliform and carangiform swimmers, however, poses major challenges to experimental studies (Tytell, 2007). One such challenge is the lack of experimental control over live fish, which thwarts the control over the governing parameters. Another challenge stems from the difficulties in obtaining the 3-D flow and pressure fields around the live fish needed to accurately measure and decompose the forces for self-propelled bodies. This is particularly hard to achieve because most particle image velocimetry (PIV) measurements provide flow fields in 2-D planes (Tytell, 2007). Therefore, the efficiency in the recent experimental studies has been obtained using hydromechanical models such as Lighthill's elongated body theory (EBT) (Lighthill, 1969). Based on EBT the thrust is produced only at the tail and the rest of the body undulations create wasted energy (Lighthill, 1969). Therefore, according to this theory carangiform swimmers are more efficient than anguilliform swimmers because most of their body undulations are restricted in the caudal fin region. This theory, however, incorporates many simplifying assumptions and neglects viscous forces that can produce thrust along the body of anguilliform swimmers (Shen et al., 2003; Taneda and Tomonari, 1974). As discussed in Tytell and Lauder, conclusions based on EBT could lead to results that are contrary to what is observed in nature (Tytell and Lauder, 2004). For example, many eels migrate thousands of miles (van Ginneken and van den Thillart, 2000) and many sharks swim steadily in the sub-carangiform (between carangiform and anguilliform) mode (Gemballa et al., 2006).
In our previous work (Borazjani and Sotiropoulos, 2009a) we helped to reconcile a number of such inconsistencies obtained from simpler theories by performing 3-D numerical simulations of tethered carangiform (Borazjani and Sotiropoulos, 2008) and anguilliform (Borazjani and Sotiropoulos, 2009a) virtual swimmers under similar conditions (Reynolds number). We found that: (1) the carangiform swimmer's Froude efficiency (ηf) is maximized as the Re* approaches ∞ while the anguilliform swimmer's efficiency is maximized somewhere in the transitional regime; (2) the temporal evolution of the net force acting on the anguilliform swimmer varies more smoothly than for the carangiform swimmer case; (3) the swimming power required for a self-propelled anguilliform swimmer is lower than that for the carangiform swimmer at the same Re*; and (4) for both swimmers the swimming power is higher than the power required to tow the respective rigid body at the same Re* (Borazjani and Sotiropoulos, 2009a).
Our previous work (Borazjani and Sotiropoulos, 2009a) pointed to a number of similarities as well as to a number of differences between the two modes of swimming. The similarities are due to similar BCF undulatory propulsion and similar governing flow parameters such as Re* and St*. The differences, however, could be due to either the difference in the respective body shapes (form) and/or the type of body undulations (kinematics). The purpose of the present study is to systematically investigate and quantify the effects of form and kinematics on the hydrodynamic performance of undulatory swimming. The main difficulty for accomplishing such an undertaking is finding a rational approach for isolating the effects of form and kinematics for different swimmers. To isolate such effects, in this work we construct virtual, self-propelled swimmers of a given body shape (fixed form) and make them swim with different kinematics. For example, to quantify the effects of kinematics on a mackerel body we can compare the performance of a mackerel swimming like a mackerel (with carangiform kinematics) with that of a mackerel swimming like an eel (anguilliform kinematics). Similarly, to isolate the effects of form we compare the performance of a mackerel body and an eel body both swimming with the same kinematics (carangiform or anguilliform). More specifically, we employ the following virtual swimmers (Fig. 1): (1) a mackerel body swimming like a mackerel (denoted as MM); (2) mackerel body swimming like a lamprey (ML); (3) a lamprey body swimming like a lamprey (LL); and (4) a lamprey body swimming like a mackerel (LM). It is of course not possible in nature to make a live mackerel swim like an eel or an eel swim like a mackerel. Thus, using virtual (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a; Kern and Koumoutsakos, 2006; Liu and Kawachi, 1999; Liu and Wassersug, 1997) or biorobotic (Barrett et al., 1999; Hultmark et al., 2007; Lauder et al., 2007) swimmers are the only feasible alternatives for such an undertaking.
To compare the performance of the four virtual swimmers, we perform self-propelled, fluid–structure interaction (FSI) simulations under the same conditions. All virtual swimmers are placed in the same, initially stagnant, fluid with viscosity ν and start undulating their bodies with the desired kinematics with a tail-beat frequency f. In the present self-propelled simulations, there is no tether to absorb the excess force F imparted by the flow [as was the case in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a)] and the fish can accelerate or decelerate depending on the sign of F. Body undulations impart a thrust-type force F on the body, which initially accelerates (on average) the virtual swimmer to higher velocities. The swimming speed continues to increase until the mean force F during one tail-beat cycle becomes zero (constant speed, inline swimming limit). Thereafter, the virtual swimmer has reached a quasi-stationary state characterized by constant values of the mean swimming speed U, Re* and St*. The swimming performance of the various swimmers is quantified at this quasi-stationary, constant-mean-speed state in terms of various performance metrics, such as power loss, Froude efficiency ηf, velocity over power (similar to mile per gallon), etc. The wake structure of various swimmers is also examined at the quasi-stationary, constant-mean-speed limit. A small portion of the results presented in this paper has recently appeared in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009b).
The rest of the paper is organized as follows. First, we briefly describe the numerical method and present the details of the fish model and prescribed kinematics. Second, we provide the results of our virtual swimmers and discuss the effects of body shape and kinematics on their performance. Finally we summarize our findings, present the conclusions of this work and outline the areas for future research.
MATERIALS AND METHODS
We have adopted the following notations for the sake of clarity throughout this manuscript: (1) boldface variables are vectors; (2) italic variables are scalars; (3) parameters with superscript * are non-dimensionalized using the self-propelled swimming speed U (m s–1) at the quasi-steady, constant-speed limit; (4) parameters with super/subscript o are non-dimensionalized using a characteristic reference velocity Uo (m s–1) (to be defined later). There is one exception to this last rule for the non-dimensional swimming speed, which is defined as U*=U/Uo.
The governing equations and boundary conditions
The equations governing the motion of an incompressible Newtonian viscous fluid are the 3-D, unsteady, incompressible Navier–Stokes equations. In this work, and in order to facilitate the numerical simulation of self-propelled swimming, the governing equations are formulated in the non-inertial frame of reference that moves with the fish.
Non-dimensionalized by a characteristic reference velocity Uo, and the fish length L, the governing equations in a non-inertial reference frame, which moves with velocity uc and angular velocity Ω, relative to the inertial frame, read in conservative form as follows (Borazjani, 2008): (3) where Reo=UoL/ν is the Reynolds number, t is the time and u, v and w are defined as: (4) where ua is the non-dimensional Cartesian absolute velocity vector of the fluid in the inertial reference frame, ur is the non-dimensional Cartesian relative velocity vector of the fluid in the non-inertial reference frame, rr is the non-dimensional position vector relative to the non-inertial frame, QT is the rotation matrix of the non-inertial frame relative to the inertial frame, and p is the non-dimensional pressure. The time derivative relative to the non-inertial frame is defined as: (5) This conservative formulation has been developed and used successfully to simulate the flow in non-inertial reference frames (Beddhu et al., 1996; Kim and Choi, 2006). Note that, as shown by Kim and Choi (Kim and Choi, 2006), the non-conservative formulation is not adopted herein because it can lead to severe instability in the numerical method. The inviscid (Euler) equations, which are also solved in this work, are recovered from Eqn 3 by letting Reo→∞.
We are interested in solving the governing flow equations in a domain containing an arbitrarily complex 3-D flexible body moving with prescribed kinematics. Therefore, the governing equations (Eqn 3) need to be supplemented with appropriate boundary conditions at the boundaries of the flow domain, which could be either fluid boundaries or the solid boundary of a moving immersed body.
In this work the non-inertial reference frame is attached to the virtual swimmer's center of mass. The motion of the center of mass is obtained by solving the Newton's Second law of motion (momentum) equations for the fish in non-dimensional form: (6) where U* is the fish swimming speed non-dimensionalized by Uo, Mred=m/ρL3 is the reduced mass (where m is the mass of the virtual swimmer and ρ is the fluid density) and is the force coefficient (where F is the force vector exerted on the virtual swimmer's body by the fluid). The position of the non-inertial frame is obtained by solving the following equation: (7) where xc is the position vector of center of mass non-dimensionalized by L, i.e. the position of the origin of the non-inertial frame relative to the inertial frame. In this work we have restricted the virtual swimmer to swim only along the streamwise direction. We are not considering the motion in other directions because we have removed all of the median and pectoral fins in our virtual swimmers that the fish use for stability during swimming.
Let the boundary of the fish body be defined by the dynamically evolving surface Γ(t). Γ(t) is discretized with K material points, which lie on it at all times and can be tracked with their global Lagrangian position vectors rkr (t) (k=1, K) in the non-inertial frame: (8) In this work the motion relative to the non-inertial frame is prescribed, i.e. position of the markers relative to the non-inertial frame is known at all times. Therefore, the shape of Γ(t) can be obtained at any time from the known relative position vectors, Eqn 8 and the calculated position of the non-inertial frame from Eqn 7 as follows (for k=1, K): (9) where is the global position vector of material point k at time t in the inertial frame of reference. With the shape of Γ(t) known at time t, boundary conditions for the Eulerian fluid velocity vector u(r,t) must be prescribed at all points of Γ(t). For viscous flow, the no-slip and no-flux boundary conditions need to be satisfied as follows: (10) This boundary condition enforces the link between the Eulerian description of the fluid motion and the Lagrangian description of the moving immersed body (for details, see Gilmanov and Sotiropoulos, 2005).
For inviscid flow only the no-flux condition is satisfied on the body, i.e. the fluid velocity normal to the body is set equal to the normal velocity of the body while the fluid velocity components tangent to the body need to be prescribed by interpolation from the interior fluid nodes. The mathematical formulation of these boundary conditions reads as follows: (11) (12) where un is the fluid velocity normal to the body, ut is the fluid velocity vector tangential to the body and n is the normal vector to the surface.
The numerical method
The flexible fish body is handled as a sharp interface immersed in the background Cartesian grid using the hybrid Cartesian immersed-boundary (HCIB) method (Gilmanov and Sotiropoulos, 2005). This method has been described in detail in Gilmanov and Sotiropoulos (Gilmanov and Sotiropoulos, 2005) and only a very brief description of the technique is given herein. The HCIB approach employs an unstructured, triangular mesh to discretize and track the position of a fish body. Boundary conditions for the velocity field at the Cartesian grid nodes that are exterior to, but in the immediate vicinity of, the immersed boundary, i.e. immersed boundary (IB) nodes, are reconstructed by linear or quadratic interpolation along the local normal to the boundary. No explicit boundary conditions are required for the pressure field at the IB nodes due to the hybrid staggered/non-staggered formulation of Gilmanov and Sotiropoulos (Gilmanov and Sotiropoulos, 2005). The HCIB reconstruction method has been shown to be second-order accurate on Cartesian grids with moving immersed boundaries (Gilmanov and Sotiropoulos, 2005). The IB nodes at each time step are recognized using an efficient ray-tracing algorithm (Borazjani et al., 2008).
This method has been validated extensively by Gilmanov and Sotiropoulos (Gilmanov and Sotiropoulos, 2005) and Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008) for flows with moving boundaries and has also been applied to simulate fish-like swimming for tethered carangiform swimmers (Borazjani and Sotiropoulos, 2008; Gilmanov and Sotiropoulos, 2005) and tethered anguilliform swimmers (Borazjani and Sotiropoulos, 2009a). As in our recent work with tethered swimmers (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a), we solve the governing equations using the fractional step method of Ge and Sotiropoulos (Ge and Sotiropoulos, 2007). The Poisson equation is solved with flexible GMRES (Saad, 2003) and multigrid as a preconditioner using parallel libraries of PETSc (Balay et al., 2004). For more details the reader is referred to Ge and Sotiropoulos (Ge and Sotiropoulos, 2007) and Borazjani et al. (Borazjani et al., 2008).
We use the partitioned approach to solve the coupled fluid–structure system of equations (Borazjani et al., 2008; Piperno and Farhat, 2001). The system of Eqns 3, 6 and 10 is partitioned into a fluid and a structure domain. Each domain is treated computationally as an isolated entity and is separately advanced in time and the interaction effects are accounted for through boundary conditions, i.e. the fluid equation (Eqn 3) and structure equation (Eqn 6) are coupled together through the no-slip boundary condition equation (Eqn 10) and the dynamic boundary condition F on the right-hand side of Eqn 6. Based on the treatment of boundary conditions two types of FSI coupling are possible: loose coupling (LC–FSI) or strong coupling (SC–FSI) [the readers are referred to the work of Borazjani (Borazjani, 2008) for the details of the different coupling methods]. LC–FSI is explicit in time while SC–FSI is implicit in time and is implemented by performing several sub-iterations at each time step (Borazjani et al., 2008). As discussed in Borazjani et al. (Borazjani et al., 2008) even SC–FSI can be unstable if the added mass is considerable relative to the actual mass of the system. To stabilize the SC–FSI we have used the following under-relaxation scheme (Borazjani et al., 2008) for all simulations reported in this work when solving Eqn 6: (13) where α is the under-relaxation coefficient, U*l+1 and U*l+1 indicate the solution after and before under-relaxation has been applied at the l+1 SC–FSI iteration. The under-relaxation coefficient is obtained dynamically during each SC–FSI iteration using the Aitken acceleration technique (Borazjani et al., 2008). The SC–FSI algorithm with Aitken acceleration typically converges within 4–5 iterations at each time step. The FSI coupling method formulated in the non-inertial reference frame is validated in Appendix B by simulating vortex-induced vibrations (VIV) of an elastically mounted cylinder in the non-inertial frame of reference.
Fish body kinematics and non-dimensional parameters
The bodies of the virtual swimmers used in this study are exactly the same as those used in our previous tethered simulations. The carangiform body was modeled after the actual anatomy of a mackerel (Borazjani and Sotiropoulos, 2008) whereas the anguilliform body was created from a lamprey computed tomography (CT) scan by Professor Frank Fish, provided to us by Professor Lex Smits from Princeton University (Borazjani and Sotiropoulos, 2009a). The lamprey was chosen to represent the anguilliform swimmer in this work because our tethered anguilliform swimmer (Borazjani and Sotiropoulos, 2009a) was a lamprey, which simplifies the comparison of tethered and self-propelled swimmers. Furthermore, the lamprey kinematics is similar to other anguilliform swimmers such as the eel (Hultmark et al., 2007). Only the caudal fin was retained for the mackerel whereas all other fins were neglected for the lamprey due to a lack of detailed experimental data (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). The bodies are meshed with triangular elements as needed by the HCIB method (Fig. 2).
The kinematics for BCF locomotion is generally in the form of a backward traveling wave with the largest wave amplitude at the fish tail (Gray, 1933), which can be described by the following equations for the lateral undulations of the fish body (all lengths are non-dimensionalized with the fish length L): (14) where z is the axial (swimming) direction measured along the fish axis from the tip of the fish's head; h(z,t) is the lateral excursion of the body at time t; a(z) is the amplitude envelope of lateral motion as a function of z; k is the wave number of the body undulations that corresponds to a wavelength λ; and ω is the angular frequency. Both modes of BCF propulsion studied herein, i.e. anguilliform and carangiform, are described by the above traveling wave equation (Eqn 14) by choosing an amplitude envelope and a wavelength (or wave number), referred to hereafter as shape parameters, that match that mode of swimming (Tytell and Lauder, 2004; Videler and Hess, 1984). In this work, we define the fish length L as the total length from the tip of the nose to the end of the tail of the straight virtual swimmer. For the mackerel and lamprey bodies the fork length is equal to 0.94L and L, respectively.
The amplitude envelope a(z) for the anguilliform kinematics was approximated by an exponential function (Borazjani and Sotiropoulos, 2009a; Tytell and Lauder, 2004): (15) where amax is the tail-beat amplitude.
For carangiform kinematics the amplitude envelope was approximated by a quadratic curve of the form (Borazjani and Sotiropoulos, 2009a): (16) For a typical anguilliform fish the coefficient amax is set equal to amax=0.1 (Hultmark et al., 2007). The following values are used for the coefficients a0=0.02, a1=–0.08 and a2=0.16 to match the experimental curve of Videler and Hess (Videler and Hess, 1984) for typical carangiform kinematics. Both kinematics have the maximum displacement at the tail amax=0.1, i.e. hmax=0.1L. Fig. 3 shows the amplitude envelopes of anguilliform and carangiform kinematics according to Eqns 15 and 16, respectively. The wave number in all simulations is based on the non-dimensional wavelength λ/L=0.642 for anguilliform (Borazjani and Sotiropoulos, 2009a; Hultmark et al., 2007) and λ/L=0.95 for carangiform (Borazjani and Sotiropoulos, 2008; Videler and Hess, 1984) swimmers. Fig. 4 shows the midlines of the fish calculated for one tail-beat cycle using Eqn 14 with the anguilliform and carangiform shape parameters and the amplitude envelopes calculated by Eqns 15 and 16, respectively. Fig. 1 shows one instant in time of the body shapes resulting when such undulations of the midline are imposed on the lamprey and mackerel bodies.
The four important non-dimensional similarity parameters in fishlike swimming are: (1) the swimming Re* defined by Eqn 1; (2) the St* defined by Eqn 2 based on the maximum lateral excursion of the tail A=2amax and the tail-beat frequency f: St*=2famax/U; (3) the non-dimensional wavelength λ/L; and (4) the non-dimensional amplitude envelope a(z/L)/L. Sometimes the so-called slip velocity or slip ratio, defined as slip=U/V=U/(ω/k), is used instead of non-dimensional wavelength. Using either parameter is correct. However, the slip velocity changes if the tail-beat frequency is changed, while the wavelength and the tail-beat frequency are independent.
Setup of numerical experiments and computational details
The four virtual swimmers are released in a given fluid with a prescribed tail-beat frequency and the swimming velocity is calculated based on the forces on the swimmers body. The simulations are continued until the quasi-steady state is reached, i.e. the mean swimming speed stays constant. The swimming speed U in the quasi-steady state depends on the tail-beat frequency. Therefore, because Re* (Eqn 1) and St* (Eqn 2) depend on the swimming speed U, both of these parameters also depend on the tail-beat frequency. Consequently, in order to obtain a quasi-steady value for Re* in a desired range, the tail-beat frequency should be chosen carefully. To accomplish this we use the information from our previous simulations with tethered swimmers (Borazjani and Sotiropoulos, 2008) as guidance to select the tail-beat frequency. Before explaining how the tail-beat frequency is chosen here, let us summarize again how we performed the tethered simulations (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a) using the notation we employ in this paper. We fixed the Reynolds number Reo by moving the tether at a constant speed Uo in a fluid of a given viscosity. Subsequently, for the prescribed Reo value, we systematically varied the Strouhal number Sto (based on Uo) by varying the tail-beat frequency until we found the critical Stouhal number at which the mean net force on the tether is zero, i.e. the virtual swimmer could self-propel itself with constant velocity equal to the corresponding velocity of the tether Uo. The self-propelled simulations can be considered as the inverse of the tethered simulations because we now fix the tail-beat frequency and find the resulting swimming velocity U in the self-propelled quasi-steady state. Therefore, if we prescribe in a self-propelled swimmer the tail-beat frequency corresponding to the value of , obtained from the tethered simulations of the same swimmer, we should obtain a self-propelled swimming speed U∼Uo (i.e. equal to the speed of the tether for ) in the quasi-steady state, assuming of course that the fluid is the same in the tethered and self-propelled simulations.
To obtain the same hydrodynamic environment in the self-propelled simulations as in the tethered simulations, we non-dimensionalize the fluid governing equation (Eqn 3) by setting the characteristic velocity Uo to be equal to the tether velocity and selecting the fluid to have the same kinematic viscosity as in the tethered simulations. We have chosen the Reo to exactly match the Reynolds number in the tethered simulations, i.e. (1) a viscous environment (Reo=300); (2) a moderately viscous environment (Reo=4000); (3) and an inviscid environment (Reo=∞), which are referred to as cases R1, R2 and R3, respectively, hereafter.
To match the tail-beat frequency f with the tethered simulations, we have non-dimensionalized the body traveling wave equation (Eqn 14) with the velocity Uo and fish length L and use the non-dimensional angular frequency ω as follows: (17) The above non-dimensional angular frequency is used along with the non-dimensional time tUo/L in Eqn 14. fo=fL/Uo is the non-dimensional frequency and Sto=f2amax/Uo is the Strouhal number based on the tail-beat frequency f, lateral excursion of the tail 2amax and the tether velocity Uo. Now for the desired value of Reo in a self-propelled simulation, we find the critical (zero force on the tether) Strouhal number from the corresponding tethered carangiform swimmer simulation (Borazjani and Sotiropoulos, 2008) and choose a value for Sto in the self-propelled simulation to be close to . Namely, we have chosen Sto=1.1, 0.6 and 0.3 for Reo=300, 4000 and ∞, respectively. We then use this Sto to find the non-dimensional angular frequency ω (Eqn 17) and prescribe it in the body traveling wave equation (Eqn 14).
The selection of the tethered carangiform swimmer (MM) (Borazjani and Sotiropoulos, 2008) as the base line for computing the non-dimensional angular frequency ω is of course arbitrary but it does not in any way limit the generality of our results. Ultimately all swimmers in each case race each other in the same hydrodynamic environment (fixed Reo) and beat their tail at the same angular frequency ω (fixed Sto). Naturally, however, this treatment will yield for every swimmer a steady, swimming speed U different than Uo, i.e. U*=U/Uo≠1 due to differences in kinematics, body shape and the different numerical approaches used to simulate the tethered and self-propelled swimmers. Because we have selected the tethered MM swimmer as the base-line case for selecting the tail-beat frequency we should expect that the non-dimensional swimming speeds we obtain for the self-propelled MM swimmer are, in general, closer to unity and this is indeed the case as discussed in the Discussion section.
The computational domain and time step for the mackerel body are exactly the same as the tethered mackerel simulations (Borazjani and Sotiropoulos, 2008). Similarly, for the lamprey body the domain and time step are the same as in the tethered lamprey simulations (Borazjani and Sotiropoulos, 2009a). The computational domain is a cuboid with dimensions 2L×L×7L, which is discretized with 5.5 million grid nodes. The domain width 2L and height L are more than 15 times the lamprey width 0.067L and height 0.066L and 10 times the mackerel width 0.2L and height 0.1L, respectively. As in the tethered simulations, the fish is placed 1.5L from the inlet plane in the axial direction and centered in the transverse and the vertical directions. The virtual swimmers start to undulate in an initially stagnant fluid and the boundary conditions on the domain outer boundaries are far-field (Neumann) boundary conditions. The reduced mass Mred in Eqn 6 is set equal to 0.01 for all of virtual swimmers.
Calculation of swimming forces and efficiency
The procedure we employ to calculate the hydrodynamic forces and efficiency has been discussed extensively in our previous publications (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). Therefore, only a brief description is given below for the sake of completeness.
In our simulations, we consider inline swimming along the x3 direction. The component of the instantaneous hydrodynamic force in the x3 direction (which for simplicity will be denoted as F) can be readily computed by integrating the pressure and viscous forces acting on the body as follows (where repeated indices imply summation): (18) where nj is the jth component of the unit normal vector on dA and τij is the viscous stress tensor.
Depending on whether F(t) is negative or positive, it could contribute to either hydrodynamic drag D(t) or thrust T(t). To separate the two contributions we adopt the force decomposition approach proposed by Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008): (19) (20) In the above equations the subscripts p and v refer to force contributions from pressure and viscous terms, respectively.
The numerical details for calculating the various surface integrals involved in the above equations in the context of the HCIB numerical method can be found in Borazjani (Borazjani, 2008). A detailed validation study demonstrating the accuracy of our numerical approach for calculating the viscous and pressure components of the hydrodynamic force in the non-inertial reference frame is provided in Appendix A.
The power loss due to lateral undulations of the fish body (Pside) is calculated as follows: (21) where is the time derivative of the lateral displacement (i=2 direction), i.e. the velocity of the lateral undulations.
The mean quantities of force, thrust, drag and power are obtained by averaging the instantaneous values over several swimming cycles. We non-dimensionalize each mean quantity in two ways: (1) non-dimensional quantities with superscript * are scaled with the calculated swimming speed U, which is the typical definition for various coefficients, such as force or thrust coefficients, found in the literature; and (2) non-dimensional quantities with superscript o are scaled with the characteristic velocity Uo, which is the velocity scale we used to non-dimensionalize the governing equations. The latter scaling is useful because using a common reference velocity for all swimmers preserves the relative proportion of a given swimming performance measure X (where X can be power, force, thrust, etc.) between different swimmers in a given environment, e.g. (where i,j=1 to 4 for the MM, ML, LM and LL swimmers, respectively). We refer to quantities with * as performance measure X coefficient (e.g. power coefficient) and to quantities with o as non-dimensional performance measure X (e.g. non-dimensional power). The various non-dimensional quantities are defined as follows.
Force coefficient and non-dimensional force in the axial direction: (22) where ρ is the density of the fluid, F is the mean force, U is the mean swimming speed and Uo is the characteristic velocity.
Thrust coefficient and the non-dimensional thrust: (23) where T is the mean thrust.
Power coefficient and non-dimensional power: (24) where P is the mean side power.
The Froude propulsive efficiency ηf based on the thrust force for constant speed inline swimming is defined as follows (see Tytell and Lauder, 2004): (25) where P is the mean power loss due to lateral undulations. We also use an efficiency measurement (called mean efficiency hereafter), (26) which is similar to the miles per gallon (MPG) fuel efficiency for cars as suggested by Schultz and Webb (Schultz and Webb, 2002).
It is important to note that the Froude efficiency equation (Eqn 25) can only be applied under inline, constant-speed swimming when the thrust force is balanced exactly by the drag force and the net mean force acting on the fish body is zero (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). Therefore, Eqn 25 is used to compute the efficiency only when the virtual swimmer has reached the quasi-steady state of constant-mean velocity.
We visualize the 3-D wake structure of different virtual swimmers for various cases using an iso-surface of the q variable (Hunt et al., 1988). The quantity q is defined as , where S and Ω denote the symmetrical and asymmetrical parts of the velocity gradient, respectively, and ∥•∥ is the Euclidean matrix norm. According to Hunt et al. (Hunt et al., 1988), regions where q>0, i.e. regions where the rotation rate dominates the strain rate, are occupied by vortical structures.
As discussed above, for a given hydrodynamic environment (fixed Reo) all four virtual swimmers are released with the same tail-beat frequency (fixed Sto), and the self-propelled, FSI simulations are continued until the swimmers reach quasi-steady state. The calculated time series of swimming speeds for the three hydrodynamic environments (cases R1, R2 and R3) and for all four swimmers are shown in Fig. 5. It can be observed from this figure that the swimmers with the lamprey body (LM, LL) consistently require more cycles to reach the quasi-steady state relative to the swimmers with the mackerel body (MM, ML). For example, in Fig. 5C the swimmers with lamprey body take about 80 cycles to reach quasi-steady state while swimmers with mackerel body reach there only after 35 cycles.
Fig. 5A shows the calculated swimming speed for the viscous case R1 (Reo=300, Sto=1.1). It can be observed that the ML swimmer reaches the highest velocity of all swimmers. Comparing swimmers with the same body shape, those with anguilliform kinematics (ML and LL) always reach velocities higher than those with carangiform kinematics (MM, LM). The swimming speed of the ML swimmer in the steady state is about 11% faster than the MM swimmer but about 29% and 42% faster than the LM and LL swimmers, respectively.
Fig. 5B shows the results for the moderately viscous case R2 (Reo=4000, Sto=0.6). As in the viscous case, the ML swimmer is the one that reaches the highest velocity for this case as well. Similar to the previous case, among swimmers with the same body those with anguilliform kinematics (ML and LL) always win. It is worth noting, however, that in this case the anguilliform kinematics win the race by a significantly narrower margin than was the case for the viscous case R1. Namely, in this case the ML swimmer is only about 3% faster than the MM swimmer and about 16% and 10% faster than the LM and LL swimmers, respectively.
Strikingly different trends are observed in Fig. 5C, which shows the results for the inviscid case R3 (Reo=∞, Sto=0.3). It is seen that for this case the MM swimmer reaches the highest velocity among all swimmers. Comparing swimmers with the same body, those with carangiform kinematics (MM, LM) now win the race against their anguilliform counterparts (ML, LL) by a wide margin. The MM swimmer is about 23% faster than the ML swimmer and about 38% and 46% faster than the LM and LL swimmers, respectively. A surprising finding in this regard is that during the initial acceleration phase, the swimmers with anguilliform kinematics (ML, LL) are seen to accelerate faster than their carangiform counterparts (MM, LM) (see Fig. 5C,D). Ultimately, however, the swimmers with carangiform kinematics overtake the swimmers with anguilliform kinematics and win the race. More specifically, it takes the MM swimmer three cycles to catch up with the ML swimmer (Fig. 5D) while it takes the LM swimmer about 28 cycles to catch up with the LL swimmer (Fig. 5C).
In the above discussion we focused on relative swimming speed without considering how much power is required by each swimmer to attain a given speed. In Fig. 6 we plot the time history of instantaneous non-dimensional power (Eqn 24) for all virtual swimmers in different cases. Recall that, by definition, the non-dimensional power preserves the relative magnitude of the power between the different virtual swimmers. Fig. 6A shows the non-dimensional power results for the case R1 (Reo=300, Sto=1.1). In this case, the swimmers with the mackerel body (ML, MM) have similar power time records, which exhibit higher fluctuations than the swimmers with the lamprey body (LM, LL). Fig. 6B shows the non-dimensional power results for the case R2 (Reo=4000, Sto=0.6). The relative values between different swimmers in this case are similar to the previous case R1 but the absolute values of non-dimensional power are one-order of magnitude lower. Fig. 6C shows the non-dimensional power the results for the inviscid case R3 (Reo=∞, Sto=0.3). The values of non-dimensional power are again one-order of magnitude smaller than the previous case R2. Therefore, we can see that the non-dimensional power is decreasing as Reo increases for all swimmers. Similar to the previous cases, R1 and R2, we see that the swimmers with the mackerel body (MM, ML) have higher power values than those with lamprey body (LM, LL). Comparing swimmers with the same body we see that the ones with carangiform kinematics (MM, LM) have higher fluctuations in power than the ones with anguilliform kinematics (ML, LL).
An interesting observation in Fig. 6C is that for the R3 case the power becomes negative at some intervals during the swimming cycle, which indicates that the swimmer is able to extract energy from the flow. Note that for this inviscid case (R3) the viscous forces do not contribute to the swimming power and only pressure forces are involved. To explore the physical reasons for this interesting finding, we visualize in Fig. 7 the instantaneous pressure field (non-dimensionalized by ) at the mid-plane of the swimmer at two instances in time when the power is positive and negative, respectively. The anterior portion of the MM swimmer body neither moves (small ) nor contributes to the power significantly. Therefore, most of the power is transmitted to the flow in the posterior region, where the arrows in the figure indicate the direction of the tail motion. It can be observed in Fig. 7A that, due to the tail motion, two pockets of high and low pressure are present on the two sides of the tail such that the resulting pressure force is against the tail motion. Therefore, in this case the swimmer has to expend power to advance its tail against the pressure-induced force (positive power). In Fig. 7B, however, when the tail starts to return, the resulting pressure force is in the same direction as the tail motion and the swimmer could use this favorable pressure to its advantage (negative power) to extract power from the flow. It is interesting to note that all of the other swimmers also exhibit power minima at the same time during the cycle but the magnitude of the power is clearly not as small as for the MM swimmer (Fig. 6C). This can be explained by examining the body shape and kinematics of the MM swimmer. The large surface area and large amplitude of motion in the posterior region relative to the other swimmers allows the MM swimmer to use the pressure pockets in the flow most effectively.
Another important observation that emerges collectively from Figs 5 and 6 is that for all three hydrodynamic environments, swimmers with the mackerel body (MM, ML) not only consistently reach higher velocities than swimmers with the lamprey body (LM, LL) but use more power as well. Reaching higher velocities does not necessarily indicate higher efficiency if more power is consumed, and more quantitative analysis is required to further probe this important point. To accomplish this, we report in Table 1 various mean quantities and efficiencies in the quasi-steady-state. In addition to the computed mean value of U* and the resulting Re* and St*, we also provide results for: (1) the power coefficient and mean power equation (Eqn 24); (2) the thrust coefficient and mean thrust equation (Eqn 23); (3) the Froude efficiency ηf equation (Eqn 25); (4) mean efficiency equation (Eqn 26); and (5) the root mean squared (r.m.s.) of the swimming speed fluctuations, which are computed once the quasi-steady, constant-mean-speed state has been reached. It can be observed from the table that in the viscous case R1, the ML swimmer is most efficient, has the smallest power coefficient, produces the most amount of non-dimensional thrust and reaches the highest speed. However, the LM swimmer uses the smallest power . In the moderately viscous case R2, the ML swimmer produces the highest non-dimensional thrust and reaches the highest swimming speed but the LL swimmer is most efficient and uses less power. Finally, in the inviscid case R3, the MM swimmer produces the highest non-dimensional thrust and reaches the highest speed and is most efficient but the LL swimmer uses the least power. In all cases the MM swimmer always exhibits velocity fluctuations with the highest r.m.s. From this table the effects of body shape and kinematics cannot be readily deduced. We will discuss these issues in more detail in the Discussion section below.
Figs 8, 9, 10 show the wake structures visualized by the q-criterion of different swimmers for the cases R1, R2 and R3, respectively. It can be observed from the Figs 8 and 9 that a double row structure is present in the viscous and transitional cases for all the swimmers. However, a single row wake structure emerges in the inertial case R3 for all swimmers (Fig. 10).
Apart from large-scale characteristics of the wake (single vs double row), the body shape can affect more subtle wake features, such as the shape of the vortex rings. For example, in Fig. 9 the rings shed by the lamprey body (LL and LM swimmers) are more circular and simpler in shape whereas the rings shed by the mackerel body (MM and ML swimmers) are more complex and less organized. Moreover, in Fig. 10, where all swimmers have a single row wake, the swimmers with the mackerel body show a remarkable vortex-within-a-vortex structure in the vicinity of the tail, which is not observed in the swimmer with lamprey body. Fig. 11 examines the vortex-within-a-vortex structure more closely for the MM and ML swimmers, which shows two vortex rings: vortex ring 1 is generated by the leading edge of the tail; while vortex ring 2 is generated by the trailing edge of the tail. The swimmers with the lamprey body do not show such structure probably due to the lack of a homocercal large aspect-ratio tail. This important wake feature of the mackerel body swimmers will be discussed further in a subsequent section. Comparing the wake structure of the swimmers with similar body shapes but different kinematics in Figs 8, 9, 10, we find that the kinematics appear to have a weaker effect on the small features of the wake than the body shape. In fact, for all cases the wake of swimmers with the same body shape but different kinematics looks fairly similar.
Self-propelled vs tethered simulations and experiments
As mentioned earlier to achieve swimming speeds close to one, which are desirable from a numerical standpoint, the tail-beat frequency of each virtual swimmer at Reo is selected close to the (Stouhal number at which net mean force on the tether is zero for Reynolds number Reo) found from tethered simulations of the MM swimmer (Borazjani and Sotiropoulos, 2008). Inversely, we can pose the question whether the swimming speed that results from the self-propelled, FSI simulations yields values for St* and Re* that are comparable with those obtained from the tethered model. The Re* and St* values reported in Table 1, i.e. St*=1.10, 0.61 and 0.25 for Re*=299, 3910 and ∞, respectively, are in excellent agreement with the values found in the tethered MM swimmer simulations reported in our previous publication (Borazjani and Sotiropoulos, 2008), i.e. , 0.6 and 0.26 for Reo=300, 4000 and ∞, respectively. Similarly, the values obtained from the LL swimmer in Table 1, i.e. St*=1.39, 0.67 and 0.47 for Re*=238, 3598.8 and ∞, respectively, are also in agreement with the values from the tethered LL swimmer with amax=0.1L (Borazjani and Sotiropoulos, 2009a), i.e. , 0.63 and 0.46 for Reo=300, 4000 and ∞, respectively. It can be observed that for the LL swimmer, the difference between St* of self-propelled and tethered simulations is somewhat higher relative to the MM swimmer because there is a larger difference between the Re* of self-propelled and tethered simulations for the LL swimmer. Note that this is because of the fact that, as already discussed above, the Sto and tail-beat frequency were chosen based on the of the tethered MM swimmer.
In spite of the differences in kinematics and shape, most fish in nature swim in a range of Strouhal numbers 0.25–0.35 (Triantafyllou et al., 1993; Triantafyllou and Triantafyllou, 1995). Our virtual swimmers with their differences in shape and kinematics get closer to this range as Re increases (inviscid simulations). Some fish, e.g. pacific salmon, however, have been observed to swim at high Strouhal number at low swimming velocities (Lauder and Tytell, 2006). Our results underscore the conclusion of Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008) that, for a given body shape and kinematics, for each Re* there is unique St* at which self-propelled swimming is possible, i.e. the pacific salmon has to swim at high St* because this is the only St* that can make it swim steadily at low Re* (swimming speed). According to all of the values given in Table 1, St* is a decreasing function of Re* for all swimmers, which is also in agreement with our pervious findings from the tethered simulations (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). Furthermore, comparing the self-propelled and tethered MM and LL swimmers Re* and St*, we observe that the self-propelled swimmers have lower Re* but higher St* than the tethered ones. For example, the self-propelled LL swimmer in the R2 case has Re*=3599, which is lower than the corresponding value for the tethered LL swimmer Reo=4000 but has St*=0.67, which is higher than the tethered LL swimmer . These results are also in agreement with our finding that St* is a decreasing function of Re*.
The Froude efficiency and power coefficient values reported in Table 1 for MM and LL are also consistent with the values previously found for the tethered swimmers (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a) and show the same trend. The efficiency of the MM swimmer increases as Re* increases but the efficiency of the LL swimmer is maximized in the transitional regime (case R2). The non-dimensional power is decreased for all swimmers as Re* is increased as in our tethered simulations. The MM swimmer has higher values than the LL swimmer in the same case, as was also the case in the tethered simulations. This is in agreement with the observations of Tytell in that trout have larger estimated wasted power than eels (Tytell, 2007). However, in Tytell (Tytell, 2007), the power loss estimates for trout decreased as swimming speed increased, similar to the trend observed here, but the estimates for eel did not vary much at least for a few swimming speeds that were reported. Note that Tytell cautioned about the fact of deriving firm conclusions from relatively few experimental measurements that are not adequate to distinguish real differences among species from random variability among individuals in a statistically meaningful manner. As discussed in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009a), non-dimensionalizing the power in case R3 as in Tytel (Tytell, 2007) (i.e. using , where S is 0.18L2 for eel and 0.54L2 for trout and assuming that the same S is good for the LL and MM swimmers, respectively) results in a power coefficient of about 0.0026 for the LL swimmer, which is remarkably close to 0.004 reported by Tytell for the eel, and 0.0013 for the MM swimmer, which is comparable with the 0.007 value reported for the trout. As discussed in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009a), the surface area of a lamprey is much closer to that of an eel than that of a mackerel to a trout. Nevertheless, in Tytell (Tytell, 2004a) the cost of wake production of steadily swimming eels was found to increase with swimming speed with power 1.7 (less than 2), i.e. the power coefficient, which is scaled by U3, decreases as swimming speed (Re*) increases.
Comparing the results in Table 1 of the present study and table 2 of Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009a), we found that the power required for undulatory swimming is higher than towing the rigid fish at the same Re*. This is in agreement with our previous results (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). However, this is in contrast with the results of Barrett et al. who found that in a narrow range the power required for swimming was lower than that of towing (Barrett et al., 1999). As discussed in our previous publications (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a) the current simulations are at much lower Re relative to the experiments of Barrett et al. (Barrett et al., 1999) and cannot provide conclusive evidence for or against those results.
It can also be observed from Table 1 that the r.m.s. of the swimming speed fluctuations (r.m.s. of U) increases as Re* decreases. This finding is in agreement with experimental observations. For instance, fish larvae swimming at low Re* show greater fluctuations in the swimming speed than the corresponding adult fish. More specifically, the swimming speed of zebra fish larvae at Re*≈300 fluctuates between 14 and 24L–1 (Muller et al., 2008) whereas for an adult eel swimming at Re*≈10,000 the swimming speed fluctuates within a much narrower range between 0.9 and 1.1U (Muller et al., 2001; Tytell and Lauder, 2004).
Self-propelled simulations vs experiments and wake structure
A single row wake structure has been observed in experiments for carangiform swimmers (Muller et al., 1997; Nauen and Lauder, 2002; Wolfgang et al., 1999) while a double row structure has been observed for anguilliform swimmers (Hultmark et al., 2007; Muller et al., 2001; Tytell and Lauder, 2004). Comparing the wake of the carangiform swimmers with the MM swimmer in case R3 (Reo=∞, Sto=0.3), which best corresponds to carangiform swimming in nature, we observe the same single row structure. The single row structure of the MM swimmer from the self-propelled simulations is similar to the single row structure of the tethered swimmer as well (Borazjani and Sotiropoulos, 2008), which is in agreement with the wake of a mackerel (Nauen and Lauder, 2002) [see Buchholz and Smits (Buchholz and Smits, 2006) for the description of the skeleton of this wake structure]. However, the wake of the LL swimmer in case R2 (Reo=4000, Sto=0.6), which best corresponds to anguilliform swimming in nature, we observe a double row structure. This is also similar to the wake of the tethered LL swimmer from our previous work (Borazjani and Sotiropoulos, 2009a), the 3-D self-propelled simulations of Kern and Koumoutsakos (Kern and Koumoutsakos, 2006) and the wake observed for anguilliform swimmers in nature (Hultmark et al., 2007; Muller et al., 2001; Tytell and Lauder, 2004) [see Buchholz and Smits (Buchholz and Smits, 2008) for the description of the skeleton of this wake structure].
The ML swimmer (mackerel body and anguilliform kinematics) is somewhat similar to sharks because most sharks have fusiform body shapes and swim in the sub-carangiform mode [except the members of family of Lamnidae that swim in thunniform mode (Gemballa et al., 2006)]. Wilga and Lauder observed a ring-within-a-ring vortex structure shed from the tail of a dogfish shark (Wilga and Lauder, 2004). They hypothesized this complex vortical structure was due to the inclined axis of rotation of the shark's tail relative to the horizontal axis of locomotion, and argued that the underlying vorticity generation mechanism in this case is similar to that at work when a vortex ring is created by a flow pulse through an inclined exit orifice. In the present study we have observed the same vortex-within-a-vortex structure for both MM and ML swimmers (Fig. 11) and its formation appears to be clearly associated with the homocercal large-aspect ratio tail: the leading edge creates the outer vortex loop while the trailing edge creates the inner vortex loop. The effect of heterocercal shape of the tail and non-symmetrical motion of the tail, i.e. similar to the tail motion of the shark (Wilga and Lauder, 2002; Wilga and Lauder, 2004), on the near-wake vorticity dynamics has yet to be understood.
Another unanswered question is the relative importance of form (shape) and kinematics on the wake structure. For example, Tytell and Lauder discuss the importance of body shape on the wake structure stating that: ‘even if a mackerel, for example, swam using the same kinematics as an eel, its wake would probably differ from an eel's due to the differences in body shape’ (Tytell and Lauder, 2004). In the experiments with the pitching and heaving hydrofoil (Koochesfahani, 1989; Triantafyllou et al., 1991) and in our previous work with virtual swimmers (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a), the St* has been found to be the primary parameter in determining the wake structure. For both tethered anguilliform and carangiform swimmers we found at a given Re* a single row of vortices at low St* while a double row of vortices was found at high St* (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a). The Re* was found to affect the St* at which the wake transitioned from single to double row structure. Note that in the tethered simulations the Re* and St* were independent parameters but in self-propelled simulations the Re* and St* depend on each other. As discussed previously, for a fixed body shape and kinematics, for each Re* there is unique St* at which self-propulsion is possible and this St* is a decreasing function of Re*. Therefore, here we cannot study the effects of Re* and St* independently, as we did with the tethered model, but we can study the effects of body shape and kinematics independently.
The results in Figs 8, 9, 10 clearly show that the conclusions we derived from our earlier simulations with the tethered MM and LL swimmers regarding the dominant role of the St* still hold for the various swimmers we considered in this work. Namely, regardless of body shape and kinematics, a double row structure is present in the viscous and transitional cases (R1 and R2), for which the St* is relatively high (∼1.1 and 0.6, respectively), while a single row wake structure emerges in the inertial case R3 for which the St* is low (∼0.3). Moreover, we also observe that when the hydrodynamic environment is fixed, i.e. when the calculated Re* and St* values are similar for all swimmers: all four swimmers have the same wake structure (single or double row). Consequently, we can conclude that the large-scale features of the wake structure (single vs double row) do not depend on body shape or kinematics but primarily depend on flow parameters Re* and St* (which in this work depend on each other and can be used interchangeably). This is, indeed, in agreement with a wide range of wake visualization results from a variety of applications that have produced similar wake structures despite the differences in the geometry. For example, high-aspect-ratio fins (Koochesfahani, 1989; Triantafyllou et al., 1991), low-aspect-ratio square fins (Buchholz and Smits, 2006; Buchholz and Smits, 2008), a vibrating sphere (Govardhan and Williamson, 2005) and spherical flapping foils (Dong et al., 2006). Our previous conclusion based on tethered simulations (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a) that the St* is the primary parameter determining the single vs double row wake structure is still in agreement with the current results, because in our self-propelled simulations the single row is observed in the case with low St* while the double row wake emerges in cases with higher St*. As mentioned in the Results section, the body shape affects the small features of the wake, such as smoothness and complexity of the loops, while the kinematics does not appear to have much effect on such subtle wake features.
We should also comment on the wake of different swimmers in the viscous case R1 (Re*∼300). In this case the wake is laminar and does not break up into small-scale structures but rather becomes weak and dissipates quickly as can be observed in Fig. 10. This is in agreement with the experimental results of Muller et al. who observed that the wake of zebra fish larvae tends to die off quickly due to low Re* effects (Muller et al., 2008). Finally, in agreement with our previous tethered simulations (Borazjani and Sotiropoulos, 2008; Borazjani and Sotiropoulos, 2009a), we also see that as Re* increases (from case R1, R2 and R3) the boundary-layer thickness decreases and a thinner layer of the plotted q iso-surface is visible on the fish surface in Figs 8, 9, 10. This finding is in agreement with a similar observation by Muller et al. (Muller et al., 2008).
Effects of swimming kinematics
To study the effects of kinematics we compare the swimmers with the same body but with different kinematics, i.e. MM is compared with ML and LL is compared with LM. From Table 1 it can be observed that the swimmers with carangiform kinematics (MM and LM) achieve in the inviscid case R3 swimming speeds 23% and 13% higher than swimmers with anguilliform kinematics and the same body shape (ML and LL), respectively. However, in the viscous and transitional regimes cases R1 and R2, the swimmers with anguilliform kinematics consistently achieve higher velocities relative to those with carangiform kinematics, e.g. ML and LL swimmers achieve 10% and 20% higher swimming speed in case R1 and 3% and 7% in case R2 than the MM and LM swimmers, respectively. In fact, even for the case R3 the swimmers with anguilliform kinematics (ML and LL) initially are ahead of the swimmers with carangiform kinematics (MM and ML) but the swimmers with carangiform kinematics ultimately catch up and reach higher velocity in the quasi-steady state. This is explained by the fact that the carangiform kinematics produces thrust via a lift-based mechanism similar to the pitching and heaving foils while anguilliform kinematics produces thrust via an undulatory pump mechanism (Blake, 2004; Muller et al., 2001; Webb, 1975), i.e. each part of the body generates thrust by accelerating the adjacent fluid by the body undulations [see also the related discussion and figs 8 and 9 of Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009a)]. If there is no flow, the lift-based mechanism (carangiform kinematics) cannot produce thrust and, as discussed in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2009a), its performance naturally improves as swimming speed increases. However, the undulatory pump mechanism (anguilliform kinematics) can produce more thrust at low swimming speeds because there is larger difference between the swimming speed U and the body wave speed V, which enables the body to push the fluid backwards more effectively. Therefore, the anguilliform kinematics initially outperforms the carangiform kinematics when the swimming speed is low but as the speed increases the carangiform kinematics works more effectively and takes over. Furthermore, these results are consistent with the observations on linear acceleration of eels, in which the eels with higher acceleration showed larger amplitude of motion throughout their body, i.e. they swam with more pronounced anguilliform-type kinematics at higher accelerations (Tytell, 2004b). This suggests that the different swimmers might change their kinematics initially during the start of swimming to achieve better performance, a trend that has also been documented in a number of experiments (Blake and Domenici, 2000; Weihs, 1973; Weihs, 1974).
From the power consumption standpoint, we can see from the Table 1 that the swimmers with anguilliform kinematics tend to use less power in the transitional (case R2) and inertial (case R3) regimes. More specifically, in case R3 the swimmers LL and ML need only half the power while in case R2 they need about 13% less power than the LM and MM swimmers, respectively. In the R1 case, the MM and ML swimmer use practically the same amount of power while the LL swimmer use about 10% more power than the LM swimmer. Concerning thrust production, we observed from the results for in table 1 the same trend as the swimming speed. Namely, in the viscous (case R1) and transitional (R2) regimes the anguilliform kinematics produces more thrust while in the inertial regime (case R3) carangiform kinematics produces more thrust. More specifically, in case R1 the ML swimmer produces 17% more thrust than the MM swimmer while the LL swimmer produces 22% more thrust than the LM swimmer. In the R2 case, the difference is smaller than case R1 with the ML swimmer producing only 5% more thrust than the MM swimmer while the LL swimmer produces only 7% more than the LM swimmer. In the R3 case, the ML swimmer produces almost half of the thrust of the MM swimmer and the LL swimmer produces about 34% less thrust than the LM swimmer.
Taking into account the efficiency values of Table 1, we find that: (1) the MM and LM swimmers are 16% and 2% more efficient than the ML and LL swimmers, respectively, in the inertial regime (case R3); and (2) the ML and LL swimmers are 20% and 27% more efficient in the viscous (case R1) and 16% and 18% more efficient in the transitional (case R2) regime than MM and LM swimmers, respectively. Using the mean efficiency (η*) as a performance metric, we find the same trends as using the Froude efficiency (ηf). Therefore, it can be concluded that the anguilliform kinematics not only reaches higher velocities but is also more efficient in the viscous (case R1) and transitional (case R2) regimes. However, in the inertial (case R3) regime the carangiform kinematics both reaches higher velocities and is more efficient than anguilliform kinematics. This is in agreement with our heuristic argument in our previous publication (Borazjani and Sotiropoulos, 2009a) where we showed that the efficiency of the anguilliform kinematics, which works with the undulating pump mechanism, decreases in the inviscid environment while the efficiency of carangiform kinematics, which works similar to a pitching and heaving airfoil, increases in the inertial regime.
Finally, it is worth noting that based on the r.m.s. of the swimming speed it can be observed that the anguilliform kinematics are characterized by lower r.m.s. relative to the carangiform kinematics for all bodies and in all cases. In our previous work we showed that the anguilliform swimmers produce thrust more smoothly than the carangiform swimmers (Borazjani and Sotiropoulos, 2009a), which in fact is in agreement with experimental observations. Observations of swimming eels (anguilliform) have revealed about 10% (Muller et al., 2001) or 4% (Tytell and Lauder, 2004) velocity fluctuations about the mean velocity U while swimming mullets (carangiform) have been found to exhibit velocity fluctuations more than 20% of the mean (Muller et al., 1997; Nauen and Lauder, 2002). Therefore, according to our results the smaller velocity fluctuations and the smoother force variation of anguilliform swimmers are at least partly due to the kinematics. The effect of kinematics on velocity and force fluctuations can be explained by the fact that the lift-based mechanism (carangiform kinematics) inherently generates much larger pressure variations than the undulatory pump mechanism (anguilliform kinematics). In the next section we will discuss the effect of form (body shape) on the velocity fluctuations as well.
Body shape effects
To elucidate the effects of body shape on the hydrodynamics of undulatory swimming, we compare the swimmers with the same kinematics but with different bodies, e.g. MM is compared with LM and LL is compared with ML. From Table 1 it can be observed that in all three hydrodynamic environments (cases R1, R2 and R3) the MM swimmer reaches a higher velocity than the LM swimmer. Similarly, in all cases the ML swimmer reaches a higher velocity than the LL swimmer. Therefore, in all cases the mackerel body achieves higher swimming speeds than the lamprey body. These swimmers, however, use more power and produce more thrust than the swimmers with the lamprey body as evident from the mean non-dimensional power and thrust values in Table 1. For example, the MM swimmer uses about 2.5 and 6 times more power than the LM swimmer in cases R2 and R3, respectively. The higher power consumption and thrust production in swimmers with the mackerel body could be explained by the fact that the mackerel body has a higher aspect ratio. That is, among virtual swimmers having the same length as the mackerel body has larger side area than the lamprey body, i.e. analogous to considering two flat plates flapping in exactly the same way: the plate with higher side area will require more power and produce more thrust.
Comparing swimming efficiency, the MM and ML swimmers have higher efficiency than the LM and LL swimmers, respectively, in the viscous (case R1) and inertial regimes (case R3) while the reverse is true in the transitional regime (case R2). Therefore, the mackerel body is more efficient in the viscous and inertial regimes while the lamprey body is more efficient in the transitional regime. We argue that this trend could be partly due to the combined effects of St* and aspect ratio at least in the Re* of order 104 range. Experiments with high (infinite) and low aspect ratio oscillating foils (at Re* of order 104 and higher) have shown that for the former efficiency is maximized at around St*=0.25 (Anderson et al., 1998; Triantafyllou et al., 1993) while for the later the maximum efficiency occurs at higher St* (Buchholz and Smits, 2008). The mackerel body is best represented by a high aspect ratio foil and in the inertial regime propels itself at St* close to that at which high-aspect ratio foils maximize efficiency. The lamprey body however is analogous to a low aspect ratio foil and should be more efficient at higher St*, which occur in the transitional regime. These experiments, however, does not explain the higher efficiency of the mackerel body in the viscous regime and further research is needed. A possible explanation might be the Re* effect in the viscous regime because the Re* difference between the swimmers with the mackerel body and the lamprey body in this regime is larger than any other regime. In fact, in the viscous regime the Re* of the MM and ML swimmers are 57% and 40% higher than the LM and LL swimmers, respectively, compared with 16% and 12% in the transitional regime. The large difference between the Re* of the swimmers with the mackerel body and the lamprey body in case R1 might prevent us from reaching a conclusion with respect to the effect of body shape in the viscous regime as the higher efficiency of the swimmers with the mackerel body might be due to the higher Re* and not necessarily the body shape effect.
Regarding the r.m.s. of the swimming speed fluctuations, it can be observed that the mackerel body always has higher r.m.s. relative to the lamprey body for all kinematics and in all cases. The larger aspect ratio and side area of the mackerel body results in larger force on the body with larger fluctuations. Consequently, the body shape is also an important parameter in the observed in nature lower velocity fluctuations and smoother force record of anguilliform swimmers. As discussed in the previous section, the anguilliform kinematics also lowers the velocity fluctuations. However, the r.m.s. values in Table 1 show that the effect of body shape is much stronger than the kinematics as the swimmers with the mackerel body and anguilliform kinematics (ML) always have higher r.m.s. than the swimmers with the lamprey body and carangiform kinematics (LM), e.g. in the R1 case the r.m.s. U* for ML swimmer is 0.013 while for the LM swimmer is 0.0049.
In this work we constructed self-propelled virtual swimmers and used them to explore for the first time hypothetical scenarios that are extremely challenging to be explored experimentally. One such scenario raised previously by Tytell and Lauder (Tytell and Lauder, 2004), when discussing the importance of body shape on the hydrodynamics and comparing anguilliform and carangiform swimmers, was ‘what if’ a mackerel swam like an eel. With our approach we carried out a systematic investigation of the hydrodynamics of a mackerel swimming like an eel and vice versa and conclusively clarified previous hypotheses regarding the effects of body shape and kinematics on the relative performance of anguilliform and carangiform swimming. Our findings make a strong case in support of the notion that hydrodynamic considerations might have played an important role in the evolution of different fish species and their modes of swimming.
It is important to mention that in this work we did not study the effects of dorsal/pectoral fins as they were removed from our virtual swimmers. Furthermore, we did not consider the dorsoventral asymmetry of tail motion of real carangiform swimmers (Gibb et al., 1999). In addition, in our study the body wave passes through the tail as well, i.e. there is no phase angle difference between the body and tail wave. Yet, it is already known that the phase angle between the body and tail waves can affect the performance of fish-like swimming (Zhu et al., 2002). We have also not studied the effect of body amplitude on swimming whereas it has been observed that eels increase their body wave amplitude as the Re* is increased (Tytell, 2004a). Similarly, the chub mackerel increases its tail amplitude with swimming speed but the kawakawa tuna does not change its tail amplitude (Donley and Dickson, 2000). Finally, we have only considered straight, inline swimming in a uniform ambient flow environment and thus neglected any effects ambient vortical structures and/or turbulence might have in swimming performance, which previous work has shown to be important (Beal et al., 2006; Liao et al., 2003). In their work the reduced mass of all virtual swimmers were equal while the reduced mass of different swimmers in nature can play a role in velocity fluctuations during steady swimming. This issue, however, was not explored in the present study so that we can focus with certainty on the effects of body shape and kinematics. In conclusion, in spite of several novel insights into the hydrodynamics of undulatory swimming we have contributed in this and our previous papers many important issues remain to be explored. It is our hope that the computational algorithms we have developed and our overall approach to the problem will provide the biological community with the tools required to tackle these problems in the future.
Additional materials and methods
Our numerical method in the non-inertial reference frame has been fully validated in Borazjani (Borazjani, 2008). In this section we report results from Borazjani (Borazjani, 2008) to demonstrate the ability of our numerical method to accurately predict the forces acting on an immersed body in the non-inertial reference frame and perform accurate FSI simulations.
(A) Validation of the method for calculating hydrodynamic forces: forced inline oscillations of a cylinder in a fluid initially at rest
To validate the ability of the method to predict the hydrodynamic force and its pressure and viscous contributions in the non-inertial reference frame, we consider the case of a circular cylinder starting to oscillate in the horizontal direction in a fluid initially at rest. This case is exactly the same as the one used in Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008) to validate the version of our method in the inertial reference frame, i.e. the background fluid mesh was fixed in space and the cylinder oscillated over the fixed mesh. The simulation in the present study is performed in the non-inertial reference frame, i.e. the background fluid mesh is attached to the center of the cylinder and the flow equation (Eqn 3) is solved in the non-inertial frame that oscillates with the cylinder.
The results from both inertial and non-inertial reference frames are reported and validated against the benchmark experimental and computational results of Dutsch et al. (Dutsch et al., 1998). The forced motion of the cylinder is described by a harmonic oscillation: (A1) where xc is the location of the center of the cylinder, fc is the oscillation frequency and Am is the oscillation amplitude. The flow induced by such oscillations is governed by two non-dimensional parameters: (1) the Reynolds number Re=UmD/ν based on the maximum oscillation velocity Um, cylinder diameter D, and the fluid kinematics viscosity ν; and (2) the Kuelegan–Carpenter number KC=Um/fcD. According to Eqn A1, the KC number is equal to 2πAm/D. The computations are performed at Re=100 and KC=5 for which both experimental and numerical results have been reported by Dutsch et al. (Dutsch et al., 1998). The domain, as in our previous paper (Borazjani and Sotiropoulos, 2008), is discretized with a mesh consisting of 721×481 nodes in the inline (oscillatory) and transverse directions, respectively. 300×100 nodes are distributed uniformly in a 3D×D box, which contains the cylinder during the oscillations. The domain outer boundaries are placed 50D from the initial position of the cylinder, and Neumann boundary condition (∂ui/∂nj=0, where ui is the ith component of the velocity and nj is the normal to the outer boundary surface) has been used.
Fig. A1 compares the computed hydrodynamic forces in the inertial frame (in red) and the non-inertial frame (in blue) with the computational results of Dutsch et al. (Dutsch et al., 1998) (in black). The solid, dotted and the dashed lines represent the total horizontal hydrodynamic force, and its pressure and viscous components, respectively. It is clear that the calculated forces in both reference frames are in excellent agreement with the results of Dutsch et al. (Dutsch et al., 1998).
Fig. A2 compares the horizontal velocity profiles at x1=–0.6D calculated in the inertial (solid lines) and non-inertial (dashed lines) reference frame with the measurements of Dutsch et al. (Dutsch et al., 1998) (square symbols) for three different phase angles (ϕ2πfct). Similar to the force comparisons, the velocity profiles are in excellent agreement with the measurements.
Finally, Fig. A3 shows the calculated instantaneous vorticity field at four different phase angles in the non-inertial reference frame, which are dominated by two counter-rotating vortices. The computed results in the non-inertial frame is identical to the computational results of Dutsch et al. (Dutsch et al., 1998) and Borazjani and Sotiropoulos (Borazjani and Sotiropoulos, 2008), which are not shown herein but are reported in the same format in their papers.
(B) Validation of the FSI approach: VIV of an elastically mounted cylinder
The flow past a single, elastically mounted two-dimensional cylinder has served as the generic VIV model problem and has been widely studied both numerically and experimentally. A 2-D circular cylinder with mass m is elastically mounted (with a spring of stiffness k and damping factor c) in a uniform flow of velocity U, and is free to vibrate in the direction perpendicular to the flow direction. The cylinder diameter is D. The Reynolds number of the flow is defined as Re=UD/ν. The dynamic equation governing the motion of the cylinder in non-dimensional form reads as follows: (B1) where xi are the components of the cylinder center of mass position vector non-dimensionalized by the diameter D, Cxi=2Fi/ρU2D is the force coefficient (Fi is ith component of force vector F), and other non-dimensional parameters are defined as follows: the damping coefficient: (B2) the reduced velocity: (B3) the reduced mass: (B4) These are the three important non-dimensional parameters that along with Re* govern the behavior of the system.
It is well known that when the natural frequency of the cylinder falls within the so-called ‘lock-in’ region, large amplitude vibrations are excited (Blevins, 1990). Within this region, the vortex shedding frequency changes to match the frequency of the structure's motion, resulting in the observed large vibration amplitudes. The synchronization frequency is not necessarily the natural frequency of the structure and has often been observed to exceed it significantly (Govardhan and Williamson, 2005; Sarpkaya, 2004).
The FSI solver in the non-inertial frame is validated by investigating the ‘lock-in’ phenomenon for a case for which benchmark numerical results are available in the literature. One set of simulations was performed in the inertial reference frame, i.e. the background fluid mesh is fixed in space and the cylinder oscillates over the fixed mesh (Borazjani, 2008; Borazjani et al., 2008). Another set of simulations are performed in the non-inertial reference frame i.e. the background fluid mesh is attached to the center of the cylinder and the flow equation (Eqn 3) is solved in the non-inertial frame that oscillates with the cylinder (Borazjani, 2008).
The computational domain is a rectangular box with dimensions 32D×16D. The cylinder is initially located on the horizontal axis of symmetry of this box 8D from the inlet. The computational domain is discretized with 281×241 grid nodes, which are clustered toward the four sides of a square box centered on the cylinder. This square box is discretized with a uniform 50×50 mesh, which corresponds to a near-cylinder grid spacing of Δh=0.02D. A time step of Δt=0.02 is used.
To investigate the ‘lock-in’ phenomenon, the Re* and reduced mass of the system are fixed (Re=150 and Mred=2) and the natural frequency of the system is systematically varied. This is accomplished by varying the reduced velocity Ured with increments of 1 within the range of 3≤Ured≤8. The resulting variation of the maximum displacement of the cylinder with Ured conditions is plotted in Fig. B1. As seen in the figure, large amplitude vibration (exceeding 10% of the cylinder diameter) is observed within 4≤Ured≤7 while outside this region the cylinder vibration amplitude is drastically reduced. Based on these results the lock-in region for this system is 4≤Ured≤7 with maximum vibration amplitude of about 0.5D attained for Ured=4. To quantitatively validate our method, we also include in this figure the recent results of Ahn and Callinderis (Ahn and Kallinderis, 2006) who employed an unstructured, finite-element ALE approach. As seen in the figure the three results are in excellent agreement with each other.
This work was supported by NSF Grants 0625976 and EAR-0120914 (as part of the National Center for Earth Surface Dynamics) and the Minnesota Supercomputing Institute. We are grateful to Professor Smits at Princeton University for providing the lamprey morphology data from the CT data of Professor Fish.
LIST OF SYMBOLS AND ABBREVIATIONS
- tail-beat amplitude
- amplitude envelope
- width of the wake
- oscillation amplitude
- body/caudal fin
- damping factor
- critical damping factor
- mean force coefficient
- mean non-dimensional force in the axial direction
- mean power coefficient
- mean non-dimensional power
- mean thrust coefficient
- mean non-dimensional thrust
- computed tomography
- cylinder diameter
- elongated body theory
- tail-beat frequency
- oscillation frequency
- natural frequency
- non-dimensional tail-beat frequency
- force vector exerted on the virtual swimmer's body by the fluid
- fluid–structure interaction
- lateral displacement of fish body
- lateral velocity of fish body
- maximum displacement of the tail
- lateral excursion of the body at time t and location z
- hybrid Cartesian immersed-boundary
- immersed boundary
- wave number
- Kuelegan–Carpenter number
- total body length
- loose coupling fluid–structure interaction
- lamprey swimming like a lamprey
- lamprey swimming like mackerel
- mass of the virtual swimmer
- mackerel swimming like a lamprey
- mackerel swimming like a mackerel
- miles per gallon
- reduced mass
- normal vector the surface
- jth component of the normal vector
- non-dimensional pressure
- particle image velocimetry
- power losss due to lateral undulations of the fish body
- rotation matrix of the non-inertial frame relative to the inertial frame
- non-dimensional position in the inertial frame
- non-dimensional position relative to the non-inertial frame
- root mean squared
- mean Reynolds number based on swimming speed U*
- non-dimensional viscosity or Reynolds number based on characteristic velocity Uo
- symmetrical parts of the velocity gradient
- strong coupling fluid–structure interaction
- mean Strouhal number based on swimming speed U*
- Strouhal number based on the characteristic velocity Uo
- mean thrust
- period of a tail-beat cycle
- non-dimensional Cartesian absolute velocity vector of the fluid in the inertial reference frame
- velocity of the non-inertial frame or the center of mass
- maximum oscillation velocity
- fluid velocity normal to the body
- non-dimensional Cartesian relative velocity vector of the fluid in the non-inertial reference frame
- fluid velocity vector tangential to the body
- mean swimming speed
- mean non-dimensional swimming speed
- characteristic (tether) velocity
- reduced velocity
- vortex-induced vibration
- position vector of center of mass non-dimensionalized by L
- location of the center of the cylinder
- axial direction measured along the fish axis from the tip of the fish's head
- under-relaxation coefficient
- a dynamically evolving surface
- mean efficiency (velocity over power)
- Froude efficiency
- kinematic viscosity
- fluid density
- viscous stress tensor
- angular frequency
- damping coefficient
- angular velocity
- © 2010.