# Aerodynamics of flapping wings with passive and active deformation

Florian Bouard, Thierry Jardin, Laurent David

► To cite this version:

Florian Bouard, Thierry Jardin, Laurent David. Aerodynamics of flapping wings with passive and active deformation. Journal of Fluids and Structures, 2024, 128, pp.104139. 10.1016/j.jfluidstructs.2024.104139. hal-04757808

HAL Id: hal-04757808

https://isae-ensma.hal.science/hal-04757808v1

Submitted on 4 Jan 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Graphical Abstract

Aerodynamics of flapping wings with passive and active deformation

Florian BOUARD, Thierry JARDIN, Laurent DAVID

<figure><img src="image_2.png" /></figure>
# Highlights

## Aerodynamics of flapping wings with passive and active deformation

Florian BOUARD, Thierry JARDIN, Laurent DAVID

* Passive and active deformations are found to promote aerodynamic performance of three-dimensional hovering flapping wings

* It is shown that performance enhancement results, to leading order, from quasi-steady effects associated with the revolving motion

* The systematic study on actively deformed wings can help determine optimal structural properties of passively deforming wings
# Aerodynamics of flapping wings with passive and active deformation

Florian BOUARD<sup>a,b</sup>, Thierry JARDIN<sup>a</sup>, Laurent DAVID<sup>b</sup>

## Abstract

This paper reports direct numerical simulations of the flow past rigid and flexible flapping wings under hovering flight conditions. Both passive and active deformations are considered. It is shown that passive deformation can help increase aerodynamic performance through significant wing bending. Bending occurs at the frequency of the prescribed flapping motion and is, in this case, characterized by moderate amplitude and phase lag with respect to the prescribed flapping motion. Bending is then actively prescribed (rather than being a result of passive deformation) with varying phase lag. This allows to decouple the role of bending amplitude and phase lag on aerodynamic performance of the flapping wing. It is shown that both lift and efficiency can be significantly enhanced for phase lags around $3\pi/2$ but this enhancement reduces with increasing pitch angle. The influence of morphing on aerodynamic performance can be explained by the concomitant role of quasi-steady and unsteady effects. These results hence demonstrate that morphing can be beneficial to the aerodynamics of flapping wings. Furthermore, they can help define structural properties that promote aerodynamic performance of flapping wings through passive deformations (with relevant amplitude and phase).

**Keywords:** flapping wings, deformable wings, morphing, hovering flight, FSI, DNS

<footnote><sup>a</sup>ISAE-SUPAERO, Université de Toulouse, 10 Av. Edouard Belin, Toulouse cedex 4, 31055, France</footnote>

<footnote><sup>b</sup>Université de Poitiers, ISAE-ENSMA, CNRS UPR 3346, Institut PPrime, 11 Bld. M. et P. Curie, TSA 41123, Poitiers cedex 9, 86073, France</footnote>
<table><thead><tr><th colspan="2">Nomenclature</th><th></th><th></th></tr></thead><tbody><tr><td></td><td>C<sub>L</sub></td><td>Lift coefficient [-]</td><td></td></tr><tr><td>α</td><td>C<sub>P</sub></td><td>Power coefficient [-]</td><td></td></tr><tr><td>α̇</td><td>C<sub>V</sub></td><td>Vertical force coefficient [-]</td><td></td></tr><tr><td>φ̇</td><td>D</td><td>Drag [N]</td><td></td></tr><tr><td>φ<sub>0</sub></td><td>E</td><td>Young modulus [Pa]</td><td></td></tr><tr><td>φ<sub>eff</sub></td><td>f</td><td>Flapping motion frequency [s<sup>-1</sup>]</td><td></td></tr><tr><td>β</td><td>f<sub>n</sub></td><td>First natural solid frequency [s<sup>-1</sup>]</td><td></td></tr><tr><td>ε</td><td>h</td><td>Wing thickness [m]</td><td></td></tr><tr><td>γ</td><td>I</td><td>Cross-section second moment of area [m<sup>4</sup>]</td><td></td></tr><tr><td>ν<sub>s</sub></td><td>L</td><td>Lift [N]</td><td></td></tr><tr><td>ω<sub>z</sub></td><td>Q<sub>P</sub></td><td>Pitching torque [N.m]</td><td></td></tr><tr><td>U&#773;<sub>ref</sub></td><td>Q<sub>R</sub></td><td>Revolution torque [N.m]</td><td></td></tr><tr><td>ρ</td><td>R</td><td>Wing radius [m]</td><td></td></tr><tr><td>ρ<sub>s</sub></td><td>r</td><td>Radial position [m]</td><td></td></tr><tr><td>ξ</td><td>Re</td><td>Reynolds number [-]</td><td></td></tr><tr><td>A</td><td>T</td><td>Flapping period [s]</td><td></td></tr><tr><td>a<sub>1</sub></td><td>t</td><td>Time [s]</td><td></td></tr><tr><td>b<sub>1</sub></td><td>U<sub>eff</sub></td><td>Effective wing tip velocity [m.s<sup>-1</sup>]</td><td></td></tr><tr><td>c</td><td>V<sub>d</sub></td><td>Displacement velocity [m.s<sup>-1</sup>]</td><td></td></tr><tr><td>c<sub>1</sub></td><td></td><td>Phase lag [-]</td><td></td></tr></tbody></table>

# 1. Introduction

Flapping wings, as extensively represented in nature, could help enhance performance of NAVs (Nano Air Vehicles) as an alternative to conventional rotatory and fixed wings. The RoboBee from Harvard Microrobotics Lab (Wood, 2008), with its very small dimensions of 3cm for the wing span and a weigth of 60mg, is one such fascinating example. Over the last three decades, researchers have thus conducted extensive studies on the aerodynamics of flapping wings aimed at understanding insect flight. Most of this research has focused on
rigid wings to characterize mechanisms responsible for lift production (e.g. Dickinson et al., 1999; Sane and Dickinson, 2001, 2002).

Conversely, their flexible counterparts have received much less attention despite their representation in natural species as shown by Wootton (1981) and later by Lucas et al. (2014). In fact, the flight of insects relies on passive deformation of their wings through natural stiffness while birds and bats can use their feathers or their muscles to actively deform their wings (Altshuler et al., 2015).

Studies on flexible wings have focused on the passive deformations as a response to aero-dynamic, inertial and elastic forces. These researches highlighted a lift augmentation for moderate flexibility, whereas too much flexibility causes a decrease in the lift production (Shyy et al., 2010). Moreover, both chordwise and spanwise flexibilities have been investi-gated as well as their combination. Heathcote et al. (2004) have investigated the effect of the chordwise flexibility on a pure plunging airfoil started from rest. They showed a higher thrust production for airfoil with intermediate stiffness. Later, Heathcote et al. (2008) stud-ied the spanwise flexibility of a rectangular wing oscillating in pure heave. They suggested better aerodynamic performance through higher effective angles of attack during the motion, again for an intermediate flexibility. Dai et al. (2012) also studied the effect of the flapping-to-natural frequency ratio of an hovering wing using 3D simulations. They showed that, when the non-dimensional stiffness (i.e. frequency ratio) is lower than 0.3, the deformation can improve lift production as well as aerodynamic efficiency. More recently, Diaz-Arriba et al. (2022) performed numerical simulations on flexible wings undergoing both revolving and pitching motions in hovering conditions. They also concluded that a moderate flexibility can improve aerodynamic performance compared to rigid wings and that a too high degree of flexibility is detrimental.

