Numerical simulation of unsteady aerodynamic characteristics of the three-dimensional composite motion of a flapping wing based on overlapping nested grids

Numerical simulations of the unsteady aerodynamic characteristics of the flapping wing composite motion are performed. To avoid negative grid sizes arising with the use of a dynamic grid and leading to divergences in the simulation and to errors in the results, an overlapping nested grid is used for the flow field background, wing, and fuselage structure. The analysis is based on the Navier–Stokes equations (N-S) and the pressure–velocity coupling method, while spatial dispersion is handled using the second-order finite volume and the adaptive step size solving strategy. The lift and resistance generated by the wing for different combinations of flow velocity, flutter frequency and amplitude, and torsion angle are determined, and the aerodynamic efficiency and flow fields are compared to find the flapping parameters that give the best aerodynamic efficiency. The simulation results show that the aerodynamic lift of a flapping wing can be greatly increased by increasing the flapping frequency, while, for a fixed frequency, the lift can be further increased by increasing the flapping amplitude, although by only a small amount. Increasing the torsion angle in the flapping of the wing can also increase the lift, but the aerodynamic efficiency will be reduced if this angle is too large. Thus, an appropriate selection of flapping wing motion parameters can effectively increase the flight lift and improve the aerodynamic efficiency.


I. INTRODUCTION
Flapping wing aircraft are inspired by the flapping of the wings of birds and insects and perform, in some respects, better than traditional fixed and rotary wing aircraft, with high efficiency, good maneuverability, low energy consumption, and low noise. The aerodynamic properties of such vehicles are, therefore, of great current interest, in particular, with regard to the lift, drag, and thrust forces that arise from the flapping wing compound motion and the corresponding aerodynamic efficiency.
Efficient flight of birds and insects relies on maximizing lift and minimizing drag while saving energy and reducing weight. [1][2][3] Lightweight feathers, hollow bones, and membranous air sacs reduce weight, while streamlined bodies and small wingtips reduce parasitic and induced drag. During flight, updraft and shear wind speeds can be used to capture the return airflow and increase the flutter interval, which can effectively reduce energy loss and allow longdistance flights. Chordwise torsion and spanwise contraction during asymmetric flapping of the wing can obviously increase the lift and thrust generated by the wing, reducing airflow resistance during flapping. 4,5 The lift force can be increased by increasing the proportion of the downward movement during the flapping cycle.
There has been much research into the mechanisms of the flapping wing flight in birds and insects, and a number of important ARTICLE scitation.org/journal/adv results have been obtained. There are basically two approaches to the investigation of flapping wing flight. The first of these approaches involves observation and measurement of wing motion and analysis of the data thus obtained. From such analyses, empirical formulas can be derived for the relationships among the parameters of the flight attitude and wing motion. [6][7][8] Weis-Fogh and Jensen et al. studied the motion of insect wings and proposed that the flapping mechanism allows insects to maintain a high lift and achieve rapid mobility. 9 From an analysis of moth flight, Ellington and Sane et al. proposed a "delay stall" mechanism in which the leading edge vortex generated by the flapping of the wing persists throughout the motion, thus providing continuous high lift. [10][11][12] Zeng and Sun took high-speed photographs of the wings of bees during flight and conducted an aerodynamic simulation based on the motion parameters thereby obtained. They calculated the changes in the aerodynamic force and torque during the flight overshoot of the bees and to analyze their efficient maneuvering mechanism. 13,14 Von Kármán, Bratt, and Jones et al. found from studies of the generation mechanism of wing thrust that the thrust on a wing is determined by the direction and position of the tail vortex. [15][16][17][18] In bird flight, a row of tail vortices is generated at the trailing edge of each wing, and the wings will be pushed if the direction of the tail vortices is the same as the direction of the incoming flow.
In the second approach, flapping wings are studied by numerical simulation, experiment, and theory. Ansari et al. studied the aerodynamic and torque characteristics of the flapping wing under the conditions of different sizes and motion parameters of the single wing by means of numerical calculation and found that the flapping frequency and amplitude of the flapping wing could effectively improve the aerodynamic efficiency of the flapping wing. 19 Tuncer et al. performed numerical simulations to analyze the dynamic stall (flutter), flow reattachment on the surface of a flapping wing, flow at the trailing edge of a flapping wing, and wake behavior. 20 The effects of flexible wing deformation on the aerodynamic characteristics of a flapping wing have also been analyzed. 21,22 It should be noted that empirical formulas obtained from observations of flapping wing flight are not always able to give accurate values of the flight parameters. 23 In addition, most numerical analyses of flapping wings have been confined to the aerodynamic characteristics of two-dimensional single-motion modes, with much less work having been done on the analysis of unsteady aerodynamic characteristics of the three-dimensional composite motion of flapping wings. 24,25 However, some simulations of the unsteady aerodynamic characteristics of three-dimensional composite motion of wings of medium and large birds have been performed, and the results obtained have helped lay the foundation for the design and control of bionic flapping wing aircraft.
In this paper, using the wingspan parameters of medium and large birds, an aerodynamic model of a flapping wing under threedimensional rigid constraints is established with the aim of analyzing the unsteady aerodynamic characteristics of the composite motion of the wing during flapping wing flight. To avoid negative grid sizes arising with the use of a dynamic grid and leading to divergences in the simulation and to errors in the results, this paper adopts an overlapping nested grid technique for the flow field background, wing, and fuselage structure. With an unstructured grid of varying size, a very fine grid can be used around the wings to ensure high precision and accuracy of the simulation results. The method of solution uses the Navier-Stokes (N-S) equations combined with the Boussinesq hypothesis and the pressure-velocity coupling method, while for the spatial dispersion a second-order finite volume and an adaptive step size solving strategy is used. 26,28 The aerodynamic behavior of the flapping wing is analyzed for different values of flow velocity, flapping frequency and amplitude, and torsion angle, and the maximum aerodynamic efficiency of the aerodynamic model is determined in terms of a size parameter. The simulation results show that the aerodynamic lift of a flapping wing can be greatly increased by increasing the flapping frequency, and once the flapping frequency has been determined, the aerodynamic lift can be increased further by increasing the flapping amplitude, although the lifting range is not large. Increasing the torsion angle during wing flapping can increase the flight lift, but excessive torsion angles cause a decrease in the aerodynamic efficiency.

