Axial propulsion with flapping and rotating wings, a comparison of potential efficiency

Interest in biological locomotion and what advantages the principles governing it might offer in the design of manmade vehicles prompts one to consider the power requirements of flapping relative to rotary propulsion. The amount of work performed on the fluid surrounding a thrusting surface (wing or blade) is reflected in the kinetic energy of the wake. Consideration of the energy in the wake is sufficient to define absolute minimum limitations on the power requirement to generate a particular thrust. This work applies wake solutions to compare the minimum inviscid propulsive power requirement of wings flapping and in rotation at wing loading conditions reflective of hover through a state of lightly-loaded cruise. It is demonstrated that hovering flapping flight is less efficient than rotary wing propulsion except for the most extreme flap amplitude strokes (  >  160°) if operating at large wake wavelength. In cruise, a larger range of flap amplitude kinematics (  >  140°) can be aerodynamically more energy efficient for wake wavelengths reflective of biological propulsion. These results imply, based on the observed wing kinematics of continuous steady flight, that flapping propulsion in animals is unlikely to be more efficient than rotary propulsion.


Introduction
Animals as varied as birds, bats, insects, turtles and some stingrays propel themselves by flapping. The impressive feats of endurance exhibited by some species, seemingly unmatched by manmade vehicles of similar scale, might suggest flapping requires less power than rotary propulsion. Whether flapping or rotating wings, the role of a propulsor is thrust generation. The associated work of the propulsor is reflected in the energy of fluid accelerated away from the thrusting plane (the plane swept out by the flapping or rotating wings) generally in the opposite direction of the produced thrust. A classic wake model, actuator disk theory (ADT) [1][2][3], reflects a conception of a propulsive wake formed purely of axially directed momentum, but this defies real fluid behavior. Along with the formation of an axial jet, there is a tendency for some fluid to flow upstream around wing tips from the pressure (downstream) surface to the suction (upstream) surface, figure 1(D). This upstream flow constitutes wasted energy described as tip loss [1]. Tip flow and an analogous unsteady flow resulting from the temporal variation of thrust in flapping propulsion, figure 1(E) [4], are unavoidable. Perhaps contrary to our expectations from observing biology, ADT suggests that the flapping strategy is at best equal in required power to rotating wings but a detailed consideration of induced flow losses yields a more nuanced result.
Because the work of a propulsor is reflected in the energy imparted to the wake, the problem of propulsive power can be explored by analyzing the wake without the need to consider most of the details of the wings' design or operation, e.g. airfoil, chord length, stall condition, etc. The structure of the propulsive wake can be modeled as sheets swept out during the motion of the wings [5], figures 1(C)-(F). In a viscous fluid like air or water, the effects of the wing pressure are preserved in the wake as regions of circulating flow which were generated on the surfaces of the wings and were shed into the wake. An analogous inviscid wake, where the details of viscosity are ignored, is an immaterial, impermeable surface of zero thickness which preserves the nature of this circulating flow as a quanti ty called circulation. The circulation in the wake is a history of that present around the wings when any particular region of the wake was generated [5]. Circulation is related to the wings' surface pressure and thrust while its gradients, along the wakes' span and in time, induce flow losses, figure 1(F). The wake moves in accordance with its own induced velocity field and the required propulsive power corresponds to the rate of work done by moving the wake sheets through space [6].
A wake which does not deform while moving aft with an apparent displacement velocity, v , has been described by Betz [7,8], to represent the minimum power solution. Various authors have developed solutions for the necessary circulation distribution to achieve undistorting wakes including for propellers [7,10], cruising birds and bats [6,11,12] and helicopters in forward flight [13]. All of these solutions focus on cruising flight. In this work, the method of Hall et al [6], which provides a numerical procedure for determining circulation distributions that satisfy the Betz criteria [7,8] for arbitrary wakes, was extended to illustrate how minimum propulsive power varies for flapping and rotation at wing loading conditions reflective of hover through cruise.
In the sections below, a brief explanation will be presented of how the approach of Hall et al [6] was applied for this study along with the development of the mathematical models. A key element of the developed analysis is the ability to apply inviscid wake results to lightly and heavily loaded flight conditions. Results will be presented for extreme cases representing hover and lightly loaded cruise with a discussion of what these extremes of performance suggest about the relative propulsive efficiency of flapping and rotary wings. Throughout the paper, attempts will be made to illustrate the changing nature of the mathematical form of the inviscid solutions with varying wing loading condition and to show these solutions converge to ADT as flapping and rotating frequencies become large.

