## Abstract

Insects are a prime source of inspiration towards the development of small-scale, engineered, flapping wing flight systems. To help interpret the possible energy transformation strategies observed in Diptera as inspiration for mechanical flapping flight systems, we revisit the perspective of the dipteran wing motor as a bistable click mechanism and take a new, and more flexible, outlook to the architectural composition previously considered. Using a representative structural model alongside biological insights and cues from nonlinear dynamics, our analyses and experimental results reveal that a flight mechanism able to adjust motor axial support stiffness and compression characteristics may dramatically modulate the amplitude range and type of wing stroke dynamics achievable. This corresponds to significantly more versatile aerodynamic force generation without otherwise changing flapping frequency or driving force amplitude. Whether monostable or bistable, the axial stiffness is key to enhance compressed motor load bearing ability and aerodynamic efficiency, particularly compared with uncompressed linear motors. These findings provide new foundation to guide future development of bioinspired, flapping wing mechanisms for micro air vehicle applications, and may be used to provide insight to the dipteran muscle-to-wing interface.

## 1. Introduction

The physiology, kinematics and aerodynamics of flapping flight are historically rich fields of research that have recently received complementary attention owing to interest in developing biomimetic or bioinspired wing motor mechanisms for micro air vehicle (MAV) applications, such as search-and-rescue and reconnaissance missions [1]. The rationales for scrutinizing and emulating insect flight mechanisms are multifold: natural selection gives credence to the effectiveness and robustness of the evolved solutions [2]; insects are of the same length scales as target MAV design spaces to favourably match aerodynamic influences [3]; engineering untethered sensory and control strategies is challenging and biological practices are prime sources of inspiration [4]. Such justifications highlight a striking aspect of MAV development. While the biological research elucidates critical physics of the flight mechanisms, realizing MAVs is a pursuit of flight itself, which has no guarantee of success. As a result, taking inspiration from nature's successful flyers appears to be an effective development strategy for MAVs.

Dipteran flight mechanisms are particularly well studied, and researchers have uncovered important aspects regarding muscular control over resulting flight dynamics [5,6]. However, the theories regarding the muscle-to-wing interface—the working wing motor—differ in terms of describing the motor as a hinged structure [7,8], an intricate gearbox [9–12] or a bistable ‘click’ mechanism [13,14]. Controversy aside, Dickinson & Tu [5] provided the insights that, whatever the mechanical means, (using an automotive metaphor) ‘the transfer of energy from the engine to the wheels [is] regulated by the configuration of the transmission’. In other words, it is well established that Diptera chiefly alter aerodynamic forces via wing stroke amplitude control related to muscle-induced changes in the motor mechanical advantage [5,15,16]. From this knowledge, it is deduced that the muscle–structure integration of the working flight mechanism plays an important role to significantly tailor wing stroke dynamics. Nevertheless, the composition of the interfacing mechanism remains a subject of continued study and attention.

In this light, the compressed or bistable click mechanism perspective of dipteran wing motor architecture has certain attraction given bistable oscillator modelling results compared with some biological evidence [14], and the operational bistability of the motor appears evident when illustrated in particular schematic form, such as in figure 1*a* as redrawn from Thomson & Thompson [14]. Researchers have also analysed the dipteran wing motor as a bistable mass-spring-damper and demonstrated that the excited nonlinear dynamics improve energy efficiency when compared with a linear motor using sinusoidal motions [17,18]; experimental results verified the analyses under certain operating conditions using a mesoscale air vehicle with bistable motor mechanism [19]. As pertains to MAV development, the opportunity for a biologically inspired, gear-free flight mechanism to provide significant and potentially continuous modulation of wing stroke amplitude for various flight dynamics has merit in terms of versatility, fabrication ease and durability. Yet, nature's solutions tend to be highly adaptable, and the fixed, bistable systems considered [17–19] lack the realistic, wide-ranging adjustability in diptera [5,12], whatever the precise source of the muscle-to-wing energy transformation. Consequently, there is impetus to assess how a compressed or bistable flight mechanism may usefully govern wing motions for flight dynamics control and a biologically inspired approach to the problem may prove fruitful. For the sake of clarity, unless specifically referring to a bistable configuration, the remainder of the text refers to this flight mechanism concept as a ‘compressed motor’.

The objective of this study is to investigate previously unexplored structural aspects of a representative wing motor mechanism which are indicative of (possibly collective) constituents in Diptera. The motivations are to exploit the potential interface between adaptive muscular control and resulting wing stroke dynamics and yield enlightened design strategies for new generations of MAV flight mechanisms. By adopting the compressed motor perspective, we seek to uncover means by which such a flight mechanism may near-continuously vary wing dynamics which are biologically relevant as well as advantageous for inspired MAV designs. Moreover, the experimentally validated modelling strategy devised here can inform engineering development by identification of critical design variable impacts and may introduce new opportunity to interpret analogous sensitivities of the biological system, albeit at coarse resolution.

