Aerodynamic investigation on shifted-back vertical stroke plane of flapping wing in forward flight

In order to properly understand aerodynamic characteristics in a flapping wing in forward flight, additional aerodynamic parameters apart from those in hover—an inclined stroke plane, a shifted-back stroke plane, and an advance ratio—must be comprehended in advance. This paper deals with the aerodynamic characteristics of a flapping wing in a shifted-back vertical stroke plane in freestream. A scaled-up robotic arm in a water towing tank was used to collect time-varying forces of a model flapping wing, and a semi-empirical quasi-steady aerodynamic model, which can decompose the forces into steady, quasi-steady, and unsteady components, was used to estimate the forces of the model flapping wing. It was found that the shifted-back stroke plane left a part of freestream as a non-perpendicular component, giving rise to a time-course change in the aerodynamic forces during the stroke. This also brought out two quasi-steady components (rotational and added-mass forces) apart from the steady one, even the wing moved with a constant stroke velocity. The aerodynamic model underestimated the actual forces of the model flapping wing even it can cover the increasingly distributed angle of attack of the vertical stroke plane with a blade element theory. The locations of the centers of pressure suggested a greater pressure gradient and an elongated leading-edge vortex along a wingspan than that of the estimation, which may explain the higher actual force in forward flight.


Introduction
Biological flyers voluntarily create a vortex-dominant flow with their flapping wings to gain sufficient lift and control forces [1][2][3][4]. Most flyers adept at hovering, for example, mainly use a large sweeping motion with a high angle of attack in a nearly horizontal stroke plane. This forms a stable leading-edge vortex (LEV) on their wings which augments a lift force [5]. Moderate but shorter wingspan than that of a man-made aircraft, which may have been evolved for their flight performance, assists the stable LEV attachment [6]. Wing pitching-up at dorsal and ventral stroke reversals also brings out an extra lift as a result of the additional circulation [7]. Accordingly, they have acrobatic maneuverability that has not been reached by a man-made aircraft.
The LEV attachment of their flapping wings is one key phenomenon explaining the lift augmentation of the flyers adept at hovering. In contrast to an airfoil in a high angle of attack, where the vortices are shed at each edge similar to blunt bodies, the sweeping motion of the wings pulls the flow over the leading edge and builds the LEV on the suction side of the wing. This LEV then reattaches the flow over the LEV and settles the shed vortex at the trailing edge. The Kutta condition is satisfied, and the lift is augmented. The LEV itself provides a lower pressure area near the leading edge, also enhancing the lift. Ellington [5] found the LEV on a wing of a hawkmoth (Re c ∼ 6000), and Dickinson et al [7] found it at a lower Reynolds number (Re c ∼ 120) with a scaled-up robotic model of a fruit fly. Tremendous follow-up studies revealed that the LEV of flapping wings is stabilized by radial pressure gradient along a wingspan, which is originally delivered by appropriate centripetal and Coriolis forces that are driven by the sweeping motion of the wings [8,9].  [15,16]).