Materials and methods
The inviscid approach used in this study is an analysis of wake geometry and the associated optimal circulation distribution. The operating condition under which these rigid or undeforming wakes were created are examples of pure axial propulsion, meaning any body motion is parallel to the thrust vector, figure 1(C). Pure axial propulsion describes hover as well as cruising propellers and buoyant flappers. Cruising aerial flappers must produce both lift and thrust, so are not purely axial propulsors, but their propulsive power requirement remains bound by the axial solution [6]. The rotary wakes considered in this study were those generated by two-bladed rigid propellers or rotors rotating at constant angular velocity, e.g. θ = 2πf . Wakes were also modeled for pairs of unbending wings flapping about a common root, figure 1(A), in a sinusoidal pattern at a constant frequency, i.e. θ = 2πf (Θ/2) cos (2πft), thus the motion of both propulsive styles can be compared at a common angular frequency. The analysis was performed far downstream of the wings in the 'infinite wake', figures 1(C) and (E), where the selfinduced velocities of the upstream wake have the same magnitude of influence as the downstream portion. Two parameters were varied in the geometry of these undeforming prescribed wakes: (1) the ratio of the wake's axial length for a single flap stroke or revolution relative to the wing span, b, a normalized wavelength of the wake described here as λ, and, for flapping, (2) the flap amplitude, Θ. The axial length of one period of the wake is defined by the apparent axial velocity of the wake divided by the frequency of flapping or rotation, U axial /f , figure 1(C).

Derivation of the governing equation
Exploration of propulsive efficiency requires that a propulsor's power can be described as a function of its thrust as well as other operating parameters reflective of the wake wavelength, λ, and flap amplitude, Θ. Hall et al [6], have derived such a relationship and a summary of the key equations are presented here.
The average thrust generated for two wings or blades according to equation (2.5) of Hall et al [6] is where dA is an infinitesimal element of the wake surface, k ·ê A is the scalar product of the axial unit vector with the area unit vector, and W represents the surface area of the entire wake while the air density, frequency, circulation and thrust are ρ, f , Γ and T respectively. This equation can be derived from the sectional lift equation as described by the Kutta-Joukowski theorem, l = ρvΓ, [14] integrated over the span of the wing and averaged in time. Hall et al [6] alternatively derived it based on Kelvin impulses. The bilateral symmetry of the wake geometries and associated circulation distributions ensure that forces and moment lateral to the axial direction are zero, figure 1(A). The areas in equation (1) can be non-dimensionalized by the span squared, b 2 , using the terms dĀ = dA/b 2 and W = W/b 2 , and a non-dimensional form of circulation can be conveniently defined as According to equation (2.13) of Hall et al [6], the average induced power, P i , is where the quantity w i ·ê A is the scalar product of the induced velocity and wake surface normal unit vector. This equation can be derived from the scalar product of sectional force and the induced velocity, which is again integrated over the span of the wing and averaged in time. Hall et al [6] alternatively derived it directly from the kinetic energy of the fluid. Evaluating equation (3) requires a description for the induced flow. The Biot-Savart law is used to calculate the self-induced flow at any point according to the circulation entrained in the wake [14]. The equation for the induced flow at any point on the wake is where r is a position vector, s is the direction across the area element of the local circulation gradient and the subscript 'o' designates quantities which are functions of the area of the wake surface to be evaluated by the integral. Introducing the non-dimensional terms described above as well as non-dimensionalizing the vectors r and r o with b, the induced velocity can be expressed as Substituting equations (2) and (5) into (3) yields an expression for the power with a non-dimensional integral term, (6) Additionally, the term λ/λ is included at this stage to introduce the wavelength into the equation. Many of the terms in equation (6) can be consolidated into one non-dimensional variable, the wake geometry parameter, and substituting equation (7) into (6), the result is Bioinspir. Biomim. 13 (2018) 036012 The inclusion of the term λ into equation (8) serves to make explicit the axial velocity at the thrusting plane, bf λ, a term which will be exploited below. To further generalize equation (8), non-dimensional coefficients for the thrust and power are defined as (9) Using these coefficients, equation (8) becomes The wake geometry parameter, H, is the primary result of this work. It contains all the effects of the wake geometry and, being only a function of the wake geometry, λ and Θ, this parameter will be used to discriminate the propulsive power requirement between various wakes. The induced velocity is a function of circulation, so the problem of interest is determining what circulation distribution will minimize the induced power for a prescribed undeforming wake geometry. Hall et al [6] provide a numerical formulation which solves for the distribution of circulation over a vortex lattice model of a wake subject to the constraint of producing the desired thrust. Their approach, which is applied here, yields not only the optimal circulation distribution required to calculate H but also the induced power, P i , and displacement velocity, v associated with satisfying the Betz criteria for optimal propulsion [7,8].

