A Computational Study on the Aerodynamic Influence of a Propeller on an MAV by Unstructured Overset Grid Technique and Low Mach Number Preconditioning

The influence of a propeller on the aerodynamic performance of an MAV is investigated using an unstructured overset grid technique. The flow regime of a fixed-wing MAV powered by a propeller contains both incompressible regions due to the low flight speed, as well as compressible flow areas near the propeller-tip region. In order to simulate all speed flow efficiently, a dual-time preconditioning method is employed in the present study. The methodology in this paper is verified as providing a reliable numerical simulation tool for all flow regimes, in the additional presence of moving boundaries, which is treated with an overset grid approach.


INTRODUCTION
For a propeller driven MAV, the ratio of the propeller's diameter to the dimensions of the wing or fuselage is relatively larger than for conventional aircraft, which means that most of the MAV is submersed in the slipstream of the rotating propeller.Therefore, the powered effect cannot be ignored during the simulation of such MAV configurations.
Precise modeling of a propeller is challenging, but crucial for a variety of reasons.The rotational motion of the blades with respect to the fuselage causes that the aerodynamic problem is naturally unsteady.The most basic method for modeling the influence of a propeller is by means of introducing an actuator disk [1], where the propeller is represented as an infinitely thin disk capable of sustaining a pressure discontinuity.Only a few works has been reported where such an approach was used for an MAV configuration.Choi and Ahn [2] studied the aerodynamic influence of a pusher propeller on a MAV by modeling the propeller as an actuator disk, and found that the slipstream has a significant influence on the performance of the MAV.Ageev [3] investigated the aerodynamic performance of a disc-wing MAV with propeller in a wing slot.However, the actuator disk methods used by both Choi and Ageev did not consider the rotation effects of the propeller, which limits the simulation of propeller effects.Because the actuator disk model turns the unsteady flow into a quasi-steady one, significant information is lost, which inspires a full modeling technique with regard to the propeller flow.Dynamic mesh techniques have been applied to solve flow problems involving geometric deformation or relative motion of multiple bodies.Batina [4] developed a spring analogy scheme which models the grid as linear springs and solves the equilibrium equations to determine the new location of the grid nodes.However, this method would result in grid crossing and negative volume when large motions occur in the flow-field.Farhat et al. [5] modified the spring analogy by adding non-linear torsion springs to avoid grid crossing problem, but the torsion springs model is difficult to retain the mesh quality.Liu et al. [6] proposed a fast deformation method based on the creation of a Delaunay graph of the original mesh, which is verified as a robust and efficient approach for dynamic mesh problems.However, both the spring model and graph mapping methods suffer from the problem of potential grid crossing for large motion, they cannot guarantee mesh quality.A robust approach to deal with the moving boundaries is the overset grid technique, which was first proposed by Steger [7], and subsequently extended by Nakahashi et al. [8] for application on unstructured grids.By hole-cutting and interpolation in overlapping meshes the negative volume can be completely avoided.These properties motivate the use of the overset grid technique in the present study.
In Computational Fluid Dynamics a rigid division exists between numerical methods for the computation of compressible and incompressible flows.For an incompressible flow, with Mach number say below 0.3, incompressible governing equations are employed.When the Mach number increases above 0.3, compressible flow governing equations need to be employed.However, for predominantly low speed flow problems with embedded high-speed regions, a purely incompressible method would not be appropriate.For compressible flow, whenever the Mach number is very small, there is a huge disparity between the convective and acoustic eigenvalues in the equations system, which will lead to poor convergence properties as well as poor solution accuracy.There are two main approaches advocated in the development of algorithms for the computation of low Mach number flows: 1) Extending incompressible solvers towards the compressible regime and 2) Modification of compressible solvers to allow for low Mach number.Several methods have been reported in literatures.The pressure-correction method which was introduced by Choin [9], has been employed effectively within several implementations like SIMPLE [10].Karki and Patankar [11] developed the SIMPLER method for compressible flows, applicable for a wide speed region.A similar approach has been adopted by Bijl and Wesseling [12], Mary et al. [13].Bijl [14] developed a computational method for all-speed flow with a staggered scheme based on pioneering work of Harlow and Amsden et al. [15,16].In the region of density based methods, two distinct techniques have been proposed to capture solution convergence for low-Mach number regimes, namely asymptotic and preconditioning.The asymptotic scheme introduces a perturbed form of the equations.Preconditioning is to pre-multiply time-derivatives by a suitable preconditioning matrix.Due to the high efficiency of preconditioning, it has becomes the method of choice for solving flow fields involving multiple scales.By modifying the governing equations without altering the steady state solution, the eigenvalues will be in the same order.Such work can be found in refs.[17][18][19][20].For unsteady problems, the preconditioning is introduced within a dual-time framework wherein the physical time-derivatives are used to march the unsteady equations and the preconditioned pseudo time-derivatives are used for numerical discretization and iterative solution.In this paper, the dual-time preconditioning framework will be employed in the simulation.
This paper is organized as follows.In section 2 the methodologies employed in this paper are presented, including the overset grid and dual-time preconditioning approaches.Subsequently, several validation cases of the methods are presented to assess the feasibility and reliability of the solver in section 3. Numerical simulation of the propeller-wing interaction of an MAV is performed in section 4, followed by conclusions in section 5.