The primary guiding cues for this investigation are found at the intersection of dipteran physiological and anatomical understanding with nonlinear structural dynamics. Using the compressed motor perspective as shown in the schematic in figure 1*a*, collective motions of motor constituents supported between the scutum and pleural wing process (PWP) are primarily ‘transverse’, in other words, up and down, similar in dynamic shape to the fundamental mode response of a simply supported beam. Thus, we refer to motor support dynamics and characteristics as ‘axial’ to distinguish them from the transverse motions of the parascutal shelf (PSS) and first and second axillary sclerites (Ax1 and Ax2) which govern wing motions. Wing supporting element elasticities and the tensioned muscles connected to them are recognized to govern wing stroke amplitude control, regardless of the interpreted mechanism architecture [5,7,13,15]. On the other hand, the precise translation from changing axial support stiffness and compression characteristics to wing stroke amplitude adjustment is still unclear. In nonlinear structural dynamics, axial support flexibility of compressed beams is known to significantly affect both static and dynamic energy storage in transverse bending, which has critical implications when the beam is compressed near the threshold of elastic stability [20,21]. Considering the potential structural analogy between the dipteran wing motor and a compressed, simply supported beam, this study hypothesizes the importance of wing motor *axial support stiffness* and *compression characteristics* in stroke amplitude control and suggests that active tuning of the characteristics may empower an exceptionally large and near-continuous range of wing motions while flapping frequency and exciting forces are unchanged. Such axial support variables as means for dipteran flight control were not included in the previous investigations adopting the compressed (and specifically bistable) motor perspective [14,17–19]. Thus, the potential influence of axial support characteristics on flapping flight is unclear and the opportunity to intelligently harness axial stiffness and compression for substantial and versatile adjustment of wing stroke in MAV designs remains unrealized.

To explore and scrutinize the importance of these previously unexplored motor structure aspects, the following sections introduce the motor architecture, describe the representative experimental system and modelling formulation, and report our investigations, discoveries and conclusions.

## 2. Representative wing motor architecture

The model wing motor architecture considered here is selected to reasonably bridge biological and engineering perspectives and to distribute the benefit of research findings. The dipteran wing motor schematic is shown in figure 1*a*; the representative structural model is shown in figure 1*b*. Rows (1), (2), (3) indicate comparable anatomical and structural configurations of the systems over the course of one-half flapping/excitation cycle; rows (1) and (3) are not to be necessarily interpreted as orientations of static equilibria but rather as dynamic states. We study a base excited simply supported beam transversely vibrating in bending indicative of the combined stiffness and inertia of the PSS, Ax1 and Ax2 which relatively translate and rotate owing to oscillatory scutellar lever excitations that are driven by deformation of the notum via longitudinal and dorsoventral muscles [7,14]. As shown later in our investigation, the base excitation of the structural wing motor is mathematically and mechanically analogous to *force* excitation. The motor dynamics determine wing motions through coupling to the radial vein and a representative location of this component as applied to the structural system is provided in figure 1*b*. In other words, assuming a rigid wing, the wing tip motion is directly related to the wing motor displacements through an expression of the moment/lever arm relation between the two motions. For the remainder of this study, wing and aerodynamic influences are accounted for via rotational inertia at the simple supports and nonlinear damping [2] which will be described in greater detail in the following sections. Beam axial boundary conditions include one end fixed in translation, whereas the other end is supported by an axial spring of stiffness *k*_{d} which statically compresses the beam by distance Δ. The combined boundary influences of the variably stiffened axial spring are representative of the scutum and PWP, where the restoring elasticities/stiffnesses are tuned by several muscles, including the pleurosternal (PSM), dorsoventral and first basalar [7,13,15]; the comparable influence of a modulated compression distance represents adjustments such as scutum positioning. Depending on axial support stiffness and compression distance relative to motor structure (beam) parameters, the system may be linear, monostable nonlinear or bistable; these are dynamic dependencies not explored in previous studies [13,14,17–19].

## 3. Experimental system, methods and preliminary evaluation

### 3.1. Experimental system and methods

The system fabricated to realize our representative structural wing motor is shown in figure 2. The right-most translational bearing (Parker 4301) is locked in place, whereas the left-most translational bearing may move to accommodate static compression distance Δ. Side and top schematics, figure 2*a,b*, respectively, indicate the location of a micrometer (Starrett 261) that adjusts the compression distance Δ, and show the aluminium cantilever support springs each having equivalent one-dimensional axial stiffness *k*_{d} [22] with respect to transverse bending motions of the spring steel beam (McMaster 9074K27, tempered 1074/5). Axial stiffness is adjusted by varying cantilever lengths. Rotational inertias at the simple supports are the consequence of the simple support axles; the inertias are modelled as cylinders of mass *m*_{a} and radius *r*_{a}, and are here representative of wing inertia. The nonlinear dissipation in the system is due to the simple support axles that oscillate in lubricated bearings (Sid Harvey F2-40 oil; McMaster 8600N7 bearings). In this study, nonlinear dissipation (drag) of the experimental system is identified by iteratively matching model results to the measurements once linear, viscous damping (beam elastic loss) is accurately determined. Therefore, the drag induced via simple support bearing oscillation [23,24] realizes a phenomenologically comparable energy dissipation (force amplitude proportional to the square of velocity) as that force which resists the vibrations of a rigid wing in air flow. Alongside the oscillating bearing inertia, the experimental system may emulate the salient inertial and damping forces encountered by a flapping, rigid wing without the need for explicit inclusion of a wing attachment during experimentation. The system is attached to an electronically controlled electrodynamic shaker (APS Dynamics 400) and excited in the axis of beam transverse vibrations *w*(*x*, *t*), as indicated in figure 2*b*,*c*. The beam is oriented so as to eliminate appreciable gravitational body forces upon transverse beam motions that could otherwise induce asymmetries in the static loading.