A. Numerical calculation method
Owing to the low flying speed of birds, the viscosity of air has a significant effect, and so the problem considered here belongs to the category of low-Reynolds-number flow. The N-S equations for three-dimensional incompressible flow are discretized by the finite volume method. In the flight of a bird, the wings generate lift force by their continuous movement. The resulting flow field is unsteady and governed by the dimensionless N-S equations where Eq. (1) is the continuity equation and Eq. (2) is the momentum equation; u, v, and w are the three components of dimensionless velocity; and p is dimensionless pressure. In dimensionless, the reference speed, length, and time are U, b, and b/U, and U is the average speed of the flap flapping for one cycle; b is the chord length. Re is the Reynolds number. The Reynolds number is given by where ρ is the air density, U is the average speed of the flap flapping for one cycle, b is the chord length, and μ is the kinematic viscosity coefficient of air. If the coordinates x, y, and z are written as xi (i = 1, 2, 3, respectively) and the corresponding components of the velocity u, namely,

ARTICLE
scitation.org/journal/adv u, v, and w, are written as Φi, then Eq. (2) can be put into the more concise form where u is micro-velocity vector, Φi represents u, v, and w, respectively, and the gradient operator is given by with i, j, and k being the unit vectors in the x, y, and z directions, respectively. Equation (4) is the general form of the N-S governing equation for the three-dimensional incompressible unsteady viscous flow field of a flapping wing in a Cartesian coordinate system. To ensure the accuracy of the time calculation for the unsteady flow field, a two-time-step implicit iterative method is adopted. Pressure-velocity coupling is used, with a Green-Gauss cellbased method being adopted for gradient interpolation in spatial dispersion. 29 Because wing flapping will produce thrust turbulence, the numerical calculation process uses the shear stress cloud (SST) model for turbulence modeling. 27 Additionally, to improve the airflow characteristics at the wing surface during wing movement, the meshes of the wing edge and surface are encrypted.

