SUMMARY
The time courses of wing and body kinematics of three freely hovering droneflies (Eristalis tenax) were measured using 3D highspeed video, and the morphological parameters of the wings and body of the insects were also measured. The measured wing kinematics was used in a Navier–Stokes solver to compute the aerodynamic forces and moments acting on the insects. The time courses of the geometrical angle of attack and the deviation angle of the wing are considerably different from that of fruit flies recently measured using the same approach. The angle of attack is approximately constant in the mid portions of a halfstroke (a downstroke or upstroke) and varies rapidly during the stroke reversal. The deviation angle is relatively small and is higher at the beginning and the end of a halfstroke and lower at the middle of the halfstroke, giving a shallow Ushaped wingtip trajectory. For all three insects considered, the computed vertical force is approximately equal to the insect weight (the difference is less than 6% of the weight) and the computed horizontal force and pitching moment about the center of mass of the insect are approximately zero. The computed results satisfying the equilibrium flight conditions, especially the moment balance condition, validate the computation model. The lift principle is mainly used to produce the weightsupporting vertical force, unlike the fruit flies who use both lift and drag principles to generate the vertical force; the vertical force is mainly due to the delayed stall mechanism. The magnitude of the inertia power is larger than that of the aerodynamic power, and the largest possible effect of elastic storage amounts to a reduction of flight power by around 40%, much larger than in the case of the fruit fly.
INTRODUCTION
To understand the aerodynamics, energetics and control of insect flight, it is necessary to know the time history of the aerodynamic forces and moments produced by the flapping wings. It is difficult, even impossible, to directly measure the forces and moments on the wings of a freely flying insect. Existing means of circumventing this limitation are to measure experimentally or to compute numerically the forces and moments on model insect wings (e.g. Dickinson et al., 1999; Usherwood and Ellington, 2002a; Usherwood and Ellington, 2002b; Sun and Tang, 2002a).
In order to use the experimental and computational methods to obtain the aerodynamic forces and moments and to study insect flight, measurements of wing kinematics and some morphological parameters are required. Other researchers have measured wing kinematics of many insects in free flight, using highspeed cine or video; and also measured morphological data of these insects (Ellington, 1984a; Ellington, 1984b; Dudley and Ellington, 1990; Willmott and Ellington, 1997). But since these reported studies used only one camera, the continuous time variation of wing orientation (geometrical angle of attack, wing rotation rate, etc.) could not be obtained. Recently, the time course of threedimensional (3D) wing motion of freely flying fruit flies was measured using three orthogonally aligned, highspeed cameras (Fry et al., 2005). Measurements of 3D wing motion of other insects are of great interest, but some limitations to Fry et al.'s work meant that morphological parameters such as weight and position of center of mass could not be measured. If these data were also measured, one could use them to test the experimental and computational models (a reasonable test of the experimental and computational models is that the measured or computed vertical force approximately balances the insect weight and, in hovering flight, the horizontal force and the pitching moment about the centre of mass of the insect are approximately zero).
In the present study, we measure the time course of 3D wing motion of hovering droneflies using three orthogonally alined highspeed cameras and also measure the required morphological data, and then employ the method of computational fluid dynamics (CFD) to compute the aerodynamic forces and moments of the wings flapping with the measured kinematics. The reasons we chose droneflies are (1) droneflies might have different wing kinematics from that of fruit flies used in Fry et al.'s study, (2) they are capable of motionless hovering, and (3) they fly even in strongly lit laboratory conditions. Comparison between computed results and the conditions of vertical force being equal to weight, horizontal force being equal to zero, and moment about the center of mass being equal to zero, can provide a test of the computational model. Analyzing the time courses of the wing motion and aerodynamic forces and moments could provide insights into how the aerodynamic forces and moments are produced and how the flight balance is achieved. Using the computed results, the power requirements can be readily estimated and the aerodynamic efficiency of flight examined.
MATERIALS AND METHODS
Animals
We netted droneflies of the species Eristalis tenax Linnaeus in a suburb of Beijing in September 2006. All experiments were conducted in the same day of capture. Only the most vigorous individuals were selected for filming.
3D highspeed film
We filmed the hovering flight of the droneflies in an enclosed flight chamber using three orthogonally aligned synchronized highspeed cameras (MotionXtra HGLE, 5000 frames s^{–1}, shutter speed 50 μs, resolution 512×320 pixels). The flight chamber, a cube of 15×15×15 cm^{3}, was built from translucent glass. We backlit each camera view using densely packed arrays of light emitting diodes (LED) covered with diffusion paper. The reason we used LED arrays as the light source was that they produced much less heat than cine lights at the elevated light levels required for highspeed filming. The filming was started only when the insect hovered in the cubic zone of approximately 5×5×5cm^{3} at the central region of the chamber, which represented the intersecting field of views of the the three cameras.
Measurement of wing and body kinematics
Using a flat calibration panel to calibrate the cameras, one can obtain the transform matrices between the world coordinate system and three image coordinate systems, which are called projection matrices. Based on the basic principles of stereo vision, the projections of the scene point in the world coordinate system onto any image coordinate system can be calculated, provided that the point's coordinate in the world coordinate system and projection matrices are known. The general structure of a trinocular stereo vision system is shown in Fig. 1.
Let (X_{w} Y_{w} Z_{w}) be the world coordinate system; let (X_{C1} Y_{C1} Z_{C1}), (X_{C2} Y_{C2} Z_{C2}) and (X_{C3} Y_{C3} Z_{C3}) be the camera coordinate systems of camera 1, 2 and 3, respectively; let (u_{1} v_{1}), (u_{2} v_{2}) and (u_{3} v_{3}) be the image coordinate systems of cameras 1, 2 and 3, respectively. Let P be an arbitrary point, whose coordinates are X_{w}, Y_{w} and Z_{w} in the world frame, X_{C1}, Y_{C1} and Z_{C1} in the frame of camera 1, X_{C2}, Y_{C2} and Z_{C2} in the frame of camera 2, X_{C3}, Y_{C3} and Z_{C3} in the frame of camera 3. Let p_{1}, p_{2} and p_{3} be the projective points of P on the three cameras; p_{1}'s coordinates in the image frame of camera 1 are u_{1} and v_{1}; p_{2}'s coordinates in the image frame of camera 2 are u_{2} and v_{2}; p_{3}'s coordinates in the image frame of camera 3 are u_{3} and v_{3}. The coordinate transformation between the world coordinate system and the three image coordinate systems are: (1) where M_{1}, M_{2} and M_{3} are projection matrices of camera 1, camera 2 and camera 3, respectively. Let m1_{i,j}, m2_{i,j} and m3_{i,j} (i, j=1,2,3) denote the elements of M_{1}, M_{2} and M_{3}, respectively. Eliminating Z_{C1}, Z_{C2} and Z_{C3} from Eqn 1 gives: (2) (3) (4) The coordinates of point P's projection onto all three image coordinate systems can be determined via Eqn 2, 3, 4.
We take a line segment, whose length is equal to the body length of the insect, to represent the body of the insect, and used the outline of the wing obtained from scanned picture of the wing (see Fig. 2A) to represent the wing. The line segment of the body and the outline of a wing are designated as the models of the body and the wing, respectively. They can be represented as a set of points. Using the method described above, we can easily computed the projection of the point set on all the three image coordinate systems.
We developed a toolbox for Matlab (The MathWorks, Inc., Natick, MA, USA) to extract the position and attitude of the body and of the wings from images captured by all three cameras. We put the models of the body and the wings into the world coordinate system, then changed the positions and attitudes of those models until the best overlap between a model's projection and displayed frame was achieved in all three views. At this point, the positions and attitudes of those models would be taken as the positions and attitudes of the body and wings of the insect. Generally, several readjustments of each model's position and attitude were required to obtain a satisfactory overlap.
Measurement of morphological parameters
The present method of measuring the morphological parameters follows, for the most part, the detailed description of the method given by Ellington (Ellington, 1984a).
The insect was killed with ethyl acetate vapor after filming. The total mass m was measured to an accuracy of ±0.01 mg. The wings were then cut from the body and the mass of the wingless body measured. The wing mass m_{wg} was determined from the difference between the total mass and the mass of wingless body.
Immediately after cutting the wings from the body, the shape of one of them was scanned using a scanner (HP scanjet 4370; resolution 3600×3600 d.p.i.). A sample of the scanned picture of a wing is shown in Fig. 2A. Using the scanned picture, wing length R (the distance between the wing base and the wing tip) (see Ellington, 1984a) and local wing chord length were measured to an accuracy better than ±0.5%. Parameters including wing area, mean chord length, radius of second moment of wing area, etc., were computed using the measured wing shape.
The wingless body (with legs removed) was scanned from two perpendicular directions (the dorsoventral and lateral views; Fig. 3). Following Ellington (Ellington, 1984a), the cross section of the body was taken as ellipse and a uniform density was assumed for the body. With these assumptions and the measured body shape, the center of mass of body could be determined. Body length (l_{b}) and distance between the wing roots (l_{r}) were measured from the dorsal views; distance between the wingbase axis and the center of mass (l_{1}) and distance between the wingbase axis and the long axis of the body (h_{1}) were measured from the lateral view (the pitch moment of inertia about the center of mass could also be determined).
Computation of aerodynamic forces and moments and power requirements
The aerodynamic forces and moments were computed using the CFD method. On the basis of studies on wing–wing interactions (Lehmann et al., 2005; Sun and Yu, 2006), aerodynamic interactions between the left and right wings could be neglected. During hovering flight, the body did not move, and it was assumed that the aerodynamic interaction between the body and the wings was negligible. Therefore in the present CFD model, the body was neglected and only the flows around one wing were computed (the aerodynamic force and moment produced by the other wing were derived from the results of the computed wing). The wing planform of the insect was obtained from the present measured data. The wing section was assumed to be a flat plate with rounded leading and trailing edges, the thickness of which was 3% of the mean chord length of the wing.
The flow equations and the solution method used were the same as those described in Sun and Tang (Sun and Tang, 2002a). The computational grid had dimensions 100×99×105 in the normal direction, around the wing section and in the spanwise direction, respectively (portions of the grid used for one of the insects in the present study are shown in Fig. 4). The normal grid spacing at the wall was 0.0015c (where c is the mean chord length of wing). The outer boundary was set at 20 chord lengths from the wing. The nondimensional time step was 0.02 (nondimensionalized by c/U, where U is the mean speed of wing at the radius of second moment of wing area; there are about 400 time steps in a stroke cycle). A detailed study of the numerical variables such as grid size, domain size, time step, etc., was conducted and it was shown that the above values for the numerical variables were appropriate for the calculations.
RESULTS
Three droneflies hovering at the central region of the flight chamber (the zone of the interesting field of view of the three cameras) were filmed. They were denoted as DF1, DF2 and DF3, respectively. For each of the droneflies, film of around 12 wing strokes were digitized. Samples of the original video sequences for DF1 are presented as Movie 1 in supplementary material.
Morphological parameters
Morphological parameters of insects filmed in free hovering are given in Table 1. Parameters in the table include the total mass of an insect (m), the mass of a wing (m_{wg}), the wing length (R), mean chord length of wing (c), area of a wing (S), radius of second moment of wing area (r_{2}), body length (l_{b}) and distance between the two wing roots (l_{r}), distance between the wingbase axis and the center of mass (l_{1}), distance between the wingbase axis and the long axis of the body (h_{1}), and the pitch moment of inertia about the center of mass (I_{y}).
Wing kinematics
We determine the stroke plane in the same way as that of Ellington (Ellington, 1984b). The stroke angle (ϕ), the stroke deviation angle (θ) and the angle of attack (α) of a wing are defined as in Fig. 5. The measured data of these angles as functions of time for the left and the right wings of DF1 are shown in Fig. 6 (for the three insects in hovering, about seven wellrepeated wing strokes in which the left and right wings moving symmetrically are captured). As seen in Fig. 6, the motion of the right wing is approximately the same as that of the left wing, as expected for hovering flight. The stroke positional angle varies with time approximately as a sinusoidal function. The angle of attack does not vary significantly in the midposition of the down or upstroke, but varies sharply during the stroke reversal. The stroke deviation angle is relatively small; it is higher at the begining and the end of a downstroke or upstroke, and lower at the middle of the down or upstroke.
Data of seven cycles for DF1 are superposed and plotted in Fig. 7. Curves of simple functions that approximately fit the data are also plotted in the figure. For a clear description of the data, we express time during a cycle as a nondimensional parameter, t̂, such that t̂=0 at the start of a downstroke and t̂=1 at the end of the subsequent upstroke. Data of the positional angle (Fig. 7A) can be well approximated by a simple harmonic function: (5) where n is the wingbeat frequency, the mean stroke angle and Φ the stroke amplitude [ and Φ are defined as following: =(ϕ_{max}+ϕ_{min})/2,Φ =ϕ_{max}–ϕ_{min}, where ϕ_{max} and ϕ_{min} are the maximum and minimum values of ϕ, respectively (see Ellington, 1984b)].
The time variation of the angle of attack (Fig. 7B) can be approximated by the following functions. It is a constant at the midposition of the down or upstroke and varies with time at the stroke reversal; the rate of its variation can be fitted by a cosine function. The constant of α at the midposition of the downstroke is denoted as α_{d}, that of the upstroke as 180°–α_{u} (α_{d} andα _{u} are the geometrical angles of attack of the downstroke and upstroke, respectively). The rotation time (nondimensionalized by stroke period) is denoted asΔ t̂_{s} for supination (α change from α_{d} to 180°–α_{u}) andΔ t̂_{p} for pronation (α change from 180°–α_{u} toα _{d}). Wing rotation may not be exactly symmetrical (wing rotation is symmetrical when the rotation occurring is symmetrical with respect to stroke reversal, i.e. half the rotation is conducted in the later part of a halfstroke and the other half in the early part of the following halfstroke); we useΔ t̂_{s,t} to denote the time (nondimensionalized by stroke period) the supination is advanced (a negative value means the rotation is delayed), andΔ t̂_{p,t} to denote the nondimensional time the pronation is advanced. The function representing the time variation of α during supination is: (6a) where a is a constant: (6b) t_{1} is the time when the wingrotation starts: (6c) The expression of the pronation can be written in the same way.
Finally, the time variation of the deviation angle (Fig. 7C) can be approximately represented by the first two terms of a Fourier series: (7)
Results for the other two insects are similar. For all the three droneflies, parameters specifying the wing motion are determined from the data. They include: and Φ, specifying the positional angle; α_{u}, α_{d} andΔ t̂_{s},Δ t̂_{p},Δ t̂_{s,t} andΔ t̂_{p,t}, specifying the angle of attack; θ_{0}, θ_{1},θ _{2}, η_{1} and η_{2}, specifying the deviation angle. The results are given in Table 2. The wing beat frequency (n), stroke plane angle (β) and body angle (χ) are also given in Table 2 [the stroke plane angle is the angle between stroke plane and the horizontal, and body angle is the angle between the long axis of the body and the horizontal (see Ellington, 1984b)].
The filmed insects are only in approximate hovering flight; i.e. some of them move at very small velocity, and the nondimensional velocity of motion, denoted by advance ratio (J) was measured (J is the velocity of motion divided by the mean wingtip speed 2ΦnR). The values of J are also included in Table 2 and are very small for the three insects.
The computed aerodynamic forces and flows
With measured wing kinematics and using the method described in Materials and methods, aerodynamic forces and moments produced by the flapping wings were computed. Since the wing motions of the three insects are similar, only the computed results for one insect, DF1, are discussed in detail here.
Let V and H be the computed vertical and horizontal forces of a wing and M the pitching moment about the center of mass. Let L and D be the lift and drag of a wing, respectively (wing lift is the force component perpendicular to the stroke plane and wing drag is the force component in the stroke plane and perpendicular to the wing span). The force and moment coefficients are defined as: (8a) (8b) (8c) (8d) (8e) where U is the mean velocity of the wing at r_{2} (U=Φnr_{2}) and ρ the fluid density. Fig. 8 gives the time courses of C_{V}, C_{H}, C_{M}, C_{L} and C_{D} in one cycle.
The downstroke produces a little more wing lift and drag than the upstroke (see Fig. 8D,E), because the angle of attack of the downstroke is a little larger than that of the upstroke (α_{d}=33.8° and α_{u}=33.0°). The vertical force is approximately the same as the wing lift (comparing Fig. 8A with D), because the stroke plane is almost horizontal (β is only 4°). From the data in Fig. 8A, the azimuthal position of the line of action of the vertical force of a downstroke or an upstroke can be calculated. The calculation shows that the line of action is near the position of the mean stroke angle, about 3° before and 5° after the position of the mean stroke angle for the downstroke and the upstroke, respectively.
As seen from the time courses of force coefficients (Fig. 8), the major part of the mean vertical force, or the mean lift, come from the midportions of the down and upstrokes (the large force peaks at t̂=0.12–0.38 and t̂=0.62–0.88) and the contribution from wing rotation around stroke reversal is relatively small (the relatively small force peaks at t̂=0.38–0.5 and t̂=0.88–1). As discussed by Sun and Tang, the large C_{L} in the midportion of a down or upstroke (during which period the wing is in pure translational motion) is due to the delayed stall mechanism, and the relatively small C_{L} peak near the end of the downstroke or upstroke is due to the fast pitchingup rotation of wing (Sun and Tang, 2002a). Using the data in Fig. 8, it is estimated that around 60% of the mean vertical force is contributed by the pure translational motion, or is due to the delayed stall mechanism.
The flowfield data can provide further evidence for the above statement. The contours of the nondimensional spanwise component of vorticity at midspan location are given in Fig. 9. The leadingedge vortex does not shed in an entire down or upstroke, showing that the large C_{L} and C_{D} in the midportion of the halfstroke are due to the delayed stall mechanism.
Flight power
With the flows computed by the CFD method, the mechanical power of a wing can be easily calculated. The aerodynamic power of a wing (P_{a}) is: (9) where M_{a} is the aerodynamic moment around the wing root andΩ the angular velocity vector of the wing. M_{a} is readily calculated using the force distribution obtained from the flow computation and Ω is known from the measured data. The inertial power of the wing (P_{i}) is: (10) where M_{i} is the inertial moment of the wing. The way to determine M_{i} has been described previously (e.g. Sun and Tang, 2002b). The component of M_{i} due to the flip rotation of wing cannot be calculated because the moment of inertia of the wingmass about the axis of flip rotation is not available. Since most of the wingmass is located near the axis of fliprotation, it is expected that this part of M_{i} is much smaller than the other parts, and is therefore neglected in the present study. The moments of inertia for the other two axes are approximately the same (denoted as I_{wg}). , where r_{2,m} is the radius of gyration of I_{wg} (Ellington, 1984a). Ellington's measured data, r_{2,m}/R=0.4, is used here, and with data of m_{wg} in Table 1, values of I_{wg} for DF1, DF2 and DF3 are computed as 1.12×10^{–11} kg m^{2}, 1.30×10^{–11} kg m^{2} and 1.72×10^{–11} kg m^{2}, respectively. With I_{wg} known and wing acceleration computed from Ω, M_{i} can be calculated.
The instantaneous nondimensional power (C_{p}), nondimensional aerodynamic (C_{p,a}) and inertial power (C_{p,i}) of a wing (nondimensionalized by 0.5ρU^{2}Sc) for DF1 are given in Fig. 10 (those for DF2 and DF3 are similar). It is interesting to note that the time course of C_{p} is similar to that of C_{p,i}; this is because the inertial power is large than the aerodynamic power. This means that elastic energy storage could be important for the droneflies.
Integrating C_{p} over the part of a wingbeat cycle where it is positive gives the nondimensional positive work (C^{+}_{w}); C^{+}_{w}=20.95. Integrating C_{p} over the part of a wingbeat cycle where it is negative gives the nondimensional negative work (C^{–}_{w}); C^{–}_{w}=–9.08. The mass specific power (P^{*}) is determined as the mean mechanical power over a wingbeat cycle divided by the mass of the insect and it can be written as: (11) where τ_{c}=t_{c}U/c (t_{c} is the period of wingbeat) (Sun and Tang, 2002b). In the case of zero elastic energy storage, with the added assumption that the cost for dissipating the negative work is negligible, and P^{*}=60 W kg^{–1} for DF1 (95 W kg^{–1} and 60 W kg^{–1} for DF2 and DF3, respectively). In the case of 100% elastic energy storage, and P^{*}=34.4 W kg^{–1} for DF1 (it is 50.2 W kg^{–1} and 38.1 W kg^{–1} for DF2 and DF3, respectively). Thus for DF1, the largest possible effect of elastic energy storage amounts to reducing the power by about 40% (for DF2 and DF3, the value is 47% and 36%, respectively).
DISCUSSION
Wing kinematics and comparisons with available dronefly data and with fruit fly data
From the measured data of wing kinematics (Table 2), we have the following main observations. The wing beat frequency (n) of the three hovering droneflies ranges from 142 Hz to 209 Hz. The stroke amplitude (Φ) is around 110°. The downstroke angle of attack (α_{d}) is around 35° [the upstroke angle of attack (α_{u}) varies with the stroke plane angle; see below]. The stroke plane angle (β) ranges from 4° to 16.4°. When β is small, α_{u} is about the same as α_{d} (e.g. for DF1, is 4°, α_{u} is 33.0°, approximately the same as α_{d}=33.8°); whenβ is relatively large (for DF2, is about 16°), α_{u} is smaller than α_{d} by about 10°. Wing rotation time for spination is a little longer than that of pronation (e.g. for DF1, the former is 32% of the wing stroke period and the latter is 26%). Wing rotation is approximately symmetrical. The mean stroke angle () also varies with β; is small when β is small (for example, DF1 β=4°, =7.1°). The deviation angle is small, around 10°. As shown in Fig. 11A, the variation ofθ gives a very shallow Ushaped wingtip trajectory.
Comparison with available dronefly data
Using one cine camera, Ellington measured the positional angle and deviation angle of many insects, including droneflies (Ellington, 1984b). His data show that the ϕ varies with time approximately as a simple harmonic function, and that θ is small and the wingtip trajectory is a shallow Ushaped curve. Our data of ϕ, θ and the wingtip trajectory agree with Ellington's [comparing Fig. 7A and Fig. 11A of the present paper with fig. 13 of Ellington (Ellington, 1984b)].
Because only one camera was used, Ellington could not measure the angle of attack. The present study has provided the time course of the angle of attack, giving data on angles of attack during downstroke and upstroke translations and rotation duration timing.
Comparison with the results of fruitflies
In a previous study, Fry et al. measured the wing kinematics of fruit flies using the same method as we use here. It is of great interest to compare the kinematics of the droneflies measured here with that of the fruit flies. For comparison, we replotted the data of the fruit flies in Fig. 12 [the data are taken from fig. 2 of Fry et al. (Fry et al., 2005); note that in Fry et al.'s paper, α is zero when the wing surface is perpendicular to the stroke plane, while in the present study α is 90° when wing is in this position, and we have converted the α values of Fry et al. to fit our definition]; data from DF1 are also given in the figure.
For both insects, ϕ varies with time approximately according to a simple harmonic function, except that for the fruit fly the duration of downstroke is slightly longer than that of the upstroke (comparing Fig. 12A with Fig. 12D). However, large differences in α and θ exist between the two insects.
For the fruit fly, after the supination (downstroke/upstroke) rotation and at the beginning of the following upstroke (Fig. 12B, t̂≈0.6), the α value is nearly 170°, i.e. the wing chord is almost parallel to the horizontal plane, and then α decreases to about 130° near the end of the upstroke; α varies throughout the halfstroke. But for the dronefly, after the supination,α approximately keeps a constant value (around 145°) in the upstroke (Fig. 12E).
For the fruit fly, stroke deviation angle is larger, peaktopeak amplitude is around 30°, whilst for the dronefly, this value is around 10°. As a result, the wingtip trajectory of the fruit fly is a deep Ushaped curve and that of the dronefly a very shallow one (compare Fig. 11A, B). Moreover,θ of the fruit fly has a fast and large decrease at the beginning of the upstroke (Fig. 12C, at t̂≈0.6), shifting the minimum point of the Ushaped curve of the upstroke forward (see Fig. 11B); but for dronefly, the minimum point of the Ushaped curve is approximately at the middle (see Fig. 11A).
As discussed below, the differences in α and θ between the two insects have significant effects on how the flight is balanced.
Aerodynamic forces and flight power and comparison with the results of fruitflies
Aerodynamic forces
In their study of hovering fruit fly flight, Fry et al. `replayed' the measured wing kinematics on a robotic model wingpair and measured the aerodynamic forces of the wingpair. In the present study, we used the measured wing kinematics of the droneflies in a CFD model to compute the aerodynamic forces and their corresponding moments. It is of great interest to compare the aerodynamic forces of the fruit flies and that of the droneflies.
When hovering with an approximately horizontal stroke plane, for the dronefly, the aerodynamic forces of the upstroke are a little smaller than that of the downstroke. For example, for DF1 (β=4.0°), the vertical force produced during the upstroke is approximately 11% less than that during the downstroke. By contrast, for fruit fly, the aerodynamic forces of the upstroke are larger than those of downstroke. For comparison, we have replotted the vertical force data of the fruit flies in Fig. 13 [using data taken from fig. 3A of Fry et al. (Fry et al., 2005)]. From the data, we estimate that the vertical force produced during the upstroke is approximately 58% larger than that during the downstroke.
For the dronefly, the vertical force that balances the insect weight is mainly produced in the mid positions of the down and upstrokes (see Fig. 8A), during which the wing is in translational motion with constant angle of attack. Since the stroke plane is approximately horizontal and the variation of the deviation angle is small in the mid portions of the downstroke and upstroke, it can be said that the weightsupporting vertical force is mainly contributed by the lift of the wings. However, for the fruit fly, a large portion of the vertical force is produced in the early upstroke (Fig. 13). In this period of the upstroke, as discussed in Fry et al.'s paper, the wing surface is almost horizontal and the wing moves downward (Fry et al., 2005). Therefore, the vertical force produced in this period is contributed by the drag of the wings. Thus we see that when hovering with an approximately horizontal stroke plane, the dronefly mainly uses the lift principle for the weightsupporting force, whilst the fruit fly uses both lift and drag principles.
Flight power
As shown above, for droneflies, the time course of mechanical power is similar to that of the inertial power (see Fig. 10), because the magnitude of inertial power is larger than the aerodynamic power. However, for the fruitflies, the inertial power is smaller than the aerodynamic power and the time course of mechanical power is similar to that of aerodynamic power [see Fig. 14 in which we have plotted the power data of the fruit flies using data from Fry et al. (Fry et al., 2005)]. As a result, droneflies could benefit from elastic energy storage much more than fruit flies: for the droneflies, the largest possible effect of elastic energy storage amounts a reduction of power by around 40%, whilst for fruit flies, the value is around 16% [calculated using data in table 1 of Fry et al. (Fry et al., 2005): (115–97)/115=0.16].
Let us explain why the inertial power is relatively large for the droneflies, but aerodynamic power is relatively large for the fruit flies. The inertial torque of a wing is proportional to: (12) where r̂_{2,m} is the radius of the second moment of wing mass, normalized by R. The aerodynamic torque is proportional to: (13) where C̄_{d} is the mean drag coefficient of wing and r̂_{d} the radius of the first moment of wing drag, normalized by R. Then, the ratio of the inertial and aerodynamic torques is proportional to: (14a) It is assumed that r̂_{d} is approximately the same for the dronefly and fruit fly wings; measured data and computations show that droneflies and fruit flies do not have large differences in r̂_{2,m}, Φ and C̄_{v}/C̄_{d}; r̂_{2,m} is around 0.4 for the droneflies (Ellington, 1984a), and around 0.48 for the fruit flies [estimated based on data of Vogel (Vogel, 1965)]; Φ is around 110° for the droneflies (Table 2), and around 140° for the fruit flies [table 1 of Fry et al. (Fry et al., 2005)]; C̄_{v}/C̄_{d} is around 1.2 for the droneflies (calculated using the computed force data), and around 0.93 for the fruit flies [calculated using data in table 1 of Fry et al. (Fry et al., 2005)]. Then the above ratio can be written as: (14b)
For the fruit flies, although n is a little larger than that of the droneflies [around 218 Hz (Fry et al., 2005); n for the droneflies is around 172 Hz], R and m_{wg}/m are very small [R=2.39 mm (Fry et al., 2005); m_{wg}/m=0.12%, estimated using data of Vogel (Vogel, 1965); note that the m_{wg} is the mass of one wing]; these values are only about one fifth and one fourth of those of the droneflies, respectively], thus the inertial torque is small compared to the aerodynamic torque, resulting the small inertial power compared to the aerodynamic power. The large inertial power (compared to the aerodynamic power) for the droneflies can be explained by similar reasoning.
Flight force and moment balance and validation of the computation model
In our previous studies on flight force requirements and dynamic stability (e.g. Sun and Tang, 2002b; Sun and Wang, 2007), wing kinematic parameters for equilibrium flight, in which the forces and moments acting on the insect are balanced, must be given. Wing angles of attack (a_{d} and a_{u}) and mean stroke angle () were not available and they were determined using the equilibrium flight conditions, which require the vertical force being equal to the weight and for hovering flight, the horizontal force and the pitching moment about the center of mass being zero.
Now that the angles of attack and the mean stroke angle have been measured in the present study (along with other wing kinematic parameters and the morphological parameters of wing and body), it is of interest to input all these parameters to the computation model and see whether or not an equilibrium flight is realized. This could provide a test of the computation model.
The nondimensional weight of an insect, denoted as C_{G}, is defined as C_{G}=mg/0.5ρU^{2}(2S). Using data in Tables 1 and 2, values of C_{G} for DF1, DF2 and DF3 are computed and given in Table 3. The results of the mean vertical force, horizontal force and pitching moment coefficients (C̄_{V}, C̄_{H} and C̄_{M}), calculated using the computation model, are also shown in Table 3. It is seen that for all three droneflies considered, the equilibrium flight conditions are approximately met: C̄_{V} is close to C_{G} (C̄_{V} is different from C_{G} by less than 6%) and C̄_{H} and C̄_{M} are close to zero. These results show that the computation model is reasonably accurate.
Implications of the present results
Wing kinematic patterns and aerodynamic mechanisms
The present measurements on droneflies show that the wing moves at a constant angle of attack during the translational phase of a halfstroke and has a shallow Ushaped wingtip trajectory. Based on the wing kinematic pattern, analysis using the CFD model shows that the droneflies primarily use lift to produce weightsupporting force during the translational phase (via the delayed stall mechanism).
Fry et al.'s measurements on fruit flies (Fry et al., 2005) showed a very different wing kinematic pattern: the wing had a large downward plunge at the start of a halfstroke, resulting in a deep Ushaped wingtip trajectory. Based on this wing kinematics, their experimental study using a robotic wing showed that the fruit flies rely heavily on drag mechanism during stroke reversal in producing vertical force.
These results demonstrate that insects having different wing kinematic patterns may employ different aerodynamic mechanisms for flight, and that in order to reveal the aerodynamic mechanisms an insect uses, detailed wing kinematics measurements should be conducted first. It is suggested that the studies on hovering of fruit flies (Fry et al., 2005) and droneflies (present study) are extended to other flight modes, such as forward flight and maneuver, and later on, other representative insects should be studied.
CFD method, combined with detailed kinematics measurements: a tool of great promise
An advantage of the CFD method in the study of insect flight aerodynamics and dynamics is that when the computational code is validated for an insect in a certain flight mode (for example, hovering), it can be readily used for other flight modes and for other insects. When changing from hovering to forward flight or other flight modes, only the input boundary conditions, determined by the wing motion, need to be changed, and this can be easily accomplished. When the code is used for other insects, the wing shape, hence the computational grid, and the Reynolds number in the Navier–Stokes equations need to be changed (Wu and Sun, 2004); both can be accomplished without much difficulty.
A second advantage is that the CFD method can provide any physical quantities that are needed for flow analysis. For example, aerodynamic force distribution on the wing is available, thus total aerodynamic force, moment about the center of mass of body (for flight balance study), and moment about wing root (for power calculation) can be readily obtained. Another example is that streamline patterns and vortices in the flow can be easily visualized.
The CFD model in the present study has been successfully validated for hovering droneflies. As discussed above, it can be readily used for study of other flight modes of droneflies and flight of other insects, provided that detailed wing kinematic patterns are measured. Combined with the method of 3D highspeed videography, the CFD method could be a tool of great promise in the future study of insect flight aerodynamics and dynamics.
LIST OF ABBREVIATIONS AND SYMBOLS
 c̄
 mean chord length
 C_{D}
 wing drag coefficient
 C_{G}
 weight coefficient
 C_{H}
 horizontal force coefficient
 C̄ _{H}
 mean horizontal force coefficient
 C_{L}
 wing lift coefficient
 C_{M}
 pitching moment coefficient
 C̄ _{M}
 mean pitching moment coefficient
 C_{P}
 nondimensional power
 C_{V}
 vertical force coefficient
 C̄ _{V}
 mean vertical force coefficient
 C_{w}
 nondimensional work per cycle
 nondimensional positive work per cycle
 nondimensional negative work per cycle
 CFD
 computational fluid dynamics
 D
 drag of a wing
 g
 the gravitational acceleration
 h_{1}
 distance from wingroot axis to longaxis of body
 H
 horizontal force
 I_{wg}
 moment of inertia of the wing about wing root
 I_{y}
 moment of inertia of the body about center of mass
 J
 advance ratio
 l_{b}
 body length
 l_{r}
 distance between two wing roots
 l_{1}
 distance from wingroot axis to body center of mass
 L
 lift of a wing
 LED
 light emitting diode
 m
 mass of an insect
 m_{wg}
 mass of one wing
 M_{a}
 aerodynamic moment of wing around wing root
 M_{i}
 inertial moment of wing around wing root
 M_{y}
 pitching moment about center of mass
 n
 stroke frequency
 P_{a}
 aerodynamic power
 P_{i}
 inertial power
 P^{*}
 bodymassspecific power
 r_{2}
 radius of the second moment of wing area
 r_{2,m}
 radius of the second moment of wing mass
 r_{d}
 radius of the first moment of wing drag
 R
 wing length
 S
 area of one wing
 t
 time
 t_{c}
 period of wing beat
 t̂
 nondimensional time (t̂=0 and 1 at the start and end of a cycle, respectively)
 U
 reference velocity (mean flapping velocity at r^{2})
 V
 vertical force of a wing
 α
 angle of attack
 αd
 constant value of α at the midposition of the downstroke
 αu
 constant value of α at the midposition of the upstroke
 β
 stroke plane angle
 Δt̂_{p}
 nondimensional time duration for pronation
 Δt̂_{p,t}
 nondimensional advance time for pronation
 Δt̂_{s}
 nondimensional time duration for supination
 Δt̂_{s,t}
 nondimensional advance time for spination
 η1
 phase angle of the first harmonic of Fourier series of deviation angle
 η2
 phase angle of the second harmonic of Fourier series of deviation angle
 θ
 deviation angle
 θ0
 constant value of Fourier series of deviation angle
 θ1
 amplitude of the first harmonic of Fourier series of deviation angle
 θ2
 amplitude of the second harmonic term of deviation angle
 ρ
 density of fluid
 ν
 kinematic viscosity of fluid
 ϕ
 positional angle
 mean positional angle
 Φ
 stroke amplitude
 χ
 body angle
 χ0
 free body angle
 Ω
 wing rotation velocity vector
ACKNOWLEDGEMENTS
This research was supported by the National Natural Science Foundation of China (10732030).
FOOTNOTES

Supplementary material available online at http://jeb.biologists.org/cgi/content/full/211/13/2014/DC1
 © The Company of Biologists Limited 2008