## ABSTRACT

Tadpoles are unusual among vertebrates in having a globose body with a laterally compressed tail abruptly appended to it. Compared with most teleost fishes, tadpoles swim awkwardly, with waves of relatively high amplitude at both the snout and tail tip. In the present study, we analyze tadpole propulsion using a three-dimensional (3D) computational fluid dynamic (CFD) model of undulatory locomotion that simulates viscous and unsteady flow around an oscillating body of arbitrary 3D geometry. We first confirm results from a previous two-dimensional (2D) study, which suggested that the characteristic shape of tadpoles was closely matched to their unusual kinematics. Specifically, our 3D results reveal that the shape and kinematics of tadpoles collectively produce a small ‘dead water’ zone between the head-body and tail during swimming precisely where tadpoles can and do grow hind limbs – without those limbs obstructing flow. We next use our CFD model to show that 3D hydrodynamic effects (cross flows) are largely constrained to a small region along the edge of the tail fin.

Although this 3D study confirms most of the results of the 2D study, it shows that propulsive (Froude) efficiency for tadpoles is overall lower than predicted from a 2D analysis. This low efficiency is not, however, a result of the high-amplitude undulations of the tadpole. This was demonstrated by forcing our ‘virtual’ tadpole to swim with fish-like kinematics, i.e. with lower-amplitude propulsive waves. That particular simulation yielded a much lower Froude efficiency, confirming that the large-amplitude lateral oscillations of the tadpole do, indeed, provide positive thrust.

This, we believe, is the first time that the unsteady flow generated by an undulating vertebrate has been realistically modelled in three dimensions. Our study demonstrates the feasibility of using 3D CFD methods to model the locomotion of other undulatory organisms.

## Introduction

Undulatory propulsion for most aquatic vertebrates is accomplished by sending alternating waves down the body towards the tip of the tail. These waves create a jet in the wake and hence a forward force, namely thrust. Undulatory swimming is a highly effective means of continuous locomotion by aquatic vertebrates. Given that most vertebrates are aquatic and have had hundreds of millions of years to adapt to that environment, it is reasonable to suppose that natural selection has closely, if not optimally, matched the kinematics of these animals to their geometry. Indeed, the sleek shape and graceful kinematics of many pelagic fishes support this view. Tadpoles, however, are characterized by an awkward shape and equally awkward kinematics. They have a globose head-body, with a laterally flattened tail abruptly appended to it. Tadpoles swim with large lateral deflections compared with those of normal teleost fishes (see Wassersug and Hoff, 1985;

Wassersug, 1989, for data on tadpoles; Videler, 1993; Rowe *et al.* 1993, for fishes). We previously used computational fluid dynamics (CFD) methods to explore the unsteady hydrodynamics of undulatory swimming by a 2D tadpole-shaped body (Liu *et al.* 1996). In that 2D study, we explored the effects of the strange shape of tadpoles on both their hydrodynamics and overall biology. Results from that study suggested that the shape of tadpoles was well matched to their kinematics. Furthermore, the 2D study suggested that the morphology and locomotion of premetamorphic anurans were directly linked to the obligatory requirement of these organisms to develop hindlimbs and ultimately to transform into frogs.

Here, we extend our 2D modelling of tadpole swimming to the more realistic 3D situation. Essentially we ask how the 3D effects of unsteady undulatory hydrodynamics by swimming tadpoles affect their locomotor performance. We approach this problem using a robust Navier–Stokes (N-S) equation solver (Liu *et al.* 1996; H. Liu and K. Kawachi, in preparation) that can simulate highly unsteady flows around an undulating locomotor that has an arbitrary 3D geometry. This N-S solver is capable of accurately mimicking the real movements of undulating vertebrates, such as tadpoles, in the aquatic medium. The validity of the 3D steady simulation was verified by comparison of flow pattern, pressure distribution and integrated drag with reliable experimental data collected on ship models (Liu, 1995). The time-accurate reliability was further confirmed by an extensive study involving grid refinement and by comparison with limited experimental results related to an oscillating hydrofoil (H. Liu and K. Kawachi, in preparation).

