Aerodynamic efficiency of a bioinspired flapping wing rotor at low Reynolds number

This study investigates the aerodynamic efficiency of a bioinspired flapping wing rotor kinematics which combines an active vertical flapping motion and a passive horizontal rotation induced by aerodynamic thrust. The aerodynamic efficiencies for producing both vertical lift and horizontal thrust of the wing are obtained using a quasi-steady aerodynamic model and two-dimensional (2D) CFD analysis at Reynolds number of 2500. The calculated efficiency data show that both efficiencies (propulsive efficiency-ηp, and efficiency for producing lift-Pf) of the wing are optimized at Strouhal number (St) between 0.1 and 0.5 for a range of wing pitch angles (upstroke angle of attack αu less than 45°); the St for high Pf (St = 0.1 ∼ 0.3) is generally lower than for high ηp (St = 0.2 ∼ 0.5), while the St for equilibrium rotation states lies between the two. Further systematic calculations show that the natural equilibrium of the passive rotating wing automatically converges to high-efficiency states: above 85% of maximum Pf can be obtained for a wide range of prescribed wing kinematics. This study provides insight into the aerodynamic efficiency of biological flyers in cruising flight, as well as practical applications for micro air vehicle design.


Introduction
Flapping wing based propulsion has several aerodynamic benefits for the application of micro air vehicles (MAVs), as demonstrated by the superior flight skills of natural flyers such as insects or birds. In particular, flapping wing produces higher lift force than conventional aerofoil at an angle of attack (AoA) above the stall angle, due to the delayed stall of the leading edge vortex (LEV) [1][2][3]. On the other hand, flying insects and birds have shown the extraordinary capability of vertical take-off, landing, hovering and manoeuvrability. These features of the flapping wing have brought a strong interest in developing MAVs that mimic the wing motion of insects or birds [4,5].   For aeronautical vehicles and flying animals, the wings produce significant lift to support weight in the air, the aerodynamic efficiency is often defined by the cost of energy to stay aloft or to travel a certain distance, which is associated with lift production. For example, the current studies on insect flight primarily focus on a 'normal hovering' state, where the aerodynamic force of the flapping wings averages to a net lift primarily to support the weight in the air. The aerodynamic efficiency of the wing is therefore mainly concerned with the production of lift for a given power. On the other hand, moving animals such as swimming fish and cruising flyers use their wings or fins to produce thrust against the drag of the fluid, for which the propulsive efficiency associated with thrust production is often used to measure the efficiency of their movement of wings or fins. For flying and swimming animals, a dimensionless parameter which describes the kinematics of their wings or fins is the Strouhal number St [6][7][8]. This dimensionless number is known to govern a well-defined structure of the wake shed from an oscillating aerofoil in the free stream, and is closely related to the efficiency of flying and swimming animals for propulsion [6,8]. In particular, optimum propulsive efficiency is found for wing motions with St lying between 0.2 and 0.4, corresponding to a stably formed wake structure and average velocity profile equivalent to a jet [6][7][8][9][10]. Further extensive literature reviews on data from flying and swimming animals showed that many animals cruise at this interval of St [6,8].
The existing studies on flapping wing efficiencies have been primarily focused on one of the above two aspects. For the normal hovering kinematics of insects, because the body assumes no forward speed, the flapping wing produces no net thrust in a flapping cycle. Thus, the aerodynamic efficiency for producing lift is of primary concern for such wing motions. On the other hand, most study on the undulatory swimming fish tails have been devoted to propulsion, while the force perpendicular to the free stream in the form of lift is of secondary factor. However, natural flyers in cruising flight use their flapping wings to generate both lift and thrust. The aerodynamic efficiency of such wing motion in terms of both propulsive efficiency and efficiency of lift production has received less study.
One particularly interesting kinematics is a wing that flaps and rotates at the same time. Vandenberghe et al. [11] studied these kinematics experimentally with a rectangular wing mounted on a vertical shaft and free to rotate horizontally. As the wing flaps up and down above a threshold frequency, it starts to rotate and finally reaches a stable speed with St ∼ 0.26 for a range of input frequencies. Similar wing kinematics has been proposed for the design of new helicopters without the reaction torque [12]. Guo et al. [13] investigated a new type of bioinspired flapping wing vehicle, namely the flapping wing rotor (FWR), for which a pair of antisymmetric configured wings is rotated by the aerodynamic thrust produced by vertical flapping motion, as illustrated in figure 1. The positive average aerodynamic lift can be obtained by varying the wing pitch angle in an asymmetric manner. Further numerical investigation has shown that the flow on the flapping and rotating wing forms compactly attached three-dimensional (3D) vortex ring structure which connects the LEV, trailing edge vortex and wingtip vortex that enhances lift production [14]. Li et al. [15] and Wu et al. [16] used quasi-steady (QS) aerodynamic analysis and CFD models, respectively, and showed that the FWR kinematics in hovering flight can produce higher lift coefficient than the conventional insect-like flapping wings and rotary wing. The aerodynamic efficiency of the FWR for lift production (defined by the dimensionless power factor P f [17,18]) is standing between the other two conventional types.
For animals in cruising flight, the equilibrium between the thrust produced by the flapping wing and the aerodynamic drag results in a wing kinematics pattern that closely resembles the wing motion of the FWR. In this study, investigations have been made into the aerodynamic efficiency for producing both vertical lift and horizontal thrust of the FWR kinematics. A 3D QS aerodynamic model which adapts the empirical coefficients obtained from high-fidelity CFD simulation is used to calculate the  aerodynamic force and power of the wing. Additionally, analysis using a two-dimensional (2D) CFD model is carried out to capture the transient status of the flow field and unsteady forces. A typical wing model of elliptical shape and aspect ratio λ = 3.6 and semi-span R = 0.098 m is chosen for the FWR model. The wing flapping motion follows a simple harmonic function (SHF) at a constant rotating speed. The wing flapping frequency is f = 10 Hz, flapping amplitude Φ = 70°corresponding to a Reynolds number of 2500.
The results of this study show that both the propulsive efficiency and efficiency for producing lift of the FWR wing peaks within a narrow interval of St: between 0.1 and 0.5 for certain range of wing pitch angles, which agrees closely with reported data of natural flyers in cruising flight [8]; however, the St for high efficiency of lift and high efficiency of propulsion in general differs. The higher propulsive efficiency corresponds to St between 0.2 and 0.5 and the higher efficiency of lift corresponds to St between 0.1 and 0.3; in particular, the St for the rotational equilibrium state of the wing lies between the maximum propulsive efficiency state and the maximum efficiency of lift state. Furthermore, systematic calculations show that high efficiency of lift (above 85% with respect to maximum P f ) can be obtained at the natural equilibrium state of the wing for wide range of prescribed wing kinematics. Insights of the results for biological flyers in cruising flight as well as for MAVs design are provided.