Besides, studies on flexible flapping wings have shed some light on the main physical reasons for the performance enhancement due to moderate flexibility. For example, Addo-Akoto et al. (2021) have shown experimentally that flexible wings produce more lift than rigid wings at the end of a stroke because of the induced cambering and twisting of the wing during this phase. These deformations result in an upward tilting of the instantaneous net force vector, thus increasing the lift production of flexible wings compared to rigid ones. Moreover, the linear twist experienced by the wing causes the LEV to be robustly attached to the wing surface, hence increasing the radial limit of delayed stall. Other works include numerical simulations on 2D chordwise flexible flapping wings under hovering conditions, as performed by Eldredge et al. (2010) and Vanella et al. (2009), where the wing consisted of two rigid segments linked by a hinge. Both concluded that chordwise flexibility has a major influence on vortex dynamics around the wing. The wake capture mechanism is fundamentally altered by flexibility with the LEV from one stroke having almost no influence on the flow of the following stroke. This permits to generate an all new LEV, increasing lift production with respect to that obtained with rigid wings. Another mechanism was highlighted by Kodali et al. (2017) who have investigated the effect of resonance of the wing structure. They considered passerine and goose wing models and observed that for both flyers, a flapping frequency that matches the dominant bending mode, is conducive to thrust generation. This suggests that birds may take advantage of resonance to fly.
On the contrary, the controlled deformation of the wing is only poorly explored. Kang et al. (2023) and Wang et al. (2022) performed numerical simulations of morphing wings in forward flight and suggested a potential increase of the overall performance through dynamic chordwise retracting and stretching and local deformation of a two-jointed wing, respectively. Cong et al. (2023) numerically studied the effect of active wing-chord adjustment of a rectangular plate under forward flight conditions with a combined heaving-pitching motion. The wing chord is adjusted during the motion to investigate the amplitude and phase of the chord variation. They showed better aerodynamic performance (for both thrust and lift) if the wing stretching is in phase with the heaving motion by delaying the detachment of the LEV. On the other hand, out of phase stretching appears to be detrimental because of the enhanced detachment of the LEV. Dong et al. (2022) have numerically explored the effect of dynamic morphing on a hummingbird-inspired flapping wing under hovering flight conditions. They showed the effect of both twisting and bending on aerodynamic performance. They concluded that bending has almost no influence on lift and efficiency while twist and camber have a large impact. We can note other studies such as Soto and Bhattacharya (2023), who have demonstrated experimentally the effect of dynamic twisting on the flow field around a heaving flat plate. When the plate is twisted in an opposite direction to the heave, the LEV is greatly increased and so the circulation during the downstroke, causing the lift and drag coefficient to be relatively higher compared to rigid wings.

To further investigate on the potential benefits of morphing in flapping wings, we here focus on the controlled deformation of a flapping wing under hovering flight conditions. Toward that end, we use Direct Numerical Simulation (DNS) to solve the flow past a revolving-pitching wing. We first conduct simulations on rigid wings, as a reference case. Then, we account for fluid-structure interactions (FSI) to assess the effect of flexibility on aerodynamic performance. The induced wing bending and twisting are quantified and correlated with performance enhancement. Specifically, the FSI-induced bending phase lag (with respect to the prescribed revolving motion) and amplitude are found to drive changes in aerodynamic performance. Finally, the measured FSI-induced bending is prescribed to a flapping (morphing) wing and the phase lag is varied to further understand its effect on aerodynamic performance. A quasi-steady lift model is further used to help understand the relative contributions of translational, rotational and added mass forces on the aerodynamic performance of rigid and morphing flapping wings. Unsteady effects are in turn inferred from unsteady results and quasi-steady analysis.

## 2. Numerical Methods

The geometry considered here is a three-dimensional wing undergoing both revolving and pitching motions, with angular speeds $\dot{\phi}$ and $\dot{\alpha}$ respectively, as shown in Figure 1 (bottom right). A sine function is used for both revolving and pitching motions, such that during half a flapping period $T/2$ (with complete flapping period $T \approx 0.098$s), the wing covers a revolving amplitude of 120° and the minimum angle of attack $\alpha$ is reached at mid-stroke. Hence, the motion is symmetric, i.e. the same angle of attack is reached at mid-downstroke and mid-upstroke as depicted in Figure 1, top right. The pitching axis is located one quarter
<figure><img src="image_3.png" /><figcaption>Figure 1: Sketch of the computational domain (left), revolving and pitching speeds for both downstroke (white area) and upstroke (gray area) for a reference case (top right) and kinematics of the flapping wing for the downstroke phase (bottom right).</figcaption></figure>

chord away from the leading edge. We stress that these conditions are similar for all cases studied in this article. In other words, revolving amplitude $\phi_0$ and flapping period $T$ remain constant across all rigid and deformed wing cases and the motion is always symmetric.

Based on the chord $c$ and the mean velocity over a stroke $\bar{U}_{ref} = 1.221\text{m.s}^{-1}$ at the $r = 0.72R$ radial position, the Reynolds number is set to 840 for the rigid wing cases. The wing is a NACA0012 profile extruded in the radial direction with a constant chord $c = 0.01\text{m}$ and no twist. The aspect ratio $AR$ is $R/c = 4$.

The StarCCM+ ver.13.04 software is used to compute the flow around the flapping wing by directly solving the incompressible Navier-Stokes equations using a cell-centered finite volume method. The wing is placed in a small cylindrical domain that moves in a larger, fixed background cylindrical domain using an overset mesh technique. Figure 1 (left) displays the computational domain considered here. The wing is treated as a non-slip wall and flaps along an horizontal axis parallel to the top and bottom boundaries of the background cylindrical domain. These two boundaries are modeled as stagnation inlet and pressure outlet, respectively. Lateral boundaries are treated as slip walls. Both temporal and spatial discretizations are achieved through second order numerical schemes. A predictor-corrector approach is used to solve the momentum and continuity equations in an uncoupled way. Specifically, SIMPLE-type (see e.g. Patankar and Spalding, 1972; Patankar, 1980) and Rhie-Chow-type algorithms are used for pressure-velocity coupling. The reader is referred to the work of Muzaferija and co-workers (Muzaferija, 1994; Demirdžić and Muzaferija, 1995) for further details about these numerical methods.

The numerical simulations are first conducted on a reference, rigid wing case with
<figure><img src="image_4.png" /><figcaption>Figure 2: Period-averaged lift coefficient $\bar{C}_L$ as a function of the flapping period $T$ for different spatial and temporal resolutions. (a) shows results obtained for different spatial resolutions $\Delta x_s = 0.02c, 0.01c$ and $0.005c$ for $\Delta t = T/1000$. (b) presents $\bar{C}_L$ for different temporal resolutions $\Delta t = T/1000, T/500$ and $T/250$ for $\Delta x_s = 0.02c$. Horizontal lines represent $\pm 3\%$ of the $10^{th}$ period $\bar{C}_L$ value for the $\Delta x_s = 0.005c$ and $\Delta t = T/1000$ case.</figcaption></figure>

$\alpha(0.25T) = \alpha(0.75T) = 45^\circ$. Three spatial resolutions are tested using cell dimensions of $\Delta x_s = 0.02c, 0.01c$ and $0.005c$ on the surface of the wing and $\Delta x_v = 0.04c, \Delta x_v = 0.02c$ and $\Delta x_v = 0.01c$ in the small moving domain. Similarly, three temporal resolutions are tested, with timesteps $\Delta t = T/250, T/500$ and $T/1000$. These tests demonstrated that, with respect to the case with $\Delta x_s = 0.02c$ and $\Delta t = T/250$, an increase in spatial and temporal resolutions do not yield significant changes in period-averaged lift coefficient $\bar{C}_L$, as shown in Figure 2. Note that horizontal lines in Figure 2 depict $\pm 3\%$ of the $\bar{C}_L$ value obtained for the $\Delta x_s = 0.005c$ and $\Delta t = T/1000$ case (highest resolutions) during the $10^{th}$ period of the motion. Furthermore, it is shown that initial transients have sufficiently decayed after two flapping cycles. Hence, results reported in this paper are those obtained during the third flapping cycle with $\Delta x_s = 0.02c$ and $\Delta t = T/250$. This is in line with previous work conducted at similar Reynolds numbers, where validation of the solver on various moving rigid body cases are reported (e.g. Jardin et al., 2012; Jardin, 2017; Jardin and Doué, 2019).