In the present study, we use this CFD approach to elucidate the unsteady 3D hydrodynamics and the instantaneous force-related quantities produced during undulating (stroking) cycles by tadpoles. We first extend our previous 2D analyses by investigating the unsteady hydrodynamics of an undulating 3D ‘model’ bullfrog tadpole. This analysis confirms that during swimming there is a ‘dead water’ region at the base of the tadpole’s tail precisely where the hindlimbs can, and do, develop. Since the CFD models quantitatively reproduce the 3D effect of fluid flowing around undulating tadpoles, they allow us to examine the cost of high-amplitude lateral oscillations of the tadpole’s head (rostral wobble) during swimming. The 3D modelling suggests that those rostral movements by the tadpole are, in fact, nicely coupled to their propulsive wave in the tail effectively to create thrust. We demonstrate how effective the tadpole kinematics is by forcing our ‘virtual’ tadpole to swim with the kinematics of a fish; i.e. with less rostral wobble and lower-amplitude waves at the tip of the tail.

Our 3D CFD analysis of tadpole locomotion provides an introduction and overview of the unsteady hydrodynamics of 3D undulatory swimming by vertebrates. Our study specifically confirms that tadpoles swim efficiently, with an elegant coupling between their specific kinematics and their unique morphology.

## Materials and methods

### Navier–Stokes equations resolution

A robust and effective solver of the incompressible N-S equations (H. Liu and K. Kawachi, in preparation) is used in the present CFD modelling of 3D undulatory swimming. This 3D model is designed to simulate viscous and unsteady flows induced by an undulating ‘virtual’ locomotor swimming in a stream of constant velocity which corresponds to the given Reynolds number. The virtual swimming model can have the geometry and kinematics of a realistic aquatic swimmer.

The governing equations are discretized using the finite volume method (FVM) and are solved in a time-marching manner using the pseudo-compressibility technique. A third-order upwind differencing scheme, using the flux-splitting method, is implemented for the convective terms in a MUSCL fashion (ultimate conservative difference scheme; Van Leer, 1977). The viscous terms are evaluated using a second-order central differencing method based on a Gaussian integration in the FVM manner. An implicit approximate-factorization method, based on the Euler implicit scheme, is employed for the discretization of the time derivatives. In deriving the time-marching (time-dependent) scheme, the time derivatives in the momentum equations are differenced using a first-order, two-point, backward-difference implicit formula and are solved iteratively such that they approach the new velocity as the divergence of velocity approaches zero.

To analyze the fluid interactions between the unsteady flow and the undulating object, a moving grid system was introduced. We made the grids fit the deforming boundaries at each physical time step. We also made sure that there was sufficient grid density to resolve the viscous and unsteady flows both around the body surface and in the wake.

In developing boundary conditions on the body surface, a no-slip condition was used for the velocity components, which is to say that the velocity of a fluid particle was equal to the velocity of the swimmer at that point. This implies that our ‘virtual’ swimmer had, at each instant, a smooth but rigid surface. To incorporate the dynamic effect due to the rapid acceleration/deceleration of the undulating swimmer, pressure divergence at the surface stenciles was used for the pressure condition, derived from the local momentum equation. Details of the algorithm can be found in H. Liu and K. Kawachi (in preparation).

All the computations were carried out using an HP Apollo 9000 (755 model) series workstation; it took approximately 120–130 h of computer time for each case.

### Definition of geometry and grid system

The geometry of our ‘virtual’ swimmer, illustrated in Fig. 1A, was based on the same *Rana catesbeiana* larva (Wassersug and Hoff, 1985) used in our previous 2D analysis (Liu *et al.* 1996). Its total length was 4.7 cm.

We assumed that our ‘virtual’ tadpole had no protruding limbs, which is a realistic assumption for a bullfrog tadpole of this size and stage of development. Tadpoles at this stage have their forelimbs covered by an opercular flap. Furthermore, from the results of our 2D study (Liu *et al.* 1996), we can be fairly confident that small hindlimbs, which do not project laterally, would have a negligible effect on the flow.

We ignored the oral disc and the denticles (keratinized teeth) on the tadpole. *Rana* tadpoles have a small mouth surrounded by an oral disc that is, indeed, small. The maximum surface area of the oral disc in *Rana* is commonly less than 2 % of the total surface area of the tadpole. We also ignored the spiracle, which is the opening on the left side of the tadpole through which water that irrigates the gills is expelled. This aperture is smaller than the oral disc and flush with the abdominal wall in bullfrog tadpoles. There is no evidence that *Rana* tadpoles eject water from this asymmetric opening during normal bouts of swimming.