Wing kinematics definition
The coordinate system to define the FWR wing motion is shown in figure 2. The kinematics of the wing is defined by three elementary motions: rotation, flapping and pitching. The wing rotates about the vertical y-axis, and the rotation speed is fixed constant and indicated byψ. The dimensionless wing rotation speed is defined as where Φ is the flapping angle amplitude and f is the flapping frequency. The wing flaps vertically with the flapping angular velocityφ described by the SHF: where τ = ft is the dimensionless time ranging from 0 ∼ 1 in a flapping period.
The wing pitches at the same frequency f, with a phase shift of π /2. The pitch motion of the wing is confined to stroke reversals. At the mid-upstroke and mid-downstroke, the wing has angle of attack (AoA) denoted by α u and α d , respectively. In the reversal phase at the end of each stroke, the pitch angular velocity of the wingα is described by the following equation: where τ r = 0 ∼ 1 indicates the dimensionless wing pitch time with respect to the flapping period, and the bracket notation [·] indicates the floor function giving the greatest bounding integer. In this study, the wing pitch time takes half of the flapping period, corresponding to τ r = 0.5.

Quasi-steady model and dimensionless parameters
The wing planform is in ellipse as illustrated in figure 3. The pitching axis of the wing (z w ) is located near the leading edge; c is the local chord length at a 2D wing strip dr; h is the vertical distance between the mid-chord axis and the pitching axis of the wing. The pitching axis of the wing is taken to be located at 0.25 chord [19,20], corresponding to h = 0.25 c.
The aerodynamic force of the wing is solved by employing a QS aerodynamic model. The details of the model are provided in the previous study [15]. On a 2D wing strip, the aerodynamic forces and pitch moment in the local wing-attached frame (x w , y w and z w , as shown in figure 2) are obtained by the equations: where ρ is the air density; U is the translational velocity; ω z andω z are the wing pitch rate and pitch acceleration; u x , u y andu y are the translational velocity components in the x w -and y w -axes and the acceleration; λ y and λ yω are the added mass force coefficients, which are given by λ y = π /4ρc 2 and λ yω = π /4ρhc 2 ; λ ω is added mass moment coefficient, which is given as: λ ω = π /4ρh 2 c 2 + π /128ρc 4 [21]. The QS aerodynamic coefficients C H , C V and C ROT are empirical coefficients for the translational force and rotational force [20,22]; C RD is the rotational damping moment coefficient;x CP is the dimensionless centre of pressure (CP) of the translational force;x RD is the dimensionless location of the rotational damping force; The relation of the empirical coefficients with the effective AoA α e are provided in our previous study [15].
The 2D aerodynamic forces and pitch moment for each wing strip are integrated along the wing-span to obtain the 3D forces and moments. The vertical lift force and rotational moment of the wing can then be obtained by projecting the 3D force and moment vectors onto the global y-axis of the coordinate, as shown in figure 2. The lift and rotational moment coefficients are defined as and where l is the vertical lift force and m is the rotational moment, i.e. moment along the global y-axis,c is the mean chord length, S is the wing area. The reference velocity U 2 is defined by where the characteristic width A is taken to be the stroke amplitude at the wingtip, and the forward speed U is taken to be the rotation speed of the wing at the wingtip [8]. For a 2D aerofoil, the reduced frequency k c is defined by where U c is the rotation speed at the specific chord.

