# Introduction FD study of airfoils to predict its lift and drag characteristics, visualisation and surveillance of flow field pattern around the body, before the endeavour of the experimental study is almost patent. In the present study aerodynamic characteristics of a welldocumented airfoil, NACA 0012, equipped with plain flap is investigated. Wing with flap is usually known as high lift device. This ancillary device is fundamentally a movable element that supports the pilot to change the geometry and aerodynamic characteristics of the wing sections to control the motion of the airplane or to improve the performance in some anticipated way. CFD facilitates to envisage the behavior of geometry subjected to any sort of fluid flow field. This fast progression of computational fluid dynamics (CFD) has been driven by the necessity for more rapid and more exact methods for the calculations of flow fields around very complicated structural configurations of practical attention. CFD has been demonstrated as an economically viable method of preference in the field of numerous aerospace, automotive and industrial components and processes in which a major role is played by fluid or gas flows. In the fluid dynamics, For modelling flow in or around objects many commercial and open source CFD packages are available. The computer simulations can model features and details that are tough, expensive or impossible to measure or visualize experimentally. These devices are primarily used to improve the maximum lift coefficients of wings with changing the characteristics for the cruising and high-speed flight conditions. As a result, it is very important to understand the characteristics of the wing having different flap angles (?) at different Mach number (M). Operating the aircraft at optimum flap angle at optimum velocity may result significant amount of fuel saving. The B-17 Flying Fortress, Cessna 152 and the helicopter Sikorsky S-61 SH-3 Sea King as well as horizontal and vertical axis wind turbines use NACA 0012 airfoil which place this specific airfoil under extensive research and study. This study does not provide any experimental data for the flow over the flapped airfoil. Therefore, to reduce the scepticism associated the results obtained, the simulation process for the study is validated instead. In the validation course the results for flow over no flapped NACA 0012 is compared with published standard data by NASA [1], as nearly same computational method is used to study flapped NACA 0012 airfoil. Many researchers have studied aerodynamic characteristics of NACA 0012 using different methods and operating conditions. The Abbott and von Doenhoff data [7] were not tripped. The Gregory and O'Reilly data [10] were tripped, but were at a lower R e of 3 million. Lift data are not affected too significantly between 3 million and 6 million, but drag data are [11]. Selecting a proper turbulence model, the structure and use of a model to forecast the effects of turbulence, is a crucial undertaking to study any sorts of fluid flow. It should model the whole flow condition very accurately to get satisfactory results. Selection of wrong turbulence model often results worthless outcomes, as wrong model may not represent the actual physics of the flow. Turbulent flow dictates most flows of pragmatic engineering interest. Turbulence acts a key part in the determination of many relevant engineering parameters, for instance frictional drag, heat transfer, flow separation, transition from laminar to turbulent flow, thickness of boundary layers and wakes. Turbulence usually dominates all other flow phenomena and results in increasing energy dissipation, mixing, heat transfer, and drag. In present study flow is fully developed turbulent and Reynolds number (R e ) is set to 6×10 6 . Spallart-Allmaras, k-? realizable, k-? standard and k-? Shear Stress Transport (SST) are primarily used to model viscous turbulent model. However, these specific models are suitable for specific flow cases. Douvi C. Eleni [2] studied variation of lift and drag coefficients for different viscous turbulent model. His study shows that for flow around NACA 0012 airfoil k-? Shear Stress Transport (SST) model is the most accurate. # II. # Theoretical Background Lowest flight velocities are encountered by an airplane at takeoff or landing, two phases that are most perilous for aircraft safety. The stalling speed V stall is defined asthe slowest speed at which an airplane can fly in straight and level flight. Therefore, the calculation of V stall , as well as aerodynamic methods of making V stall as small as possible, is of vital importance. The stalling velocity is readily obtained in terms of the maximum lift coefficient, as follows. From the definition of C L , L= q ? SC L = 1 2 ??? ? 2 SC L Thus, V ? =? 2?? ?? ? ???? ?? (1) [In case of steady level flight, L=W] Examining Eq. ( 1), we find that the only alternative to minimize V ? is by maximizing C L for an airplane of given weight and size at a given altitude. Therefore, stalling speed resembles to the angle of attack that yieldsC L,max : V stall =? 2?? ?? ? ???? ??,??????(2) Figure 2 : When a plain flap is deflected, the increase in lift is due to an effective increase in camber and a virtual increase in angle of attack In order to decrease V stall , C L,max must be increased. However, for a wing with a given airfoil shape, C L,max is fixed by nature, that is, the lift properties of an airfoil, including maximum lift, depend on the physics of the flow over the airfoil. To assist nature, the lifting properties of a given airfoil can be greatly enhanced by the use of "artificial" high-lift devices. The most common of these devices is the flap at the trailing edge of the wing, as sketched in Fig. 2. When the lap is deflected downward through the angle ?, as sketched in Fig. 2b, the lift coefficient is increased for the following reasons: ? The camber of the airfoil section is effectively increased, as sketched in Fig. 2c. The more camber an airfoil shape has at a given angle of attack, the higher the lift coefficient. J trailing edge of the lap, points A and B, respectively, in Fig. 2d. Line AB constitutes a virtual chord line, rotated clockwise relative to the actual chord line of the airfoil, making the airfoil section with the deflected lap see a "virtual" increase in angle of attack. Hence, the lift coefficient is increased. For these reasons, when the flap is deflected downward through the flap deflection angle ?, the value of C L,max is increased and the zero-lift angle of attack is shifted to a more negative value, as shown in Fig. 3. In Fig. 3, the lift curves for a wing with and without laps are compared. Note that when the flaps are deflected, the lift curve shifts to the left, the value of C L,max increases, and the stalling angle of attack at which C L,max is achieved is decreased. However, the lift slope remains unchanged; trailing-edge laps do not change the value of Range (R) is characterized by the maximum distance that an aircraft can travel with a full tank of fuel. Range is technically defined as the total distance (measured with respect to the ground) traversed by the airplane on a tank of fuel.All the way through20thcentury aviation, range has been avital design factor, especially for transcontinental and transoceanic conveyors and for tactical bombers for the army. The range formula for jet airplane which gives a quick, practical estimate for range and which is generally accurate to within 10 to 20 percent is given by R = 2 ? 2 ?? ? ?? 1 ?? ?? ?? ?? 1/2 ?? ?? (?? 0 1/2 ? ?? 1 1/2 ) (3) From Eq. ( 3) that to obtain maximum range for a jet airplane, we want the following: ? Minimum thrust-specific fuel consumption c t . ? Maximum fuel weight W f . ? Flight at maximum ?? ?? 1/2 /C D . ? Flight at high altitudes, that is, low ? ? . Endurance (E) is defined as the entire time that an airplane stays in the air on a tank of fuel. In different applications, it may be desirable to maximize one or the other of these characteristics. The parameters that maximize range are different from those that maximize endurance. The formula for endurance is given by E= 1 ?? ?? ?? ?? ?? ?? ln ?? 0 ?? 1 (4) From Eq. (4) that for maximum endurance for a jet airplane, we want: ? Minimum thrust-specific fuel consumption c t . ? Maximum fuel weight W f . ? Flight at maximum C L /C D . # b) Mathematical Formulation of Turbulance Model Equations for mass and momentum are solver by the solver for all flows. In case of turbulent flow transport equations are also solved additionally. The equation representing the conservation of mass or continuity equation, can be written as follows: ???? ???? + ? ? ??? u ? = S m (5) Eq. 5 is valid for incompressible as well as compressible flows which is the general form of the mass conservation equation. S m is the source of the mass added to the continuous phase from the dispersed second phase (for instance, due to vaporization of liquid droplets) and any user-defined sources. Momentum conservation in an inertial reference frame can be described by Eq. J ?? ???? (?? u )+? ?(?? u u )=-???+? ? ( ? )+?? g + F (6) where p is the static pressure, ? is the stress tensor (expressed below) and ?? g and F are the gravitational body force and external body forces , respectively. F contains additional model-dependent source terms like porous-media and user-defined sources as well. The stress tensor ? is given by: ? =?? ??? u + ? u ? ? ? 2 3 ? ? ? u ?? (7) Where, µ is the molecular viscosity, I is the unit tensor, and the second term on the right hand side is the consequence of volume dilation. The continuity equation for 2-D, steady and incompressible flow is given by: ???? ???? + ???? ???? = 0 (8) Fluent facilitates with various turbulent model having various characteristics suitable for various specific field of study. As stated earlier, no single turbulence model is generally recognized as being superior for all courses of problems. Choice of turbulence model depends on contemplations such as the physics incorporated in the flow, the conventional practice for a definite sort of problem, the level of exactness required, the obtainable computational resources, and the amount of time offered for the simulation. To make the most apposite choice of model for required work, one requires to comprehend the competencies and limitations of the various options. However, Douvi C. Eleni [2] shows in his study that the most accurate among Spalart-Allmaras Model, k-?? realizable Model and k-? Shear Stress Transport (SST) Model, is k-? SST Model for 2D NACA 0012 airfoil simulation process. Therefore, for this study k-? SST Model is employed. # c) The k-? shear stress transport (SST) model The k -? shear-stress transport (SST) model was proposed and developed by Menter [9] to effectively blend the vigorous and precise formulation of the k -? standard model in the near-wall region with the free-stream liberation of the k -? standard model in the far field. This is achieved by the conversion of the k -? model into a k -? formulation. The k -? SST model is comparable to the standard k -? model, but following enhancements are included: ? A blending function was multiplied to both of the standard k -? model and the transformed k -? model and then added together. In the near-wall region the blending function is one activating the standard k-? model. Away from the surface it is zero, which activates the transformed k -? model. ? A damped cross-diffusion derivative term is incorporated in the ? equation of SST model. ? The modified definition of the turbulent viscosity is used to account for the transport of the turbulent shear stress. ? The constants of modeling are made different. These features make the SST k -? model more accurate and reliable for a wider class of flows (e.g., adverse pressure gradient flows, airfoils, transonic shock waves) than the standard k -? model. The SST k? model has a similar form to the standard k -? model: ??(???? ) ???? + ??(?? ?? ?? ??) ?? (?? ?? ) = P -? * ??k + ?? ???? ?? [(µ + ? k µ t ) ???? ???? ?? ] (9)??(???? ) ???? + ??(?? ?? ?? ??) ??(?? ?? ) = ?? ?? ?? P -??? 2 + ?? ?? ?? ?? [(µ + ? ? µ t ) ???? ???? ?? ] +2(1-F 1 ) ?? ?? ?? 2 ?? ???? ???? ?? ???? ???? ??(10) Where, P = ?? ij ???? ?? ???? ?? ?? ij = µt?2?? ???? ? 2 3 ?? ?? ?? ?? ?? ?? ?? ???? ? ? 2 3 ?????? ???? ?? ???? = 1 2 ? ???? ?? ???? ?? + ???? ?? ???? ?? ? Turbulent eddy viscosity can be expressed as: ?? ?? = ?? ?? 1?? max (?? 1 ?? ,??? 2 ) Blending of inner (1) and outer (2) constant for each of the constants are done by: ?? = ?? 1 ?? 1 + (1 ? ?? 1 )?? 2 Other functions are given by: F 1 = tanh(?????? 1 4 ) arg 1 = min??????? ? ??? ?? * ??ð??"ð??" , 500 ?? ð??"ð??" 2 ?? ? , 4????2?????? ???? ð??"ð??" 2 ? CD k? = max?2???? ??2 1 ?? ???? ???? ?? ???? ???? ?? , 10 ?20 ? F 2 = tanh(??????2 2 ) arg 2 = ?????? ?2 ??? ?? * ??ð??"ð??" , 500 ?? ð??"ð??" 2 ?? ? The model constants are: ? 1 = ? 1 ? * ? ? ? 1 ? 2 ?? * ? 2 = ? 2 ? * ? ? ? 2 ? 2 ?? * III. # Computational Method The well documented airfoil, NACA 0012, is utilized in this study. As NACA 0012 airfoil is symmetrical, theoretical lift at zero angle of attack, AoA (?) is zero. In order to validate the present simulation process, the operating conditions are mimicked to match the operating conditions of NASA Langley Research Center validation cases [1]. Reynolds number for the simulations is R e =6x10 6 , the free stream temperature is 300 K, which is the same as the ambient temperature. The density of the air at the given temperature is ? = 1.225kg/m 3 and the viscosity is ?=1.7894×10 -5 kg/ms. Flow for this Reynolds number can be labelled as incompressible. This is a supposition close to reality and there is no necessity to resolve the energy equation. A segregated, implicit solver, ANSYS Fluent 12, is utilized to simulate the problem. The airfoil profile is engendered in the Design Modeler and boundary conditions, meshes are created in the pre-processor ICEM-CFD. Pre-processor is a computer program that can be employed to generate 2D and 3D models, structured or unstructured meshes consisting of quadrilateral, triangular or tetrahedral elements. The resolution and density of the mesh is greater in regions where superior computational accuracy is needed, such as the near wall region of the airfoil. As the first step of accomplishing a CFD simulation the influence of the mesh size on the solution results should be investigated. Mostly, more accurate numerical solution is obtained as more nodes are used, then again using added nodes also escalates the requisite computer memory and computational time. The determination of the proper number of nodes can be done by increasing the number of nodes until the mesh is satisfactorily fine so that further refinement does not change the results. Fig. 5 depicts the variation of coefficient of lift with number of grid cells at stall angle of attack (16°). 120000 quadrilateral cells with C-type grid topology is applied to establish a grid independent solution (Fig 6). From fig. 4 it is evident that 120000 cells are quite sufficient to get a stable and accurate result. Moreover, Douvi C. Eleni [2] was able to generate accurate results using only 80000 cells. The domain height and length is set to approximately 25 chord lengths. This computational model is very small compared to that of NASA's validation cases (fig. 7). Tominimize problems concomitant with the effect of far-field boundary (which can particularly influence drag and lift levels at high lift conditions), the far-field boundary in the grids provided have been located almost 500 chords away from the airfoil. But then again, simulation of NASA's specification of the large computational domain requires very high computer memory. Furthermore, far-field boundary contributes very little on the result. : Mesh of the computational domain around NACA 0012 airfoil (top) and closed detail to the airfoil (bottom) Ansys recommends turbulence intensities ranging from 1% to 5% as inlet boundary conditions. In this study it is assumed that inlet velocity is less turbulent that pressure outlet. Hence, for velocity inlet boundary condition turbulence intensity is considered 1% and for pressure outlet boundary5%. In addition, Ansys also recommends turbulent viscosity ratio of 10 for better approximation of the problem. For accelerating CFD solutions two methods were employed on the solver. The pressurebased coupled solver (PBCS) introduced in 2006, reduces the time to overall convergence, by as much as five times, by solving momentum and pressure -based continuity equations in a coupled manner. In addition, hybrid solution initialization (fig. 8 (a) and 8 (b)), a collection of recipes and boundary interpolation methods to efficiently initialize the solution bas ed purely on simulation setup, is employed -so the user does not need to provide additional inputs for initialization. The method can be applied to flows ranging from subsonic to supersonic. It is the recommended method when using PBCS and DBNS (density -based coupled solver) for steady -state cases in ANSYS Fluent 13.0. This initialization may improve the convergence robustness for many cases [6]. IV. # Validation of the Simulation Process To validate the computational method stated earlier, results obtained by the 2D simulation of NACA 0012 for zero flap angle (?) is compared with NASA's result. The lift curve, drag polar, pressure coefficient (C P ) curve (AoA 0, 10 and 15 degree) for present study is obtained and overlapped on the standard curves provided in NASA's website [1] to observe the fit of current study data. As NASA recommended the definition of the NACA 0012 airfoil is slightly altered so that the airfoil closes at chord = 1 with a sharp trailing edge. To do this, the exact NACA 0012 formula is used, then the airfoil is scaled down by 1.008930411365. The scaled formula can be written: From -16 degree AoA to 16 degree AoA the lift curve is almost linear. Throughout this regime no separation occurs and flow remains attached to the airfoil. At stall AoA lift coefficient is reduced drastically due to intense flow separation generation. Slight deviation from Abbott and Von Doenhoff's [7] unstripped experimental results occurs (almost 3%), as the computational domain of the current study is nearly 1/20 th of the original computational domain under experiment of NASA. Fig. 10 depicts the conformation of drag polar of current validation study with NASA's validation cases. Present study results tie reasonably well with Abbott and Von Doenhoff's [7] unstripped experimental results. At zero angle of attack (AoA) surface pressure coefficients matches with all experimental data particularly well having slender deviations at the trailing edge of the airfoil (fig. 11). However, surface pressure coefficients for flow having AoA 10 degree and 15 degree appear to (fig. 12 and fig. 13) conform to data of experiment conducted by Ladson et al [8]. Leading edge upper surface pressure peak do not appear to resolve well in both cases. Additionally, present study depicts higher pressure than Ladson study [8] on the lower surface on the leading edge of the airfoil primarily due to assuming zero surface roughness of the wall. 15 (see appendix). As NACA 0012 is a symmetric airfoil, for zero AoA it can be observed that velocity profile, pressure profile and streamlines are also same on both upper surface and lower surface of the airfoil. As a consequence, lift generation is also zero for this case (fig. 9). However, with changing AoA the position of stagnation point also changes (fig. 14). J At stagnation point pressure is maximum and velocity is zero which is characterized by distinct red point on the velocity contour plots. It is also apparent that with positive AoA stagnation point moves toward trailing edge on the lower surface of the airfoil. This pressure deviation on the upper and lower surface of the airfoil principally creates significant amount of positive lift. Moreover, separation of flow is also evident at high angle of attack (?) from fig. 15 (see appendix). In turn this flow separation phenomenon creates another source of aerodynamic drag, called pressure drag due to separation. That is why high lift usually associates with high drag. Two major significances of separated low over the airfoil can be noted. The first is the loss of lift. The aerodynamic lift is derived from the net component of a pressure distribution in the vertical direction. When the flow is separated higher pressure is created on the top surface pushing the airfoil downward, thus creating less lift. V. # Results and Discussion NACA 0012 airfoil having different flap angles (?) was subjected to flow of varying Mach number (M). Flow having Mach number greater than 0. Lift coefficient (C L ) escalates with increasing Mach number (M) but a dramatic downslope is obtained at free-stream velocity (V ? ) approaching to sonic velocity. When Mach number (M) is in between 0.8 to 1.2, the flow is said to be transonic which is characterized by some very complex effects. This problem of drastic increasing in drag coefficient (C D ) (fig. 17 and 18) and decreasing in lift coefficient (C L ) can be dealt by using thin airfoil or supercritical airfoil. A rise in critical Mach number (M cr ) usually means an upsurge in the dragdivergence Mach number. Hence, before encountering drag divergence a transonic airplane having a thinner airfoil can fly at a higher Mach number if everything else being equal. The drag coefficient (C D ) remains somewhat constant at low Mach number (M). But, very sudden and dramatic escalation is observed when Mach number (M) approaches to unity (fig. 17). This phenomenon can be also observed in fig. 18 which depicts variation of drag coefficient (C D ) with Mach number (M). However, fig. 18 is attained for zero angle of attack (AoA). As in this study flap angle (?) is varied, the virtual AoA is also changed (fig. 2). As a result, fig 17 gives dissimilar outcomes form fig. 18 to some extent. 19c where presence of shockwave is depicted. Shock waves themselves are dissipative occurrences, which results in an escalation in drag on the airfoil. Moreover, sharply increase of pressure across the shock waves creates an adverse pressure gradient, causing the flow to separate from the surface. This flow separation also contributes to the drag substantially. However, with high flap angles (?) (i.e. 40 and 50) this trend occurs somewhere at Mach 0.5 (fig. 17). This is mainly due to increasing flap angle (?) associates with increasing frontal area of the airfoil. 1![Figure 1 : Typical high lift devices Some high lift devices are illustrated in fig. 1. These devices are primarily used to improve the maximum lift coefficients of wings with changing the](image-2.png "Figure 1 :") 2![When the flap is deflected, we can visualize a line connecting the leading edge of the airfoil and the Year 013](image-3.png "? 2 Global") 3![Figure 3 : Illustration of the effect of flaps on the lift curve. The numbers shown are typical of a modern medium-range jet transport](image-4.png "Figure 3 :") 4![Figure 4 : Typical values of airfoil maximum lift coefficient for various types of high-lift devices](image-5.png "Figure 4 :") 2![](image-6.png "6 2 Global") 5![Figure 5 : Variation of lift coefficient with number of grid cells [2]](image-7.png "Figure 5 :") 67![Figure 6](image-8.png "Figure 6 Figure 7 :") 8![Figure 8 : (a) Standard Initialization, 279 Iterations (b)Hybrid Initialization, 102 Iterations[6] ](image-9.png "Figure 8 :") ![y = 0.594689181[0.298222773??? -0.127125232x -0.357907906x 2 + 0.291984971x 3 -0.105174696x 4 ] developed turbulent flow is simulated in Fluent to match NASA's criteria. Variation of lift coefficient (C L ) with angle of attack (?) for the simulation can be observed from fig. 9.](image-10.png "") 9![Figure 9 : Lift curve of NACA 0012 airfoil](image-11.png "Figure 9 :") 11![Figure 11 : Variation of pressure coefficient (C P ) for 0 degree AoA](image-12.png "Figure 11 :") 13![Figure 13 : Variation of pressure coefficient (C P ) for 15 degree AoA Pressure and velocity contours along with streamlines for different AoA (?) are presented in a tabular form in fig. 15 (see appendix).As NACA 0012 is a symmetric airfoil, for zero AoA it can be observed that velocity profile, pressure profile and streamlines are also same on both upper surface and lower surface of the airfoil. As a consequence, lift generation is also zero for this case (fig.9). However, with changing AoA the position of stagnation point also changes (fig.14).](image-13.png "Figure 13 :") 14![Figure 14 : Variation of stagnation point position with AoA](image-14.png "Figure 14 :") ![3 is considered compressible. Density based solver utilizing k-? Shear Stress Transport (SST) modelin Fluent facilitates mimicking compressible flow over the body under experimentation very accurately. The resultant forces are typically resolved into two forces and moments. The component of the net force acting normal of the airfoil is lift force (F L ) and acting horizontal to the airfoil is drag force (F D ). The curves showing variation of lift coefficient (C L ) and drag coefficient (C D ) with different flap angles (? ) and Mach number (M) are analysed to realize the aerodynamic behaviour of plain flapped NACA 0012 airfoil. Curves of C L /C D and ?C L /C D are analysed further, as these are crucial factors affecting range (R) and endurance (E) of aircrafts. However, flow at high flap angles (?) (i.e. 30, 40 and 50 degrees) are very unstable and it remains unconverged even after 5000 iteration in Ansys Fluent flow solver. Hence, flow for flap angle (?) 30, 40 and 50 degrees are slightly erratic. At fig. 16 variation of lift coefficient (C L ) with Mach number (M) can be observed. Higher lift coefficient (C L ) is obtained for higher flap angles (?) at any Mach number (M). However, atypical behavior of the curves are evident at increasing Mach number (M).](image-15.png "") 1716![Figure 17 : Variation of drag coefficient (C D ) with Mach number (M) for different flap angle (?) The airfoil subjected to the flow passes through three distinct phases with Mach number (M) represented by point a, b and c in the fig. 18. At point a, free stream Mach no is characterised by M ?