B. Flapping-wing aerodynamic calculation method
The calculation of the aerodynamic lift is based on a number of micro-elements, as shown in Fig. 1. It is assumed that the aerodynamic direction is always perpendicular to the wing airfoil during flapping flight. During flutter, the instantaneous velocity of a microelement is ωr, and the instantaneous aerodynamic lift generated by the micro-element is where ρ is the air density, ω is the flapping angular velocity, r is the distance from the airfoil element to the X-axis, CL(α) is the lift coefficient at angle of attack of micro-element, and c(r) is the chord length when the distance between the wing element and the X-axis is r. The instantaneous aerodynamic force of the entire flapping wing can be calculated by integrating over all micro-elements of the flapping wing element as follows: r is the ratio of the distance between the wing element and the X-axis to the wingspan,ĉ(r) is the chord length when the airfoil microelement isr,c is the average chord length of the flapping wing, R is the wing length, and Swing is the airfoil area. Similarly, the aerodynamic drag FD, thrust F T , power Pw, and propulsion efficiency η during flapping can be obtained as follows: where CD(α) is the drag coefficient at an angle of attack of microelement, C T is the thrust coefficient, C T is the flutter power factor, M is the flapping wing torque, H(t) is the expression for the time dependence of flapping motion, and θ(t) is the expression for the time dependence of torsional motion.

III. NUMERICAL SIMULATION MODEL
A. Three-dimensional composite motion model of the flapping wing The simulation model adopted in this paper is shown in Fig. 2. It comprises of three parts: the left and right wings and the body (the fuselage, if this were an aircraft). The left and right wings flutter around the axis of the body and at the same time can achieve synchronous twisting. The coordinate system O-XYZ is established with the center of mass of the body as the coordinate origin. The X axis coincides with the axis of the body, with the positive direction being taken as that pointing from head to the tail; the Z axis is perpendicular to the X axis, pointing toward the left wing, and the Y axis is then determined according to a right-hand coordinate system. U∞ is the vector representing the wind speed and direction. During wing flapping, if the wing is located above the OXZ plane, the flapping amplitude angle α is positive, and vice versa. The torsion angle β is viewed from the positive direction of the Z axis, with clockwise being positive, and counterclockwise negative. Wing_Right, Wing_Left, and Wall_Body represent the right wing, the left wing, and the  Table I. Flapping of the wing is a composite motion consisting of span flutter and chord torsion. The lift and thrust are increased by chord torsion during flapping, the air resistance is reduced, and the aerodynamic efficiency is improved. The wing movement is governed by the following equations for flutter and chord torsion, respectively: where H(t) is the expression for the time dependence of flapping motion and θ(t) is the expression for the time dependence of torsional motion. A and f are the flutter amplitude and frequency, B is the torsion amplitude, the deflections ξ 0 and ψ 0 represent the initial phases, andα 0 and β 0 are the initial flutter and torsion angles. In this paper, we take α 0 = 0 and β 0 = 0. According to the flapping wing aircraft prototype designed by us, its maximum flapping frequency is 5 Hz, maximum flapping amplitude is 45 ○ , maximum flying speed is 5 m/s, and the maximum wing chord torsion angle is 10 ○ . Therefore, in the process of numerical simulation, the influence of the air inlet velocity (1-5 m/s), flapping frequency (1 Hz-5 Hz), flapping amplitude (15 ○ -45 ○ ), and torsion angle (0-10 ○ ) on the overall aerodynamic performance of the flapping wing is analyzed, which is the basis of flapping wing flight dynamics modeling and flapping wing flight attitude control methods.

B. Overlapping meshing and boundary conditions
The size of the calculation region is taken to be sufficiently large (15-20 times the size of the body and the spanwise direction) that the influence of the boundary conditions on the flow field can be ignored when solving the equations for the fluttering motion to calculate the flow speed. The flow field dimensions thereby determined are 12.5 m × 9.3 m × 6.1 m, as shown in Fig. 3(a). To ensure an accurate calculation of the aerodynamic characteristics at the wing surface, the mesh size there needs to be encrypted. Since the overall flow field is large, if a uniform grid size was adopted throughout, stringent requirements would need to be imposed on calculation time and computer hardware. Therefore, to ensure the continuity of the flow field and accuracy of the calculation at the airfoil, and to increase the calculation speed, an O-shaped mesh is used for the outer flow field and adopts an unstructured tetrahedral mesh for the wing airfoil. The mesh scales at the interface between the airfoil and the flow field are the same to ensure continuity of the flow field, as shown in Fig. 3(b). The basic parameters of each part of the grid are shown in Table II.
Thus, the flow field and each part of the flapping wing model are covered by overlapping grids, so it is necessary to assemble the grid for each part in an appropriate manner when constructing the overall simulation model. At the same time, an overlapping interface between each grid needs to be established.
As the boundary condition on the flow velocity at the wall surface, the component in the positive X direction is taken to be equal to the oncoming flow velocity, while the components in the Y and Z directions are taken to be zero, i.e., UX = U∞, UY = 0, and U Z = 0. The pressure far-field boundary condition is assumed at the outlet. The no-slip boundary condition is applied at the surface of the flapping calculation model, i.e., the oncoming flow velocity near the wing surface is taken to be equal to the wing wall motion velocity.