To define realistically the shape of our virtual tadpole, we based its geometry on previously published, high-quality side- and top-view illustrations of *Rana* tadpoles (Fig. 1B,C). We note that the body outlines (cross sections) of a large number of aquatic undulating swimmers can be well approximated by a continuous analytical function, if one ignores all fins and other appendages. We based such a function on the 2D shape of multiple cross sections along the body length of tadpoles taken from these illustrations. The shape of each cross section was determined by two curves; one for the dorsal line on the sagittal plane (Fig. 1B) and one for the lateral body line (Fig. 1C) on the horizontal plane. In other words, we used a side view to digitize the dorsal and ventral margins of the tadpole, and a top view to digitize the lateral margins (i.e. side line) of the head-body and tail.

Generalized pond tadpoles, such as our model *Rana* larva, have a density close to that of water (Campeny and Casinos, 1989) and normally swim horizontally without any apparent inclination of the body (Wassersug and Hoff, 1985; Wassersug, 1989). These facts suggest that the slight dorso-ventral asymmetry of real tadpoles, even if it leads to some discrepancy of pressure distribution across the horizontal plane (*xz* plane), would not produce a significant vertical force. Thus, to simplify our problem, we considered our 3D virtual tadpole to have a symmetrical horizontal plane (*xz* plane) as illustrated in Fig. 1A, which reduced our computing time by half. We used a single elliptical function to describe the head-body because it fits a typical body of revolution (Fig. 1A). However, the tadpole tail (Fig. 1A,B) is composed of dorsal and ventral fins that are of different height. We therefore introduced a smoothed thin plate, with a cross section composed of two elliptical curves, to define the dorsal and ventral fin margins. The heights of these two curves were taken from the digitized dorsal and ventral lines along the longitudinal axis of the tail (Fig. 1A) with an appropriate spline fit applied. As a result, we were able to construct a complete 3D configuration for the tadpole based largely on realistic images and partially derived analytically (Fig. 1C). (It should be noted that more complicated functions could be introduced to match realistically the shape of 3D organisms with more complicated geometries.)

Considering the complicated geometry of any undulating swimmer, we elected to use the so-called C-O type grid topology (pattern) for the FVM. This grid topology is capable of closely matching the girthwise lines of the body (i.e. the O-type at each section), particularly at the snout and the tail tip. It is also effective in resolving the boundary layer on the surface of the posterior tail and in solving the vortex in the wake (i.e. the C-type). To ensure grid generation with high quality, a multi-block grid system was introduced, with each block corresponding to each portion (i.e. the head-body, the tail and the wake) of our undulating swimmer. A grid system of approximately 70 000 grids was used in the present simulation. In addition, to ensure grid quality (i.e. to avoid overlapping of grid lines when regridding during swimming), the configuration at the tail tip was slightly modified, as illustrated in Fig. 1B. Specifically, for our simulations, the tail is truncated with a non-zero height, instead of tapering to a point. To minimize the influence of this modification, we kept the height at the truncation quite small – one-fifth of the maximum height of the dorsal fin – and kept the body length unchanged. This led to a negligible increase in surface area, i.e. less than 1 % of the total surface area of the tadpole. As discussed below in detail, this slight modification at the tail tip caused insignificant changes in flow and force production. Furthermore, it departs little from reality, since *Rana* tadpoles in nature often have the tips of their tails bitten off by predators and can survive even greater caudal amputations with little increased risk of predation (Wilbur and Semlitsch, 1990; for *Hyla*, see also Figiel and Semlitsch, 1991).

### Kinematics of undulatory swimming

The kinematics of undulatory swimming used in the present 3D analysis was based on the forward locomotion in a straight line of a real *Rana catesbeiana* larva (Wassersug and Hoff, 1985), as was used in our previous 2D analyses (Liu *et al.* 1996). As pointed out by Wassersug and Hoff (1985), the common head-bodies of tadpoles are of a spheroidal shape that does not noticeably deform with locomotion. Thus, we assumed that the head-body of our tadpole-shaped model could oscillate, but not deform. In contrast, the thin tail fins, which lack any cartilaginous or osseous support, are extremely flexible (see Long and Nipper, 1997) and could, in principle, be passively deflected out of the vertical plane by water pressure during swimming. Fortunately, high-speed films of bullfrog tadpoles swimming in a flow tank show little passive bending of the fins (Wassersug and Hoff, 1985). Consequently, we have assumed that a two-dimensional cross section through a tadpole changes little in shape as one moves more dorsally or ventrally from the midplane of the body. Similarly, we assume that the dorsal and ventral fins are held vertically and are sufficiently stiff that they show no strain during undulatory swimming, other than the lateral undulation driven by the caudal muscle itself. This is tantamount to treating our virtual tadpole as an oscillating plate during swimming and requires two basic assumptions: (1) that the lateral travelling wave is one-dimensional at each cross section, and (2) that the body elongates slightly on alternate sides during undulatory swimming.

