Abstract
Numerical simulations are carried out to study the fluid dynamics of a complex-shaped low-aspect-ratio pectoral fin that performs the labriform swimming. Simulations of flow around the fin are achieved by a developed immersed boundary (IB) method, in which we have proposed an efficient local flow reconstruction algorithm with enough robustness and a new numerical strategy with excellent adaptability to deal with complex moving boundaries involved in bionic flow simulations. The prescribed fin kinematics in each period consists of the power stroke and the recovery stroke, and the simulations indicate that the former is mainly used to provide the thrust while the latter is mainly used to provide the lift. The fin wake is dominated by a three-dimensional dual-ring vortex wake structure where the partial power-stroke vortex ring is linked to the recovery-stroke ring vertically. Moreover, the connection of force production with the fin kinematics and vortex dynamics is discussed in detail to explore the propulsion mechanism. We also conduct a parametric study to understand how the vortex topology and hydrodynamic characteristics change with key parameters. The results show that there is an optimal phase angle and Strouhal number for this complicated fin. Furthermore, the implications for the design of a bioinspired pectoral fin are discussed based on the quantitative hydrodynamic analysis.
1. Introduction
In recent decades, the bionic propulsion systems that employ mechanisms obtained from fish swimming have been increasingly applied in the propulsion of underwater vehicles [1–4]. The swimming categories of fish are typically named after the body and/or caudal fin (BCF) swimming modes and median and/or paired fin (MPF) swimming modes [5, 6]. Due to the relatively simple design and high propulsive performance, oscillating caudal fin has become the most popular biomimetic propulsion system [7–11]. However in nature, pectoral fins are used mostly in thrust production, maneuvering, reversing, and rapid stop, and unsuspected diversity has been revealed in the musculoskeletal morphology of the pectoral fin structure and this diversity has clear functional implications [12]. The rapid development of underwater vehicle industry has motivated some investigations on the propulsion mechanism of the pectoral fin. For instance, the flapping foil, as a simplified model of the pectoral fin, has been the focus of considerable theoretical, experimental, and numerical works. It has been shown that, for a two-dimensional foil, optimal thrust condition coincides with the formation of a well-organized inverse Kármán vortex street [13–17]. The studies on finite-aspect-ratio flapping foils [18–22] indicate that the wake of a three-dimensional foil is dominated by two sets of vortex rings that convect downstream at oblique angles to the wake centerline.
The above studies on flapping foils clearly show that observations drawn from two-dimensional foils do not simply carry over to lower aspect-ratio three-dimensional ones. As we know, the shape of a pectoral fin is more complex than that of a flapping foil, and according to the studies by Gibb et al. [23] and Westneat and Walker [24], pectoral fins usually perform a compound rotational motion. Since a pectoral fin represents a significantly more complex situation, it is expected that the wake evolution can change dramatically between an ideal flapping foil and a more realistic pectoral fin. The experimental investigations of Lauder et al. [25, 26] on the bluegill sunfish pectoral fin have shown the presence of a distinct leading-edge vortex on the fin dorsal edge during abduction, and their particle image velocimetry measurements at selected planes also reveal that the fin wake has a highly complex structure. The study of Ramamurti et al. [27] on the digitized pectoral fin of bird wrasse also indicates the emergence of a large leading-edge vortex and the shedding of a pair of counterrotating vortices at the end of the upstroke. The simulations on the bluegill sunfish pectoral fin [28, 29] have found that a number of distinct vortex structures are produced by the fin stroke, and they are subject to mutual induction effects, leading to deformation of the vortex filaments and creation of the highly complex conglomeration of vortices.
On the other hand, in order to develop a biomimetic propulsive system which provides performance comparable to the fish fin, some investigations on pectoral fins have been devoted to the hydrodynamic performance of the labriform swimming mode [2, 30–32] and the development of the mechanical pectoral fin. These studies are focused on the pectoral fin with a combination of locomotion formulated by several rotational motions in a constant velocity free stream. The fin performance within certain range of parameters is tested experimentally or computed numerically and the effect of kinematic parameters on hydrodynamic characteristics is discussed.
At present, the number of studies which have systematically investigated the three-dimensional vortex wake dynamics and tried to carefully explore the propulsion mechanism underlying the generation of forces during the pectoral fin stroke is limited. Since kinematic studies indicate that pectoral fin motions depend on the fish and its travel speed [33, 34], the wake structure and propulsion mechanism of pectoral fins employed by different kinds of fishes probably vary significantly. Are the fluid dynamics obtained from the aforementioned pectoral fins of bluegill sunfish or bird wrasse applicable to the pectoral fins of other kinds of fishes? Comparative analysis across species of fishes is the next logical step towards understanding the generality of research results of fish pectoral fins and this will require studying hydrodynamic characteristics and wake evolution of pectoral fins used by many kinds of fishes, elaborately. The current work is centered on the fluid dynamics of the pectoral fin of black bass (Micropterus salmoides) and provides new insight in this particular field. Black bass is selected for the present study because its pectoral fin is large enough to observe by a video camera and this kind of fish is easy to breed in an aquarium. Furthermore, the bass has been described as a generalist blessed with the functions of accelerators, cruisers, and maneuvers in the study of Webb [35, 36], and it is one of the typical fish using the labriform swimming mode which is the focus of this paper and a number of previous researches [2, 30, 37].
For simulating flows with complex moving boundaries, the choice of an accurate, efficient numerical method with high fidelity is quite significant. The IB method [38, 39] discretizes the flow governing equations on a fixed Eulerian grid, while the IB is explicitly tracked as a curved surface in a Lagrangian approach. This method retains the advantage of explicit interface information from the Lagrangian method and the simplicity of the fluid motion description in Cartesian coordinates and makes grid generation greatly simplified. Thus, the IB method has obvious advantages over the conventional body-fitted method in biomimetic flow simulations with highly mobile boundaries and complex geometries, as shown by the study of Kazakidi et al. [40] where a result-based comparison between IB and body-fitted methods is provided. In our work, based on the basic idea of the IB method and advanced CFD techniques, an improved ghost-cell IB method is proposed. Our current numerical method possesses the following characteristics:(1)Although various reconstruction algorithms have been proposed to enforce the desired boundary conditions, their robustness and efficiency have not been adequately addressed in most IB studies. Herein a simple, fast, and stable reconstruction scheme is developed and it is robust enough to deal with a variety of situations.(2)For moving boundary problems, one case encountered is that a node which is in the solid emerges into the fluid with the time advancement due to boundary motion. We have proposed a new numerical strategy based on the modified Shepard interpolation of Thacker et al. [41] to treat this.
In this paper, we set up a fin computational model approximating the pectoral fin kinematics of black bass. Based on the developed IB method, firstly we focus on exploring the propulsion mechanism of the complex-shaped low-aspect-ratio fin. The connection between the fin kinematics, wake dynamics, and force production is analyzed in great detail so as to determine the complex flow physics that is key for the propulsive performance. Furthermore, the wake flow field of the fin is also compared with available results on three-dimensional flapping foils to reveal the wake feature of the fin with superior hydrodynamic performance. Following this, since the fin performance strongly depends on kinematic parameters, our interest is centered on the effects of phase angle and Strouhal number on the hydrodynamic performance and wake structure. This allows us to address the practical question of how the fin performance varies with key parameters and to understand the fundamental mechanism of the variation. Finally, based on our simulations, we comment on the overall fluid dynamics of the fin and corresponding implications for engineering designs of bionic pectoral fin propulsors.
2. Computational Model
Based on the observation for the pectoral fin of black bass [32, 37], the fin model employed in the current computations is depicted in Figure 1 with a maximum chord length of 0.155 m, which is chosen as the length scale for the current flow, and an aspect ratio of about 1.21. The fin model is 6.2 times larger than the pectoral fin of the observed black bass [32, 37]. Kinematically, labriform swimming is a combination of an anteroposterior (forth-and-back) rowing motion, a dorsoventral (up-and-down) flapping motion, and a feathering motion that denotes a twisting motion of the fin pitching [42]. Although, in certain situations, a fish may exhibit pure rowing motion (which is also called drag-based swimming) or pure flapping motion (lift-based swimming) [43, 44], in general fish rarely perform pure flapping or rowing movements, and instead, in most cases a coupled rotational motion is used [23, 24]. Blake [45] has shown that fish which are propelled by drag-based labriform swimming mode are more probably to have triangular pectoral fins than rectangular or square ones. Theoretical models have predicted that fish which employ drag-based labriform propulsion should have wedge-shaped pectoral fins with relatively blunt distal edges [46]. Referring to above two studies on the effect of pectoral fin shape, black bass can be considered as one of the typical fish using drag-based labriform propulsion in terms of its pectoral fin shape (Figure 1), and the experimental analysis of the pectoral fin of black bass by Kato and Furushima [47] has also demonstrated this. Therefore, we focus on drag-based labriform swimming mode of the black bass pectoral fin, which consists here of a rowing motion and a feathering motion.
The rowing angle is defined as the rotation angle around the -axis, as shown in Figure 2(a), whereby the coordinate system is transformed into the coordinate system and the uniform inflow velocity , the velocity scale for the current flow, is oriented with the arrow along the negative -direction. The feathering angle is defined as the rotation angle around the -axis, as displayed in Figure 2(b), whereby the coordinate system is transformed into the coordinate system. Numerical simulations are carried out with cosinusoidal kinematic mode being applied for rowing and feathering motions, in accordance with the experimental observation and analysis of the bass pectoral fin [47]. The rowing motion of the fin is defined aswhere is the average angle of the rowing motion, is the rowing amplitude, and is the angular frequency of the fin. If we use to represent the motion frequency of the fin which is chosen as the time scale for the current flow, the relationship between and is .
(a)
(b)
The feathering motion of the pectoral fin is then defined aswhere is the average angle of the feathering motion, is the feathering amplitude, and is the phase angle between rowing and feathering (hereinafter referred to as phase angle).
Although the pectoral fin kinematics have not been measured in detail by relevant equipment (e.g., high-speed cameras), the present study on labriform swimming mode retains the complex shape of the black bass pectoral fin and adopts the formulated kinematical description with rationality and representativeness shown by previous investigations on pectoral fins [32, 37, 47]. Therefore, it is expected that some significative insights would be gained into the wake vortices and hydrodynamics of the fin, thereby guiding us towards biomimetic pectoral fin propulsors with excellent performance. The study on the fin digitized using high-speed, high-resolution videos from a black bass pectoral fin will be presented in a separate paper.
Referring to the coordinate system, in the current work, the thrust and lift coefficients are defined aswhere is the fluid density, is the planform projected area of the fin, is the thrust, and is the lift (the -axis is downward). Correspondingly, the time-averaged thrust and lift coefficients are expressed aswhere is the motion period of the fin and is the number of cycles over which the mean is computed.
The propulsive efficiency of the fin is given bywhere is the average power output during one cycle and is the average power input over one period. We compute them as follows [29]:where is the cycle-averaged thrust, denotes the fin surface, is the local surface force, and is the local velocity of the fin surface.
In the current study, the angular velocity in the coordinate system is expressed using the angular velocity components of the rowing and feathering motions, and , as
The motion of the fin is determined by the resultant angular velocity. Eight phases during one period (every ) are shown in Figure 3 via plots using three views of the fin stroke cycle. The left and middle columns present the side and back view, respectively. The right column shows the three-dimensional view of the fin and an ideal fish body in fixed position. Clearly to see is that (1) this fin performs the labriform swimming that consists here of two rotational motions around moving axes, rather than pure rowing motion or feathering motion, (2) at early stage, the fin flaps backward quickly, and (3) the area normal to the incoming flow is adjusted by the feathering motion. The connection between these fin characteristics, vortex dynamics, and force production will be presented in detail below.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
In terms of the length scale , velocity scale , and time scale for the current flow, two nondimensional parameters for the fin, Strouhal number () and Reynolds number (), are defined as and , respectively.
3. Numerical Methods
3.1. Governing Equations and Numerical Procedure
The fluid motion is governed by three-dimensional unsteady, incompressible Navier-Stokes equations aswhere is the fluid velocity, is the pressure, is the fluid density, and is the dynamic viscosity. The above equations are discretized using a cell-centered arrangement of the primitive variables (). CV and CS denote the control volume and control surface, respectively, and is the unit vector normal to the control surface.
The proposed IB method will be presented in the framework of a finite-volume, fractional-step, Navier-Stokes solver for incompressible flow. Both convective and diffusive terms are advanced implicitly in time using a Crank-Nicolson scheme [48] which eliminates the numerical stability constraint, and all spatial derivatives are discretized using central differences [48]. In the following section we will focus on the treatment of complex moving boundaries.
3.2. Immersed Boundary Treatment
A multidimensional ghost-cell methodology is employed to incorporate the effect of the IB on the flow. Referring to the review paper by Mittal and Iaccarino [39], this method falls into the category of “discrete forcing” IB methods. Most researches on bioinspired flow configurations involve complex geometries and kinematics, so the developed method is aimed at simulating flows around arbitrarily complex moving boundaries. Hence, an unstructured mesh with triangular elements is used to represent the surface of the IB, and the choice also makes interpolation along the normal direction readily computed, which is beneficial to the subsequent reconstruction algorithm of boundary condition. The surface of the pectoral fin represented by surface triangulation is shown in Figure 4(a). The unstructured surface grid of the fin is then “immersed” into the Cartesian volume mesh of the flow field, as depicted in Figure 4(b).
(a)
(b)
3.2.1. Reconstruction Scheme of Boundary Condition
This numerical method proceeds by identifying so-called “ghost-cells” of which nodes are inside the solid but have at least one neighbor in the fluid. We express the local flow variable in terms of a polynomial and employ it to derive the value at the ghost cell. Although polynomials of higher degree are expected to be more accurate, the use of them leads to the easier occurrence of numerical instability and boundedness problems [49]. Higher-order approximations also involve an excessively large interpolation stencil that increases the complexity of numerical algorithms [48]. Additionally, in the light of previous studies [50, 51], employing an approximation of one order lower accuracy at the boundary does not reduce the overall accuracy of the whole numerical method. Moreover in the IB methods of Fadlun [52], Yang and Balaras [53], and Liao et al. [54], the linear flow reconstruction has been used and therein a number of validation tests are provided. In order to minimize the probability of numerical instability and save computational time, a linear reconstruction scheme is adopted.
As shown in Figure 5, a normal line segment is extended from the node of ghost cells into the fluid to an image point such that the interface intersection point is midway between the ghost node and the image point. The point is the point at which the boundary condition is to be satisfied. The local flow variable around the point is approximated by the following interpolation polynomial:where , denote the polynomial coefficients. The interpolated variable value at the point has the following form:where is one of the four data points and is the corresponding weight associated with the interpolation polynomial. The interface intersection point and three of the eight nodes surrounding the image point are used as the interpolation data points to determine these weights.
Following this, along the normal line segment, by a linear interpolation and a central-difference approximation which incorporate the Dirichlet condition (for the velocity) and the Neumann condition (for the pressure), respectively, the following formulae are obtained:where is the length of the normal line segment. Equations (10)-(11) complete the numerical description for the ghost cell and they are solved in a fully coupled manner with the governing equations (8) for the neighboring fluid cells.
Having accomplished the simple, fast, and stable reconstruction procedure, we now move on to discuss the capability for our IB method to deal with special situations probably encountered. For instance, a case may be encountered where one of the eight nodes surrounding the image point is the node of the ghost cell itself (see the left ghost cell in Figure 5). This case does not cause any additional problems since three of other seven nodes along with the interface intersection point can be chosen as the interpolation stencil.
It may also be the situation that the eight nodes surrounding the image point for a given ghost cell contain nodes of other ghost cells (see the right ghost cell in Figure 5). In this situation, if the number of the fluid nodes among the eight nodes surrounding the image point is greater than or equal to three, three of the fluid nodes along with the interface intersection point can be chosen as the interpolation stencil. If not, the interpolation stencil consists of the fluid nodes, the interface intersection point, and other ghost nodes. Although this does imply that some of the ghost node values are coupled to each other, it does not lead to any consistency issues since the equation for the ghost cell is solved in a fully coupled manner with the governing equations for the surrounding fluid cells along with the equations for other ghost cells. Based on the above analysis, the proposed reconstruction scheme is robust enough to tackle a variety of situations.
3.2.2. Boundary Motion
In the moving boundary case, we need to move the IB from its current location to the new location at every time step. This is achieved by moving the nodes of the surface triangles with a known velocity. Thus we use the following formula to update the position of a surface element vertex:where is the position vector of the vertex and is the vertex velocity.
One problem associated with moving boundaries is the so-called “fresh cell” issue [53, 55]. This refers to the case where a cell that is in the solid emerges into the fluid at the next time step (or, vice versa) due to boundary motion. Plotted in Figure 6 is the emergence of two fresh cells caused by boundary motion from time level to . When the time advancement is carried out for a fresh cell, some of the required values or derivatives from the previous time step are not physical because of the fact that the boundary changes locations. For the explicit or semi-implicit time-stepping schemes, due to the CFL restriction, the boundary cannot move by more than one cell deep in each time step. In this situation, Yang and Balaras [53] have employed a field-extension methodology, and a procedure of constructing the interpolation stencil by use of the normal probe has been proposed by Mittal et al. [55]. However, when the implicit time-stepping scheme is employed or the IB is moving with large amplitude, the layer of fresh cells may be more than one grid cell. For biomimetic flow simulations with complex moving boundaries, these cases may be involved, and therefore we need to develop a methodology suitable for implicit large time-stepping schemes and large amplitude movements in terms of the “fresh-cell” problem. In the current solver, based on the modified Shepard interpolation [41], a new numerical strategy is proposed.
Once the flow computation for the current time step is performed, the following numerical treatment is carried out to facilitate the time advancement for fresh cells at the next time step. As shown in Figure 7, firstly, two layers of data points are employed to construct the Shepard interpolator. The first layer consists of interface intersection points on the IB where the velocities are determined by the prescribed labriform swimming. The second layer is made up of ghost cells inside the body of which the velocities are computed by the flow reconstruction scheme as presented in Section 3.2.1. Then we check the cells inside the body and identify those which would become fresh cells at the next time step. Finally, the constructed Shepard interpolator is used repeatedly for all these cells identified to obtain the required values or derivatives to resolve the time-advancement difficulty for fresh cells at the next time step. For the mathematical details of the modified Shepard interpolation, the reader is referred to the work of Thacker et al. [41]. Since the Shepard algorithm is independent of the mesh topology and very fit for scattered data interpolation, the proposed numerical strategy is robust enough to deal with a variety of situations, even if adopting large time step size and the IB moving with large amplitude.
4. Results and Discussion
4.1. Validation and Assessment of the Numerical Method
4.1.1. Grid and Domain Independence Studies
Firstly, we present the nominal grid adopted in the current simulations. 17996 elements are used for the surface triangulation of the fin geometry and 8.83 million elements are employed for the background fixed Cartesian grid. The dimensions of the domain are × × . The element spacing sizes at the inner region are , , and , respectively, and the region embodies a × × (5.87 million elements) grid. The grid used at the inner region is uniform, and beyond the region, the mesh is rapidly stretched in the - and -directions. In the -direction, the grid stretching is rapid in the upstream region of the fin. However in the near wake region, the stretching ratio of the mesh is kept below 5% so as to keep relatively high streamwise resolution.
We have carried out comprehensive studies to assess the sensitivity of the numerical solution to the element spacing and domain size. Grid refinement studies are conducted by setting the element number as one and half times the original one in all three coordinate directions in the inner region, with the overall background element number of about 26.5 million and the doubled fin surface grid. Domain independence studies are conducted by doubling the domain size in all three coordinate directions with the whole element number of about 17.1 million.
Simulation results of the temporal variation of the thrust coefficient during the last three periods () are shown in Figure 8(a), where a stable periodic state has been achieved and we observe a close result among the nominal grid, the finer grid, and the domain enlargement grid. Thus the nominal grid is employed for the whole following simulations of flow past the fin to ensure accuracy and minimize the computational time and required memory.
(a)
(b)
(c)
(d)
4.1.2. Accuracy Tests and Verification of Capability to Deal with “Fresh-Cell” Problem on Well-Known Examples
We use the developed IB method to simulate flow past a sphere for . This is a well-known example that allows us to test the accuracy of the IB approach for three-dimensional flows. We adopt two grid densities: a nominal grid with 1003 grid elements and a finer grid with 2003. The drag coefficient for different Reynolds numbers is calculated, as presented in Figure 8(b). We compare the obtained with the studies of Fornberg [56] and Fadlun et al. [52] and find that the sensitivity of the numerical solution to the grid has been removed and the numerical results agree very well with these previous investigations.
The second validation study is conducted on a transversely oscillating cylinder in a free stream. The choice of kinematic parameters in the canonical example reproduces the conditions in the numerical simulations by Yang and Balaras [53], and the specified motion of the cylinder is formulated bywhere is the cross-stream location of the cylinder, is the oscillation amplitude, and is the excitation frequency. The oscillation period is divided into 350 time steps. Under the same flow conditions we have quantitatively compared the time-averaged drag coefficient with the corresponding values from the study by Yang and Balaras [53], as displayed in Figure 8(c), for (—natural shedding frequency). It can be seen that the present IB calculations match the corresponding results from [53], thereby providing further validation in the accuracy. To verify the capability to cope with “fresh-cell” problem especially in the case of the layer of fresh cells being more than one gird cell, we divide the oscillation period into 150 time steps. At this time the layer of fresh cells can be about two cells deep in each time step due to the oscillation amplitude of (—the diameter of the cylinder) and element spacing around the cylinder of . The simulation results obtained by 150 time steps/period are also shown in Figure 8(c), where our numerical results are still in good agreement with Yang and Balaras’s data.
In order to check the ability to deal with the IB moving with large amplitude, the third validation test is carried out on a flapping foil with large heave amplitude. The foil employed in the experiment of Techet et al. [57] features a NACA0012 cross section, with an average chord () of 5.5 cm. We compare our simulation results with the experimental data of the foil with dimensionless heave amplitude of 1.5 which is three times the value used in a number of previous studies [21, 58, 59] on flapping foils. The grid size used in the present simulation has a minimum spacing around the foil of . The flapping period is divided into 500 time steps. At this time the layer of fresh cells can be up to eight cells deep in each time step based on these above choices. Figure 8(d) shows the comparison between numerical and experimental results of the time-averaged thrust coefficient under the same flow conditions. The simulation results of are found to match well with those in the study of Techet et al. [57].
4.1.3. Validation Study on the Pectoral Fin with the Same Geometry
In order to further validate the IB method used in this study, we compare the computed thrust coefficient during one period with published experimental and numerical studies [2, 37, 60], as shown in Figure 9, where the flow has attained a steady state through seven cycles in the simulation. Since the experimental and numerical studies used for comparison have used the same fin geometry as that in the current work, the validation study is quite convictive. Due to the fact that the fluid viscosity and the thickness of the pectoral fin are considered, the IB method and CFD simulation of Wang et al. [60] can calculate the hydrodynamic performance of the pectoral fin in the real environment more accurately than the potential theory method of Su et al. [2]. Moreover in comparison with numerical results of Wang et al. [60], a smaller fluctuation can be observed in our computed results, and the present simulation matches well with experimental results of Kato [37].
4.2. Wake Structure and Propulsion Mechanism
The objective in this section is to provide a comprehensive description of the wake structure and to discuss the propulsion mechanism of the pectoral fin. Vortices in three-dimensional simulations are identified by the criterion [61]. A positive value of is a measure for any excess of rotation rate with respect to the strain. Therefore, the flow exhibits a swirling motion within a region where > 0 as shown by Chakraborty et al. [62]. Earlier researches on flapping foils [63, 64] have indicated the significance of the leading-edge vortices in the fin performance. The kinematics and geometry of the pectoral fin are, however, more complex than those of flapping foils, and the flow around the fin should be dominated by more than one distinct and strong vortex structure. Therefore, we need a detailed analysis to obtain insight into the propulsion mechanism underlying the force generation. For this analysis, we focus on the case with , , , , , , and , and similar values have been used in previous studies on the hydrodynamic characteristics of pectoral fins [30, 65]. All the results presented here and in the following sections have been obtained by simulating the fluid flow over seven motion periods of the fin. When mean quantities are computed, we have discarded the first four cycles and all instantaneous quantities plots correspond to the fifth period, by which time the flow field has achieved a time-periodicity state since the variation in mean hydrodynamic force coefficients in the following strokes is less than 1%.
Figure 10 shows the temporal variation of the thrust and lift coefficients. The thrust peaks once in each period, and so does the lift. Approximately, the pectoral fin mainly provides the thrust in the first half cycle (power stroke), and the lift in the second half cycle (recovery stroke). It should be noted that the time-averaged value of (7.87) is larger than the time-averaged value of (4.91). For comparison, In Kato’s experiment [37] on a mechanical pectoral fin, the mean thrust coefficient was about 1.41 times the mean lift coefficient under the nominal experimental condition. For a pectoral-fin propulsive system developed by Wang et al. [32], the mean value was from 0.18 to 1.75 times in compared to that of for different combinations of kinematic parameters. In a reality of underwater swimmers, the produced lift force may be counteracted with the paired pectoral fin, other fish fins, or body undulations. In a study on modulating the behavior of the robot Madeleine propelled by multiple appendages, four-flippered and two-flippered gaits were used to achieve the starting, cruising, and stopping by the fore-aft and port-starboard interactions [66]. Sfakiotakis et al. [67], who investigated the swimming of a multiarm robot, found that forward propulsion along an approximately straight line might be generated by synchronized sculling movements of the arms, using various patterns of arm coordination. Following this, in order to explore the propulsion mechanism, we analyze the connection between the fin kinematics, vortex dynamics, and force production, which is based on Figure 11. This figure presents the pectoral fin at eight phases in one period (every ). For each phase, three plots are shown in Figure 11. The left-hand-side plot presents a perspective view of the wake vortex. The middle and right-hand-side plots show contours on the downstream and upstream surface of the fin that correspond to the pressure, respectively. Also plotted are vectors on the fin that represent the magnitude and direction of the local surface force.
(a) Thrust coefficient
(b) Lift coefficient
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
Figure 11(a) presents the vortex structures, pressure distribution, and surface force at that is the beginning of a motion period. The dorsal edge vortex (DEV1), ventral edge vortex (VEV1), and tip vortex (TV1), which were formed during the recovery stroke of the previous cycle, have shed from the fin. Moreover, we identify ring shed from the fin during the power stroke of the previous cycle and shed during the period before last. On the downstream surface, a local high pressure region emerges at the fin tip due to the shedding of the tip vortex. Vectors of the local surface force show that the thrust is mostly produced near the fin tip.
In Figure 11(b), at the phase of , after the vortices formed at the previous motion period have completely shed from the fin, new attached vortices are formed at the edge of the fin due to the rapid backward motion. The new dorsal edge vortex (DEV), ventral edge vortex (VEV), and tip vortex (TV) are connected together to constitute a vortex structure which is reminiscent of a horseshoe-shape vortex, and this leads to a large low pressure region in the outer half of the upstream surface. Thus, there is a large pressure difference between the downstream and upstream surface at the phase. In addition, the fin surface is almost vertical to the inflow so that the surface force is tilted in the thrust direction, as shown by the force vectors. Hence, the thrust almost attains the maximum at this phase with small negative lift (referring to Figure 10).
From the phase of to , with backward motion of the fin, the attached vorticity is continually accumulating. The shedding process of the dorsal edge vortex (DEV) and tip vortex (TV) begins at about the phase of , as shown in Figure 11(c), and rings and are convecting downstream with some dissipation. Different from the phase of , on the upstream surface, the obvious low pressure region is only created near the ventral edge. This is because the shedding process of the dorsal edge vortex (DEV) and tip vortex (TV) has started, while the ventral edge vortex (VEV) remains attached to the ventral edge. However, the fin surface still remains good verticality to the inflow at this phase and the shedding process of DEV and TV has only just started, so considerable thrust is still produced at this phase, while the lift is nearly zero (referring to Figure 10).
The next phase of is shown in Figure 11(d), and the tip vortex (TV), partial dorsal edge vortex (DEV), and partial ventral edge vortex (VEV) have shed from the fin and this leads to a smaller low pressure region than that of the previous phase. Moreover, the fin surface is at an oblique angle to the inflow at this phase, which makes the proportion of the pressure difference force assigned in the thrust direction decrease. These are the reasons that the thrust of is much smaller than that of (referring to Figure 10(a)). Furthermore, ring has evolved into a vortex structure which is reminiscent of a horseshoe-shape vortex due to the viscous dissipation.
In Figure 11(e), at the phase of , the vortex structure that consists of the dorsal edge vortex (DEV), ventral edge vortex (VEV), and tip vortex (TV) has completely shed from the fin, and the new tip vortex (TVW) is being created at the fin tip.
Figure 11(f) corresponds to , and TVW is connecting to the VEV + TV + DEV vortex structure shed at the previous instant to form a vortex ring TVW + VEV + TV + DEV. The new dorsal edge vortex (DEVW) and new ventral edge vortex (VEVW) are produced, and they have opposite direction against that of DEV and VEV. A low pressure region appears near the dorsal edge on the downstream surface under the influence of DEVW. According to the spatial angle of the fin at this phase, the pressure difference is assigned to generate lift force and negative thrust force. Since the included angle between the fin surface and the inflow decreases as compared to the previous phase, the proportion of the pressure difference force assigned in the lift direction becomes larger, leading to an increase in the lift.
The next phase shown in Figure 11(g) is at , and TVW is being stretched due to the axial induced velocity of ring RW. With the forward motion of the fin, the dorsal edge vorticity and ventral edge vorticity are increasing continuously. Ring has become two separate columnar vortices with the viscous dissipation. An obvious low pressure region appears in the vicinity of the dorsal edge on the downstream surface, which is clearly associated with strengthening DEVW, and the local surface force vectors show that there exists large pressure difference force. Due to the spatial angle of the fin at this phase, the pressure difference force is reoriented in the direction of lift, and the lift close to the maximum is generated (referring to Figure 10(b)).
Figure 11(h) corresponds to , and DEVW, VEVW, and TVW are being shed from the fin and these vortices correspond to DEV1, VEV1, and TV1 in Figure 11(a) at the next cycle, respectively. Moreover, under the impact of the axial induced velocity of ring RW that corresponds to ring in Figure 11(a) at the following period, VEVW and TVW are being stretched.
Previous studies [21, 68] have shown that three-dimensional flapping foils create wide divergent wakes characterized by two oblique jets because of mutual vortex induction mechanisms. In contrast to the intense, narrow jets which are characteristic of superior hydrodynamic performance, highly divergent foil wakes imply low levels of thrust and propulsive efficiency of finite aspect-ratio foils. However, above vortex dynamics analysis reveals a three-dimensional dual-ring vortex structure of the pectoral fin where the partial power-stroke vortex ring (/RW) is linked to the recovery-stroke ring (DEV1/DEVW + VEV1/VEVW + TV1/TVW) vertically. Thus, it is interesting to investigate the wake field produced by the fin, and Figure 12(a) shows one isosurface of the streamwise velocity . This jet does not seem to show any significant transverse expansion and keeps rather compact. Figure 12(b) presents contours of the streamwise velocity on one streamwise plane (plane shown in Figure 12(a)). We can see that, in comparison with the wake of a three-dimensional flapping foil, the wake of the fin has only one region (as shown by dashed line) of concentrated streamwise momentum where streamwise velocity in the wake is about 49% larger than that in the free stream. Therefore, the wake of the fin is more compact and this can be considered as a sign of the higher propulsive performance of the fin.
(a)
(b)
4.3. Effect of Key Parameters on Wake Topology and Hydrodynamic Performance
4.3.1. Phase Angle
Previously, the effect of phase angle between heaving/rolling and pitching for a flapping foil has been well examined, and some interesting phenomena were reported. von Ellenrieder et al. [19] have indicated that the shedding timing of the leading-edge and trailing-edge vortices is strongly affected by the phase angle. In the studies of Polidoro [69] and Guglielmini and Blondeaux [17], it was shown that the phase angle with the optimal foil performance is . Read [70] has found that the variation in phase angle can change the width of the foil wake. Since the phase angle is an important parameter associated with the vortex structures and hydrodynamic performance, in this section we address the effect of phase angle on the complex-shaped low-aspect-ratio fin.
Figure 13 shows the vortex structures for different phase angles at the end of the motion cycle, when the fin stroke has completed and a multitude of distinct vortex structures are clearly seen. Herein our simulations cover a wide range of phase angles around the nominal value of . In particular, we consider two cases with phase angles of and and these represent a −50% and +50% variation over the nominal value of the phase angle. Similarities can be observed in the vortex structures for all cases, although the tip vortex (TVW) for is stronger than that for the and cases. Interestingly, the wake structure for the lower-phase-angle case is wider and shorter as compared to other two cases with phase angles of and and this means that a wide divergent jet is created for the lower-phase-angle case, which indicates low levels of thrust and propulsive efficiency as discussed in Section 4.2. It is also noted that for a columnar vortex (HV) can be seen to be formed of helical vortex filaments that are reminiscent of the tip vortex structure of a lifting wing [71], and in fact 21% more mean lift coefficient is produced when the phase angle is increased from to based on the current simulations. In the wake of a bluegill sunfish pectoral fin, the helical structure of the abduction tip vortex was also observed in a study by Dong et al. [29].
(a)
(b)
(c)
The temporal variation of the thrust and lift coefficients for different phase angles is presented in Figure 14. Interestingly, although the and cases produce slightly higher peak thrust, they do so with the 76% and 131% lower trough value in comparison with the nominal case, respectively. Furthermore, a smaller fluctuation in thrust is observed for . These suggest a larger mean thrust for the nominal case. Note also that the time at which the trough is found for is much earlier than for two higher-phase-angle cases. This implies that the fin has undergone a rapid transition from thrust to drag producing behavior for the lower-phase-angle case, leading to the fin being subjected to more drag effect in the stroke, and further causing low level of mean thrust. We now move on to the plot of the lift coefficient. It is of interest to observe an approximate peak value of the lift for two lower-phase-angle cases, and however the case has a remarkably higher peak lift. This is associated with the aforementioned helical vortex filaments similar to the tip vortex structure of a lifting wing in the higher-phase-angle case. In fact 54% more peak lift is reached as the phase angle is increased from to . As compared to two higher-phase-angle cases, the case shows an obviously lower trough value, which implies low level of mean lift. It should be noted that the nominal case of has a flatter peak and this suggests the steady output and long duration of relatively high level of lift force. Therefore, only 21% more mean lift is generated, as the phase angle is raised from to , although 54% more peak lift is obtained. There is also a significant reduction in from 7.87 to 4.58 as is decreased from the nominal value to . From the viewpoint of biomimetic fin propulsor, this kind of relatively flat peak is beneficial to the mechanical design since the fluctuant load which is withstood by the structure reduces.
(a) Thrust coefficient
(b) Lift coefficient
Figure 15 plots the variation of the time-averaged thrust coefficient and propulsive efficiency with the phase angle, where a noticeable result is revealed. The simulations show that the highest thrust case with a thrust coefficient of 4.91, the case, is also the most efficient case with a propulsive efficiency of 44%. As the phase angle is reduced to , there is a significant decrease in thrust coefficient to 2.23 and the propulsive efficiency is greatly reduced to 15%. The lower-phase-angle case produces low levels of thrust and propulsive efficiency, and this is remarkably linked to the wide divergent jet at this phase angle. There is also some reduction in thrust coefficient to 3.57 and propulsive efficiency to 24% as the phase angle is increased to , but this decrease is not as much as that observed at the lower phase angle. Consequently, the present simulations indicate that there is an optimal phase angle of for this complex-shaped pectoral fin.
(a) Time-averaged thrust coefficient
(b) Propulsive efficiency
4.3.2. Strouhal Number
In flapping-foil fluid dynamics and aquatic swimming, the Strouhal number is usually taken as a key parameter governing well-defined series of vortex evolution and propulsive performance [21, 59, 72]. Moreover for many kinds of fishes, the increase in swimming speed with their pectoral fins is achieved primarily by the variation in the motion frequency of the fins [28, 73]. On the other hand, the presence of an optimal range of the Strouhal number has been well established for flapping foil and caudal fin propulsion [68, 74, 75]. Herein we deal with the practical question of how the hydrodynamic performance and wake structure of the fin are expected to vary with change in Strouhal number and try to establish optimal for the current black bass pectoral fin.
In addition to simulations at which have been described so far, two additional cases are simulated with Strouhal numbers of 0.25 and 0.85 around the nominal value of 0.55. The wake vortices at the end of the motion period for different Strouhal numbers are displayed in Figure 16. It can be seen that the vortices become stronger as the Strouhal number is increased, especially the recovery-stroke vortex ring (DEVW + VEVW + TVW) near the fin. This makes sense because an increase in Strouhal number implies an increase in the fin stroke velocity relative to the surrounding fluid velocity, thereby leading to the increase in the vortex strength. It is interesting to note that the power-stroke vortex ring (RW) detaches from the recovery-stroke ring for the lower-Strouhal-number case. Considering as a ratio of the time scale of the convection to the flapping time scale of the fin, a decreasing suggests an increasing flapping time scale and this allows the power-stroke ring to convect downstream farther during the labriform swimming of the fin and detach from the recovery-stroke ring. Another important observation is that, with the increase in Strouhal number, the vortex rings in the wake (RW and previously shed ring) become more oriented in the streamwise direction with decreasing included angle between the -axis and ring axis, which increases the streamwise momentum of the jet.
(a)
(b)
(c)
Figure 17(a) presents the variation in the thrust coefficient with time for different Strouhal numbers. Interestingly, all the plots show the same phase in time; for example, the location of the peak is reached at about 14% of the cycle and the location of the trough is found at about 84% of the stroke. It should be noted that the peak value becomes remarkably larger with increasing Strouhal number. There is a 97% increase in peak thrust, as the Strouhal number is increased from the nominal value to 0.85, in contrast with 66% lower trough value. This implies the increase in the mean thrust with Strouhal number. The variation in the lift coefficient with time for different Strouhal numbers is shown in Figure 17(b). Similar to the thrust behavior, for all cases temporal variation of shows the same phase; for example, the time at which the trough occurs is about 6% of the period and the time at which the peak appears is about 79% of the cycle. It can be seen that the peak lift increases as the Strouhal number is raised. It is however interesting to note that the 108% lower trough value is found, as is increased from 0.55 to 0.85, in spite of a 69% increase in peak lift. This suggests that the increase in the mean lift with Strouhal number is not so obvious as the increase in with .
(a) Thrust coefficient
(b) Lift coefficient
The variation of the time-averaged thrust coefficient and propulsive efficiency with the Strouhal number is shown in Figure 18. A finding is that, for the fin, with the increase in Strouhal number, the thrust increases monotonically. This is associated with the strengthening of the wake vortices and the orientation change of the vortex rings (referring to Figure 16), leading to the increase in the jet streamwise momentum which is directly proportional to the thrust produced by the fin. Such a phenomenon has been well reported for two-dimensional and three-dimensional flapping foils [22, 76–78], and the present study clearly reveals that the complex-shaped fin with labriform swimming mode also shows similar behavior. This is in consistency with the fact that many kinds of fishes increase their swimming speed using larger thrust produced by the increase in the fin motion frequency. However, the propulsive efficiency in Figure 18 attains a peak value at a Strouhal number of about 0.55 first and then decreases after that point. In particular, as is increased from 0.25 to 0.55, there is a significant increase in efficiency from 22% to 44%, which is remarkably linked to the increasing thrust. Also a not precipitous reduction in efficiency to 32% is found, as the Strouhal number is further increased to 0.85, due to the more viscous cancellation of opposite signed vorticity with the emergence of complex interconnections among vortex rings (referring to Figure 16). Thus, the present computations indicate that there is an optimal Strouhal number range for the low-aspect-ratio pectoral fin. The optimal value of about 0.55 is highly close to the 0.54 value predicted in the previous study on the bluegill sunfish pectoral fin by Bozkurttas et al. [28]. Also Hover et al. [74] have found that, depending on the choice of kinematic parameters, the optimal Strouhal numbers of a flapping foil can change from 0.3 to 0.6.
(a) Time-averaged thrust coefficient
(b) Propulsive efficiency
5. Conclusions
Numerical simulations are carried out to study the propulsion mechanism and hydrodynamic performance for the pectoral fin of a black bass during the labriform swimming. Simulations of flow around the complex-shaped low-aspect-ratio fin are achieved by the developed IB method, in which we have proposed an efficient local flow reconstruction scheme with enough robustness for tackling special cases probably encountered for the fin cutting through Cartesian mesh of the flow field. Moreover, a new numerical treatment with excellent adaptability based on modified Shepard interpolation is developed to deal with the role variation of the Eulerian grid points near the IB caused by the fin motion. The good agreement with reference results in the literatures shows that the developed IB method can simulate the unsteady motion of the complicated pectoral fin quite well.
In each cycle the prescribed fin kinematics consist of the power stroke and the recovery stroke, and the computational results indicate that the former is primarily employed to provide the thrust while the latter is primarily employed to provide the lift beneficial to the maneuverability. The vortex dynamics analysis interestingly reveals that the fin wake is dominated by a three-dimensional dual-ring vortex wake structure where the partial power-stroke vortex ring is linked to the recovery-stroke ring vertically. In particular, the power-stroke ring is oriented obviously more downstream and is responsible primarily for the thrust production of the black bass pectoral fin. However, the recovery-stroke vortex ring is oriented remarkably more ventrally and serves to generate lift force. Previous studies [21, 58] on finite-aspect-ratio flapping foils have shown that the foil creates a vortex ring every half-period and these rings are released in alternating order, and finally the divergent foil wake is dominated by two sets of inclined vortex rings. Consequently, the dual-ring vortex wake structure linked vertically with distinct directions of ring axes means that the pectoral fin of black bass has a highly compact wake (referring to Figure 12 and corresponding discussion) and produces a jet oriented more posteriorly and ventrally, leading to the generation of more thrust and lift as compared to an engineered flapping foil. Drucker and Lauder [79] have found that the pectoral fin of bluegill sunfish creates a single vortex ring every fin stroke at low forward speeds and a pair of linked vortex rings (with one ring attached to the fish body and only partially complete) at maximal swimming speeds. The fully three-dimensional, volumetric imaging [80] has revealed that the wake of the shark tail is made up of one set of dual-linked vortex rings formed per half-cycle. Thus despite the absence of the flexibility and adjacent fish body, it seems that the predicted wake structure does share some rudimentary features with real fish fins. For the development of a biologically inspired robotic fin, although we cannot fully replicate the exceedingly complex structure and the whole kinematical minimal features of actual pectoral fins, the above key features of vortex structures should be embodied. Analysis on the connection between the fin kinematics, vortex dynamics, and force production shows that the high thrust generated during the power stroke is driven by several distinct vortex structures attached to the fin edge and the large area normal to the incoming flow. However, the produced lift during the recovery stroke is motivated by the increasingly strong dorsal edge vortex and pressure difference force tilted in the lift direction. Furthermore in contrast to a three-dimensional flapping foil, the study on the wake field indicates that the fin wake remains fairly compact and has only one area of intense streamwise momentum, and this implies high levels of thrust and propulsive efficiency.
We also conduct a parametric study to understand how the propulsive characteristics of the pectoral fin change with key governing parameters and to explore the vortex dynamics mechanism underlying the variation. The simulations show that the case with a phase angle of is of the highest thrust and is the most efficient case. The lower-phase-angle case creates wider and shorter wake that is responsible for low levels of thrust and propulsive efficiency. However, the columnar vortex with helical vortex filaments is formed under the higher-phase-angle case and this supports more produced lift. Therefore, the current investigation suggests that the phase angle should be adjusted close to when a mechanical pectoral fin mimicking the labriform swimming is developed, in order to obtain best propulsive performance during steady forward locomotion, whereas for making use of larger lift to maneuver, the fin motion should be actuated with further increased phase angle.
With the increase in Strouhal number, the strength of the vortex rings in the wake increases and these rings become more oriented in the streamwise direction, leading to increased thrust and further causing rapidly increasing propulsive efficiency. Interestingly, the present results suggest that a Strouhal number of 0.55 which is the nominal value in this study makes the pectoral fin operated under the highest efficiency. After that point a relatively slowly decreasing efficiency is observed, which is associated with the more viscous cancellation due to complex vortex interaction, and the fin still has a considerable efficiency at relatively large Strouhal numbers. To our knowledge, the presence of an optimal interval of has been well documented for flapping foil and caudal fin propulsion, and the simulations show that such a range does exist for the black bass pectoral fin. Thus, the fin should be operated within the optimal range of this parameter when the biomimetic underwater vehicle propelled by it is cruising, so as to maintain optimal propulsive efficiency and increase endurance time and mileage. However for a starting maneuver, the fin can swim at a relatively high Strouhal number to provide larger thrust while keeping considerable propulsive efficiency.
The work contributes to the understanding of the fluid dynamics of biomimetic pectoral fin propulsion. However, the flow associated with actual fish locomotion is more diverse and complex than the flow over the pectoral fin undergoing the labriform swimming considered here. Ideally, simulations of flow around fish body with flexible fins are desirable. This will be a subject of our future study.
Competing Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This paper has been supported by NSFC (National Natural Science Foundation of China) funds (Grant no. 51479039).