Simulations on a flexible (passive deformation) flapping wing are then conducted using the same approach as that reported by Diaz-Arriba et al. (2022). The wing flexibility is defined by two parameters, the mass ratio $\gamma = \rho_s h/(\rho c)$, which relates inertial to aerodynamic forces. Here, $h$ is the thickness of the wing and $\xi = f/f_n$ is the frequency ratio, where $f = 1/T$ is the frequency of the flapping motion and $f_n$ is the first natural frequency of the wing. $f_n$ is directly related to the Young modulus $E$ and density of the solid $\rho_s$ using the classical Euler-Bernoulli beam theory. Specifically, $2\pi f_n = (1.875/2\pi R^2)\sqrt{EI/\rho_s A}$ where $A$ and $I$ are the area and the second moment of area of the wing cross-section. Based on these parameters, the wing deformation is computed at each time step following the principle of virtual work. The total Lagrangian displacement finite element formulation is used. Strain
is assumed to be small and the material is modeled as isotropic, linear-elastic with the linear stress-strain relationship given by Hooke's law. The resolution of wing deformation is strongly coupled to the resolution of the fluid motion (i.e. two-way coupled approach). Accordingly, the mesh within the moving cylinder is deformed at each time step. The approach has been previously validated against experiments in Diaz-Arriba et al. (2022).

Finally, simulations on a morphing (active deformation) flapping wing are conducted. As will be shown thereafter, torsion induced by fluid-structure interactions is negligible with respect to flexion. Hence only spanwise bending is prescribed in the morphing wing case. The latter is applied normal to the chord. We stress that the approach is similar to that used for passive deformation except that deformation is prescribed (i.e. fluid-structure interactions are not solved for). That is, the mesh within the moving cylinder is deformed at each time step following the prescribed deformation. Similar to the rigid wing case, a mesh sensitivity study has been conducted on morphing wings. This study (not reported here for the sake of conciseness) showed that, with respect to the case with $\Delta x_s = 0.02c$ and $\Delta t = T/250$, an increase in spatial and temporal resolutions do not significantly impact the period-averaged lift coefficient $\bar{C}_L$. Same conclusions as for rigid wings can be drawn and results reported for this configuration are obtained during the third flapping cycle.

In what follows, results are analyzed in terms of the mean vertical force coefficient $\bar{C}_V$ and the mean vertical force-to-power coefficient $\bar{C}_V/\bar{C}_P$:

$$C_V = 2(L \cos \beta + D \sin \beta)/(\rho S \bar{U}_{ref}^2), \quad (1)$$

$$C_P = -2(\dot{\phi} Q_R + \dot{\alpha} Q_P)/(\rho S \bar{U}_{ref}^3). \quad (2)$$

where $L$, $D$, $Q_R$, $Q_P$, $\dot{\phi}$ and $\dot{\alpha}$ are the lift, drag, revolution and pitching torques, and revolving and pitching angular speeds, respectively. $\beta^1$ is the stroke plane angle and is computed as: $\beta = \arctan(\bar{D}/\bar{L})$. The mean value of each coefficient is period-averaged over the 3rd period of the flapping kinematic.

# 3. Results

## 3.1. Rigid wings

Figure 3 shows the mean vertical force $\bar{C}_V$ as a function of the mean vertical force-to-power ratio coefficient $\bar{C}_V/\bar{C}_P$ obtained for the rigid wing for different angles of attack $\alpha$. Results are reported for mid-stroke angles of attack varying from 0° to 75° with a step of 5°. First, both $\bar{C}_V$ and $\bar{C}_V/\bar{C}_P$ increase with pitch angle, from 0° to 15°. Maximum $\bar{C}_V/\bar{C}_P$ is reached at 15° after which it decreases while $\bar{C}_V$ continues to increase. Maximum in $\bar{C}_V$ is reached at 50°. Beyond 50°, $\bar{C}_V$ decreases with increasing $\alpha$. The lift-optimal pitch angle of approximately 50° has already been reported for other motion laws (see Diaz-Arriba et al. (2021); Oyama et al. (2009); Sane and Dickinson (2001)), which suggests that it is relatively

<footnote n="1">¹For a symmetric movement, $\beta$ is close or equal to 0. Thus, the drag has almost no contribution to the vertical force.</footnote>
<figure><img src="image_5.png" /><figcaption>Figure 3: $\bar{C}_V$ as a function of $\bar{C}_V/\bar{C}_P$ for different $\alpha$ (with a step of 5°). Results for rigid, passively and actively deformed wings are shown. Dotted horizontal and vertical lines represent $\pm2.5\%$ of the $\bar{C}_V$ and $\bar{C}_V/\bar{C}_P$ of the passively deformed wing simulation (FSI) with $\alpha = 45°$, respectively.</figcaption></figure>

robust to revolving and pitching velocity waveforms. Note that the effective angle-of-attack of the wing is lower than 50° due to the downwash induced as a reaction to vertical force production during previous strokes. Yet, it appears that stall (in the time-averaged sense) is delayed to larger pitch angles than conventionally observed on 2D wings. This phenomenon is known to be partly caused by strong rotational effects, which promote the attachment of the leading edge vortex (LEV) that forms on the suction surface of the wing (Lentink and Dickinson, 2009; Jardin, 2017). Similarly to that observed on delta wings, this LEV supports a low pressure region that contributes to generating a strong normal, aerodynamic force.

Formation of a stable conical LEV can be observed in Figure 4 which displays a time sequence of Q-criterion isosurfaces obtained on the rigid wing case with $\alpha = 45°$. The values displayed are $Q_c^2/U_{tip}^2 = 1$ and 10 in light gray and cobalt, respectively. $U_{tip} = 2\phi_{tip}/T$ where $\phi_{tip} = \phi R$ is the amplitude covered by the wing tip. The strength of the LEV increases with $\alpha$, hence the normal force. However, the relative contribution of the normal force to the vertical force decreases with $\alpha$ and hence $\bar{C}_V$ decreases at some point. It can also be observed from Figure 4 that, in the outboard region of the wing, the LEV bursts into small scale structures. This phenomenon is known to result from the interaction of LEV core flow (directed outboard) with the tip vortex flow (directed inboard) (Medina and Jones, 2016). The region of LEV burst expands inboard as the wing both decelerates and pitches up, which here occurs after $t/T = 2.25$ for the downstroke phase (and 2.75 for the upstroke phase). Hence, when compared to Diaz-Arriba et al. (2021) where revolving and
pitching (polynomial) laws with a constant revolving speed phase were used on similar wings, significant differences in the LEV structure are observed. Specifically, the flow exhibits small scale structures over a broader spanwise region in the present case, which seemingly results in stronger wing-wake interactions at stroke reversal.

Finally, Figure 5 shows non-dimensional spanwise vorticity contours obtained in the spanwise cross-section $r = 0.72R$ for the same $\alpha = 45^\circ$ reference case. The first row displays vorticity contours for the early downstroke phase, while the second row shows the early upstroke phase. These snapshots help appreciate the complex wing-wake interactions (i.e. wake capture) at stroke reversals. In the first row, the wing revolves from left to right while pitching down. The bursted LEV that formed during the previous half-stroke appears as clockwise rotating (blue) small scale structures ($t/T = 2.02$) that then sweep on the right side of the wing, both upwards towards the leading edge, and downwards towards the trailing edge. Vorticity spots that sweep upwards strongly interact with the new counter-clockwise rotating (red) LEV. That is, they locally increase the inflow angle and velocity and thus contributes to the rapid formation of the new LEV that tends to lift off from the wing surface ($t/T = 2.06$ and $2.08$). Similar dynamics is observed on the second row. Finally, the third row displays vorticity contours obtained during the upstroke. Lift off of the LEV from the wing surface is visible at $t/T = 2.60$. However, the LEV then moves back closer to the surface ($t/T = 2.70$) and rapidly burst under the influence of deceleration and pitch up motions.