*a*

_{i}(

*x*) represents amplitude at point

*i*in the centre plane, λ is wavelength,

*T*is period,

*h*

_{i}(

*x*,

*t*) is the elevation of the centre plane,

*t*is time and

*x*is the coordinate in the

*x*direction, as illustrated in Fig. 1A, corresponding to the body length. Note that this ‘virtual’ tadpole swims like a waving plate along a single horizontal plane, and the lateral deviation of its centre plane

*h*

_{i}(

*x*,

*t*) is a function of longitudinal position

*x*and time

*t*only, and has nothing to do with movement in the

*y*and

*z*planes. As an aside, the present definition of the travelling wave can in principle be extended, by taking the

*y*coordinate in equation 1 into consideration, to mimic full 3D movements of an undulatory swimming animal capable of moving in all three (

*x, y*and

*z*) planes. Following Liu

*et al.*(1996), the amplitudes of

*a*

_{i}(

*x*), as illustrated in Fig. 2, are again determined by using the spline interpolation from five original maximum amplitudes,

*a*

_{i}, along the length,

*L*, of the tadpole, taken directly from Fig. 4 in Wassersug and Hoff (1985). These values are: at the snout (

*x*=0,

*a*=0.05

*L*), at the otic capsule (

*x*=0.19

*L, a*=0.005

*L*), at the base of the tail (

*x*=0.384

*L, a*=0.04

*L*), at the midtail (

*x*=0.692

*L, a*=0.1

*L*) and at the tail tip (

*x*=1.0

*L, a*=0.2

*L*).

*L*, which represents the length of the tadpole, and a reference velocity

*U*that is its forward velocity, the Reynolds number (

*Re*) is defined by: where 𝒱 is kinematic viscosity, with a value of 1.533×10

^{−6}m

^{2}s

^{−1}. Given a body length of 4.7 cm for our virtual bullfrog tadpole, and a reasonable and realistic forward swimming speed of 5

*L*s

^{−1}, we calculate an

*Re*of 7200, as in our previous 2D analyses (Liu

*et al.*1996). Furthermore, by introducing a reduced frequency (

*k*) of 2π

*fL*/2

*U*, where

*f*is frequency, we can restate equation 2 in a simplified form: Subscript i denotes the grid points of the centre line on the centre plane. The reduced frequency

*k*is evaluated to be 5.843, corresponding to the forward swimming speed of 5

*L*s

^{−1}, from a plot of forward velocity against tail-beat frequency in Wassersug and Hoff (1985). The overall propulsive wavelength is taken as 0.87

*L*, on the basis of empirical data of 0.8±0.1

*L*from the same study.

Note that our virtual tadpole could be made to swim more like other aquatic vertebrates, including various teleost fishes, by appropriate modifications to the parameters of wavelength, reduced frequency and the amplitudes of the propulsive wave at five or more points along the length of the animal. We compare the hydrodynamics of a tadpole swimming like a fish with that of one swimming like a tadpole by forcing our ‘virtual’ tadpole to swim like a saithe, *Pollachius virens* (see Fig. 2). The kinematic data for the saithe are taken directly from Videler (1993). The fish’s key parameters were a wavelength of 0.86*L*, a body length of 0.35 m, a forward swimming speed of 1.05 m s^{−1} and, hence, a calculated reduced frequency of 3.653. This yielded a *Re* of 2.4×10^{5}.

## Results and discussion

### CFD-visualized flows of tadpole locomotion

The flow past a non-undulating simulated 3D tadpole, with its tail kept straight back at *Re*=7200, was computed first to obtain an overview of the steady flow around a tadpole-shaped body and to assess dead drag. Iso-speed contours at three cross sections (the middle of the head-body, the base of the tail and the midtail) are illustrated in Fig. 3. The figure illustrates the boundary layer development, as well as the velocity profiles on the horizontal and vertical planes.