Unstructured Overset Grid Technique
Unstructured overset grid methodology is a mature flow simulation technique that has been used for over 20 years to simplify grid generation for complex geometries as well as an enabling technique for the simulation of bodies in relative motion [21].As described in ref. [8], the overlapping technique consists of two main procedures: Chimera-hole cutting and the intergrid-boundary definition.First, the distance from each node in the volume grid to the wall boundaries (named as wall-distance in literature) is measured.Secondly, the nodal activity can be determined for all nodes by comparing the wall-distances.Thirdly, the cell activity can be determined by the nodal properties.The unstructured overset grid scheme used in the present paper is based on referenced literatures [22][23][24].

Governing Equations
The time-dependent Reynolds-Averaged Navier-Stokes (N-S) equations for moving bodies can be written in integral form as follow: Here, W=( , u, v w e) T , F(W) is the convective fluxes, F v is the viscous fluxes, vg is the motion of the control volume V.In order to close the N-S equations, the Spalart-Allmaras (S-A) turbulence model is employed.

Preconditioning Dual Time-Stepping Algorithm
In order to obtain reliable physical solutions for unsteady flow with low Mach preconditioning, it is necessary to use a pseudo time approach [20].Here, a preconditioning pseudo time derivative term is added in Eqn(1) as follows: For optimizing the performance of the pseudo iterations, the pseudo time derivative should be written in terms of the primitive variable vector Q=(P, u, v, w, T) T by introducing a preconditioning matrix as below: According to Venkatesearan et al. [25], and Philip et al. [26], the low Mach preconditioning matrix form of is defined by: where p ' = 1 ' 2 , and is the speed of sound.Here , where is the preconditioning parameter.
The definition of the preconditioning matrix is determined by the parameter of , in the present case, this parameter is defined as follows: = Min(Max(M max 2 , min ),1.0)… (5) here, M max is the local maximal Mach number defined as: where V and V g are the local fluid and grid moving velocities, respectively.Subscripts cell and neighbors represent the current control volume and its adjoined cells.min = 3MAX(Ma 2 , Ma g,mean 2) … (7) where Ma is the Mach number of free stream, while Ma g,mean is the mean Mach number of the moving boundaries.
Hence, the eigenvalues of the matrix can be obtained as follows: where is the flow velocity normal to the surface of the control volume: = V in = un x + vn y + wn z .
When =1 no preconditioning is applied.
Xiao [27] gives a detailed description spatial and temporal discretization schemes.

VALIDATION
In order to apply the overset grid method and preconditioning approach described above, a solver was developed at the MAV Research Center in Nanjing University of Aeronautics and Astronautics (NUAA).To assess the accuracy of the solver, two cases are studied: 1) Flow around a 2-D NACA0012 airfoil is simulated to verify the accuracy of the preconditioning method in steady flow.
2) A plunging NACA0012 in free stream is employed to assess the feasibilities of the overset grid scheme as well as the unsteady dual-time strategy.