## 3.2. Passive deformation

The Fluid-Structure Interaction (FSI) simulation is then conducted for a wing with mass ratio $\gamma = 0.5$ and frequency ratio $\xi = 0.5$. This corresponds to the following dimensional wing properties: Young modulus $E = 4.02$ MPa, Poisson ratio $\nu_s = 0.48$ and density $\rho_s = 2.542$ kg/m$^3$. Diaz-Arriba et al. (2022) showed that a similar wing (NACA0012 profile and $AR = 4$ with constant $c = 0.01m$) with these properties enhances lift compared to rigid cases.

The aerodynamic performance of the flexible wing is reported in Figure 3 using the upward pointing triangle symbol (‘Passive Deform.’) for three $\alpha$ values corresponding to maximum $\bar{C}_V$, maximum $\bar{C}_V/\bar{C}_P$ and an intermediate case at $\alpha = 30^\circ$. For maximum $\bar{C}_V$ ($\alpha = 45^\circ$), it is shown that a moderate flexibility enhances the aerodynamic performance of a flapping wing, supporting previous results from the literature. In the present case, flexibility has virtually no influence on the $\bar{C}_V/\bar{C}_P$ ratio but it has a non-negligible impact on the $\bar{C}_V$ value, which is increased by about 7% with respect to the rigid wing case. For maximum $\bar{C}_V/\bar{C}_P$ ($\alpha = 15^\circ$), the opposite trend is observed. Flexibility has no influence on the $\bar{C}_V$ value while $\bar{C}_V/\bar{C}_P$ decreases by about 10%. Finally, for the $\alpha = 30^\circ$ case, flexibility has a positive impact on $\bar{C}_V$, which is increased by approximately 5.5% with respect to the rigid wing case. On the contrary, the $\bar{C}_V/\bar{C}_P$ ratio is very similar to that obtained in the rigid wing case, with a difference of less than 1%. In what follows, we will focus on the $\alpha = 45^\circ$ case where, conversely to other cases, flexibility is found to have a non-negligible positive impact on aerodynamic performance.
<figure><img src="image_6.png" /><figcaption>Figure 4: Q-criterion isosurfaces at different instants $t/T$ for the $\alpha = 45^\circ$ reference rigid wing case. $Q_c^2/U_{tip}^2 = 1$ and 10 isosurfaces are displayed in light gray and cobalt respectively.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Figure 5: Non-dimensional spanwise vorticity contours in the $r = 0.72R$ spanwise cross-section for the $\alpha = 45^\circ$ reference rigid wing case. Early downstroke, early upstroke and upstroke phases are displayed.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 6: Non-dimensional wing tip displacement (a) and twist angle (b) with respect to the rigid wing position as a function of non-dimensional time. Dashed-line in (a) is a sinusoidal fit of the measured displacement. (b) displays the twist angle at three spanwise locations. White regions depict the downstroke phase while gray ones represent the upstroke phase.</figcaption></figure>

The increase in lift for the $\alpha = 45^{\circ}$ case results from passive deformation, which can be decomposed into two principal displacements: the bending and twisting of the wing. Bending is quantified by monitoring the displacement of a point located on the pitch axis, at the tip of the wing. The displacement is measured in a direction normal to the chord and with respect to the rigid wing case. Figure 6(a) shows the displacement normalized by the chord $c$ as a function of non-dimensional time (non-dimensionalized by the period of the flapping motion $T$). It is observed that the amplitude of bending reaches almost 40% of the wing chord and that the deformation slightly lags the prescribed revolving motion. The bending deformation over ten consecutive flapping cycles is found to be reasonably well approximated by a sine function of the form:

$$ Pos(t) = a_1 \sin(b_1 t + c_1). \qquad (3) $$

where $a_1 = 0.00388m$ is the amplitude, $b_1 = 63.95s^{-1}$ is the angular frequency<fnref n="2" /> and $c_1 = -3.006$ is the phase lag of the wing tip position compared to that of the wing root. Note that the fit is performed from the sixth to the last period such that the transient observed on the tip location in the time interval $t/T = [2, 4]$ is not taken into account. However, it was found that the change in tip location during this transient does not significantly affect the period-averaged lift. This is shown in Appendix A where $\bar{C}_L$ value of the FSI simulation is plotted for ten flapping cycles.

Wing twisting is quantified by monitoring displacement of points located at the leading and trailing edges of the wing, for three spanwise locations ($r = 0.5R, r = 0.75R$ and $r = R$). The line connecting points at the leading and trailing edges makes an angle $\epsilon$ with respect to that obtained in the rigid wing case. Figure 6 (b) shows $\epsilon$ obtained for the three spanwise

<footnote n="2">²Actually, $b_1$ is the angular frequency related to the sinusoidal, revolving motion of the wing: $b_1 = 2\pi/T$.</footnote>
locations during the seventh period of the flapping motion. First, the same trend is observed for every wing location. Moreover, the maximum twist angle is found to be at the tip of the wing, which is expected since deformation increases with wing loading (which increases from root to tip). However, the twist angle never exceeds 1°. From Figure 3 it can be inferred that such a small deformation will have a weak impact on aerodynamic performance. As such, bending is here found to be the predominant deformation mechanism.

### 3.3. Prescribed deformation

In this section, bending is prescribed to the flapping wing (i.e. active deformation or morphing), rather than being induced by fluid-structure interactions. Note that only bending is applied since twisting was found to be weak compared to bending in the passive wing case. We first prescribe the same bending as that obtained passively in the previous section, albeit with some simplifications. The corresponding displacement velocity is directly computed from the fitted position in equation 3:

$$V_d(r,t) = \frac{r}{0.04}a_1b_1\cos(b_1t + c_1). \quad (4)$$

Note that $V_d$ is here taken as a linear function of the radial (spanwise) coordinate $r$. However, because the aerodynamic force is a quadratic function of $r$, the wing actually takes an approximately parabolic shape. Nonetheless, simulations with both deformations (linear and quadratic) were found to provide similar results (see Appendix A) for this type of wing. From Figure 3, one can notice the relative similarity between the passive and morphing wings’ aerodynamic performance. That is, both active and passive deformations yield similar aerodynamic performance despite that only bending is prescribed in the morphing wing case and that it is defined as a linear function of $r$. Therefore, insight into the role of FSI-induced deformation phase lag and amplitude can be provided by dynamically morphing the wing for different values of parameters $a_1$ and $c_1$. In what follows, we focus on the role of morphing phase lag on the aerodynamic performance of the flapping wing.

#### 3.3.1. Overall analysis

Figure 7 depicts the influence of $c_1$ on wing tip displacement. The morphing wing is shown in lightgray while the rigid wing is shown in gray, for reference. It can be observed that a phase lag of $\pi$ induces backward deformation of the wing throughout the stroke while being in the rigid wing position at stroke reversal. On the contrary, a phase lag of $2\pi$ causes the wing to bend forward during the stroke while remaining at the rigid wing position at stroke reversal. Finally, for an intermediate value of $3\pi/2$, the wing undergoes maximum deformation at stroke reversal and moves back to the rigid wing position at mid-stroke. Hence, the wing tip lags inboard sections during the first half of downstroke (resp. upstroke) and leads inboard sections during the second half of downstroke (resp. upstroke). Because of the phase lag, the morphing wing has a lower or larger speed than the rigid wing, depending on the value of $c_1$. For example, if $c_1 = \pi$, the morphing wing has a lower speed during the first half stroke and a larger speed during the second half stroke compared to the rigid wing case. Conversely, if $c_1 = 2\pi$, the morphing wing has a larger speed during the
<figure><img src="image_9.png" /><figcaption>Figure 7: Illustration of the wing during downstroke ($\alpha$ = 45°) for three $c_1$ values. Gray wing depicts the rigid wing while morphing wings are shown in lightgray.</figcaption></figure>