Wake displacement velocity and definition of the wake wavelength
As described by the Betz criteria, a minimum inviscid power wake appears to move as a rigid structure with a displacement velocity [7,8], v . The magnitude of the displacement velocity can be described as where the product of the thrust and half the displacement velocity yields equation (8). It should be noted that equation (11) is not an approximation but is a result of the wake solution. While not explicitly stated in Hall et al [6], equation (4.8) can be manipulated and inserted into equation (4.3) to show that this is the case. This displacement velocity corresponds to a wake which is infinite in length, but at the thrusting plane, where the wake is semi-infinite (extending aft from the thrusting plane), approximately half of the displacement velocity is present [8].
The wake axial length, as illustrated in figure 1(C), is defined to be U axial /f , so the normalized wake wavelength is defined as where U is the flight speed. Applying equations (11) and (12) as well as the definition for the thrust coefficient, equation (9), yields the expression where Ū = U/bf . For lightly-loaded propulsion, where wing loading is low and the induced velocity and thrust coefficient can be assumed to be zero, this equation reduces to simply λ =Ū , which is traditionally termed the advance ratio [9]. Exploring the alternative extreme of hover, where the free stream velocity is zero, yields The changing nature of the relationship between λ and C T as a function of the wing loading conditions will be important in interpreting results below. This relationship varies smoothly between λ ∝ C 0.0 T in lightly loaded cruise to λ ∝ C 0.5 T in hover. Terms like reduced frequency and advance ratio are often used to define the operational space of flapping and rotary propulsors. Both are related to the wavelength parameter λ but the geometric definition of the wavelength allows for it to be applied more generally across the space of varying propulsor wing loading.

Reduction to the ADT solution
An important benchmark for comparison with the wake results is ADT. The value of H corresponding to ADT can be found through a limiting process which reduces the wake wavelength, λ, to zero. Equation (8) can be shown to have the same form as the ADT solution as hover is approached. The term bf λ, employing equation (12) with the forward velocity set to zero, can be replaced by 0.5v yielding The use of 0.5v is defining the wake length, which is generated by the induced velocity present at the thrusting plane, i.e. generated by a semi-infinite wake, and thus should be half the displacement velocity of an infinite wake. To achieve the wake of an actuator disk in hover, the propulsive frequency must tend toward infinity (and the wavelength, λ, toward zero) in order to collapse the wake into a disk, and the wake geometry parameter, H, must be evaluated as λ → 0. Additionally, since the power is simply the product of the thrust and half the displacement velocity [6], v can be replaced by the term 2P i /T thus yielding which can be further transformed to where A is the area swept out by the wings and Θ is the stoke amplitude, which is π for rotary propulsion. The ADT power derived from momentum theory is [1,3] Setting equation (17) where the amplitude for rotary flight is π. The induced power coefficient for ADT can be defined by substituting equations (19) into (10), Readers familiar with hovering rotary propulsion may perceive equations (10) and (20) to be of the incorrect form, comparing it to equation (2.16) in Leishman [1] for example. For the loading condition of hover in particular, as demonstrated by equation (14), the wake wavelength is a dependent variable of the thrust coefficient and through substitution, equation (20) transforms into the more expected ADT form:

Grid resolution
The solution approach for the circulation distribution is based on an undeforming vortex lattice method requiring that the wake be discretized into smaller elements, small vortex 'rings' of varying circulation strength (see figure 1) and it is the collective strength of these vortex rings which form the circulation distributions.
For the current study, each wake period was discretized in the spanwise and temporal directions for a sufficient number of periods of wake to yield the behavior of an infinitely long wake. Figure 2 demonstrates the resulting variation with refined resolutions for rotary kinematics. Results are plotted for the term H/λ from equation (10). It is seen that the solutions converge exponentially with increasing spanwise and temporal resolution across the entire range of wavelength, figure 2(A). For the smallest wavelength, the inter-period resolution influences the extent of convergence to the ADT solution (black line), figure 2(B). As wavelength becomes smaller, the influence of one period of the wake on the next becomes greater because wake periods are spaced more closely together. Therefore, more periods of wake are needed to ensure convergence to the infinite wake solution.
For this study, a grid resolution of 36 equally spaced spanwise and temporal elements per period were selected. Each wake also included 20 periods up and down stream of the wake period of interest, 41 wake periods in total. Based on the rate of exponential convergence, it is estimated that there is less than a 0.1% discrepancy with the fully converged solution for all wake geometries due to the intra-period resolution.
There is greater divergence at the extreme of smallest wavelength, e.g. λ = 0.05, and this can be seen clearly relative to the ADT solution. Given the tendency of the wake solutions to converge to the ADT solution in this extreme region, solutions with sufficient inter-period resolution to ensure convergence were not deemed necessary.
Wake solutions were found for wake wavelength λ = 10 n π, with n spanning from −3 to 2, varying in resolution over 50 locations in the domain. Only values of the domain highlighted by the gray shaded region will be featured in later figures. Flap amplitudes are presented from 60° to 180° in 10° increments.
While the value of the wake geometry parameter, H, is derived from determining the optimal wake circulation distributions, an explicit discussion of these distributions is saved for the supplementary material (stacks.iop.org/BB/13/036012/mmedia).

Wake geometry
To begin exploring the results, we consider the influence of wake geometry on H, while ignoring the actual operating condition. The axial length of a wake results from the combination of the propulsor's forward speed and self-induced displacement velocity, therefore, the same wake geometry can represent both hover and various axial speeds so long as the wake wavelength are equal. All the propulsors in figure 3 approach their absolute minimum value of H as λ → 0, reflecting the ADT prediction. The largest amplitude of flapping, Θ = 180°, converges toward the same ADT solution as rotation making these effectively equivalently efficient at small wavelength. Contrary to ADT's predictions, which ignore many of the nuances of the wake geometry, as wavelength grows, some flapping amplitudes require less power relative to rotary propulsion with smaller amplitudes progressively eclipsing the rotary H with increasing wavelength.
Viewing the components of the wake geometry parameter, H, associated with work performed on the wake surface in the axial and swirl directions makes clear that it is primarily swirl energy contributing to divergence from ADT, figure 3. The axial component reflects the energy predicted by ADT with additional self-induced axial flow while the swirl component is all in addition to ADT. The propulsive power is proportional to the product of the circulation (which is proportional to lift) and the component of induced velocity perpendicular to the local wake surface (see equation (7)). For smaller wavelength, wakes are axially compact meaning circulation is oriented on the wake surface such that it does not induce much swirling flow. Induced power grows with increasing wavelength because of increases in both the induced velocity and wake surface projected in the swirl direction, and thus resulting in divergence from ADT. For flapping at a given wavelength, the slopes of the surface normals change throughout the flap stroke with some regions where surface normals are oriented more closely with the flight direction, i.e. in the midst of each flap stroke rather than at stroke reversals where the surface normals have large swirl oriented components. This results in a slower rise in the wake geometry parameter with wavelength relative to rotary kinematics. Note there is never a radial component to the wake surface normals considered in this study which precludes induced velocity from performing work in the radial direction.