Aerodynamic power and efficiency measures
The mean aerodynamic power over a flapping cycle T can be obtained by summing the aerodynamic power of each independent axis of rotation:P is the mean aerodynamic power of the ith axis. The mean aerodynamic power coefficient can be defined asC The mean aerodynamic powerP is the total power required for the wing to overcome the fluid forces. Therefore, when measuring the efficiency of lift production,C P can be used directly to define the dimensionless power factor P f [17,18], which measures the power efficiency of flying animals and vehicles for sustaining a specific weight: (2.14) However, when measuring the propulsive efficiency, the mean aerodynamic power of each independent axisP i needs to be treated differently. The propulsive efficiency of the wing η p is defined by the ratio of aerodynamic power output (for propulsion) to the power input: When the flapping wing is producing positive propelling moment, this definition agrees with the usual definition of propulsive efficiency for oscillating foils in the free stream [6,7].

Comparison of 3D quasi-steady and 2D unsteady forces
This study employs a 3D QS aerodynamic model for modelling the forces and power of the FWR wing. Additionally, a 2D CFD analysis is carried out using the commercial CFD solver ANSYS Fluent, which solves the 2D unsteady, incompressible Navier-Stokes equations based on a finite volume method. The accuracy of this solver has been extensively validated against several experimental and numerical studies in flapping wing aerodynamics [23]. The 2D aerofoil is chosen as flat plate of 2% thickness. The choice of the simplified models is due to the efficiency for computing the various wing kinematic cases.
In this investigation, the QS model and 2D CFD model are compared with the high-fidelity 3D CFD results from the previous study. A particular kinematic case from Wu et al. [16] is chosen with the kinematic parameters of the wing specified by: Φ = 70°, α u = 60°, α d = −20°and the rotation speed η = 1.10. In this case, the wing model is of rectangular shape and with aspect ratio 5.8 (λ = R/c). The wing semi-span R and the flapping frequency f are given by R = 0.098 m and f = 10 Hz, corresponding to Re ∼ 1600 (Re = U 2c /υ, where U 2 is the mean flapping velocity at the radius of the second moment of wing area r 2 , andc is mean chord length). For the 2D CFD analysis, a series of wing chords located along the wing-span (ranging between 0.2 and 0.7 wing-spans) is taken for investigation. To compare the 2D model with 3D results, the thrust coefficient of the 2D results (C T = T/0.5ρŪ 2 c, where T is the thrust force andŪ is the local mean flapping velocity) is converted to the 3D rotational moment coefficient C m using a scale factor of λ(I 3 /I 2 ) obtained from a standard blade element analysis, λ is the wing aspect ratio and I k is the kth dimensionless moment of wing area defined by the equation: The time courses of the forces and rotational moments by different models are shown in figure 4a, and the flow field for the 2D wing chord is shown in figure 4b. All the cases have the same St = 0.45; the 2D CFD result presented here is taken at 0.35 wing-span, with reduced frequency k c = 1.15. The full spectrum of 2D results at different span-wise locations ranging between 0.2R and 0.7R is provided in the electronic supplementary material.
From the results in figure 4a, it is clear that the time courses of the 3D QS forces and moments agree very well with the 3D CFD results (a more comprehensive validation of the QS model for both transient and mean aerodynamic forces is provided in [15]); while for the given wing chord at 0.35 wing-span, the 2D CFD results appear larger due to the well-known downwash of the wingtip vortex for a 3D wing. However, the variations of the 2D transient forces agree qualitatively with the 3D transient forces.
In figure 4b, the flow over the 2D wing chord shows a dynamic formation and shedding of vortices that closely resembles the dynamic stall of conventional aerofoil. In the downstroke, a strong LEV first forms on the upper surface until a reverse flow emerges from the trailing edge and the LEV starts to shed from the surface. The shedding of the LEV is near the end of the downstroke, corresponding to a decrease in lift coefficient.
The particular choice of wing chord at 0.35R to represent the flow field follows from the observation that the shedding of LEV is at the end of the downstroke. In previous studies on 2D flapping wing, the frequency of LEV shedding was shown to play a significant role in determining the transient forces.
Lewin & Haj-Hariri [24] and Wang [25] used numerical simulation to analyse an aerofoil in 2D flow undergoing heaving oscillating. They show that the dimensionless reduced frequency k c and the Strouhal number St serve to govern the time scale associated with the growth and shedding of the vortices on the wing. The k c is the primary factor governing the LEV shedding and St the secondary factor related to the growth of LEV. For a smaller value of k c , the LEV tends to separate and advect along the freestream, leading to an early separation of the LEV; while for larger k c , the LEV tends to separate later and stays longer on the wing in each flapping stroke. Similar flow phenomenon is observed in this study. The LEV on the outer wing chords (with smaller k c ) tends to separate early while it stays longer on the inner wing chords with larger k c (see the electronic supplementary material).
In this study, in order to obtain the qualitative characteristics of the flow, the 2D CFD analysis is taken at different span-wise locations of the FWR wing, which has the same St but different k c (ranging between 0.6 and 2.0 for different St cases). The particular cases of the 2D results with LEV shedding frequency that matches with the flapping frequency of the wing (i.e. LEV shed at the end of the downstroke) is then chosen to represent the flow field, the resulting transient forces of the 2D calculations are found in qualitative agreement with 3D models.

