Induced airflow in flying insects I. A theoretical model of the induced flow

Flapping birds and insects are often likened to revolving propeller blades or rotors because their wings generate lift by steadily pushing air downward. Two influential aerodynamic models of flight in insects(Ellington, 1984c) and birds(Rayner, 1979) drew much inspiration from the extensive theoretical work on rotor aerodynamics. These models focused primarily on the far-field wake of flapping wings using the conceptual abstraction of a pulsed actuator disk, albeit with appropriate modifications to account for the changes in stroke angles and other correction factors appropriate for birds or insects(Rayner, 1979; Ellington, 1984c). Although successful in capturing some interesting far-field features due to flapping wings (e.g. vortex gaits) and stroke-averaged parameters such as power,efficiency etc., subsequent researchers were unable to synthesize these far-field theories with near-field aerodynamics of flapping wings, primarily due to a lack of understanding of the major components of aerodynamic force generation in high-angle of attack flapping flight.

In the last decade, several groups focused their attention on detailed aerodynamic measurements (Dickinson et al.,1999; Sane and Dickinson,2001; Usherwood and Ellington,2002; Dickson and Dickinson,2004; Maybury and Lehmann,2004) and flow visualization on scaled flappers(VandenBerg and Ellington,1997b; Birch and Dickinson,2001) and flapping insects(Dickinson and Gotz, 1996; Willmott et al., 1997; Srygley and Thomas, 2001; Bomphrey et al.,2005). In addition to these experiments, there are also several new Computational Fluid Dynamic (or CFD) simulations of flow around flapping insect wings (Liu et al.,1998; Ramamurti and Sandberg,2002; Sun and Tang,2002; Wang et al.,2004; Miller and Peskin,2005) and a few analytic models(Minotti, 2002; Zbikowski, 2002). Together,these researches have vastly improved our understanding of the near-field mechanisms of force production by flapping wings. However, unlike earlier descriptions of wakes behind insects and birds, these studies focused almost exclusively on the proximate mechanisms of aerodynamic force generation (for reviews, see Sane, 2003; Wang, 2005).

Due to these advances on the experimental and numerical front, it is now possible to address how far-field flow is related to near-field aerodynamic mechanisms. In addition, because the flow environment greatly influences many biological processes such as olfaction or mass flow, such calculations may be especially useful to researchers who work on the interface of flight and sensory-motor physiology. For instance, a detailed understanding of induced flow can allow us to determine the odor plume dynamics over the antennae during flight and understand how moths and other insects accurately track odors. It may also be important from the physiological perspective to understand, for instance, how convective heat loss is enhanced by induced flow, thus influencing the overall flight energetics in endothermic insects.

The main aim of this two-part paper is to estimate the magnitude of the gross flows around an insect body using the near-field approach and discuss their biological importance. The first part of this study presents a derivation of induced airflow using helicopter (or rotor) theory and a blade element-momentum approach modified for application to hovering insects. This model is used to predict the magnitude of mean self-generated air flow in flying insects. In the second paper (Sane and Jacobson, 2006), we test some of the theoretical predictions via systematic measurements of the magnitude of self-generated airflow along the insect body using hot-wire anemometry and show that in addition to the mean induced flow predicted by the theory outlined here, there are additional higher frequency components due to flapping wing motion. We will discuss the relevance of these mean flows and higher frequency flow fluctuations to animal flight studies.

Many recent studies show that at Reynolds number over ca. 100, the flows and forces due to flapping wings are well approximated by Euler equations (Wang, 2000; Ramamurti and Sandberg, 2002; Sane and Dickinson, 2002; Sun and Tang, 2002). Therefore, the model presented in this paper assumes the fluid to be essentially inviscid and incompressible. The effects of vortices, which have their origins in viscous interactions at solid-fluid boundaries, are incorporated into the model through use of empirical coefficients. For example, the influence of leading edge vorticity and tip vorticity is implicitly included within the main derivation of circulation through the use of lift coefficients measured on the finite wings. Although this procedure considerably simplifies the model in the initial stages, there is evidence from both propeller (Leishman,2000) and insect flight literature(Birch et al., 2004) that tip losses contribute directly and significantly to the overall flow structure. Several approximate empirical formulae have been previously proposed but none offer exact solutions to account for tip losses. The earliest of these formulations, developed by Prandtl (for reviews, see Johnson, 1980; Seddon, 1990; Leishman, 2000; see also Jones, 1990 for a review on finite wing theory) for a two dimensional (2D) far-wake resulting from to a finite number of rotor blades, modeled the effect of tip loss in terms of a decrease in the effective actuator disk radius. It predicted a reduction in effective blade radius as the number of blades in the propeller system decreased. As the number of blades tends to infinity, the propeller blades approximate an actuator disk and the effect of tip losses is minimal. Given the sparse data on sectional lift characteristics of flapping insect wings,Prandtl's approach is arguably the simplest to incorporate within a blade-element model and will be tentatively incorporated in the theory described here (see Discussion). Finally, because the current model is based on rigid propeller theory, the aerodynamic effects of wing rotation at either end of the stroke (Dickinson,1994; Dickinson et al.,1999; Sane and Dickinson,2002) are ignored but can be incorporated in future extensions of this model.