Unlike the results we obtained in our 2D analyses – in which a pair of large-scale vortices as well as strong back flow were observed over the whole tail – a separated region, namely the ‘dead water’ zone, is detected in the 3D analysis only at the base of the tail. This dead water zone is quite small. However, it is precisely where the hindlimbs develop, endorsing our view that the crease between the head-body and the tail of tadpoles provides a space where the hindlimbs can grow without impeding flow.

The shape of the 3D simulated tadpole in side view (as in Fig. 1B) is quite close to a ‘streamlined’ body and, thus, produces a very thin boundary layer at the dorsal fin. This 3D effect of the strong crosswise (girthwise) velocity gradient, and hence pressure gradient, appears to create strong secondary flow (cross flow) and hence significantly reduces the separated region behind the head-body over the whole tail (Fig. 3). As a result, the 3D tadpole-shaped object has relatively low overall dead drag (a drag coefficient *C*_{d} of 0.0287) with similar and equally low contributions from pressure (0.0169) and from friction (0.0118), respectively. In this respect, a non-undulating 3D simulated tadpole performs better than predicted from a 2D representation (Liu *et al.* 1996). This low dead drag for a non-undulating tadpole may be ecologically significant for tadpoles resting on the bottom in streams (see also Dudley *et al.* 1991).

Unsteady flow at the same Reynolds number was elucidated by ‘swimming’ our tadpole model (Fig. 1A) using the travelling wave given in equation 3. Flows around this virtual tadpole are illustrated in Fig. 4, visualized at the two moments when the tail is in a C-shaped curve (Fig. 4A) and an S-shaped curve (Fig. 4B). The figure shows iso-pressure contours on the body surface together with pressure contours and velocity profiles on the horizontal symmetrical midplane. Comparing these figures with similar illustrations from our previous 2D study (Fig. 6A,B in Liu *et al.* 1996) reveals very similar flow patterns at the same instants in the tailbeat cycle. Strong jets (accelerated fluid) are again detected in the wake, and the distributions of both positive and negative pressures very closely match our 2D results (Liu *et al.* 1996).

To capture an overall 3D image of the flow around the 3D undulating simulated tadpole, perspective views are given in Fig. 5A,B for the tadpole in both the C-curve and the S-curve positions. In this figure, instantaneous streamlines were released over the body surface while the surface was concurrently covered with iso-pressure contours in a fringe manner. The instantaneous streamlines show smooth flow over the whole undulating body and no large-scale longitudinal vortex at the dorsal fin even in the extreme posture, when the angle of attack for the posterior tail is greatest (i.e. in the C-curve position). We also see that the flow is three-dimensional in structure, with strong secondary flow (cross flow) at the snout and at the edges of the dorsal and ventral tail fin.

Focusing our attention on the flows at the tail, illustrated in Fig. 6A,B, we see that the instantaneous streamlines appear quite two-dimensional, with little vertical cross flow over most of the tail, except within a small region limited to the dorsal and ventral tail fin. Also, the pressure distributions are very two-dimensional over a large portion of the tail surface. This agrees with the results of Cheng *et al.* (1991), based on the 3D waving plate theory, that undulatory motion can reduce 3D effects. This supports the assumption for the present 3D simulation that the plate model for the 3D simulated tadpole tail is reasonable and helps to validate conclusions from our previous 2D study. Combined with the flow, illustrated in Fig. 3A,B, Fig. 6 confirms that a tadpole should be able to develop hindlimbs at the base of the tail without the limbs handicapping undulating locomotion. In other words, the unique shape of tadpoles is specially adapted to the requirement of these organisms to grow limbs and transform into frogs (see Liu *et al.* 1996; Wassersug, 1997).

However, the 3D results raise the question of why there is no large-scale vortex produced along the dorsal and ventral fin margins, given such large lateral deflections of the tail by the undulating tadpole. Note the very steep pressure gradient over the whole dorsal fin of the tail (Figs 5A,B, 6A,B). We believe that, unlike the situation of a fast-start from still water, a strong free stream around the tadpole during undulatory swimming at constant velocity is very likely to suppress flow separation at the dorsal and ventral fin margins. The lack of separation implies that, if some vortex were produced there, it could not persist in the strong longitudinal flow and would quickly collapse or be reduced to a small-scale phenomenon.