Du
Dt Equation (1) is the dimensionless Navier-Stokes equation scaled for the wingbeat motion of biological flyers, in which r and u are the position and velocity of the fluid volume, and Ω is the angular velocity of the rotating frame. J, A * , AR, Eu, and Re are the dimensionless numbers of the advance ratio, dimensionless stroke amplitude, aspect ratio, Euler number, and Reynolds number, respectively (refer to Lentink and Dickinson [8] for details). As shown in equation (1), the centripetal and Coriolis forces are the functions of AR and J, indicating that the AR and J determine the LEV characteristics. Indeed, several studies have successfully revealed the effect of AR on the LEV characteristics [10][11][12][13][14]. Kruyt et al [12] investigated the LEV at seven different ARs and revealed the existence of a radial limit of LEV stall along a wingspan, which is approximately four times of chord length. Han et al [13] examined a watertank test with a scaled-up robotic model, and showed that the wing of AR = 3.0 outperforms the others with the largest lift and reasonable efficiency at Re ∼ 10 4 . They also showed that an excessive spanwise flux at AR < 3.0 rather degrades lift production. Jardin and Colonius [14] found an optimal AR similar to the study of Han et al [13]. They further found a lifted-off LEV from the surface at approximately three times a chord. Complicated wingbeat kinematics, undulating surroundings, and unsteady wake vortex dynamics of the wings of biological flyers, however, have confined the studies to that in hover with J = 0, despite the fact that staying particular point in space is very rare for most biological flyers.
In forward flight, two additional kinematic parameters, i.e. an inclined stroke plane angle χ and a shifted-back stroke plane angle β, must be considered apart from those in hover [15][16][17]. The J, which denotes the ratio between the wingbeat speed and the freestream velocity, also should be seriously taken into account because it governs the two forces along a wingspan as aforementioned.
The χ creates because biological flyers tilt a force vector to gain a thrust. This would complicate analyses in most cases, because it inserts the angle in between the freestream and the sweeping direction thereby introducing different angles of attack along the wingspan, in contrast to that in hover. The biological flyers also shift the stroke plane backward when they move forward, similar to a sweptback wing of a man-made aircraft, creating the β as shown in figure 1. This could help to stabilize the LEV with an additional spanwise flow, thereby compensating for the loss of force production caused by an increase in J. This even could provide longitudinal flight stability in forward flight, because shifting back a stroke plane will directly send a point of action of the aerodynamic force behind the center of gravity, as to how we give stability to a fixed-wing aircraft (refer to figure 1 to see the location of the center of gravity of a hawkmoth as an example [15]). In terms of J, the J plays the same role as the AR on LEV characteristics as shown in equation (1). One previous study showed that an increase in J destabilizes an LEV similar to that on higher AR wings, but still needs further studies, as it was analyzed only on the one model wing moving in a horizontal stroke plane (Zimmerman wing planform at AR = 3.0 [18]). Here, decomposing aerodynamic forces into steady, quasi-steady, and unsteady components with an aerodynamic model will definitely help to understand the LEV, wake vortex dynamics, and the overall aerodynamic characteristics of a flapping wing in forward flight.
In this paper, we investigated the aerodynamic characteristics of a flapping wing in a shifted-back vertical stroke plane in a freestream. A scaled-up robotic arm in a water towing tank was used to collect the forces on a model flapping wing, and a semi-empirical aerodynamic model was employed to extract steady and quasi-steady force components of the model flapping wing. The findings in this paper will be a stepping stone for future studies on an effect of an AR, Rossby number, and other variables of a flapping wing in forward flight as well as of the J and/or other variables with a freestream.

Material and methods
All the apparatuses in this study had been introduced in the author's papers published elsewhere (see [19] as an example); some key statues are provided here.  The model wing was mounted on an output shaft of this arm. The wing was made of a 2 mm-thick acrylic flat plate and has a spanwise length b of 250 mm and an AR of 3.0. A Zimmerman wing planform, which was formed by two elliptical curves intersecting at a quarter chordline, was adopted to the model wing. A flapping axis was located at 0.2b from the wingroot. The distance to the wingtip R from the pivot was 1.2b and a corresponding dimensionless second-moment wing arear 2 was approximately 0.57. A six-axis force/torque sensor (Nano17 IP68, ATI) was mounted on the wingroot as shown in a red part in figure 2(a). The maximum force and torque ranges of the sensor were ±12 N and ±120 N mm, respectively. The model wing created a vertical sweeping motion about X B -axis within each constant β, which was given about Z B -axis, as shown in figure 2(b). A wing pitching motion was not applied and the χ was fixed as 90 deg for the vertical stroke plane; this is to avoid the cases with non-perpendicular χ to the freestream, which could blend the effect of the β by adding another slope in the inflow. Figure 2(c) shows a towing system on the water tank. This provided two freestream velocities of 5.3 and 9.1 cm s −1 to the model wing; the corresponding J were 0.58 and 1.0 (as discussed later).