where ϵ is infinitesimally small. This method allows the integral in Eqn 4 to attain a convergent form whose value depends on how Γ varies with r (e.g. Prandtl and Tietjens, 1957; Milne-Thomson,1973, p. 80).

The net distribution ofΓ _f(r̂) along the span of a rotating rectangular finite propeller blade is thus a combination of linear and elliptic distributions.

Here, to accommodate the quasi-steady assumption, we assume that C_L(r̂) is the steady lift coefficient, i.e. it varies only spatially along the wing span but does not explicitly depend on time or flow history.

For infinite wings or wings with high aspect ratio, the distribution of lift coefficient is usually assumed to be constant along the entire length of the wing. For flat or twisted finite rotor blades or low aspect ratio wings however, the lift coefficients may be substantially modified by a local variations in sectional inflow leading to a spanwise variation in effective angles of attack from base to tip (Birch and Dickinson, 2001; Usherwood and Ellington, 2002). For flapping insect wings, we must also include further modifications due to the geometry of the wing and twist on local angles of incidence.

where K₀ is a proportionality constant.

where x̂ is any arbitrary point on or near the moving wing.

Eqn 18 provides a general form for induced flow without specifying its explicit dependence on wing geometry. Because there is considerable interspecific variation in wing geometry in various insects, it is useful to develop a procedure for estimating non-dimensional chord length ĉ(r̂), which can be used to derive analytical estimates of induced flow speeds. This is possible with an appropriate choice of ĉ(r̂) that not only describes the wing shape, but also allows us to avoid the singularity at 1 or 0.

Eqn 22 provides a general relationship induced velocity due to a single flapping wing at a point located at distance x̂ from the base of the wing. The contribution from the second wing must be added to this with the appropriate modified value of x̂ with reference to the other wing. Because p>0 and q≥0.5 for most insect wings, the integral in Eqn 22 takes on finite values for all values of x̂<0 and x̂>1, but in the domain of r̂(1<x̂<0),this integral must be converted to a principal value integral and solved as shown by Eqn 4.

Below, I will describe how Eqn 22 applies to Drosophila and Manduca, for which there is some experimental data on induced flow. To measure the values of c(r), the outlines of a Drosophila wing and a Manducawing were traced on to a graph paper and divided respectively into 28 and 50 strips of equal thickness. The length of each strip c(r) was plotted against the position r of that strip with respect to the base(Fig. 2A,B, filled circles) for a wing of normalized length. The resulting plot depicts how c(r) varies as a function of the non-dimensional position (r).

Along the length of the wing, this integral is difficult to solve analytically because it has a Cauchy-type kernel that goes to infinity for x̂=r̂. A method to solve such integrals is outlined by Prandtl in a paper by Betz(Milne-Thomson, 1973), but not without making the derivations tedious and the conclusions difficult to apply generally.

Eqn 27 offers some simple and testable predictions for induced axial flow along the body of the insect. First, for constant stroke amplitude, the induced axial velocity is directly proportional to wing beat frequency. Second, for constant wing beat frequency, induced axial velocity is directly proportional to stroke amplitude. Third, if both stroke amplitude and wing beat frequency are variable, these two predictions can be combined to obtain a more general prediction that the induced axial velocity is directly proportional to the wing velocity.

where α′(x̂) is the sectional angle of attack and α is the morphological angle of attack.

Because the flapping wing is the only source of vorticity generation in the field, we may treat the rest of the fluid outside this region as irrotational. Assuming the fluid to be inviscid and incompressible, we can apply the Bernoulli theorem for the flow within these regions. Thus, the main results of the classical models based on standard actuator disk theory should hold in the non-vortical regions A and B (Fig. 1).

The standard actuator disk theory is based on the assumption that a rotor(or propeller) may be approximated by a 2D disk, which generates uniform and steady pressure difference above and below the disk, thereby generating sufficient forces for weight offset during hovering. This theory, when applied to insect (Ellington, 1984c)and bird flight (Rayner,1979), retained the assumption of radial uniformity of pressure distribution but replaced the steady pressure difference with pressure pulses generated over each stroke. Although the theory outlined in previous sections depends on the radial increase in velocity for a rotating propeller blade or wing, such non-uniformity is smeared out in the far field due to the presence of the slightest viscosity in actual fluids. Hence, when extending the near field theory to far-field, it is convenient to retain most assumptions of the actuator disk theory. The section below repeats some of the theory from previous publications (Rayner,1979; Ellington,1984c) but is included here for the sake of completeness.