first half stroke and a lower speed during the second half stroke. Finally, in the $c_1 = 3\pi/2$ case, the morphing wing has a larger speed during the whole stroke compared to the rigid wing case. Therefore, one could expect a larger lift production in the early downstroke and upstroke phases in the $c_1 = 2\pi$ case, compared to the rigid wing case, due to larger contributions from quasi-steady lift. Similarly, the rigid wing case would produce more lift than the $c_1 = \pi$ case during these phases. These points will be further addressed thereafter. In some cases, like in the $c_1 = 3\pi/2$ case, the deformation of the wing at stroke reversals implies that the wing covers a larger amplitude $\phi_{eff}$ than the 120° covered by the rigid wing during each stroke phase. Thus, a new reference velocity $U_{eff} = 2\phi_{eff}R/T$ is defined to normalize the aerodynamic coefficients<fnref n="3" />, namely $\bar{C}_{V,eff}$ and $\bar{C}_{P,eff}$.

Figure 8 shows the aerodynamic performance of the morphing wing for different values of phase lag $c_1$. Data are compared to those obtained for the passively deformed and rigid wings at different pitch angles. At low pitch angles, the phase lag has a drastic influence on both $\bar{C}_{V,eff}$ and $\bar{C}_{V,eff}/\bar{C}_{P,eff}$. For example, when $\alpha = 10°$, a phase lag of $14\pi/10$ increases $\bar{C}_{V,eff}$ by 29% and $\bar{C}_{V,eff}/\bar{C}_{P,eff}$ by 31% when compared to the rigid wing case. The evolution of $\bar{C}_{V,eff}$ against $\bar{C}_{V,eff}/\bar{C}_{P,eff}$ with phase lag follows a clear elliptical shape at $\alpha = 10°$, which progressively collapses as $\alpha$ increases to 70°. In other words, phase lag has a major influence at low angles of attack but only a weak impact at high $\alpha$. Yet, a phase lag value between $\pi$ and $3\pi/2$ enhances both vertical force and efficiency for all $\alpha$. As previously observed for $\alpha = 10°$, maximum efficiency is found to be around $c_1 = 3\pi/2$ for every configurations. This phase lag also generally leads to the largest vertical force

<footnote n="3">³Note that, in the rigid wing case, $U_{eff} = U_{tip}$.</footnote>
<figure><img src="image_10.png" /><figcaption>Figure 8: $\overline{C}_{V,eff}$ as a function of $\overline{C}_{V,eff}/\overline{C}_{P,eff}$ for $\alpha = 10^\circ, 20^\circ, 30^\circ, 45^\circ$ and $70^\circ$. Filled markers depict results obtained for rigid wing cases, for each angle of attack. Open markers display aerodynamic coefficients for morphing wings with different phase lag, varied by step of $\pi/10$. Results from passively deformed wing simulation (FSI) with $\alpha = 45^\circ$ case are also shown.</figcaption></figure>

coefficient for a given $\alpha$. The largest $\overline{C}_{V,eff}$ among all tested cases is obtained for $\alpha = 45^\circ$. In the next sections, we will focus on quasi-steady and unsteady analysis to help understand these trends.

### 3.3.2. Quasi-steady analysis

In this section, we use the quasi-steady, semi-empirical model from Lee et al. (2016) to help gain insight into the influence of phase lag on aerodynamic performance. The model, which derives from the model of Sane and Dickinson (2001), was developed based on data obtained for a rigid trapezoidal wing with a flat and infinitesimally thin wing profile, which somewhat slightly differs from our wing (with NACA0012 profile). Although this may contribute to discrepancies between present simulations and the model, the trends in the evolution of aerodynamic performance with phase lag are expected to be similar with both approaches.

The model separates the instantaneous aerodynamic force $F_{inst}$ acting on the wing into four components such that:

$$F_{inst} = F_a + F_{rev} + F_{rot} + F_{wc} \qquad (5)$$

where $F_a$ is the force due to added mass, $F_{rev}$ is the quasi-steady force due to the revolving motion, $F_{rot}$ is the rotational force and $F_{wc}$ is the force due to wake capture. Because the wake capture term is not straightforward to model, it will not be taken into consideration
in this investigation. Nevertheless, wake capture does not exist during the very first stroke and hence the analysis can be performed on this specific phase. This point will be further addressed in what follows.

The first term is derived from the work of Sedov (1965) who provided the added mass force experienced by a two-dimensional flat plate in an inviscid fluid. Integrating this two-dimensional force along the span and adding empirical correction factors f yields:

$$F_a = f_{\lambda,a} f_{AR,a} f_{Re,a} \left( \rho \frac{\pi}{8} R^2 c^2 (\ddot{\phi} \sin \alpha) + \ddot{\alpha} \rho \frac{\pi}{16} Rc^3 \right). \quad (6)$$

where $f_{\lambda,a} = 47.7\lambda^{-0.0019} - 46.7$ is the taper ratio $\lambda = 1$ related term. In the present case $f_{\lambda,a} = 1$. $f_{AR,a} = 1.294 - 0.590AR^{-0.662}$ is the added mass correction factor associated to the wing aspect ratio $AR = 4$. $f_{Re,a} = 0.776 + 1.911Re^{-0.687}$ is the correction factor due to Reynolds number effects. In Lee et al. (2016), Re is taken at the $R_2 = \sqrt{\frac{1}{S} \int_0^R cr^2 dr} = R/\sqrt{3}$ position, which in our case is Re = 677.

The second term is related to the revolving motion of the wing:

$$F_{rev} = f_{AR,rev} f_{Ro,rev} C_L(\alpha, Re) \frac{1}{6} \rho \dot{\phi}^2 c R^3. \quad (7)$$

where $C_L(\alpha, Re) = (1.966 - 3.94Re^{-0.429}) \sin(2\alpha)$ is the revolving lift coefficient depending on the Reynolds number and angle of attack. $f_{AR,rev} = 32.9 - 32AR^{-0.00361}$ is the aspect ratio correction factor related to the revolving motion. $f_{Ro,rev} = -0.205 \arctan((0.587(Ro - 3.105)) + 0.870)$ is the correction factor due to Rossby number effects. Here, at the $R_2$ position, $Ro = R_2/c = 2.31$.

Finally, the last term pertains to rotational effects:

$$F_{rot} = f_{r,rot} \gamma_{rot} \frac{1}{2} \rho \dot{\phi} \dot{\alpha} c^2 R^2, \quad (8)$$

where $\gamma_{rot} = 0.842 - 0.507Re^{-0.1577}$ is the rotational force coefficient and $f_{r,rot} = 1.570 - 1.239(x_{rot}/c)$ is a correction factor accounting for the rotational (or pitch) axis position. Here, $x_{rot}/c = 0.25$. Note that $F_a$ and $F_{rot}$ act normal to the chord, whereas $F_{rev}$ is normal to the wing (revolving) velocity. These terms are thus projected on the vertical axis to obtain the vertical force.

The model is straightforward to use in the rigid wing case where the motion is constrained to the horizontal plain. In such a case, $\ddot{\phi}$, $\dot{\phi}$ and $\alpha$ are defined as shown in Figure 1. Note however that, in the morphing wing case, the wing can move along any arbitrary direction. $\ddot{\phi}$, $\dot{\phi}$ and $\alpha$ should thus be replaced by their corresponding value in a coordinate system whose x-axis is collinear with the wing motion. In other words, revolving acceleration, speed and angles of attack are taken in a coordinate system that is tilted with respect to that used for the rigid wing case (where the x-axis is collinear with the horizontal direction).