A. Overall average aerodynamic characteristics of flapping wings
To analyze the overall aerodynamic performance of a flapping wing, the influence of various parameters needs to be considered. To achieve higher flight efficiency and in response to changes in wind strength and direction, birds adjust the amplitude, frequency, and torsion angle of flutter to change the lift generated by their wings and thereby reduce the resistance. Therefore, we investigate the effective lift and resistance that can be generated by different flutter amplitudes and frequencies and different torsion angles under the influence of different wing speed disturbances, and we determine the aerodynamic efficiency and the load capacity of the wings that can be achieved under these different conditions.
The aerodynamic performance of a flapping wing is affected by the flapping frequency and amplitude, the incoming velocity, and the torsion angle. Therefore, in the following analysis, the influence of each of these factors on aerodynamic performance is examined in turn, with the others held constant. Figure 4 shows the aerodynamic characteristics of the wing when the torsional angle is held constant at θ = 0 ○ , while the flapping frequency and amplitude and the incoming flow velocity vary (the data are averaged over two beat cycles after transients have stabilized).
As can be seen from the variations in the average lift in Fig. 4(a) The lift generated by the wings can thus be effectively increased by increasing the flapping frequency, although this results in an increase in energy consumption. Figure 4(b) shows the variations in the overall average drag of the flapping wing, with the negative sign indicating that the drag is opposite to the flow velocity. It can be seen that the drag increases with increasing flapping amplitude and frequency, with the frequency having the stronger influence.
From Fig. 4, it can be seen that at a frequency f = 1 Hz and amplitude φ = 15 ○ , the average lift, drag, aerodynamic efficiency, and load capacity do not vary greatly as the other parameters change. However, as the flutter amplitude and frequency increase, the aerodynamic forces generated by the wings change significantly. The flutter frequency has the most obvious effect on the wing aerodynamics. When f = 3 Hz and φ = 45 ○ , the aerodynamic efficiency reaches 36.37. However, the effective load of the wing under this parameter is only 1562.3 g (including its own weight). When the flapping frequency of the wing increases to 5 Hz, the overall aerodynamic efficiency is reduced due to the increase of air resistance during flapping. However, the average lift force generated by the wing at this time is 2.2 times of that at f = 3 Hz. Under this parameter, the payload capacity of the wing reaches 3451 g, which significantly increases the payload capacity. It can be concluded that the flapping frequency of the wing can be increased to obtain a larger payload capacity during the flapping flight. However, if the flapping frequency is too large, the system vibration will be increased, and the aerodynamic efficiency will be reduced; thus, the flapping wing flight stability will be reduced and energy consumption will be increased. Therefore, the flapping frequency should be adjusted appropriately according to different flight conditions.
During the flight of a bird, the wings are twisted through a certain chord angle during a flapping cycle. Figure 5 shows the aerodynamic characteristics of the wing as a function of frequency, amplitude, and oncoming flow velocity when the torsion angle θ = 5 ○ .
From a comparison of the average lift and resistance generated in the absence of torsion and at a torsion angle θ = 5 ○ (Figs. 4 and 5, respectively), it can be seen that the chord torsion helps to provide lift during flight when the wing flaps up and down, but at the cost of an increase in flight resistance. Figure 5

ARTICLE scitation.org/journal/adv
Although wing torsion increases the flight lift, it also leads to an increase in the aerodynamic drag. When the flapping frequency f = 1 Hz, the average resistance generated by flapping increases with increasing flapping amplitude and oncoming velocity. The average resistance generated by the wing is greatest for U = 5 m/s and φ = 45 ○ , with a value of −1.11 N (the negative sign indicates that the resistance is in the opposite direction to the oncoming flow) for a torsion angle θ = 5 ○ , while the average resistance for the same parameters but in the absence of torsion is −0.229 N. At higher flapping frequency, the average resistance decreases with increasing flow velocity. When f = 3 Hz, the average resistance for U = 1 m/s, φ = 45 ○ , and a torsion angle of the flapping frequency has a maximum negative value of −2.8 N, compared with only −0.618 N in the absence of torsion. For the parameter combination that gives the maximum average resistance of −8.23 N in the presence of θ = 5 ○ torsion. However, f = 5 Hz, U = 1 m/s, and φ = 45 ○ , the average resistance without torsion is only −1.816 N. It can be found that the chord wise torsion of the wing flapping up and down can increase the lift generated by the wing flapping, but increasing the lift also increases the drag. This is due to the increase in the projected area of the wing in the forward flight direction due to the chord wise torsional motion of the wing, resulting in much greater resistance than when there is no torsion.
From the above analysis, it can be seen that if the wing is subjected to chord torsion while flapping up and down, this