Because the flows are eventually dissipated by viscosity in the case of real fluids, the region of applicability for the formula given by Eqn 39 is somewhat limited.

To illustrate how this theory can be applied to experimental data, this section will focus on data from Drosophila and Manduca, two of the best studied examples for insect flight. These examples also offer us data for a rigorous test of the theory because of the different sizes of these insects as well as their different wing planforms. For each example, I will focus only on those predictions for which experimental data is available.

Recent calculations of circulation derived from Digital Particle Image Velocimetric (DPIV) measurements of the flow around wing sections for model Drosophila wings (Birch et al.,2004) provide a direct experimental data to test the model for induced flow.

Using p=0.75, q=0.5 and B_fly=3.2, Eqn 20 allows us to calculate chord as a function of wing position. Using Eqn 20 in the current model, I estimated the spanwise circulation as a function of non-dimensionalized distance from the base. This estimate is derived from the assumptions of quasi-steady state and a semi-elliptical distribution of forces along the wing span. The normalized functional forms of ĉ(r̂),

A similar approach for deriving spanwise circulation was proposed tentatively by Ellington(1974), who also suggested a linearly varying component of the circulation arising from spanwise changes in velocity. However, the nearly-ellipitic component in Ellington(1974) arises from a variation in chord as a function of span wise position, c(r). In contrast, in the above described framework, the nearly elliptic variation arises from a combination of c(r) and a quasi-steady elliptic distribution of circulation on any finite airfoil. Thus, the essential form of the function -its slow rise to a maximum value followed by a roll-off towards the tip -would hold even for a finite rectangular wing. This is true, for instance, in a rectangular helicopter blade (Conlisk,1997).

Fig. 3B shows a comparison of the circulation calculated using the above model and the circulation measured from sectional flow vector fields in(Birch et al., 2004). The measurements of circulation were performed on a model Drosophila wing flapping in a tank of oil filled with air bubbles as seeding particles and systematically imaging the flow fields around each span wise section using Particle Image Velocimetry (Birch et al.,2004). The functional form of circulation in the model is constructed from physical considerations and scaled using force coefficients,kinematic parameters and wing morphology, whereas the measurement and subsequent estimation of circulation were carried out without a priori knowledge of the forces or force coefficients or any assumptions about the nature of spanwise variation. Thus, in principle, these two estimates are entirely independent. Hence, the similarity in measured and predicted magnitudes validates this model. The location of the peak of Γin real wings should not necessarily match the theoretical model that neglects tip loss.

Using the values of tip velocity from Birch et al.(2004) of 0.26 m s^-1, mean chord (c̄) of 7 cm, and lift coefficient of 1.8 and the analytic function for the wing shape of a fly described above, the spanwise circulation as a function of non-dimensional distance from the base is shown in Fig. 3B. The results of these two estimates are in close agreement with respect to the maximum values of circulation achieved and the rate of loss of circulation from maximum to zero. However, unlike the experimental results the circulation measured in the model Drosophila wing falls to zero at r̂=0.85 as compared to the theoretical circulation that becomes zero only at r̂=1. Most likely, this discrepancy arises due to exclusion of wing tip losses from the theoretical model.

For any finite flapping or translating airfoil, it is necessary to take into account the effect of tip vortices when addressing the overall structure of the flow fields. Physically, ignoring these tip losses amounts to assuming that the aerodynamic effect of wing tip is small compared to the chordwise circulation and that each section of the wing continues to generate lift all along the span. As first pointed out by Prandtl(Prandtl and Tietjens, 1957)this assumption may not be valid at points closer to the tip of a finite wing,where a tip vortex enhances leakage of the fluid around the airfoil tip thus significantly reducing the ability of that region to generate any lift. In the framework described above, the effect of tip vortices is not considered explicitly, but is implicitly incorporated through the measured lift coefficients, which include the effect of tip vorticity on the finite wings and through the semi-elliptic distribution resulting from the loss of lift at wing tips.