Cruising propulsion
While valid across a spectrum of pure axial propulsion, equation (10) can be interpreted differently to explore the two operating extremes, lightly-loaded cruise and hover. The limiting case of cruise with a 'lightly-loaded' propulsor assumes that the induced flow is small enough to ignore in the formation of the wake, i.e. the length of a period of the wake is due purely to the flight speed. To compare lightlyloaded propulsors at cruise, it is desirable to redefine the coefficients C Pi and C T based on the flight speed, U, rather than wing angular velocity. In this case, C PU i = C Pi /λ 3 , C TU = C T /λ 2 C TU = C T /λ 2 and equation (10) becomes C PU i = HC 2 TU . From the behavior of H, figure 3, one would conclude that propulsors should operate at the smallest wavelength possible, but the largest amplitude when flapping, to minimize power. This is consistent with past findings [6,11] but the conclusion ignores that real  propulsors must balance the total power, including viscous and additional unsteady effects as well as other physiological or mechanical constraints.
Charting where existing flyers reside on the solution space of H reveals what trends are relevant to animals and existing vehicles. It has been observed that many birds, bats and even some insects seem to have flapping strokes at cruise which reflect a relatively limited range of Strouhal number, St = (bf /U) sin (Θ/2), and the range St = 0.2-0.6 appears to approach the minimal total power kinematics for many species [15,16,17]. This Strouhal number range is depicted in figure 4 as is a range for cruising propellers marketed for radio controlled aircraft (see supplementary tables 4.1). The propellers generally operate at wavelength smaller than the flappers where overall power is lower than the flappers. At the wavelengths where the flappers reside, a range of large amplitudes require less power than rotary propulsion, Θ = 140-180°. This range of amplitudes is larger than typically reported for steadily cruising animals. Given this, cruising animals would typically expend more effort than would be required of a rotary propulsor for the same wavelength.

Hover
In hover, the axial flight speed is zero and the wake wavelength is defined by the induced displacement velocity only, allowing again for a change in the governing equation. The displacement velocity, and therefore the wake wavelength, λ, are functions of the thrust coefficient, equation (14), thus, equation (10) becomes C Pi = √ 2HC 1.5 T . As induced velocities are wholly responsible for the wake length, a phenomenon appropriately ignored in lightly-loaded cruise must be accounted for in hover, induced swirl effects on thrust production. The induced swirl velocity is oriented in the same direction as the wing rotation and serves to diminish the relative velocity the wings experience. This induced drag reduces the effective thrust produced for a given circulation distribution and compensating for this increases circulation and power. Thus equation (10) becomes where C T eff is a design parameter reflecting the load the propulsor must carry, e.g. fuselage drag, vehicle weight, etc, while C T is reflective of the actual circulation strength of the wake. Defining the design variable C T eff allows for 'apples to apples' comparison of the varying propulsive configurations by matching the propulsive load. The set of curves in figure 5 reflects the term √ 2H(C T /C T eff ) 1.5 and since C T is a dependent variable of λ (equation (14)), the wake wavelength defines the thrust coefficient. The ratio C T /C T eff is always greater than one and serves to increase the power (see supplementary sections 1 for the analytical development of the induced swirl correction). The general trends seen in cruise, figure 4, are preserved in hover, figure 5, with the exception of 180° amplitude flapping. The induced swirl is so low for this amplitude that it more significantly out performs rotary propulsion than in cruise. As with cruise, the largest amplitudes offer the lowest power but again, these magnitudes are not typically