Figures 3(a) and (b)
show the wingbeat motion profile φ(t) and corresponding angles of attack α, where τ denotes a dimensionless time based on a wingbeat period T, i.e. τ = t/T. A near-triangular motion profile for φ(t) with a small decel-acceleration duration Δτ φ of 0.12 gave a constant stroke velocity in a translational phase (refer to [20] for the decalacceleration duration). Under conditions of a wingbeat frequency of 0.2 Hz and a stroke amplitude φ 0 of 70 deg, the stroke velocityφ in the translational phase was fixed as 30.68 deg s −1 (figure 3(a)). An advance ratio J was defined as equation (2) in this study, and two J of 0.58 and 1.0 were selected as variables. The α at the reference point (figure 2(b)) were approximately 60 deg and 45 deg at J = 0.58 and 1.0 at β = 0 deg, respectively. Because the actual α is the arctangent of a freestream U ∞ and distributed stroke velocitiesφR, i.e. α = tan −1φ R/U ∞ , the α was also distributed along the wingspan as shown in figure 3(b). The wingtip velocity-based Re had a range of 1.22 to 1.44 × 10 4 .
As a way of avoiding an effect of an underdeveloped wake, a single run was constructed of (1) a flight forward with three continuous wingbeats, (2) rewinding, and (3) an intermission of 60 s to settle down the wakes. Each case was composed of 27 runs to eliminate white noise. An ensemble average was also nearly converged within this number of runs (refer to [19] for more details).
The spatial resolutions in the stroke axis were 0.38% (90 ± 0.088 ± 0.25 deg), considering an encoders' resolution and a backlash of the gearbox. The force/torque sensor also has the resolution of 1/320 N and 1/64 N mm, respectively; the uncertainties in the forces and moments were assessed as ±1.30% and ±0.74% at the maximum load. The change in the temperature, which alters the density of the water, was negligible, and the minor uncertainty sources such as a structural vibration from the traverse system and the spanwise bending on the model wing were not counted; these were not observed during the tests.
A quasi-steady aerodynamic model that had been empirically revised with an identical wing planform at similar Re was employed (refer to [21] for example). Based on a blade element theory as shown in figure 4, translational, rotational, and added-mass force components of the aerodynamic model F trans , F rot , and F added were individually calculated in each element, as shown in equations (5)-(7).
The translational force coefficient C T was one that in hover (J = 0) provided by Han et al [21]. This C T reflects the maximum force that the wing can produce with the equally distributed α. The rotational force coefficient C R was assumed as π/4, which is a theoretical value when a pitching axis is located at half chord, reflecting no forced pitching motion in this study [22]. The coefficient of the added mass force C A was also based on a theoretical value, i.e. π/4. The term for the angular acceleration along a pitch axis (the second term in the theoretical added mass description [22]) was neglected as Dickson and Dickinson [23]  suggested. The acceleration of the inflowU inflow consisted only of a component in the X W -axis, meaning that it was considered one perpendicular to the wing surface. As discussed later, these two components of F rot and F added were completely not negligible, even the model wing did not have any rotational motions in entire wingbeat cycles. 25 blade elements in this study were sufficient to converge the estimations.

Results and discussion
For analyses, we used one force coefficient C F , which denotes the force perpendicular to the wing surface about X W -axis. This force is only one dominant component acting on the wing, simply fluctuating with a wingbeat motion. Note that this directly reflects the lift on a body as C L = C F cos φ(t) because of no pitching motion in this study. Also note that an aerodynamic model does not cover viscous forces which act parallel to the surface.