Many recent experimental results clearly document the effect of tip vortices on the overall flow field in flapping wings. Willmott et al.(1997) showed there were significant tip vortices during both upstroke and downstroke on the wings of moths actively flapping in a smoke rake, confirming the presence of bound circulation and force generation on both strokes. In CFD simulations on`virtual' hawk moths, Liu et al.(1998; Liu and Kawachi, 1998)confirmed this flow around the tip. While visualizing flow around flapping model wings, Birch and Dickinson(2001) noted a distinct spanwise change in the downwash velocity likely influencing the sectional angles of attack.

Second, because the formula was developed for a rotor with variable number of blades, it is not immediately clear if insect wings should be considered as a two-blade or a single-blade system. For a single (up or down) stroke, the net angular excursion of the two wing-system is similar to that of a single-blade rotor in full revolution i.e. N=1. For a complete wing cycle, however the two-wing system behaves like a two-blade rotor because each wing goes through an angular excursion of (approximately) 2π, i.e. N=2.

For N=1, using representative values of C̄_L=1.8, c̄=0.7 mm, R=2.5 mm, and B_flys_i(0)=0.79 (where B_fly=3.2), the ratio of the effective wing span to actual wing span k is about 0.72. For N=2, the value k is approximately 0.86. Thus, depending on the value of N chosen, the calculations estimate that about 14-28% of the wing span from the wing tip produces no lift, as compared to the experimental values of about 15%. Although Prandtl's basic approach could enable us to approximate tip losses in flapping insect wings, it is important to view the above result with some caution because the assumption of uniform inflow may not be valid in most cases.

In the few cases where flow was visualized and quantified around an insect body, the insects were placed within wind tunnels in the presence of an ambient flow typically equal to or higher than the induced axial velocity to ensure that the seeding particles streamed along the insect body(Dickinson and Gotz, 1996; Bomphrey et al., 2005). The ambient flow in Dickinson and Gotz(1996) is 20 cm s^-1as compared to the axial inflow estimate of 14.8 cm s^-1. In Bomphrey et al. (2005), the ambient flow is 3.5 m s^-1 and 1.2 m s^-1 as compared to the estimated axial inflow range of 0-0.35 m s^-1. Hence, I was unable to use these studies to test the induced flow model described above.

When the wing rotates steadily around an axis, the post-downwash lift coefficients vary between 0 to 1.25(Usherwood and Ellington,2002). Using the representative values of Φ=2π/3 radians, n=25 Hz, c̄=0.02 m., B_moths_i=0.8673 (where B_moth=1.8), the range of mean axial velocity calculated by Eqn 27 falls between 0 and 0.35 m s^-1. This range matches rather well with the observed range for the induced axial inflow in flapping Manduca sexta from the anemometric measurements (fig. 5 in the accompanying article, Sane and Jacobson,2006).

Thus, the theory presented here is capable of approximately predicting the mean induced flow speeds due to flapping wings in both the near- and far-field. It shows that this flow is directly proportional to the stroke frequency and amplitude, provided the wing shape and angle of attack remain constant or change only slowly. These equations may be useful to a broad variety of researchers studying such diverse topics as odor tracking behavior,convective heat loss or gas exchange in flying insects, all of which may be influenced by the structure of flow over their bodies.

List of symbols
 
• A

  area of actuator disk

 
• AR

  aspect ratio

 
• B_i

  chord rescale ratio

 
• c̄

  mean chord length

 
• c(r)

  chord length as a function of radial position

 
• ĉ(r̂)

  non-dimensional chord length

 
• C_L

  lift coefficient

 
• C̄ _L

  mean coefficient of lift (average of sections from base to tip of a finite wing)

 
• C_L(r̂)

  spanwise coefficient of lift

 
• \({\dot{{\phi}}}\)

  angular velocity of the wing

 
• F_L

  lift

 
• g

  gravitational acceleration

 
• k

  ratio of the effective wing span to actual wing span

 
• K₀

  proportionality constant

 
• m

  mass of insect

 
• n

  wing beat frequency

 
• P

  mean pressure above actuator disk

 
• p,q

  powers of the Beta function describing wing shape

 
• P₀

  ambient pressure

 
• r̂

  non-dimensional radial position

 
• r

  position along the wing measured from base to tip

 
• R

  wing length

 
• s_i(0)

  a shape parameter at x̂=0

 
• s_i(x̂)

  a shape parameter as a function of x̂

 
• u_∞(r)

  wing velocity

 
• v_axial

  axial flow velocity

 
• v̄ _far

  mean far field induced flow

 
• v_i

  induced velocity

 
• v_i(x̂)

  induced velocity at x̂

 
• v̄ _i

  mean near-field induced flow

 
• x, x̂

  distance from wing base of an arbitrary point on or near the moving wing

 
• Φ

  stroke amplitude

 
• Γ(r)

  circulation as a function of radial position

 
• Γ0

  spanwise circulatory component of the series expansion of circulation

 
• Γ2,Γ4...

  second, fourth etc. derivatives of the circulation with respect to r̂

 
• Γ_2D(r)

  spanwise circulation of a wing section

 
• Γ_f(r̂)

  spanwise circulation of a finite wing

 
• α

  morphological angle of attack

 
• α′(x̂)

  spanwise effective angle of attack

 
• ρ

  air density

I am grateful to Tom Daniel, Michael Dickinson and Ty Hedrick and Sudip Seal for providing useful comments. This work was supported by an NSF Inter-Disciplinary Informatics grant to Sanjay P. Sane and an ONR-MURI grant to Tom Daniel.