Propulsive efficiency versus the efficiency of lift
The FWR kinematics makes use of aerodynamic thrust produced by flapping motion to drive the wings to rotate about the vertical axis. At the same time, lift is obtained by biasing the pitch angle of the wing in the upstroke and downstroke. In this study, both the propulsive efficiency η p and the efficiency for producing lift P f are investigated for the FWR kinematics. These two efficiency measures are calculated in different kinematic conditions defined by the dimensionless St and wing pitch angles. An FWR wing model with wing shape illustrated in figure 3, wing semi-span R = 0.098 m, aspect ratio λ = 3.6, flapping amplitude Φ = 70°and flapping frequency f = 10 Hz at Re ∼ 2500 is taken for this study. The wing pitch angles are defined for four cases, varying from symmetric pitching to asymmetric pitching as: α d = −15°and α u = 15°, 30°, 45°, 60°. The S t for each of the above cases are chosen to vary between St = 0.1 ∼ 1, which determines the rotation speeds of the wing uniquely (η = 0.5 ∼ 5 for the given St range). The computed results for aerodynamic efficiencies (η p and P f ) and force (moment) coefficients (C l andC m ) against the St are shown in figure 5. The St for maximum η p and P f at different α u cases is given in table 1.
As shown in figure 5a, the variation of the propulsive efficiency η p and the efficiency of lift P f against the St follows a similar trend. As the St increases from 0.1 to 1, both η p and P f first increase rapidly to the maximum values and then decrease. It is also observed that both the maximum η p and P f occur within a narrow interval of St = 0.15 ∼ 0.5 when the wing upstroke pitch angle is within 15°∼ 45°. However, these two efficiencies appear to be complementary with respect to each other, hence do not reach maximum at the same time. In one of the extreme cases at the lower end of St = 0.19, the flapping motion of symmetric up-and downstroke pitching (α d = −15°, α u = 15°) leads to the optimal propulsive efficiency η p = 0.37 at the cost of zero mean lift coefficient (figure 5b) and efficiency P f = 0. When the flapping motion becomes asymmetric with α u = 30°∼ 45°, figure 5b shows that theC l ,C m and also P f increased dramatically, but the η p is reduced. The complementary nature of η p and P f can also be seen from figure 5a,b that the maximum P f always occur at small St with negativeC m , where the FWR wing produces net drag instead of thrust.
Previous studies have shown that the propulsive efficiency of flapping aerofoils is closely related to the evolution of the flow structure on the wing. Triantafyllou et al. [6,7] studied the propulsive efficiency of 2D oscillating aerofoil and proposed that the optimal efficiency is obtained when an aerofoil is flapped at the frequency that results in the maximum amplification of the shed vortices, and the velocity profile behind the aerofoil is in the form of an inverted von Kármán vortex street indicative of a jet. Later, by using an inviscid panel method to investigate the wake structure of a 2D oscillating aerofoil, Jones et al. [26] noted a remarkable similarity between the simulated wake and the experiment   wake structure, indicating that the formation of the well-defined wake structure is essentially an inviscid phenomenon.
At low Re and large AoA, the well-defined structure of the wake is complicated by the flow separation at the leading edge and interactions with the vortices shed from the trailing edge. Wang [25] studied the flow over an impulsively started 2D aerofoil using CFD method and observed that the thrust production is correlated with the time scale that governs the shedding of the LEV. He proposed that optimal efficiency is obtained when the duration of the flapping stroke is inside the 'thrust window' that exists before the LEV is shed.
In the view of the previous results in 2D, we have taken the optimal kinematic cases (shown in figure 5 and table 1) that result in the highest propulsive efficiency and efficiency of lift to investigate the 2D flow and forces. The 2D calculations are conducted for the cases with symmetric pitch angles: In figure 6, large thrust is produced in both the upstroke and downstroke for the two cases. In the symmetric pitching case, the production of thrust is equal in the up-and downstroke; while in the asymmetric pitching case, the thrust produced in the downstroke is larger than in the upstroke. Figure 6a shows that the LEVs of the symmetric pitching case are of equal strength on both sides of the wing in the up-and downstroke, which contributes equally large thrust in each stroke. By contrast, figure 6b shows that asymmetric LEVs are formed on the wing for the asymmetric pitching case. In the downstroke, a strong LEV is formed on the upper wing surface which contributes significant lift and thrust; while in the upstroke, a weak LEV is formed on the lower wing surface which contributes small thrust and negative lift.
It is observed that the scales of the LEVs associated with the above cases are relatively small. This is due to the small effective AoA at the given St. In general, the formation of LEV decreases the propulsive efficiency, as indicated in previous studies for oscillating 2D aerofoils [27][28][29]. However, LEVs also serve as a source of thrust production. It is noted that the St serves to balance these two effects. As shown in figure 5b, when St decreases from moderate to small value, a transition from large rotational moment to negative is observed, indicating a decrease of the effective AoA from large value to negative. Maximum propulsive efficiencies occur at medium St where small positive effective AoA forms small LEV, as shown in figure 6. Figure 7 shows the time courses of lift, rotational moment and 2D flow structure of the wing for the maximum P f case with an even larger pitch angle (α u = 45°). The forces and moments of the 2D CFD results follow qualitatively the trends of the 3D QS forces and moments. However, disturbances of the 2D forces are observed due to the shedding of LEVs. In this case, lift is produced in both upstroke and downstroke. Large anti-rotating moment is produced in the upstroke and only a small rotational moment is produced in the downstroke. It is noted that in both up-and downstroke, large LEVs are formed on the upper surface of the wing and significant vortices shedding are observed near the end of each stroke. This is most likely due to the fact that the wing has consistently large effective AoA in a whole flapping cycle at this St.
From the above analysis, it is noted that as the wing pitch angles change from symmetric to asymmetric in the up-and downstroke, a transition of the LEV structure from symmetric to asymmetric is observed. The asymmetry of the LEV structure on the wing surface results in net lift production. However, as the flow forms large LEV in the downstroke, higher power is required to overcome the vertical lift force, which leads to a diminished propulsive efficiency. The transition of the flow structure from symmetry to asymmetry thus indicates a transfer of flapping energy from propulsion to weight suspension.