Measured and identified experimental system parameters are provided in table 1. The dynamic data acquired during experimentation are shaker acceleration and beam centre velocity and displacement. Shaker acceleration in the direction of excitation is measured using an accelerometer (PCB 352 C33), whereas a laser interferometer simultaneously captures the transverse velocity and displacement of the beam along the axis of excitation (Polytec OFV 3001 S, OFV 303). With the addition of a rigid wing attachment to a simple support bearing axle (like that in figure 1*b*), one could determine the wing tip displacements using a moment/lever arm relationship between the beam centre and the wing tip; then, having the knowledge of wing tip oscillation, one could calculate the aerodynamic force generated for a given vibration of the wing motor [2]. Experiments are conducted sampling the data at 4096 Hz; the data are filtered using a sixth-order lowpass infinite impulse response filter with cut-off frequency of 100 Hz. All fundamental mode dynamics of the beam considered here (whether pre- or post-buckled) are sufficiently within this bandwidth. Apart from stationary/steady-state experiments, excitation level sweep evaluations are undertaken by programming very slow triangular wave modulations of shaker acceleration (triangular wave frequency 0.003 Hz), such that the system is excited almost quasi-statically over the course of many excitation/flapping cycles. The data are then processed over short durations of time using fast Fourier transforms. According to the linear natural frequency of the beam (approx. 48.2 Hz), the data are thus sufficient for comparison with analytical predictions that assume pure stationary excitation. Finally, the data employed towards the presentation of the following experimental results may be obtained in the electronic supplementary material.

### 3.2. Experimental structure viability

A first assessment of the experimental system efficacy and viability to reasonably represent essential wing stroke dynamic characteristics of Diptera and MAVs is needed. To reiterate from §3.1, the beam centre motions are geometrically related to wing stroke dynamics, and hence aerodynamic forces, through a scalar multiple determined by a moment arm between the beam centre and a wing attachment free tip. Figure 3 presents measurements of the beam motor structure displacements and velocities over several excitation cycles at 26 Hz when the ratio of axial compression to Euler buckling load is 0.95. Four suspension stiffnesses are considered, representing slightly less than two orders of magnitude of difference from the least to most stiff axial spring, as seen in table 1. This is more stiffness variation than has been evident in dipteran flight muscle tensioning, but comparable to that achieved by species of other orders [25].

The displacement response amplitudes, figure 3*a*, and hence aerodynamic forces, are seen to be significantly modulated by adjustment of the axial support characteristics, whereas excitations remain the same. This is consistent with established understanding that Diptera primarily change aerodynamic forces via wing stroke amplitude control [5,15,16]. Additionally, the displacement trajectories appear to be mostly sinusoidal, but close inspection of corresponding velocities indicates a reduction in peak absolute values of velocity at middle points through the half-stroke, labelled as A in figure 3*b*. This is evident in figure 3*b*,*c* by comparison of two examples of the beam structure dynamics with an oscillator undergoing simple sinusoidal motions of the same displacement amplitude (grey long-dashed curves). Comparable to the aerodynamic drag force acting on the flapping wings of Diptera which suppresses the peak absolute velocities during a cycle of flapping [26,27], the drag dissipation of the lubricated simple support bearings of the experimental system also reduces the peak velocity of oscillation. It was reported that the dipterans *Eristalis* and *Episyrphus* exhibited approximately a 1.35% reduction in root-mean-square (RMS) angular velocity when compared with simple harmonic motion of the same stroke amplitude [26]; in comparison, the experimental structure with the softest axial suspension stiffness exhibits a 1.99% reduction in RMS beam centre velocity with respect to the sinusoidal motion. As a result, Ellington [26] observed that for dipteran hovering the wing ‘accelerations and decelerations at the ends of the [stroke] cycle are consequently greater’. Indeed, this is phenomenologically emulated in the experimental system: the structural system yields a 21.2% greater peak absolute value of acceleration than the simple sinusoidal motions at the ends of the stroke cycle (91.3 versus 75.3 m s^{−2}; corresponding position labelled as B in figure 3 for the system employing the softest axial suspension stiffness). These are characteristics observed not only in the biological data [26,27], but also more generally in successful MAV platforms [28]. Therefore, the experimental system faithfully emulates the salient dynamic characteristics of wing flapping, even without a wing attachment, which underscores that the inertial and dissipation characteristics introduced by the lubricated simple support bearings and axles accounts for similar dynamic influences as a rigid, flapping wing extension. For a first assessment of the experimental system, the above-mentioned results show the motor structure demonstrates reasonable efficacy and viability to be considered in the following as an experimental platform for evaluation of the numerous roles played by adapting axial support characteristics.

## 4. Flexible motor modelling

### 4.1. Dynamic governing equation of motor structure

With the representative wing motor architecture, figure 1*b*, we can model and analytically evaluate the structure using mathematical tools that yield direct prediction of excited motor dynamic amplitudes. When considered alongside data acquired using the experimental system, the aim is to develop more rigorous conclusions as to the respective influences of compressed motor characteristics, namely axial support stiffness and compression distance. Such predictive methods thus serve as the basis for design strategies for inspired MAV platforms and introduce a new approach from which to formulate insight regarding the dipteran muscle-to-wing interface.