ARTICLE
scitation.org/journal/adv increases both lift and resistance. To understand the influence of wing torsion on the overall aerodynamic performance of a flapping wing, it is necessary to analyze the variations in aerodynamic efficiency and load capacity, as shown in Fig. 6. It can be seen that when the wing twists, the aerodynamic efficiency of the wing is much higher for f = 5 Hz and U = 5 m/s than for other parameter combinations. When the flutter angle φ = 45 ○ , the efficiency reaches a maximum of 27.6, and the lift generated by the wing at that time is capable of carrying a weight of 3966 g in flight. Without torsion, the aerodynamic efficiency of the wing for this parameter combination is 22.03, and the wing can carry 3294.4 g. Although the aerodynamic efficiency of a wing with no torsion can reach 36.4 for f = 3 Hz, U = 5 m/s, and φ = 45 ○ , it then has a load capacity of just 1519.4 g, well below that with torsion. According to the aerodynamic efficiency curve of the flapping wing under different parameters, the chord wise torsion of the wing increases the average resistance, so that the overall aerodynamic efficiency is lower than when there is no torsion. However, the wing chord torsion can effectively improve the flying load of the flapping wing. Therefore, it is necessary to determine the appropriate flapping frequency and chord torsion angle in the subsequent flapping wing flight control, so as to achieve higher flapping wing flight efficiency. When the torsion angle θ = 10 ○ , the aerodynamic characteristics of the aerodynamic force of the wing in a flapping cycle for different parameter combinations are shown in Fig. 7. It can be seen that the lift force increases with increasing flutter amplitude when the flutter frequency is constant. The larger the flutter amplitude, the greater is the lift generated. When the flutter angle φ = 45 ○ , the maximum lift generated at flutter frequencies of 1 Hz, 3 Hz, and 5 Hz is 6.977 N, 19.432 N, and 58.89 N, respectively: when the frequency increases from 1 Hz to 5 Hz, the lift generated by the fluttering process increases by a factor of 8.44. Thus, the flutter frequency is an important factor determining the lift of the flapping wing. To obtain a greater lift, the overall aerodynamic performance can be changed by increasing the flutter frequency and amplitude, which is consistent with observations of birds, which increase both the flutter frequency and wing swing angle when they need greater lift. From an analysis of the aerodynamic data in Figs. 4-7, it is found that adding chord torsion to the flapping of a wing can effectively increase the lift force generated, but at the cost of increased aerodynamic resistance and reduced overall aerodynamic efficiency. Therefore, in the design of a bionic flapping wing aircraft, to increase the lift generated by flutter, the torsion angle of the wing can be increased, but not by too much, and the effect on aerodynamic efficiency needs to be taken into account to achieve the best flight performance.
The results of the above numerical simulations of a flapping wing under different conditions have shown that the flutter frequency and amplitude and the torsion angle of the wing have a great influence on lift and resistance. Therefore, the choice of appropriate flutter parameters is key to improving lift and increasing thrust.