Figure 9 shows the effective lift coefficient predicted by the model for the $\alpha = 10^\circ$ rigid wing case (a) and for the morphing wing case with $c_1 = 3\pi/2$ (b). Total lift and contributions from added mass, revolving and rotational (or pitch) motions are displayed as
<figure><img src="image_11.png" /><figcaption>Figure 9: Effective lift coefficient $C_{L,eff}$ predicted by the model and its different components ($C_a$, $C_{rev}$ and $C_{rot}$) as a function of non-dimensional time for rigid (a) and morphing wing ($c_1 = 3\pi/2$) (b) with $\alpha = 10^\circ$. Results from numerical simulations obtained during the first and third downstroke phases of the flapping motion are also shown.</figcaption></figure>

a function of non-dimensional time during downstroke. Results obtained during the 1st and 3rd downstroke phases using direct numerical simulations are also depicted.

First, the model predicts with reasonable accuracy the global trend in lift produced during the 1<sup>rst</sup> downstroke phase of the flapping motion, for both cases. In particular, a first bump is observed near $t/T \approx 0.15$, followed by a valley near $t/T \approx 0.20$ and a peak around $t/T = 0.35$. The trend is different during the 3<sup>rd</sup> downstroke. In the rigid wing case, a wide valley is observed around $t/T = 0.15$. In the morphing wing case there is no clear evidence of both bump and valley at $t/T \approx 0.15$ and $0.20$, respectively. As previously mentioned, discrepancy between the trend from the model and that from simulations of the 3<sup>rd</sup> downstroke partly results from the presence of strong wake-capture in the simulations.

Moreover, it is observed from the model that lift production is mainly due to the revolving motion. Added mass and rotational effects have comparatively less influence, although their contributions on instantaneous lift are not negligible. As expected, it can further be noted that those effects have a very small net contribution on the time-averaged lift. While the trend in lift is relatively well predicted by the model, the latter generally overestimates lift when compared to direct numerical simulations. This discrepancy may arise from unsteady effects, which are not accounted for in the quasi-steady model. In particular, unsteady development of the wake as the wing moves from rest tends to attenuate lift, which is in line with observations from Figure 9. In classical unsteady airfoil theory, this effect is taken into account using the Wagner function. Its well-known counterpart for harmonic motions in the frequency domain is the Theodorsen function, which is a transfer function that attenuates and lags quasi-steady lift depending on the frequency of the motion.

Given the reasonable trend provided by the model, it is then used to compute $\overline{C}_{L,eff}$ as a function of $\phi_{eff}$ for different phase lag $c_1$, as shown in Figure 10(a) for $\alpha = 45^\circ$ and $\alpha = 10^\circ$. Results from numerical simulations are also shown (averaged over the 1<sup>rst</sup>
<figure><img src="image_12.png" /><figcaption>Figure 10: $\bar{C}_{L,eff}$ predicted by the model as a function of $\phi_{eff}$ for different $c_1$ values, for $\alpha = 45^\circ$ and $\alpha = 10^\circ$. Results obtained during the 1$^{rst}$ downstroke using numerical simulations are also shown (a). $\bar{C}_{L,eff}$ for the morphing wing are displayed using open markers and shown for different $c_1$ values, varied by step of $\pi/10$. Filled markers depict $\bar{C}_{L,eff}$ obtained for the rigid wing. Component of the model for $\alpha = 45^\circ$ are displayed in (b).</figcaption></figure>

downstroke phase) for comparison, where open and filled markers are used for morphing and rigid wings, respectively. Furthermore, for $\alpha = 45^\circ$, time-averaged contributions from added mass, revolving and rotational quasi-steady forces are displayed in Figure 10(b).

First, for both $\alpha = 45^\circ$ and $\alpha = 10^\circ$ cases, the evolution of $\bar{C}_{L,eff}$ with $\phi_{eff}$ (as $c_1$ varies) reveals a clear elliptical shape, as discussed from simulation results in the previous section. Thus, while $\bar{C}_{L,eff}$ values from the model do not accurately match those from simulations, the model is able to qualitatively predict the influence of $c_1$ on $C_{L,eff}$. In turn, Figure 10(b) provides insight into the contributions responsible for the tilted elliptic shape in Figure 10(a). Specifically, it appears that the elliptic shape is induced by added mass effects. Moreover, the rotational component does not influence the period averaged $C_{L,eff}$, i.e. the period-averaged value of $C_{L,rot}$ is null. Finally, quasi-steady forces due to the revolving motion are found to drive the minimum and maximum $\bar{C}_{L,eff}$ values (for a given $\alpha$), hence the inclination and position of the elliptic shape in the $(\bar{C}_{L,eff},\phi_{eff})$ space.

As previously discussed, discrepancies between the model and the simulations may partly arise from unsteady effects. For example, the Wagner effect may be responsible for the overestimation of $\bar{C}_{L,eff}$ by the model, as observed on Figure 9. We thus focus on the unsteady analysis of the flow in the next section.

### 3.3.3. Unsteady analysis

We now focus on the $\alpha = 10^\circ$ case for which the largest differences between the morphing and rigid wings' performance are observed (cf. Figure 8). Figure 11 displays the instantaneous lift coefficient obtained for the rigid wing case and for the morphing wing cases $c_1 = \pi$, $2\pi$ and $14\pi/10$ (the most efficient case) during the third period of the flapping motion. The instantaneous displacement of the point located at the wing tip and on the pitch axis for $c_1 = 14\pi/10$ is also shown (see scale on the right-side ordinate axis). Compared to the rigid
wing case, this phase lag induces a slight backward deformation at the beginning of the flapping motion. That is the morphing wing has a lower speed than the rigid wing over a small portion of the stroke. Points $v_{do}$ and $v_{up}$ mark the instants from which the morphing wing has a larger speed than its rigid counterpart (due to deformation) for the downstroke and upstroke phases, respectively.

First, it is observed that, for all cases, the time evolution of the lift coefficient during upstroke is relatively similar to that obtained during downstroke. We thus focus our analysis on the upstroke only.

In the rigid wing case, the lift coefficient slightly increases between $t/T = 2.50$ and $t/T = 2.57$ and then drops to negative values. As inferred from the previous quasi-steady analysis and Figure 9, the initial increase due to combined

<figure><img src="image_13.png" /><figcaption>Figure 11: $C_{L,eff}$ as a function of non-dimensional time obtained during the third period for the rigid and morphing wings with $\alpha = 10^{\circ}$. Results for three morphing wings with $c_1$ value of $\pi$, $14\pi/10$ and $2\pi$ are shown. The non-dimensional wing tip displacement with respect to the rigid wing is also shown for $c_1 = 14\pi/10$. The white area is the downstroke phase, the upstroke phase is represented by the gray area.</figcaption></figure>

added mass and quasi-steady forces from the revolving motion are counter-balanced by quasi-steady rotational forces and wake capture, the latter being ultimately responsible for the lift drop. The lift increases again after $t/T = 2.70$, i.e. after the effect of wake-capture has decayed and as both contributions from quasi-steady rotational force (due to pitch-up motion) and force due to revolving motion increase. Finally, as the wing decelerates and reaches very high angles of attack, $C_L$ severely drops, reaching approximately zero at the end of the stroke.