Aerodynamic efficiency of passive rotating wing
When the FWR wing is free to rotate horizontally, the rotation speed will reach an equilibrium state when the mean rotational moment over a flapping cycle is zero. The rotational equilibrium state of the FWR kinematics has been proposed for practical design of MAVs. Several previous studies have shown that this kinematics may have certain benefits in terms of system simplicity and aerodynamic efficiency compared with other conventional type (fixed wing, rotary wing and insect-like flapping wing) when applied for MAV design [13,15,16]. Apart from practical applications, the passive rotating kinematics also has a notable similarity with the cruising flight of natural flyers, where the flapping wings produce both lift and thrust, and the cruise speed is determined by the equilibrium of the body drag and the flapping propulsive thrust.
The present investigation focuses on the rotational equilibrium state of the FWR wing. As the wing produces no net propelling moment at this state, the efficiency for producing lift P f is of particular interest. In this investigation, the aerodynamic efficiency at the equilibrium state is compared with the maximum aerodynamic efficiency, which is the highest aerodynamic efficiency that can be obtained by actively tuning the rotation speed (or equivalently the St, as shown in figure 5) of the wing. The efficiency for producing lift at the rotational equilibrium state is denoted by P fe , while the maximum efficiency for the given wing pitch angles (α u and α d ) is indicated by P fm .
In this investigation, the two efficiencies P fe and P fm are calculated with the downstroke wing pitch angle α d prescribed between −45°and −10°and the upstroke wing pitch angle α u prescribed between 25°and 70°. Specifically, the individual cases when α d is fixed at −15°, −30°and −45°are taken out to analyse the efficiency variations with α u . For all the cases, the wing semi-span is taken as R = 0.098 m, wing aspect ratio λ = 3.6, the flapping frequency f = 10 Hz, flapping amplitude Φ = 70°and the Re is about 2500. The variation contours of P fe , P fm and the ratio P fe /P fm with α u and α d are shown in figure 8b. The efficiency curves with fixed α d are shown in figure 8a. In figure 8b, the negative lift regions correspond to wing kinematic cases with larger negative α d in the downstroke than positive α u in the upstroke, which consequently yields negative mean lift force in a flapping cycle.
The results in figure 8 shows that, in general, the efficiencies (P fe and P fm ) are sensitive to the variation of wing pitch angles (α u and α d ). Increasing the downstroke AoA α d generally leads to decreases in the aerodynamic efficiencies; while the optimal upstroke AoA α u for P fe and P fm always takes intermediate values between 30°and 60°, depending on α d . However, despite these variations with wing pitch angles, the efficiency at the equilibrium states P fe appears to be very close to the maximum efficiency P fm . In most of the investigated kinematic cases (i.e. α d between −45°and −15°), the ratio P fe /P fm is above 85%; furthermore, when the downstroke wing pitch angle is large (i.e. α d between −45°and − 30°), the ratio P fe /P fm reaches even higher value of above 90%.
The above results imply that for the passive rotating wing, the forces equilibrium of the flapping propulsive thrust and anti-rotating drag in the up-and downstrokes results in a wing kinematic state of high aerodynamic efficiency. It is therefore expected that the passive rotation kinematics may be favourable for bioinspired MAV design, because the rotation speed automatically converges to a highefficiency state. Furthermore, because flapping wing flyers in cruising flight are in a state of equilibrium where the production of thrust by their flapping wings balances with the drag from the body, the above results imply that cruising flyers may also benefit from this natural equilibrium to gain high aerodynamic efficiency of lift production at this state.
In order to fully characterize the kinematics of the three statuses of the wing, i.e. the maximum propulsive efficiency state, the state with maximum efficiency of lift and the equilibrium state, the St of these respective states at different wing pitch angles are further calculated. It should be noted that  because the St serves to determine the production of aerodynamic forces from propulsive thrust to anti-rotating drag, the St at equilibrium state stands in the middle of the other two states. Figure 9a shows the variations of St at the typical states with α u for fixed α d cases; the full contours of typical St for different wing pitch angles are shown in figure 9b. In figure 9, the St for maximum η p states are denoted by S tp , for maximum P f states are denoted by S tm , and for equilibrium states are denoted by S te .
The variations in figure 9a show that the St for typical states (S te , S tm and S tp ) increase monotonically with the increase of α u . Figure 9a,b shows that reducing the wing pitch angle (i.e. decrease α u and increase α d ) lowers the St for the high-efficiency states, indicating that higher efficiencies (η p and P f ) are obtained at higher rotation speed at smaller wing pitch angles (because St is inversely proportional to the rotation speed).  shows that for most of the kinematic cases with α u between 25°and 70°and α d between −45°and −10°, the St at equilibrium states S te and maximum P f states S tm fall in the interval between 0.1 and 0.5; while the St at maximum η p states appears to be higher between 0.2 and 0.9 for most of the cases. Particularly, the distribution of the data shows that when the upstroke AoA is small, i.e. α u less than 45°, nearly all the St for high-efficiency states S tm , S tp and the equilibrium state S te converge to the interval of 0.1 ∼ 0.5. The lower end of St (between 0.1 and 0.3) corresponds to higher P f states, while the higher end (between 0.2 and 0.5) corresponds to higher η p states. The current results are in close agreement with previously reported data of flying animals in cruising flight, which show that many natural flyers cruise with St between 0.2 and 0.4 [8]. The above results imply that the various wing kinematics of flapping wing flyers in cruising flight may result in high aerodynamic efficiency states for both lift production and propulsion, although these two efficiencies cannot be optimized at the same time.