Figures 5(a) and (b)
show the time-course C F in three continuous wingbeat cycles at J = 0.58 and 1.0. The C F at J = 0, which was originally obtained for the aerodynamic model [21], was also given in figure 5 as black dashed lines. Comparing to the C F at J = 0, the model wing at both J was found to slightly outperform that of J = 0 (at β = 0). This implies a favorable effect of freestream as an increase in the inflow velocity, although there was also a huge loss in the α (α at J = 0 was equally distributed as 90 deg and gradually decreased as J increased; refer to figure 3(b)).
All the curves started with a large force peak, which appeared at every end of a stroke reversal. The level of each peak in the two J cases was nearly settled from the second half stroke (τ > 0.5), implying that the wakes in each case had been nearly fully developed within the half stroke (τ < 0.5). Only the first peak showed lower than the others, indicating that the fully developed wake helped to produce aerodynamic forces in contrast to that in hover. This is in line with the previous study on the effect of J with a horizontal stroke plane at J > 0.5 [18].
With respect to β, figure 5 clearly shows timecourse changes in C F in every half wingbeat cycle. While the C F at β = −20 deg slightly decreased after the peaks, the C F at β = 45 deg gradually increased during the stroke, which was common in both J cases. These are comparable to the C F at β = 0 deg that showed no change with time. The peaks also appeared to be influenced by the β. The peaks at β = −20 deg were the highest and the peaks at β = 45 deg had the lowest values, respectively. These suggest that the effect of β appears throughout a wingbeat cycle.
In order to see where the time-course changes in C F come from, we investigated a kinematic relationship in aerodynamic variables with respect to β. If the relationship fails to explain the changes in C F , it can be seen that there is a substantial unsteady effect, which would be due to the wakes in most cases [22]. Figure 6 shows a geometrical relationship between a flapping wing and U ∞ . It was found that the U ∞ is decomposed into two vector components acting perpendicular to the wing surface with β. Geometrically, the spanwise component of U ∞ , i.e. U ∞ sin β (blue arrows), leaves U ∞ sin β sin φ (red arrows) in an entire wingbeat cycle as shown in figure 6. Thus, the inflow velocity at each blade U inflow becomes equation (9). When the wing is in a dihedral position as shown in the upper part of figure 6, this U ∞ sin β sin φ acts in a negative way, i.e. a resultant inflow is reduced by this component. As the wing continues the motion and enters an anhedral position as shown in the bottom part of figure 6, the U ∞ sin β sin φ became beneficial with its opposite direction, i.e. the inflow increases. The level of peaks after each stroke reversal also appear to originate from such a resultant inflow velocity depending on the U ∞ sin β sin φ. This suggests that the time-varying C F in this study mostly stemmed from the time-varying changes in the inflow velocity by the β, most of which the aerodynamic model may cover.  figure 7(c). The β and J here were 45 deg and 1.0, respectively. As shown in figures 7(a) and (b), both the U inflow and the α show a gradual increment during the stroke. This is due to the non-zero β and U ∞ sin β sin φ component as aforementioned. Note that the change in α was not monotonic as shown in equation (2). This always leftα andα components, and brought out the F rot and F added even in the translational phase, where the model wing had a constant stroke velocity (see the enlarged window in figure 7(c)). At β = 45 deg, the F trans , F rot and F added account for 78.1%, 8.8%, and 13.1% of the estimated aerodynamic force production in the quasi-steady model, respectively. Figure 8 compares the force components of the quasi-steady aerodynamic model to the direct measurements at four example β of −20, 0, 20, and 45 deg. The differences between the actual forces and the estimations, i.e. F measured − F trans − F rot − F added , which are shown as purple dotted lines in figure 8, are enlarged two times for better visibility.
In the translation phase, the aerodynamic model seemed to roughly estimate the time-course change in C F . For the quantities, however, the aerodynamic model underestimated the C F , even it included the F rot and F added . The difference clearly appears at J = 0.58, where the wing had a higher α, and slightly decreased at J = 1.0. This was also gradually reduced as β increased (as discussed later).
Note that the C T of the aerodynamic model was built with the same wing planform at an identical Re, and the only major difference from this study is the α along the wingspan, i.e. the same α distribution when the model was built versus an increasingly distributed α in this study [21]. This α distribution may have intensified a radial pressure gradient and strengthened the LEV on the model wing, thereby resulting in the higher actual force in the measurements. Of note is that such an α distribution always appears when a wing moves with an inclined stroke plane, which adds an angle between a freestream and a sweeping direction. This suggests that an inclination of a flapping wing in forward flight is beneficial to stabilize an LEV with such an increasingly distributed α and pressure gradient. This also can directly compensate for an effect of J in forward flight, which attenuates an LEV in general [18].
The peaks during the stroke reversal were not properly matched with those of the measurements in all cases. The actual peaks were lower than the estimations, and developed after a relatively long time from each stroke reversal (approximately 0.1τ later). The  wake of the previous stroke seemed the main source for the peaks in the measurements and the differences, because the wake in this instance would directly impinge the wing surface with high angles of attack, causing vortex shedding at all the edges of the model wing. This would be sufficient to break the assumptions in the quasi-steady concept of the aerodynamic model (e.g. Kutta condition). This is one unsteady phenomenon found in this study, and is in line with peaks in previous studies that used such a short Δτ φ [20,24]. Further investigations on flow structures at this instance are needed. Figure 9 shows mean C F in the translational phase; the range is given as orange color boxes in figure 8 (2.20 < τ < 2.40 and 2.70 < τ < 2.90). Because the measurements include other non-steady components even it moved with the constant stroke velocity as aforementioned, we extracted a translational force component in the measurements (blue triangles) by subtracting the F rot and F added , i.e. F extracted = F measured − F rot − F added . The F extracted then solely contains a steady force component with the LEV; the unsteady wake effect at the stroke reversal would also be excluded in the mean C F with the limited ranges shown in figure 8.
Irrespective of β, all the F extracted were higher than the F trans . This clearly indicates the positive effect of the vertical stroke plane and the increasingly distributed α. The effect was prominent in all β cases at J = 0.58, but at J = 1.0, it was drastically reduced at β > 0 deg. Considering that the wing in this case had a much stronger spanwise flow than those of J = 0.58 due to the higher U ∞ , it can be seen that the excessive spanwise flow, which overly convected the vorticity of an inboard LEV to the wingtip, reduced the coverage of the LEV and the benefit in the lift production, even it can more stabilize the overall LEV system (refer to [25,26] for the structure diagnosis on the LEV). Relatively lower α at J = 1.0, which provided less pressure difference on the upper surface, must have stimulated the convection along a wingspan. This is identical to the LEV behavior of the flapping plates at lower AR in hover (AR < 3) [13], and further implies that the reduction in C F could be recovered by a higher AR wing. Further studies on J, AR, and their correlations in forward flight are highly recommended. Figure 10 shows the locations of centers of pressure CP in a translational phase with respect to β (refer to the orange color boxes in figure 8 for the range of extraction). All the CPs in this study were found to be closer to the wingtip than the CP in hover [18]. This suggests an elongated LEV due to the increasingly distributed α, which could increase the radial pressure gradient on the suction side of the wing and the force. As shown in the right enlarged window, the CPs at J = 0.58 were closer to the trailing edge than that of J = 1.0. This reflects a higher α at lower J. An increase in β moved the CP to the wingtip in both J cases. This implies that the spanwise flow, which gradually increased with β as U ∞ sin β cos φ (refer to figure 6 for details), further elongated the LEV to the wingtip. Note that several studies on a lateral gust found that a spanwise flow reinforces the LEV [19,27,28]; this also supports the slight increment in C F with increasing β at J = 0.58. Such characteristics in CPs are comparable to that of the estimations near the wingtip, which rather approached the wingroot with an increase in β (left enlarged window). This seems mainly due to the limitation of a blade element theory which cannot take additional effects along a wingspan into account; a proper decomposition methodology, which can cover the distributed effect on a wingspan, is highly necessary for a flapping wing in forward flight.

Conclusion
In this paper, effect of a shifted-back vertical stroke plane on a flapping wing in forward flight has been investigated. A semi-empirical quasi-steady aerodynamic model has been used to estimate quasi-steady force components of the model flapping wing, and a scaled-up robotic arm in a water towing tank has been used to obtain an actual force on the flapping wing. It was found that the shifted-back stroke plane left a part of the freestream, resulting in the timecourse change in an inflow and an angle of attack, and gradual growth in the actual force during the stroke. The aerodynamic model showed a presence of the rotational and added-mass force components even the wing moved with a constant stroke velocity due to the shifted-back stroke plane. The vertical stroke plane in this study forms an increasingly distributed angle of attack along a wingspan with a freestream; the aerodynamic model underestimated the actual force, although a blade element theory in the model can take the increasingly distributed angle of attack into account. The centers of pressure, which were located closer to the wingtip, suggested a greater pressure gradient and an enhanced LEV along a wingspan. This may explain the higher actual forces of the model flapping wing in forward flight than that of the estimation.