NACA0012 Steady Flow
To begin with the steady flow around a NACA0012 airfoil is used to assess the accuracy of the computational approach.The C-mesh used in this case is shown in Fig. (1a).Pressure coefficient distributions for four different Mach numbers (0.6, 0.1, 0.01 and 0.001) are shown in Fig.
(2) both with preconditioning and without preconditioning, respectively.As we can see from Fig. (2a, e), for Ma=0.6, preconditioning gives essentially the same pressure distribution as when no preconditioning is applied.As expected the preconditioning method has no advantages for high Mach number flow.With decreasing Mach number as shown in the other figures, the contour lines start to display irregularities, especially notable for Ma=0.001,Fig. (2d).With the preconditioning method, on the other hand, smooth contour lines are obtained, even for the lowest Mach number value as visible in Fig. (2h).Also, the preconditioning will increase the calculation efficiency, as evidenced in Fig. (1b), which shows that the computation energy residual with preconditioning has decreased six orders while the no-preconditioning case only has decreased three orders after 5000 physical iterations.

NACA0012 Plunging in Flow
This case is targeting to test the accuracy of the overset lapping scheme and dual-time preconditioning for unsteady flow.A periodic plunging NACA0012 foil is simulated in the present solver, with the motion parameters given in Table 1.An overview of the overlapping meshes can be seen in Figs.(3,4) plots the wake vortex pattern between the present simulation and the results from study of Platzer and Jones [28].With the same vortex contours, the solver can be viewed a reliable tool to simulate unsteady flow which contains moving boundaries.

Computational Set-Up
An MAV with an inverse-Zimmerman wing layout developed by the MAV research center in Nanjing University of Aeronautics and Astronautics (MAV RCNUAA) has been studied in both numerical and experimental way.An image of the MAV configuration and its relevant dimensions can be found in Fig. (5) and Table 2.
For the numerical simulation a 3-D model of the MAV and its propeller were established in CATIA (scaled with 1:1) as shown in Fig. (6).Note that the propeller in Fig. (6) is not exactly the same as the real propeller, since that we cannot measure all the specifications of the real one.The rotation direction is opposite which does not influence the lift, drag and pitching moment in the free stream direction.The rolling moment of the MAV is not considered in the present study due to the limitation of the solver.
The calculation meshes were generated using the grid generator software Gambit in view of the simple configuration of the MAV.Fig. (7) illustrates the boundary conditions and mesh scheme.The mesh of the wing is regarded as the background mesh in the overlapping grid.The mesh of the propeller is moving with the propeller, and information translation is carried out by interpolation in the intergrid boundaries.

Results
Unsteady flow computations have been performed in two cases with propeller on and off, and with the angle of attack ranging from -4 o to 40 o with 4 o increment.The velocity of the free stream is 12m/s, while the Reynolds number is 1.15 10 5 .The propeller rotates right-handed along the +x direction with a RPM=16.000 in the computational domain.
Wind tunnel experiments on the MAV have been performed in the low-speed wind tunnel in NUAA, where timeaveraged aerodynamic forces with propeller on or off were measured during the experiment.Experimental data are used to assess the accuracy of the computational results.

Aerodynamic Coefficients
The aerodynamic coefficients with the propeller on or off, and also experimental data are plotted in Fig. (8).The time-average of the coefficients will be used in the subsequent discussion, and the comparison to the experimental data.It is easily to observe that the computational results have larger stall angles (4 o more) than  coefficient, in order to guarantee the similarity the measured point is located at the 0.3 chord length from the leading edge on the mid-plane.The influence of the separation suppression is commonly regarded to be beneficial for the performance of the aircraft, however, the lift-to-drag ratio decreases with the propeller on as plotted in Fig. (8d), that is because the propeller increases the turbulence intensity level on the wing, which in turn increases the viscous drag.At low angle of attack this may have an adverse effect on the lift-todrag ratio curve, while at higher angle of attack the delay of separation, due to the propeller, helps the overall drag to decrease which has a favorable effect as seen in Fig. (8d), for Angle of Attack >15 degree.All the aerodynamic coefficients of the simulation nicely correspond with the experimental data which proves the solver as a reliable tool.