The governing equation of motion for the motor structure, figure 1*b*, is derived using energy principles. The beam is excited by harmonic motions of the surrounding frame *ä*(*t*) = *A*_{a}cos*Ωt*, where *A*_{a} is acceleration amplitude and *Ω* is excitation frequency; this excitation is equivalent to a distributed modal force as derived using the present model formulation. Considering that beam transverse *w*(*x*, *t*) and axial *u*(*x*, *t*) displacements are functions of the beam length coordinate *x* and time *t*, the kinetic energy *T*, the potential energy *U* and the dissipation function *D* may be expressed as [20,21]4.14.24.3Operators (·) and ( )* _{x}* indicate differentiation with respect to time and to beam length coordinate, respectively. Other parameters are as follows, where subscript ‘b’ indicates the term is related to the beam: volumetric density

*ρ*

_{b}; area

*A*

_{b}; natural beam length

*L*

_{b}; rotational mass/axle radius

*r*

_{a}; Young's modulus

*E*

_{b}; moment of inertia

*I*

_{b}; viscous damping per beam length

*c*

_{b}. The drag damping force (stemming from the lubricated simple support bearings) is introduced after reduction of the continuum model using the following procedure. To better facilitate model predictive capacity by eliminating asymmetries, static body forces are neglected.

The excited, periodic response of the beam is primarily harmonic at the same frequency of the excitation. The linear vibration modes are used as spatial expansion functions4.4and4.5where are generalized coordinates of the axial and transverse responses. Note that is related to the amplitude of wing tip motions through a moment arm expression and is therefore related to aerodynamic forces. Following substitution of equations (4.4) and (4.5) into (4.1)–(4.3) and employing Lagrange's equations4.6where *L* = *T* − *U*, leads to direct determination of coefficients 4.7and4.8where the axial spring stiffness ratio is *κ* = *k*_{d}*L*_{b}/*E*_{b}*A*_{b} which is the ratio of the axial spring to beam axial stiffnesses. Non-dimensional motor displacement and time *τ* = *ω*_{o}*t* are introduced, where *ω*_{o} is the uncompressed beam fundamental natural frequency expressed in equation (A 1). Continuing application of Lagrange's equations to coordinate leads to the governing equation for the transverse beam response4.9where key terms include *P*_{cr}, critical Euler buckling load; *z*, non-dimensional excitation level; *ω*, non-dimensional excitation frequency. The operator ( )′ indicates differentiation with respect to non-dimensional time *τ*. The remaining terms and further details for the prior are given in appendix A. Note that by the model formulation, the base acceleration is equivalent to a generalized force excitation, indicating the one-to-one correspondence between the two excitation forms. The final derivation step for this study is to introduce nonlinear damping related to the damping force in the lubricated simple support bearings4.10where *η*_{a} is an aerodynamic (drag) dissipation factor identified from experimental data. The viscous loss factor *η* is initially identified by removing the beam from the fixture, clamping one end to a fixed support, and conducting ‘ring-down’ tests using an initial displacement at the other, free end of the beam; the loss factor is then determined from the measurements using the logarithmic decrement approach [29]. Then, the deviation between amplitudes predicted by the model and those as measured in experimentation (for moderate levels of excitation) is caused by the nonlinear damping related to the simple support bearing oscillations. Through iteratively updating modelled results with respect to the measurements, the equivalent aerodynamic (drag) dissipation factor *η*_{a} is approximately identified. Instead of relating this term to an additional ‘equivalent viscous’ damping factor [24], the nonlinear projection of the phenomenon is retained in the model.