To obtain an overview of the vortex structure, we drew the iso-vorticity contours at two cross sections of the tadpole and also the iso-vorticity surfaces in the wake for the tadpole with its tail in the C-curve (Fig. 7A) and S-curve (Fig. 7B) positions. The resultant wake vortices are clearly detected. Large vorticity is concentrated behind the tail tip in the wake and is released at the end of each tail stroke. The vortices are quite 3D in structure, with a width approximately equal to the greatest amplitude of the tail tip (0.2*L*) in the stroke plane and a height equivalent to the greatest height of the tail. These vortices closely match the reversed von Karman street visualized in the wake of a *Rana* tadpole swimming through a thin layer of milk, as illustrated in Wassersug and Hoff (1985).

The vortex structure in the wake is a result of shedding of vortices from both sides of the plate-like tail (Fig. 7A,B). The vortices are initiated by the shear flow over the tail fin, starting just caudal to the base of the tail. They grow down the tail and are shed at the tail tip into the wake. Note that the vorticity has a rather small magnitude on the body surface compared with that in the wake. This quantitatively confirms the result, reached above, that the longitudinal vortices over the tail fins are of quite small scale, even with such rapid lateral compressing motion. An overview of the streamlines flowing over the undulating tadpole, as shown in Fig. 5A,B, further reveals that these longitudinal vortices merely lead to the freestream-wise flow spiralling down the tail to be shed into the wake.

Minor modifications of the tail tip geometry barely affect the flow in both the steady and unsteady simulations. In the steady dead drag situation, as shown in Fig. 3, flow proceeds quite smoothly from the trailing edge without any separation. The velocity profiles and the iso-speed contours show a very thin boundary layer along the dorsal and ventral fins compared with those at the midplane, because of the ‘streamlined’ shape of the tadpole in the vertical plane (Fig. 1B). In the case of undulatory swimming, this ‘streamlined’ shape in the vertical plane is evidently important in creating the strong longitudinal (freestream-wise) flow and hence reducing the large-scale separation/vortex at the dorsal and ventral fins due to the large lateral oscillations. Thus, as seen in Figs 4A,B, 5A,B, unsteady flows also show very smooth shedding into the wake along the trailing edge.

There is one caveat here concerning the cost of truncating the tail tip. In the present simulation, our virtual tadpole had a only small portion of its tail amputated and it swam with the kinematics of a tadpole whose tail was completely intact. Wassersug and Hoff (1985) found that freely swimming *Rana* tadpoles with a large portion of their tail amputated had disproportionately large-amplitude body movements and a much wider wake – indicative of a high Strouhal number and low mechanical efficiency (see Triantafyllou *et al.* 1993) – compared with free-swimming tadpoles without damaged tails.

### Leading-edge suction at the snout

As mentioned above, the flow shows quite strong 3D structure, i.e. cross flow, over the head-body, leading to a very strong girthwise (crosswise) pressure gradient. However, in contrast to what generally occurs with steady flow, a longitudinal pressure gradient, the so-called leading-edge suction region (Wu, 1971), at the snout is barely detectable during the whole undulating cycle. This is because the idea of the leading-edge suction region was based on the case of steady flow around a ‘streamlined’ fish, without accounting for any lateral oscillation of the snout. However, our CFD analysis shows that the suction region does exist, but it alternates laterally from side to side in the undulating tadpole. We believe that the large-amplitude (maximum 0.05*L*) lateral motion of the tadpole’s snout is more likely to push or suck the incoming flow to either side of the snout rather than towards its upper or lower surfaces. This can be further confirmed by the visualized flow in the side view at the instant when the tail is in the C-curve position (Fig. 6A,B). A large suction region, with a lower-pressure ‘island’ (shown in blue on Fig. 6A) is detected on the left side (looking backwards along the model), which is quite 3D in structure. A high-pressure region is also observed on the corresponding right side. Furthermore, tracing the flows through complete tailbeat cycles reveals that this low-pressure region, as well as the positive pressure region, is readily produced at each moment, but oscillates in the vertical *xy* plane. However, the local inclination of the pressure at the snout, i.e. the suction, appears not to be large enough for this low-pressure region to contribute much forward force (i.e. to thrust).

Additionally, it is obvious that the lateral oscillation of the undulating tadpole, as seen in Fig. 6A,B, leads to a very large pressure difference across the snout. This result strongly supports the analysis by Lighthill (1993) that pure sideslip of the head of clupeid fish would generate a great pressure difference, but simultaneous yawing could enormously reduce the effective pressure difference and any associated cross-flow effects. Although we did not consider the yawing angle in our CFD analysis, we did find cross flow mostly across the snout and not over the tail, where it would be most detrimental in increasing resistance to the tadpole’s motion. Also, the cross flow that we observed changed significantly with each instant in the undulating cycle (see Fig. 5A,B). This further implies that limiting this cross flow, as pointed out by Lighthill (1993), with appropriate simultaneous yawing, should reduce the resistance during undulatory swimming.