For the $c_1 = \pi$ case, the wing has a lower speed than the rigid wing case during the first half of the upstroke phase (i.e. from $t/T = 2.50$ to $2.75$). Thus, from simple quasi-steady
considerations discussed in the previous section, $C_{L,eff}$ is expected to be lower than that obtained for the rigid wing case. Accordingly, it can be observed from Figure 11 that $C_{L,eff}$ rapidly drops after $t/T = 2.50$, reaching values below -0.5 around $t/T = 2.60$. Conversely, the morphing wing has a larger speed than the rigid wing during the second half of upstroke (i.e. from $t/T = 2.75$ to 3) and hence experiences a larger lift force. $C_{L,eff}$ peaks at $t/T = 2.87$, similar to the rigid wing case, but to a value approximately 16% larger, and then drops for similar reasons to those described above.

For the $c_1 = 2\pi$ case, the wing has a larger and lower speed than the rigid wing case during the first and second half of the upstroke phase, respectively. Again, from simple quasi-steady considerations discussed in the previous section, the morphing wing is thus expected to experience larger and smaller lift than the rigid wing between $t/T = 2.50$ and 2.75, and $t/T = 2.75$ and 3, respectively. This is, again, in line with observations from Figure 11, albeit with some differences in the precise timing at which curves from the rigid and morphing wings cross each other. This difference in timing results from effects that are not related to quasi-steady forces associated with the revolving motion (unsteady effects, wake capture, added mass). Then, $C_{L,eff}$ peaks at $t/T = 2.87$, similar to previous cases, but reaches a lower value (approximately 16% smaller than that obtained in the rigid wing case).

In the last $c_1 = 14\pi/10$ case, the morphing wing has lower and larger speed than the rigid wing before and after $v_{up}$, respectively. Therefore, one can expect $C_{L,eff}$ to be below and above that of the rigid wing case before and after (approximately) $v_{up}$, respectively. This can again be correlated with observations from Figure 8. However, it is interesting to note that in this case, $C_{L,eff}$ slightly drops after $t/T = 2.50$ but quickly increases from approximately $t/T = 2.60$ such that it never reaches negative values, conversely to other cases. The differences between the $c_1 = 14\pi/10$ case and other cases slightly after $t/T = 2.60$ and before $t/T = 2.75$ is significant. This is typically within the time interval where wake-capture phenomena are strong, and where the latter were found to be ultimately responsible for negative lift values observed in other cases (see previous section). This point will be discussed in the next paragraph. Finally, $C_{L,eff}$ peaks around $t/T = 2.87$ at a lower value than that obtained for the $c_1 = \pi$ case, despite the fact that the wing has a larger speed at that time. Again, effects that are not related to quasi-steady forces associated with the revolving motion are hence presumably responsible for this lower peak. In particular, in the $c_1 = 14\pi/10$ case, the wing experiences a stronger deceleration than in that in the $c_1 = \pi$ case at $t/T = 2.87$, and hence added mass effects have a detrimental impact on $C_{L,eff}$.

Figure 12 shows a sequence of non-dimensional spanwise vorticity contours obtained at the $r = 0.72R$ spanwise cross-section for both rigid and $c_1 = 14\pi/10$ deformed wing case with $\alpha = 10^\circ$. Snapshots are obtained during the upstroke phase of the third flapping period. First, when compared to Figure 5 which displays similar snapshots for $\alpha = 45^\circ$, the flow during the early upstroke phase exhibits a much less complex structure. At the first instant $t/T = 2.52$, it is characterized by a highly coherent LEV generated during the previous stroke. While the flow in both rigid and morphing wing cases appear very similar, slight differences exist in the position of this LEV relative to the wing. As the wing pitches down, this slight difference, together with different wing kinematics, induces different wing-wake
<figure><img src="image_14.png" /><figcaption>Figure 12: Non-dimensional spanwise vorticity contours in the $r = 0.72R$ spanwise cross-section at different instants during the upstroke phase for a rigid and morphing ($c_1 = 14\pi/10$ case) wing with $\alpha = 10^\circ$. The black outline depicts the position of the wing tip.</figcaption></figure>

interactions. At $t/T = 2.56$, a new LEV forms on the upper surface of the wing in both cases. However, because of the different wing-wake interactions, it develops differently from $t/T = 2.58$ to $t/T = 2.60$. Specifically, while at $t/T = 2.60$ the LEV loses coherency and appears to detach from the rigid wing, it remains closely attached to the morphing wing. This difference in LEV dynamics persist at $t/T = 2.65$ where large differences in $C_{L,eff}$ were observed between the $c_1 = 14\pi/10$ case and other cases (Figure 11). That is, the attached, coherent LEV in the $c_1 = 14\pi/10$ case presumably contributes to the larger lift between $t/T = 2.60$ and 2.75. In turn, this helps explain initial observations from Figure 8 where the $c_1 = 14\pi/10$ case was found to be that with the largest lifting force.

## 4. Conclusion

Extensive works have been conducted on flapping wings these last few decades. Most works have focused on rigid wings, although it has been observed in nature that wings of birds and insects typically deform during the flapping motion. Accordingly, more recent studies have analyzed the role of passive wing deformation on the flow field and resulting aerodynamic performance. These studies generally point towards a beneficial impact of flexibility on lift and/or efficiency. Conversely, the role of active deformation is poorly documented. In the present work we have therefore analyzed the influence of active bending on the aerodynamic performance of three-dimensional flapping (revolving and pitching) wings under hovering flight condition and we have compared it with that obtained from rigid and passively deforming wings.
First, it has been shown that, in the present setup, bending is the dominant deformation in the passively deforming wing case. Moreover, it has been confirmed that bending helps increase aerodynamic performance and that this increase results from a larger effective flapping amplitude and the introduction of a phase lag with respect to the rigid wing motion (hence a phase lag between the revolving and pitching motion). In addition, it has been shown that enhancement in aerodynamic performance is highly dependent on angle of attack.

Second, measured passive wing bending has been actively prescribed to the wing, which we referred to as ‘active deformation’ or ‘morphing’. Active deformation was parametrized by using a sinusoidal law with similar frequency to that of the prescribed flapping motion. The amplitude was kept equal to that measured on the passively deforming wing and the phase lag was systematically varied. It has been shown that phase lag has a strong influence on aerodynamic performance, and that this influence reduces with increasing angle of attack. Both quasi-steady and unsteady analysis demonstrated that this significant increase is correlated, to leading order, with the increase in quasi-steady forces associated with the revolving motion. Furthermore, added-mass was also found to play a non-negligible role, as well as wake-capture mechanisms, which ultimately contributed to optimal performance.

Overall, the present work has shown that active bending can drastically help increase aerodynamic performance of flapping wings. In addition, our parametric study can provide guidelines towards the use of specific structural properties that would help passively deforming wings to deform in a way that best promotes aerodynamic performance (e.g. with proper phase lag). Future work will focus on testing active deformation experimentally.

# Appendix A. Linear bending model

The difference between (active) linear and quadratic bending has been assessed on the period-averaged lift coefficient, as shown in Figure A.13. Horizontal lines depict ±2% of the last $\bar{C}_L$ value obtained from the simulation on the passively deforming wing (FSI) at $\alpha = 45^\circ$. It can be observed that, after two periods, all values are within these bounds and hence both active linear and quadratic bending accurately reproduce results from passive bending.

# Acknowledgements

This work was supported by the French government program “Investissements d’Avenir” (EUR INTREE, Reference No. ANR-18-EURE-0010). The authors are grateful to GENCI-IDRIS, GENCI-CINES (Grant A0122A07178) and CALMIP (Grant 2022- p1425) for providing HPC resources.

# References

* Addo-Akoto, R., Han, J.S., Han, J.H., 2021. Roles of wing flexibility and kinematics in flapping wing aerodynamics. Journal of Fluids and Structures 104, 103317. doi:10.1016/j.jfluidstructs.2021.103317.
<figure><img src="image_15.png" /></figure>

Figure A.13: $\bar{C}_L$ as a function of non-dimensional time for the simulation on passively deforming wing (FSI) and actively deforming wings with linear and quadratic bending. Horizontal lines depict ±2% of the last $\bar{C}_L$ value obtained for the simulation on the passively deforming wing.