Surface Flow Characterization
For a low aspect ratio and highly swept back wing as used in the present study, flow separation behaves more complex than on a conventional airplane.This section will focus on the flow separation under the influence of the propeller.Fig. (9) indicates the pressure distribution on different span-wise sections, which reveals that the "Plateaus" in the surface pressure distribution caused by the boundary layer separation (Fig. 9a) disappear when the MAV is mounted with a powered propeller, that is because the complex slipstream caused by the high speed rotating propeller changes the flow behavior on the upper surface.In order to provide more insight to the flow separation phenomenon, the topological approach is employed.Legendre [29] proposed that a pattern of skin friction (also called limited stream) could be viewed as trajectories of flow particles close to the body surface, having properties consistent with those of a continuous vector field.As a typical example, the skin friction on both upper and lower surfaces is shown in (b) lower surface, left: without propeller; right: with propeller for =8 o .It is easy to see that in the case of propeller on the skin friction pattern becomes asymmetric because of the propeller swirling effect.On a clear fuselage (b) lower surface, left: without propeller; right: with propeller (Fig. 10a, Left), four confluence lines appear on the upper surface, with a pattern of twice separation and twice reattachment.

Tip Vortex
For low-aspect-ratio wings, the lift can be divided into two parts as 2-D linear lift which is respect to the airfoil and 3-D non-linear due to the wing-tip vortex.

Propeller Wake
The flow field which contains a rotating propeller is typically an unsteady situation.In order to obtained the unsteady characteristics of the flow induced by the propeller, a simulation on a single rotating propeller is carried out, refer to Fig. (13a), in which the vortex ring behind the propeller is captured in isolation.The dimension of the vortex ring decreases with increasing distance from the propeller due to the contraction of the slip stream tube.Fig. (13b, c) give a clear picture of the vortex contour of the 3-D flow field where propeller and fuselage are both involved with two different angles of attack.Most part of the MAV is submersed in the vortex ring, which means the propeller has a significant influence on the wing.Note that when in a low angle of attack the MAV is mainly influenced by the propeller vortex ring, as the increasing angle of attack, wingtip-vortex is generated suddenly.

CONCLUSIONS AND FUTURE WORK
Dual-time preconditioning and overset grid methodologies are employed during the present study, motivated by the special flow characteristics of the propellerwing interference.The solver was developed based on the above two approaches and also is approved as a reliable and efficient tool for unsteady flow field which contains moving boundaries, especially in low Mach regime.The aerodynamic influence of a propeller on a fixed wing MAV has been studied in this paper.Conclusions are as follow: 1. Aerodynamic coefficients of the MAV increase because of the influence of the propeller.

2.
Stall angle increases 2 degrees with a powered propeller which can delay the flow separation.

3.
Lift to drag ratio is slightly decreased which is contrary to commonly known benefits of the delayed flow separation.

4.
Wing-tip vortex gives a significant influence on the aerodynamic of a MAV

5.
Vortex ring generated by the propeller can significantly influence the flow field around the wing.
Since now the limiting solver cannot obtain enough information of the flow performance, future work will engage on the development of the present methodology to provide a full understanding of the effect of flow fields which contain moving boundaries.

Fig. ( 8 )
Fig. (8).Aerodynamic Coefficients.The measured point of the moment coefficient is located at the 0.3chord of the MAV from the leading edge on the mid-plane.
Fig. (9b, c) show the pressure distribution of both sides of the powered MAV.The pressure on the upper side is almost the same while only different distribution occurs on the leading edge of the wing.
Fig. (10a, right) indicates that the first separation is delayed a small distance and the other confluence lines are disappeared because of the slipstream produced by the rotating propeller.On the lower surface as shown in Fig. (10b), the limited streams are similar to each other in most part of the surface, only a delayed separation occurs in the tail region when the MAV is operated with a rotating propeller.Obviously, flow separations retarded on both sides of the MAV caused by the propeller give a larger stall angle.
Fig. (11) shows the pattern of wing-tip vortex in different angle of attack with rotary propeller, while the vortex generated by the propeller is ignored.The vortex generated by the wing tip occupies a huge part of the upper surface of the MAV.The strength and dimensions of the vortex increase as the increasing angle of attack, that provide extra energy to submerse more air to delay the flow separation, which is the main reason why the low-aspect-ratio wing displays a large stall angle.The flow streamlines near the wing tip during one rotation cycle of the propeller on the angle of attack of 34 o are plotted in Fig. (12), as we can see the flow behavior does not change too much during one cycle.