### Thrust production

The variation of the thrust coefficient (Thrust^{*}), power coefficient (Power^{*}) and drag coefficient (Drag^{*}) due to friction are plotted against time in non-dimensionalized form in Fig. 8. The thrust coefficient (see the List of symbols) is readily produced during the whole tailbeat cycle, but unlike in our 2D tadpole model where a net thrust coefficient was observed, the sum of the thrust coefficient and drag coefficient, i.e. the net thrust coefficient, now has a very small mean magnitude of 0.01. Since our undulating fish and tadpole models are forced (restrained) to swim in a stream with constant velocity at a given *Re*, zero net thrust (i.e. the sum of thrust and drag coefficients) therefore indicates trimmed swimming. In this regard, the present 3D models far more realistically represent real tadpoles swimming steadily at a given *Re* than does our previous 2D model. The drag due to friction shows little fluctuation, unlike the distinct periodicity seen for drag in the previous 2D analyses (Liu *et al.* 1996). There, the periodicity was thought to be due to the large region of flow separation, whose magnitude changed with time, as well as to the unsteady boundary layer. Here, the reduced fluctuations in drag are probably the key reason that the ‘dead water’ zone at the base of the tail and the unsteady boundary layer over the tail are 3D in structure (Figs 5A,B, 6A,B). Additionally, the oscillating body could enhance the velocity gradient inside the boundary layer and hence increase friction-based resistance. This is also confirmed by the fact that the friction drag coefficient during undulatory swimming, calculated as 0.030, is approximately 1.8 times higher than the contribution to drag from skin friction in the dead drag situation. This means that using the dead drag in calculating the force-related quantities concerning an undulatory swimming object would seriously underpredict the real drag produced (see Webb, 1975, for a discussion of dead drag *versus* drag during swimming in fishes).

Furthermore, the time course (Fig. 8) of the thrust coefficient and the required power coefficient closely match in shape those from our previous 2D analyses (Liu *et al.* 1996). This implies that the simplification of the 2D analyses did not distort the phase of force production compared with the more realistic 3D model.

We obtained a mean Froude efficiency (*Fr*) of approximately 0.45, following a conventional method (see the List of symbols) for calculating *Fr*, by using the mean thrust multiplied by the forward swimming speed and then dividing by the required power. This *Fr* value is significantly lower than that found in our previous 2D analyses. As discussed above, we believe that thrust production, which is dominated by the pressure distribution on the tail, is reduced to some extent by the 3D flow along the edges of the dorsal/ventral fins (Fig. 6A,B). *Fr* may also be sensitive to our assumption that the tadpole swims with a rigid lateral compressing motion, without considering any elastic (passive) deformation. Such locomotion may require more power output. Clearly, a truly realistic 3D model of undulatory swimming should incorporate the girthwise (crosswise) discrepancy in kinematics as well as elastic surface effects, which remain to be investigated and may be studied by a direct extension of the present method.

### Hydrodynamics of the 3D tadpole swimming in a fish mode

In our 2D CFD analyses, we reached the conclusion that, although tadpoles swim with unusually high-amplitude travelling waves, their kinematics is closely matched to their size and peculiar shape. Specifically, we argued that there was an effective interaction of flows between the globose head-body and the plate-like tail to enable tadpoles actively to produce thrust. To confirm that this is true for a realistic 3D model of a swimming tadpole, we computed the flows around our ‘virtual’ tadpole swimming in a fish-like fashion.

As noted above, we used kinematic data for the saithe in this simulation. The amplitudes of the fish’s travelling wave are plotted in Fig. 2 in a non-dimensional form. The fish has a maximum amplitude at the snout of 0.025*L* and at the tail tip of approximately 0.1*L*. These values are approximately half those of a free-swimming bullfrog larva (Fig. 2). In this final simulation, we employed the same reduced frequency of 3.484 achieved by the saithe in order to ensure that our model tadpole undulated like the fish. However, we adjusted the body length and the swimming velocity of the saithe to maintain *Re* at 7200, as before. The wavelength for the travelling wave of the saithe is taken as 0.86*L*, slightly longer than that of the tadpole.