Altshuler, D.L., Bahlman, J.W., Dakin, R., Gaede, A.H., Goller, B., Lentink, D., Segre, P.S., Skandalis, D.A., 2015. The biophysics of bird flight: functional relationships integrate aerodynamics, morphology, kinematics, muscles, and sensors. Canadian Journal of Zoology 93, 961–975. doi:10.1139/cjz-2015-0103.

Cong, L., Teng, B., Chen, L., Bai, W., Jin, R., Chen, B., 2023. Aerodynamic performance of low aspect-ratio flapping wing with active wing-chord adjustment. Journal of Fluids and Structures 122, 103964. doi:10.1016/j.jfluidstructs.2023.103964.

Dai, H., Luo, H., Doyle, J.F., 2012. Dynamic pitching of an elastic rectangular wing in hovering motion. Journal of Fluid Mechanics 693, 473–499. doi:10.1017/jfm.2011.543.

Demirdžić, I., Muzaferija, S., 1995. Numerical method for coupled fluid flow, heat transfer and stress analysis using unstructured moving meshes with cells of arbitrary topology. Computer Methods in Applied Mechanics and Engineering 125, 235–255. doi:https://doi.org/10.1016/0045-7825(95)00800-G.

Diaz-Arriba, D., Jardin, T., Gourdain, N., Pons, F., David, L., 2021. Numerical investigation of three-dimensional asymmetric hovering flapping flight. Physics of Fluids 33, 111907. doi:10.1063/5.0069840.

Diaz-Arriba, D., Jardin, T., Gourdain, N., Pons, F., David, L., 2022. Experiments and numerical simulations on hovering three-dimensional flexible flapping wings. Bioinspiration and Biomimetics 17, 065006. doi:10.1088/1748-3190/ac8f06.

Dickinson, M.H., Lehmann, F.O., Sane, S.P., 1999. Wing rotation and the aerodynamic basis of insect flight. Science 284, 1954–1960. doi:10.1126/science.284.5422.1954.

Dong, Y., Song, B., Xue, D., Yang, W., 2022. 3D numerical simulation of a hovering hummingbird-inspired flapping wing with dynamic morphing. Journal of Applied Fluid Mechanics 15, 873–888. doi:10.47176/jafm.15.03.33365.

Eldredge, J.D., Toomey, J., Medina, A., 2010. On the roles of chord-wise flexibility in a flapping wing with hovering kinematics. Journal of Fluid Mechanics 659, 94–115. doi:10.1017/S0022112010002363.

Heathcote, S., Martin, D., Gursul, I., 2004. Flexible flapping airfoil propulsion at zero freestream velocity. AIAA Journal 42, 2196–2204. doi:10.2514/1.5299.

Heathcote, S., Wang, Z., Gursul, I., 2008. Effect of spanwise flexibility on flapping wing propulsion. Journal of Fluids and Structures 24, 183–199. doi:10.1016/j.jfluidstructs.2007.08.003.

Jardin, T., 2017. Coriolis effect and the attachment of the leading edge vortex. Journal of Fluid Mechanics 820, 312–340. doi:10.1017/jfm.2017.222.
Jardin, T., Doué, N., 2019. Influence of pitch rate on freely translating perching airfoils. Journal of Fluid Mechanics 873, 49–71. doi:doi:10.1017/jfm.2019.421.

Jardin, T., Farcy, A., David, L., 2012. Three-dimensional effects in hovering flapping flight. Journal of Fluid Mechanics 702, 102–125. doi:10.1017/jfm.2012.163.

Kang, L., Bifeng, S., Zhihe, W., Ang, C., Dong, X., 2023. Effects of dynamical chordwise retracting and stretching on the propulsive performance of a pure plunging flat plate. AIP Advances 13, 025347. doi:10.1063/5.0136580.

Kodali, D., Medina, C., Kang, C.k., Aono, H., 2017. Effects of spanwise flexibility on the performance of flapping flyers in forward flight. Journal of The Royal Society Interface 14, 20170725. doi:10.1098/rsif.2017.0725.

Lee, Y.J., Lua, K.B., Lim, T.T., Yeo, K.S., 2016. A quasi-steady aerodynamic model for flapping flight with improved adaptability. Bioinspiration and Biomimetics 11, 036005. doi:10.1088/1748-3190/11/3/036005.

Lentink, D., Dickinson, M.H., 2009. Rotational accelerations stabilize leading edge vortices on revolving fly wings. Journal of Experimental Biology 212, 2705–2719. doi:10.1242/jeb.022269.

Lucas, K.N., Johnson, N., Beaulieu, W.T., Cathcart, E., Tirrell, G., Colin, S.P., Gemmell, B.J., Dabiri, J.O., Costello, J.H., 2014. Bending rules for animal propulsion. Nature Communications 5. doi:10.1038/ncomms4293.

Medina, A., Jones, A.R., 2016. Leading-edge vortex burst on a low-aspect-ratio rotating flat plate. Physical Review Fluids 1(4), 044501. doi:10.1017/S0022112010002363.

Muzafarija, S., 1994. Adaptive finite volume method for flow prediction using unstructured meshes and multigrid approach. Ph.D. thesis. University of London.

Oyama, A., Okabe, Y., Shimoyama, K., Fujii, K., 2009. Aerodynamic multiobjective design exploration of a flapping airfoil using a navier-stokes solver. Journal of Aerospace Information Systems 6, 256–270. doi:10.2514/1.35992.

Patankar, S.V., 1980. Numerical heat transfer and fluid flow. McGraw-Hill, New York.

Patankar, S.V., Spalding, D.B., 1972. A calculation procedure for heat, mass and momentum transfer in three-dimensional parabolic flows. International Journal of Heat and Mass Transfer 15, 1787–1806. doi:https://doi.org/10.1016/0017-9310(72)90054-3.

Sane, S.P., Dickinson, M.H., 2001. The control of flight force by a flapping wing: lift and drag production. Journal of Experimental Biology 204, 2607–2626. doi:10.1242/jeb.204.15.2607.

Sane, S.P., Dickinson, M.H., 2002. The aerodynamic effects of wing rotation and a revised quasi-steady model of flapping flight. Journal of Experimental Biology 205, 1087–1096. doi:10.1242/jeb.205.8.1087.

Sedov, L., 1965. Two-Dimensional Problems in Hydrodynamics and Aerodynamics. New York: Interscience.

Shyy, W., Aono, H., Chimakurthi, S.K., Trizila, P., Kang, C.K., Cesnik, C.E.S., Liu, H., 2010. Recent progress in flapping wing aerodynamics and aeroelasticity. Progress in Aerospace Sciences doi:10.1016/j.paerosci.2010.01.001.

Soto, C., Bhattacharya, S., 2023. The effect of dynamic twisting on the flow field and the unsteady forces of a heaving flat plate. Bioinspiration and Biomimetics 18, 026010. doi:10.1088/1748-3190/acb7ba.

Vanella, M., Fitzgerald, T., Preidikman, S., Balaras, E., Balachandran, B., 2009. Influence of flexibility on the aerodynamic performance of a hovering wing. Journal of Experimental Biology 212, 95–105. doi:10.1242/jeb.016428.

Wang, C., Liu, Y., Xu, D., Wang, S., 2022. Aerodynamic performance of a bio-inspired flapping wing with local sweep morphing. Physics of Fluids 34, 051903. doi:10.1063/5.0090718.

Wood, R.J., 2008. The first takeoff of a biologically inspired at-scale robotic insect. IEEE Transactions on Robotics 24, 341–347. doi:10.1109/TRO.2008.916997.

Wootton, R.J., 1981. Support and deformability in insect wings. Journal of Zoology 193, 447–468. doi:10.1111/j.1469-7998.1981.tb01497.x.