Discussion
Among all the kinematic patterns possible, the selection of sinewave-based flapping kinematics might be considered arbitrary. One alternative wake geometry extreme would have wing half-strokes happen instantaneously with pauses between halfstrokes to maintain particular wavelengths. This approach is modelled exactly by Prandtl's tip loss correction to ADT for rotary and 180° flapping [8]. In these 'tip loss' formulations, swirl losses are completely eliminated but axial loss still grows with wavelength. Alternatively, wings might be modeled to achieve pitch reversals instantaneously and transverse half-wing stokes at constant velocity. For 180° flapping, this is equivalent to having counter-rotating wings, a problem previously explored by Theodorsen [10] and which results in lower swirl loss due to the mutually inhibitory interaction between the axisymmetric wakes of opposite circulation strength. In either extreme, realizing these kinematics is implausible yet trends of the wake geometry parameter with wavelength and amplitude would be similar and bounded by ADT. Many factors are ignored in formulating the inviscid solution and their impact on total aerodynamic power must be considered. Inviscid wake distortion [1,[7][8][9][18][19][20] (see supplementary sections 2) will tend to increase power requirements beyond the nondistorting inviscid optimum, though this unavoidable distortion will more greatly impact intrinsically unsteady flapping. Viscous and inertial effects [16,[21][22][23] will cause the total power optimum solution to diverge from that of the inviscid optimum, For example, a reduction in flap amplitude from the extreme inviscid optimum would likely be accompanied with reductions in viscous drag thus benefiting overall power consumption. Airfoil profile drag generally increases more slowly with higher average lift coefficients, at least short of stall, than it decreases due to the reduction of the wings average velocity, to which drag is proportional to the square. In total, such factors defy easy generalization. In one experimental study exploring hover and the stability of the lift enhancing leading edge vortex, rotary kinematics required lower total aerodynamic power than equivalent flappers across a broad range of size [22]. If extending interest to cruising flapping flight, with simultaneous lift and thrust production, opposite extremes of flap amplitude are preferred to generate either force [6,8]. This means the inviscid solution for continuous steady non-bouyant flapping flight is compromised [6] though an aerodynamically more efficient flapping system could derive from the lower viscous drag of having only one wing pair rather than a wing and propeller. Including more of the physical effects, Hall [12] demonstrated one example of cruising flapping propulsion with a very similar propulsive efficiency to a rotary system, but radically superior efficiency will not derive from continuous steady flapping.
The results above do not apply only to linearize aerodynamic behavior and important unsteady effects for flapping, e.g. dynamic stall, can be accounted for. The circulation distributions which spawn induced flow and determine the wake geometry parameter could result from quasi-steady airfoils operating in a linear regime or vorticity deposited in to the wake from a dynamic stall vortex. Within the identified Strouhal number regime of efficient forward flight flapping and plunging propulsors, it does appear that a particular wake structure emerges [15,17], sometimes termed a 'reverse Von Karman vortex street'. While the structure of these wakes appears much more like a series of discrete vorticies of alternating circulation strength rather than sheets of circulation, an inevitable consequence of wake instability, the patterns of induced flow are reminiscent of these vortex lattice solutions.
This work is a relative comparison of the potential to produce a particular thrust at low power given a prescribed flapping or rotary wake. The results speak to the divergence of real wakes, with tip flow and swirl, from the ADT ideal of a uniform axial jet of momentum. Thus the focus of this study are finite span wings and blades. While 2D aerodynamics studies show great promise for achieving high propulsive efficiencies [13], plunging airfoils lack many of the loss mechanisms important for evaluating the performance of finite span wings.

Conclusion
The results of this study provide a framework for comparing the performance of axial flapping and rotary propulsion seamlessly across wing loading conditions reflective of hover through a state of lightlyloaded cruise. The inviscid wake analysis ignores most details of the wing design to glean the minimum power associated with the generation of a particular thrust. Flapping power decreases with growing flap amplitude while both configurations experience trends of increasing power with wavelength, primarily reflecting the growing influence of swirl effects. Flapping offers options for more efficient operation within the wavelength range exploited by existing propulsors in both hover and cruise, but at flapping amplitudes larger than are typically observed in nature. This implies that for continuous steady propulsion, flapping offers animals minimal to no intrinsic aerodynamic energetic advantage over manmade designs; the source of apparent superior performance must reside in other aspects of the biological system.