Hysteresis of Thrust^{*}, Power^{*} and Drag^{*} due to friction are plotted in Fig. 9. A negative net thrust, i.e. drag, is obtained during the whole tailbeat cycle, with a mean value of 0.031, which is even larger than the 0.0287 dead drag value for the tadpole reported above. To understand this force production, we visualized flows for the instant when the tadpole was in the C-curve position (Fig. 10A,B). Two side views of the iso-pressure surface of the body and the instantaneous streamlines are shown. Significantly different flows and pressure distributions are observed compared with those achieved by the tadpole swimming in the tadpole mode (Fig. 6A,B). In Fig. 10, we have lost the considerable pressure distributions on the tail that created the thrust seen in Fig. 6A,B. Particularly conspicuous is the total loss of the positive pressure region at the posterior tail, which pushed fluid backwards and hence yielded the jet for forward propulsion. Furthermore, the reduction in oscillation at the snout leads to changes in the alternating regions of reduced pressure (see Wu, 1971) on the sides of the head-body. These changes detrimentally affect the separation and reattachment of flow from the dorsal tail to the posterior tail tip. All told, a tadpole swimming like a fish is seriously handicapped in thrust-production.

Finally, note that the friction-based resistance is 0.019 (Fig. 9), i.e. slightly larger than that found for the tadpole in the dead drag situation (0.0169), but reduced to approximately half that observed when our virtual tadpole swam in the tadpole mode (0.030). This reduction in drag is consistent with the fact that the saithe makes much smaller lateral excursions of its head and tail during swimming than does the tadpole. These changes in drag again confirm the fact that the friction-based resistance in undulatory swimming is dominated by the unsteady boundary layer, which can enhance the velocity gradient and hence increase the drag.

To the best of our knowledge, the present study is the first to examine unsteady viscous flows around a realistic representation of a 3D undulating vertebrate. Our 3D CFD analyses demonstrate the feasibility of CFD modelling for investigations of morphology and function in other undulating organisms. Comparisons of the integrated quantities between the present CFD analyses and other conventional hydrodynamic models should be carried out for other organisms. We are optimistic that they will reconfirm the validity of the present model. Such studies should also increase our understanding of how viscosity influences the propulsive performance of undulatory swimming.

## Acknowledgements

We thank Dr M. Matsui for kindly allowing us to reprint the illustration of a bullfrog tadpole shown in Fig. 1B. This research was supported by the Japan Scientific and Technology Corporation and the Natural Science and Engineering Research Council of Canada.

## List of symbols

*a*amplitude

*a*_{i}(*x*)amplitude at point

*i*on the centre plane (non-dimensional)- Drag*
(

*C*_{d}) drag coefficient =*D*/0.5ρ*U*^{2}*S*, where the drag is the sum of both pressure and stress forces on the swimming body (non-dimensional) *D*drag

*f*frequency (Hz)

*F*_{r}Froude efficiency =

*ThU*/*P*, the averaged hydrodynamic propulsive efficiency*h*_{i}(*x*,*t*)elevation of the centre plane of a tadpole (non-dimensional)

*k*reduced frequency

*L*body length (non-dimensional)

*P*power

- Power
^{*}power coefficient =

*P*/0.5ρ*U*^{3}*S*, where power is the hydrodynamic power *Re*Reynolds number

*S*surface area of the swimming body

*t*time

- Thrust
^{*}thrust coefficient =

*Th*/0.5ρ*U*^{2}*S*, where the thrust is the sum of the pressure force on the swimming body surface *T*period (non-dimensional)

*Th*thrust

*U*forward swimming speed (non-dimensional)

- λ
wavelength (non-dimensional)

- 𝒱
kinematic visocity

- ρ
water density

## References

*Alytes obstetricans*

*Zool. Anz.*

*J. Fluid Mech.*

*Copeia*

*Can. J. Zool.*

*ASME/FED*

*J. exp. Biol.*

*Am. Zool.*

*Frogs and Toads of Japan*

*Phil. Trans. R. Soc. Lond. B*

*J. Fluids Structures*

*J. comp. Phys.*

*Am. Zool*

*Am. Zool*

*J. exp. Biol.*

*Bull. Fish. Res. Bd Can.*

*Am. Nat.*

*Rana*tadpoles

*Copeia*

*J. Fluid Mech.*