B. Flow field analysis of flapping wing motion
The design of a bionic flapping wing aircraft requires the flapping wing to have a large load capacity and high flight efficiency. The analysis in Subsection IV A has shown that by increasing the flutter frequency and amplitude and the torsion angle, the effective lift generated during the flapping process can be increased to improve the load capacity of the wing. However, for a flapping wing aircraft with a large wingspan, an excessive torsion angle increases the flight resistance, resulting in a decrease in the overall aerodynamic efficiency. Therefore, for the following analysis, we choose the parameter values f = 5 Hz, Ux = 5 m/s, φ = 45 ○ , and θ = 5 ○ . With this combination of parameters, the payload capacity of the wing can reach 3.9659 kg, with an aerodynamic efficiency of 27.55, which meet the aerodynamic characteristics required for a bionic flapping wing aircraft of the given size.  Figure 9 shows pressure nephograms at several typical moments during one flapping cycle. The pressure nephograms of the YOZ section is on the left and the right is the pressure nephograms of the upper and lower wings. At the initial moment, the wing is in the middle position, from which the flapping starts. This was immediately preceded by the end of the last flapping period, so there is a high-pressure area on the lower surface of the wing at T 0 /12, while a low-pressure area forms on the upper surface. Then, the wing moves downward, and a starting vortex forms from the upper surface. The impact of the oncoming flow below the wing causes this starting vortex to develop into a leading-edge vortex and awake vortex. The leading edge vortex extends from the wing root to the tip, and it continues to grow while always remaining attached to the upper airfoil, as a result of which a low-pressure area forms on the upper surface of the wing. In Fig. 9, from T 0 /12 to T 0 /4, there is a clear change in the upper airfoil vorticity, and finally a low-pressure area forms on the upper airfoil and grows larger as the wing is lowered. At T 0 /4, the downward stroke of the wing is fastest, the leading-edge vortex is strongest, and the pressure difference between the upper and lower wing surfaces reaches its maximum value, with the result that a peak occurs in the lift, as shown in Fig. 8(a). During the flapping, the wing is twisted, and the fluttering also produces a thrust in the forward direction, at which point the thrust also reaches a maximum. As the leading edge vortex gradually develops and expands toward the wing tip, the leading edge vortex and the wake vortex fall off at the end of the wing downstream, reducing the upper and lower pressure difference.
During the upward stroke of the wing, the pressure on the airfoil becomes larger, and the lower wing surface has a low-pressure area from the root. With increasing flutter amplitude and flutter velocity, the high-pressure vortex on the upper airfoil separates from the latter. During the entire upward stroke, the pressure difference between the lower wing and upper wing surfaces is negative, resulting in a large resistance during this stroke. However, as a result of torsion during the flapping of the wing, the angle of attack and the projected area of the airfoil change continually during the upward stroke, so the resistance is reduced. The flapping angle is large, and the equivalent negative angle of attack of the wing relative to the incoming flow is reduced, so the adverse effect of the upward stroke is weakened.
The instantaneous flow lines of the model at different times during a flapping cycle are shown in Fig. 10. The wing starts to flutter from the middle position, and a vortex appears at the root of the wing. As the flutter angle and flutter speed increase, the vortex on the upper surface of the wing gradually increases in strength while remaining attached. At time T 0 /4, the leading edge of the wing reaches its lowest point. The vortex detaches from the tip of the wing, as can also be seen from the nephograms in Fig. 9. During the upward stroke, the pressure on the upper wing surface is large, and a starting vortex appears at the root of the lower wing surface. With increasing flutter amplitude and flutter velocity, a vortex line gradually develops and expands toward the wing tip, reaching the tip at time 3T 0 /4, which is consistent with what is observed in the flight of medium and large birds. Figure 11 shows the vortexes generated by wing flapping for one period. At the end of the flapping motion of the wing, a devortex will occur at the tip of the wing, as shown in Fig. 10 at the 3T 0 /4 moment, where there is a vortex at the tip of the wing. When the wing flaps downward from the highest position, the leading edge vortex appears at the leading edge of the wing. As the wing flaps faster, the leading edge vortex continues to develop toward the wing tip and the trailing edge of the wing and disengages from the airfoil when the wing flaps to the lowest position. As shown in Fig. 11, the leading edge vortex of the wing gradually develops toward the wing tip and trailing edge positions. The vortex is disengaged from the airfoil when the wing is moving to the lowest position.

V. CONCLUSIONS
In this paper, the aerodynamic changes occurring during the three-dimensional composite motion of flapping wings in an unsteady flow field have been investigated for different parameters of the motion. The mechanisms responsible for the generation of lift and drag have been determined, and preliminary values of the relevant parameters have been provided to aid in the design of bionic flapping wing aircraft. The following conclusions can be drawn from this study: • By increasing the flapping frequency of the wing, the overall time-averaged aerodynamic force and torque can be significantly changed. To obtain a large amount of lift in a short time, the flutter frequency can be increased, although too high frequency is not suitable for long flights. • The flapping amplitude can affect the instantaneous lift and drag on the wing during the entire flapping cycle. However, for a fixed flapping frequency, an increase in the flapping amplitude has little effect on the lift. At the same time, increasing the flapping amplitude will increase the overall vibration of the aircraft and decrease the stability of the system. • For a bionic flapping wing aircraft with large wingspan, an appropriate value of the torsion angle can be chosen, bearing in mind that although increasing the torsion angle will increase the lift force and forward thrust generated by the wing, although an angle that is too large will lead to a reduction in the overall aerodynamic efficiency.