Conclusion
The aerodynamic efficiency of a novel FWR kinematics which combines vertical flapping motion and passive horizontal rotation is investigated by using a QS aerodynamic model and 2D CFD analysis. The propulsive efficiency η p for producing horizontal thrust and the efficiency P f for producing vertical lift of the wing are investigated for a wing model of elliptical shape with wing semi-span R = 0.098 m, wing aspect ratio λ = 3.6, flapping vertically with a frequency of f = 10 Hz and amplitude Φ = 70°at the Re of 2500. The calculated data show that both the propulsive efficiency η p and efficiency of lift P f depend on the dimensionless St and wing pitch angles (α u and d u ). For small wing pitch angles (upstroke AoA α u less than 45°), both efficiencies η p and P f are found to peak at St between 0.1 and 0.5. However, these two efficiencies are complementary to each other: when maximum η p is obtained, the P f is relatively low; while maximum P f always occurs when the flapping wing produces net drag instead of thrust. In particular, high efficiency of lift production is found when St is between 0.1 and 0.3, which is, in general, lower than the St for high propulsive efficiency (between 0.2 and 0.5). Further analyses of the 2D flow of the wing in typical kinematic states show that the production of lift and thrust is closely related to the LEV structure on the wing. Maximum η p occurs when the wing forms small and symmetric LEVs in the up-and downstroke; while asymmetric LEVs large in the downstroke are associated with the production of lift and decrease in propulsive efficiency η p .
The rotational equilibrium state of the FWR kinematics is investigated. The results show that the aerodynamic efficiency at this state P fe is above 85% compared with the maximum efficiency P fm for a wide range of wing kinematics. Furthermore, systematic calculations show that most of the St for the high-efficiency states (maximum P f state and maximum η p state) and the equilibrium state of the wing are within the interval of St between 0.1 and 0.5. These results show that the natural equilibrium between thrust and drag of the flapping wings result in a state of high aerodynamic efficiency. The agreement of our results with reported data of cruising flight of animals indicates that flapping wing flyers may benefit from this equilibrium to gain high efficiency for both lift and thrust production, although these two efficiencies cannot be optimized at the same time. The above results also have implications for bioinspired MAV design, because for the passive rotating kinematics of FWR, no fine tuning of the rotation speed is needed to achieve a high-efficiency state.
Data accessibility. All supporting data are made available either in the article or the electronic supplementary material. Authors' contributions. H.L. wrote the Matlab codes for QS aerodynamic model, ran the CFD analyses, wrote most of the paper; S.G. conceived of the study, coordinated the study and helped draft the manuscript. Both authors participated in development of concepts and critical revisions and gave their approval for publication.