Note that coefficients *β*, *γ* and *κ* are non-negative. Thus, it is observed that the nonlinear governing equation for the wing motor, equation (4.10), may have negative linear stiffness when *k*_{d}Δ/[*P*_{cr}(1 + *κ*)] > 1. In this condition, the beam is compressed beyond the threshold of elastic stability and exhibits bistability. The opposite inequality indicates the structure is monostable, but still nonlinear. In the absence of compression, Δ = 0, and for sufficiently small displacements, the dynamics of the beam are approximately linear. The term *k*_{d}Δ/*P*_{cr} is the axial pre-load ratio, defined as the ratio of pre-loading force to the fundamental critical buckling load. Depending on the axial spring stiffness ratio, *κ* ∝ *k*_{d}, the effect of the axial spring may be to support an axial compressive load greater than the critical Euler load *P*_{cr}. Such are the static influences of the axial stiffness and compression. Yet, dynamic influences are what concern aerodynamic forces and the dynamic effects of axial characteristics are our interest in the following.

### 4.2. Validation and accuracy of model formulation

Prior to developing a solution strategy for the dynamic nonlinear governing equation (4.10), the model is directly validated with experimental results. To this end, we numerically integrate equation (4.10) using a fourth-order Runge–Kutta algorithm and compare model simulations alongside the measurements. The top row of figure 4 presents the numerical results, whereas the corresponding experimental data are shown in the bottom row. The system is excited with frequency 26 Hz and amplitude 30.4 mN, where the excitation force that drives the fundamental mode of the structure, *F*_{m}, is related to the accelerations of the shaker by *F*_{m} = 2*ρ*_{b}*A*_{b}*L*_{b}*A*_{a}/*π*. In all cases presented, the beam is compressed using an axial pre-load ratio of *k*_{d}Δ/*P*_{cr} = 1.10, whereas the spring stiffness *k*_{d} changes for each trajectory shown; thus, to maintain the fixed pre-load ratio, the compression distance Δ is modified from case to case to counterbalance change in stiffness *k*_{d}. Figure 4 shows the predicted time series and phase portrait trajectories are in good qualitative agreement with the data; additionally, there is overall quantitative agreement between the measurements and simulations in terms of amplitude and phase differences among the cases. The same trends regarding axial stiffness influence, for constant pre-load ratio, are evident from the modelling as has been identified from the experimental system. The more nuanced variation of velocity over time observed in the experimental data owing to minor asymmetry is not exhibited by the model which presumes a symmetric motor architecture. Nevertheless, the overall agreement between model and experiment validates the theoretical formulation as a viable foundation for the analytical solution approach to be hereafter applied.

### 4.3. Analytical solution strategy to predict wing motor dynamics

For harmonic excitation of the motor, a suitable approximation of the response *g*, governed by equation (4.10), is a Fourier-series expansion [21]. Given that wing motions are primarily harmonic with respect to the frequency of excitation (single periodic motions), a one-term Fourier expansion is employed4.11

Substitution of equation (4.11) into (4.10), eliminating higher-order harmonic terms, and assuming slowly varying coefficients leads to an equation system given by4.12*a*4.12*b*4.12*c*where *A* = 1 − *k*_{d}Δ/[*P*_{cr}(*κ* + 1)], *B* = *β* + *γ**κ*/(*κ* + 1), and reduction of the nonlinear damping follows [30]. Assuming steady-state response, the constant wing motor displacement bias *g*_{0} may be solved in terms of the dynamic amplitude *g*_{1} using equation (4.12*a*). It is then found that either or . Depending on the compression of the motor as quantified by the linear stiffness term *A*, one or three response solutions are determined from the solution to (4.12*a*). For example, if *A* is positive, then the steady-state solution to equation (4.12*a*) is physically meaningful only for : in other words, the system is monostable, because the motor is not compressed beyond the limit of elastic stability. If the linear stiffness term *A* is negative, then the motor is bistable, and all three equilibria exist , where the case of *g*_{0} = 0 is then an unstable configuration. Substituting a selection of *g*_{0} into equations (4.12*b*,*c*), the two equations are squared and summed to yield a third-order polynomial equation in terms of . The roots of the polynomial are computed, and the coefficients *g*_{0} and *ψ* are then determined.

Not all mathematically predicted response coefficient sets **a** = [*g*_{0}, *g*_{1}, *ψ*]^{T} are stable and physically meaningful (some roots may be or be complex valued). Stability may be checked by considering equation (4.12) expressed using **Pa′** = **Ga**, where **P** and **G** are matrices determined accordingly. Given a solution set **a***, stable responses are those for which the real components of the eigenvalues of the Jacobian **J** = [*∂*(**P**^{−1}**G**)/*∂***a**]|** _{a}_{=}_{a*}** are all negative [31].

The stable responses predicted from this approach, for a given system and excitation parameter set, represent the steady-state beam (wing motor) centre point displacement amplitudes. Note that by this solution method, and assuming a rigid wing is attached to the rotational support of the wing motor (as in figure 1*b*), wing tip velocities are directly related to the beam centre velocity through a moment arm relationship and steady-state flapping forces are then proportional to the square of the wing tip velocity amplitudes [2].

## 5. Flexible and compressed motor suspension enhances flapping versatility and effectiveness

The first results shown in figures 3 and 4 reveal a principal influence of adjusting axial suspension stiffness *k*_{d} and compression distance Δ, in those cases when maintaining constant pre-load ratio *k*_{d}Δ/*P*_{cr}. It is seen that axial stiffness reduction can greatly amplify motor displacement and velocity amplitudes, and hence aerodynamic force production, whereas excitation level and frequency remain unchanged. Biologically, there are limits imposed on Diptera to adjust the thorax resonance frequency (i.e. flapping frequency [2,32]) to modify wing stroke amplitude, so other means likely play important roles to tailor aerodynamic force [15]. For engineered systems, there are advantages to applying minor structural controls to yield dramatic stroke amplitude adjustments which are likewise less easily realized by modifying flapping frequency [33]. The following sections therefore investigate the influences of the flexible axial support and its compression characteristics in determining wing motor dynamic responses critical to a versatile and effective aerodynamic force generation. The parameters explored here include axial stiffness ratio *κ* ∝ *k*_{d}, axial compression pre-load ratio *k*_{d}Δ/*P*_{cr}, normalized excitation frequency *ω* and excitation level *z*. Note that as we explore influence in changing axial pre-load ratio while *κ* remains constant, this indicates change in compression distance Δ according to the definition of pre-load ratio.

### 5.1. Motor excitation level

Axial support influences on motor structure velocity amplitudes are now investigated as the generalized excitation force amplitude is varied. The experimental sweeping strategy described in §3.1 is used. The fundamental, undamped, linear natural frequency of the uncompressed beam is approximately 48.2 Hz, and, prior to buckling, the linearized resonance shifts downward in consequence to axial compression of the motor structure [21]; after buckling, the linearized resonance increases. In the following studies, the ratios of excitation frequency to the undamped, uncompressed natural frequency, *ω* = *Ω*/*ω*_{0}, are examined for the cases of *ω* = 0.33 or 0.54, although the ratios with respect to the undamped, *linearized* resonance frequencies after taking into account axial compressions considered here are close to 1. This frequency selection complies with the recognition that dipteran wing motors are operated at a frequency near to the wing–thorax system linearized resonance [2,32].

Figure 5 presents the analytically predicted (left column) and experimentally measured (right column) results of beam centre velocity amplitudes as excitation force changes in the horizontal axis for an excitation frequency of 16 Hz, representing frequency ratio *ω* = 0.33. From figure 5*a*–*d* in the left column and *e*–*h* in the right column, the suspension stiffness ratio *κ* decreases from a nearly fixed constraint (*κ* = 0.3) to a much more flexible suspension (*κ* = 0.0053). Solid curves in the left column denote symmetric dynamics predicted by the model (which may be either bistable dynamics that cross the unstable configuration or monostable oscillations), whereas the dashed curves indicate asymmetric responses specifically induced by bistability. Provided throughout figure 5 is the influence of changing pre-load ratio *k*_{d}Δ/*P*_{cr} from 0.75 to 0.95 to 1.10 indicated by the black, green (light grey) and red (medium grey) curves or data points; arrows are provided as visual aids to follow the trend of increasing pre-load ratio. In the left column, the predicted cases of zero pre-load (Δ = 0) are shown as dotted curves.

When axial compression is absent from motor design, Δ = 0, the motor dynamics are comparatively negligible across all excitation force levels. This indicates that motor compression is a critical first step to enable broad adjustability of flight mechanism dynamics. Regardless of the axial support stiffness *k*_{d}, a lack of static compression distance effectively eliminates opportunities to effect a useful ‘transformation’ of energy from the input excitation to output wing dynamics by varying stiffness. In agreement with numerical and experimental results presented in earlier sections, the reduction in axial stiffness, *κ* ∝ *k*_{d}, leads to amplification of the velocity amplitudes, as observed by comparing results for any given excitation level and pre-load ratio, and then viewing the plots in a column from top to bottom. The higher pre-load ratios consistently induce larger amplitude dynamics. Significant variation in motor velocity is evident by comparing the case of *k*_{d}Δ/*P*_{cr} = 1.10 (red/medium grey curves or circle data points) as axial support stiffness ratio decreases *κ* = 0.3–0.0053, top row and bottom row, respectively. These results suggest that the elimination of axial stiffness *k*_{d} → 0 is one means to maximize dynamic amplitudes. On the other hand, the feasibility to infinitely compress the motor, Δ → ∞, to maintain the same pre-load ratio *k*_{d}Δ/*P*_{cr} poses obvious practical challenges. Nevertheless, figure 5 shows that even an order of magnitude stiffness change, such as between the second and fourth rows of the figure, indicate flapping motions transmitted through a compressed wing motor will be significantly modulated in amplitude. For example, at 30 mN, the experimentally measured amplification is approximately 3.6 times for pre-load ratio 0.95, shown by triangle data points in figure 5*f*,*h*.

The results also indicate a distinctive characteristic of nonlinear systems, namely coexistence of multiple dynamic forms for the same excitation conditions. For example, this is evident in the left column of figure 5 where asymmetric (dashed curves) and symmetric (solid curves) motor motions are predicted for the same excitation level. Coexistent low and high amplitude symmetric dynamics are also predicted in some cases, which represent two dynamic forms that are out of phase. Some experiments for this frequency of excitation observed coexistence in the form of *nonlinear hysteresis*, for example, in figure 5*h* around excitation force 10 mN, owing to near quasi-static up and down sweeping strategy for excitation level. Thus, depending on the initial configuration (i.e. initial displacement and velocity of the system), the lower or higher amplitude dynamic can be realized. Therefore, although modulating suspension and compression of the wing motor may enable large change in resulting wing stroke dynamic amplitudes, concurrent attention to the resulting dynamic form (i.e. symmetric, asymmetric and the phase) is required to ensure that modulation of stiffness *k*_{d} or compression distance Δ induces the desired variation in aerodynamic force. Such coexistence may strictly be a matter of implementation related to engineered compressed wing motors, however, because it may be anticipated that the higher damping [2] and degree of muscular control [5] exhibited in a compressed dipteran motor architecture would sufficiently alleviate concerns of realizing one dynamic over another. A final observation relates to the multiple dynamic forms, whether coexistent or individually realized, in the context of engineered MAV flight mechanisms. Transitioning among the forms suggests means to vary wing stroke amplitude and potentially stroke symmetry for the same axial suspension characteristics. Varying the stroke amplitude influences aerodynamic force, whereas a deviation from symmetric oscillations indicates capacity for steering. The versatility of changing flapping forms and amplitudes, whereas all other design parameters and exciting conditions remain the same, has merit for MAV applications.

To more clearly explain the characteristic of coexistent dynamics and influence of changing the pre-load ratio, figure 6 presents results of numerically integrating equation (4.10) for the softest axial stiffness ratio considered here, *κ* = 0.0053, when the wing motor is excited with 5.82 mN of force amplitude at a frequency of 16 Hz. This operating state corresponds to that highlighted by the shading in figure 5*d*. Uncompressed (Δ = 0) wing motor dynamics are omitted from figure 6. The time series of the motor displacements, figure 6*a*, indicates that depending on the pre-load ratio *k*_{d}Δ/*P*_{cr}, dynamics of very different amplitudes may be activated; this was previously evident in figure 5*d*. On the other hand, the time series, figure 6*a*, and displacement–velocity trajectories, figure 6*b*, more clearly illustrate the impact of asymmetry if it should occur. For the pre-load ratio of 1.10 (dashed curves in figure 6), the motor structure is bistable, and the level of excitation is insufficient to excite large amplitude symmetric motions. Instead, two asymmetric dynamics are possible for the same excitation. It is important to note that the asymmetric motions (in this wing motor architecture) may exhibit either a positive or negative offset value; the results in figure 6 are merely one representative case when the low amplitude asymmetric dynamic has a negative displacement offset, whereas the larger amplitude response exhibits a positive value of offset. In the context of a structural motor with flapping wing extension, the offset induced by asymmetric flapping motions would induce a torque upon the hovering body, for example, a roll or pitch depending upon whether one or both, respectively, of the wings flapped asymmetrically. Although Diptera employ more intricate muscular means to contort the thorax and wing motor for purposes of flight control [5–7], the utilization of asymmetric flapping dynamics activated via the flexible and compressed wing motor structure presents one novel opportunity to simplify steering strategies for MAVs.

The results of investigations conducted using excitation/flapping frequency 26 Hz (frequency ratio *ω* = 0.54) are presented in figure 7 which applies the same plotting conventions as employed in figure 5. At this frequency, there is likewise significant amplification of the velocity responses by reduction of axial stiffness and corresponding adjustment of compression distance. More apparent, however, is the increased likelihood for coexisting dynamic forms, observed plainly in both the analytical and experimental results. When compared with excitation at 16 Hz, the range of excitation levels for which multiple dynamics may coexist is much broader when the motor structure oscillates at 26 Hz. The importance of these findings is related to the concept of domains or basins of attraction which is a means to characterize the likelihood that certain initial system configurations lead to one or other dynamic response form [31]. As one of the coexistent forms approaches a bifurcation owing to independent variable change, for example, in figure 7*d* reduction in excitation force to approximately 26 mN for the high amplitude dynamic with approximately 500 mm s^{−1}, the likelihood of realizing that dynamic reduces for the same initial condition space [34]. Consequently, the broad coexistence range observed in the analytical and experimental results of figure 7 in the third and fourth rows for higher pre-load ratios are evidence of greater challenge to consistently realize one preferred flapping motion and amplitude over the other, particularly if the system begins vibrating/flapping from rest (i.e. no excitation force). As earlier described, this result may find greater meaning in the context of engineered compressed motors, and is primarily a consequence of the low damping evaluated in the experimental motor structure which more closely represents fabricated MAV motor architectures than a dipteran motor.

Finally, figures 5 and 7 indicate an overall good agreement in the qualitative trends and quantitative values of measured and analytically predicted velocity responses of the representative wing motor structure. The non-dimensionality of parameters explored open the interpretation of results to more than engineered flight mechanisms. In addition, the findings underscore the critical importance of axial support stiffness and compression characteristics to provide broad variation in flapping types and amplitudes without otherwise modifying flapping frequency or the driving force.

### 5.2. Aerodynamic effectiveness of compressed motors

The model is then explored to investigate the effectiveness of compressed motor designs in bearing a given mass. The metrics of supported mass and mass-specific aerodynamic power are employed, to first evaluate operating conditions of motor designs which carry the same mass and second to find the power requirements to operate the motors under such conditions. For hovering flight, the supported mass and mass-specific aerodynamic power (which is the sum of induced and profile powers) are given in equations (5.1) and (5.2), respectively [2].5.1and5.2where is peak-to-peak wingbeat amplitude in radians, *n* is wingbeat frequency in hertz, *R* is wing length in metres, *C*_{L} is a mean lift coefficient, *Λ* is aspect ratio, *C*_{D,pro} = 7 *Re*^{−1/2} is steady-state profile drag and *Re* = 4*ΦnR*^{2}/*vΛ* is the Reynolds number for kinematic viscosity *v*. The leading coefficient values are due to assumption of flapping in air with conventional gravitational acceleration. In agreement with Ellington's derivation and utilization of equations (5.1) and (5.2) [2], we assume the kinetic energy involved in driving the wing motor (the shaker oscillations) is negligible owing to an elastic energy storage and release in the motor itself (strain energy in the beam); therefore, only the aerodynamic power needed to maintain the system vibrations is considered towards the net power required to hover with a given supported mass. Because the choice of the dimensioned parameters influence only the absolute values of the results, when compared with changing the dynamic amplitude trends, we arbitrarily select the dimensional parameters for the following. Thus, *v* = 15.68 × 10^{−6} m^{2} s^{−1}, *R* = 4*L*_{b}, *Λ* = 7, *C*_{L} = 2.5. In the following, the motor damping constant is selected to be closer to that reported for *Drosophila* [2], *η* = 0.24, than the value employed in the prior sections which was specific for the experimental motor structure. The remaining values required to compute equations (5.1) and (5.2) are determined from our model: excitation/flapping frequency and wingbeat amplitude which is related to motor displacement amplitude through the moment arm equation indicated above.

Figure 8 presents two sets of results, where the left column shows the flapping frequency required to support a given mass and the right column reports the aerodynamic power necessary to sustain that frequency. The top row of figure 8 presents the results for pre-load ratio *k*_{d}Δ/*P*_{cr} = 0.95 and excitation force 50 mN; the bottom row plots results for pre-load ratio *k*_{d}Δ/*P*_{cr} = 1.10 and excitation force 200 mN. For the curves, an increasing lightness of plotted curve shading (from black to light grey) indicates a two orders of magnitude increase in the axial stiffness ratio *κ* = [0.001, 0.1]. Symmetric and asymmetric dynamics are indicated by solid and dashed curves, respectively. The corresponding results for an uncompressed, linear motor are shown as the thick dotted curves.

To understand the relative trends of bearing a given mass and the corresponding power requirements, a few examples are beneficial. The circle/square symbols in figure 8*a*,*b* provide a comparison of two compressed motors; the motors are both monostable, because the pre-load does not induce a negative linear stiffness and hence bistability. It is seen in figure 8*a* that when bearing the same 16 mg mass, the motor with axial stiffness ratio *κ* = 0.001 (circle) requires an excitation/flapping frequency of 10 Hz, whereas the motor with *κ* = 0.1 (square) flaps at 30 Hz. Then, in figure 8*b*, by identifying the location on the plots corresponding to those stiffness ratios and respective flapping frequencies, it is observed that to sustain this same 16 mg mass, the motor with axial stiffness ratio *κ* = 0.1 requires more than 70% greater power than the motor with softer axial suspension. In agreement with earlier findings regarding velocity amplitude enhancement, the efficiency benefits of certain softer, compressed axial suspensions are clear. Also evident in figure 8*a* is the fact that for this excitation force, the uncompressed linear motor is not capable of bearing a 16 mg load. From this example, we find that motor compression leads to improved aerodynamic performance in terms of potentially bearing a greater mass and also relative power savings based on the axial stiffness employed.

In figure 8*c*,*d*, a greater pre-load ratio *k*_{d}Δ/*P*_{cr} = 1.10 and excitation level 200 mN are considered. For this pre-load, all of the motors are bistable with the exception of *κ* = 0.1 which has zero linear stiffness and is hence ‘essentially nonlinear’ (i.e. non-linearizable). The light-shaded circle/square symbols compare two cases of the compressed axial motors bearing a 250 mg mass, when the motors have approximately one order of magnitude difference in axial support stiffness. It is again apparent that the softer suspension (bistable) requires less power to operate in bearing the mass at approximately 33 Hz flapping frequency when compared with the stiffer axial suspension *κ* = 0.1 that flaps at approximately 57 Hz. This and the above findings indicate that axial motor supports having greater compliance may enhance aerodynamic performance whether the wing motor is mono- or bistable following the compression. Figure 8*c*,*d* also compares the results for the compressed motor having axial stiffness ratio *κ* ≈ 0.002 (dark-shaded circle) and the uncompressed linear motor operated at resonance (triangle). At resonance, the linear motor is shown to support a mass of approximately 130 mg flapping at a frequency of approximately 48 Hz, whereas the compressed motor will flap at approximately 14 Hz to support the same mass. But the disparity between aerodynamic power requirements between the motor designs is striking, figure 8*d*: the uncompressed linear motor requires almost twice the power as the compressed motor. Taking into account the collective experimental and analytical evidence, wing motor architectures that modulate axial suspension stiffness and compression characteristics may significantly enhance the range of flapping amplitudes and types that may be realized and can exploit such nonlinear dynamics for energy efficiency improvements. Whether developing flight mechanisms for MAV applications or considering the dipteran wing motor in new light, the above results show the previously unexplored aspects of axial stiffness and compression in the flapping wing motor are strongly related to the efficiency and functionality of the system.

## 6. Concluding remarks

Realizing small-scale flapping vehicles is an endeavour that is well advised to take inspired cues from nature's successful flyers. Through analyses and experimentation, a representative dipteran wing motor architecture explored in this study shows that axial stiffness and compression of the flight mechanism are key to achieving a broad range of wing stroke dynamics and improving the load-bearing effectiveness of a given design. By varying the motor-supporting axial stiffness and/or the compression thereof, a versatile collection of flapping motions is achievable to adjust aerodynamic force as well as to potentially steer. It is also found that tailoring axial suspension characteristics can reduce motor power requirements to hover with the same load. These conclusions are shown to apply whether the motor is compressed beyond the buckling threshold, representative of the bistable click mechanism, or if the compressed motor remains monostable. These results have a clear connection to the engineering context of MAV development, and in addition, the biological interpretations of the findings may enlighten new perspectives on variably stiffened and compressed dipteran wing motors as a potential explanation for the historically intriguing transmission of energy ‘from the engine to the wheels’.

## Funding statement

This research is supported by the Air Force Office of Scientific Research (FA9550-13-1-0122) under the administration of Dr David Stargel and by the University of Michigan Collegiate Professorship.

## Appendix A

Additional variables include *A*_{b} = *h*_{b}*t*_{b}, beam area; and , beam moment of inertia.A 1A 2A 3A 4A 5A 6A 7A 8

- Received December 13, 2014.
- Accepted January 2, 2015.

- © 2015 The Author(s) Published by the Royal Society. All rights reserved.