Multi element airfoils are high lift devices and provide improved aerodynamic characteristics which are beneficial for several applications such as aircraft wings, wind turbine blades. Simulations were performed using the computational panel code developed in MATLAB for the airfoils, 30P30N, GA (W)-1, RAF16, NLR 7301 configurations. The results were obtained for varying angle of attacks from -10 to 20deg to include the effects of turbulence. The lift and drag coefficients are evaluated for the airfoil configurations, 30P30N, for Reynolds number range of 1.6x10⁵-1x10⁶. The velocity and pressure contours distributions are illustrated to predict the laminar to turbulent flow characteristics on the airfoil surface. The maximum lift coefficient for single slotted flap element is found to decrease for large angle of attack beyond the stall region and for slat element the increase in lift is found to be marginal at high angle of attack. The maximum lift coefficient obtained including flap and slat elements is 3.62 and the drag coefficient is 0.368.

Keywords: airfoil, lift coefficient, drag coefficient, pressure, chord, slat, flap

High lift devices affect the aircraft performance during the landing, takeoff situations. The addition of such devices along the wing span ensures the flow characteristics are modified to generate the desired lift required for landing, takeoff conditions. Presence of flap element on the aircraft wing will help achieve rolling and pitching moment stability by controlling the desired lift produced during the cruising condition. Dong et al.¹ performed the numerical analysis for determining the flow behavior over the airfoil configurations which included the leading edge slat and trailing edge flap devices. The influence of flap gaps and the trailing edge flap angles were investigated on such configurations to quantify the change in the lift and drag coefficients. Wings of cargo aircrafts usually consist of leading edge slat and trailing edge flap configurations to provide the high aerodynamic efficiency although at the cost of structurally complex design.¹ Sorensen et al.² conducted wind tunnel experiments on thick multi element airfoil configurations at subsonic flows over airfoil elements and observed the effects of slat, circulation, boundary layer, in order to understand how such configurations provide high lift coefficients compared with conventional airfoils.² In the next sections, the numerical panel methods are discussed which require relatively less computational effort to predict the physics of flow around the airfoils and aim to study the influence of high lift devices for static and dynamic stall conditions experienced by the aircraft wings. The key parameters involved in the multi element airfoils design are deflection or orientation angle, overlap and gap distance and its variation will enable to observe its influence on overall performance of the airfoil configurations required for intended application.^3,4

Computational panel method

Traditional methods for modeling flow around slender bodies of any shape include potential flow which utilizes the superposition of source and sink on x axis and in uniform distributed flow. However, the theory does not predict accurate values for flow whose leading edge has rounded shapes. Basic panel methods were developed by Hess and Smith at Douglas aircraft in late 1950s⁵ for aircraft industry. Panel methods model the potential flow by distributing sources over the body surface. A source is point at which the fluid appears in the field at uniform rate while a sink is point which disappears at uniform rate, m3/s. Each source or sink has specific strength and vortex strength denoted by circulation, Simple 2D uniform lifting flows⁶ can be described using the following equations

$\text{φ}=\text{by}-\text{ax}\dots \text{stream line function}$

$\varphi =\text{bx}+\text{ay}$ … Velocity potential function

The resultant velocity V, at any point along the flow direction can be written as            

$\sqrt{\left[a^2+b^2\right]}$

The point source or sink is distributed uniformly in all directions of the flow field and obey the continuity equation and irrigational motion everywhere except at the point itself. The total velocity potential² for numerical panel method can be rewritten as follows

$\varphi =\text{Ux}+\text{μ}$

Where, µ is the perturbation potential away from free stream conditions. The airfoil geometry is discretized into finite number of panels over the surface. The panels are represented by the 2D shape of the surface by series of straight line segment.^2,4 The following procedure describes the calculation for 2D lifting flows.

1. Numbering of end points or nodes of the panels from 1…N
2. The center points of each panel are chosen as collocation points. The boundary condition of zero flow orthogonal to surface is applied to the points.
3. Panels are defined with unit normal and tangential vectors, $\widehat{n}$ , $\widehat{t}$ .
4. Velocity vector, denoted by V_ij are estimated by considering the two panels, i & j the source on the panel j which induce a velocity on panel i. The perpendicular and tangential velocity components to the surface at the point I, are given by scalar products of ${\text{v}}_{\text{ij}}\text{.}\widehat{n}$  and ${\text{v}}_{\text{ij}}\text{.}\widehat{t}$
5. The above quantities represent the source strength on panel j and expressed mathematically as

${\text{v}}_{\text{ij}}\text{.}\widehat{n}={\text{σ}}_{\text{j}}{\text{N}}_{\text{ij}}$
${\text{v}}_{\text{ij}}\text{.}\widehat{t}={\text{σ}}_{\text{j}}{\text{T}}_{\text{ij}}$

Where ${\text{N}}_{\text{ij}}$  and ${\text{T}}_{\text{ij}}$ are the perpendicular and tangential velocities induced at the collocation panel i and known as normal and tangential influence coefficients. The surfaces represented by the panels are solid and the following conditions are applied for the normal and tangential velocities at each of collocation points consisting of sources strengths, vortices, and oncoming velocity, U.

${\displaystyle \sum _{j=1}^N{\sigma }_j{N}_{ij}+\gamma N_{i,N+1}+\overrightarrow{U}\cdot {\widehat{n}}_i=v_{n_i}}$  (1)
${\displaystyle \sum _{j=1}^N{\sigma }_j{T}_{t,j}+\gamma T_{t,N+1}+\overrightarrow{U}\cdot {\widehat{t}}_i=}{v}_{s_i}$    (2)
${\displaystyle \sum _{j=1}^N{\sigma }_j{N}_{ij}+\gamma N_{i,N+1}+\overrightarrow{U}\cdot {\widehat{n}}_i=0}$   (3)

The above system of linear algebraic equations are solved for the N unknown source strengths, $\sigma$ i, using matrix system and expressed as

   M. a = b …  (4)

Where N is an N+1 x N+1 matrix containing the Nij and ${\text{σ}}_i$  is column matrix of N elements and A is the column matrix of N elements of unit normal velocity vectors. Matrix inversion procedures available in MATLAB are applied to solve for the source strengths using the above system of equations. The pressure acting at collocation point i is given by the Bernoulli equation as⁵

${\text{C}}_{\text{pi}}=1-{\left[\frac{{\text{v}}_{\text{Ti}}}{\text{U}}\right]}^2$ (5)

Where ${\text{v}}_{\text{Ti}}$ the tangential velocity vector is determined using the influence coefficients. The influence coefficients are important for panel method in order to determine the pressure distribution over the surface of the any given airfoil coordinates. The airfoil trailing edge presents a unique condition for the flow field parameters. Using panel method, the following criteria is applied for the stream lines around the airfoil

1. The streamlines leave the trailing edge with a direction along the bisector of the trailing edge angle.
2. The velocity magnitudes on the upper and lower surfaces near the trailing edge of airfoil approach the same limiting values.
3. The trailing edge angle is modeled as the stagnation point for finite value of trailing edge angle hence the source strength must be zero at the trailing edge.

The above assumptions are known as the Kutta condition which is essential for the successful evaluation of velocity vectors and pressure for 2D uniform flows. It can be written in the algebraic form of equation as

${\displaystyle \sum _{j=1}^N{\sigma }_j{T}_{t,j}+\gamma T_{t,N+1}+\overrightarrow{U}\cdot {\widehat{t}}_i=}$ $-\left({\displaystyle \sum _{j=1}^N{\sigma }_j{T}_{t,j}+\gamma T_{t,N+1}+\overrightarrow{U}\cdot {\widehat{t}}_i}\right)$  (6)

It must be noted that the velocity components induced at any point P due to sources on panel centered at point Q can be expressed mathematically as scalar product of vector form

${\text{v}}_{\text{PQ}}={\text{v}}_{\text{xQ}}\text{.}\widehat{t_j}+{\text{v}}_{\text{yQ}}\text{.}\widehat{n_j}$       (7)
${\text{N}}_{\text{ij}}={\text{v}}_{\text{PQ}}\text{.}\widehat{n_i}$ (8)
${\text{T}}_{\text{ij}}={\text{v}}_{\text{PQ}}\text{.}\widehat{t_i}$ (9)

Therefore, the resulting velocity components along with known source strengths and influence coefficients are added for each panel in order to obtain pressure distribution over the airfoil surface. The number of panels used (order) in the simulation ranged from 150-250 for airfoils in the MATLAB routine foil m Individual airfoil elements are discretized into several panels in usual manner, for evaluating the pressure characteristics over the integrated airfoil configuration. In Figure 1 the graphical illustration of Kutta condition over the airfoil trailing edge is shown with the limiting velocity vectors. The unit normal vector and tangential velocity are approximated with equal and opposite in directions on the trailing edge panel.

Figure 1 Graphical illustration of kutta condition.

Airfoil geometry

Slats are airfoil elements which are located in front of main airfoil and deflect the flow in the chord wise direction. They enable the flow separation to occur at higher AOA which as result provide improved lift characteristic over the wing span.² Greenwell studied the influence of trailing edge flaps which increase the profile drag as result of the flap gap and flap deflection angle however, this is countered with the reduction of induced drag. The coordinates were obtained from the UIUC airfoil database. The flap gap at the trailing edge is higher for GA (W)-1, NLR 7301 airfoils compared to 30P30N airfoil elements.^3,7 There are no leading edge slat elements for the 30P30N, NLR7301, RAE 16 elements. Hence, the total lift produced by them is relatively lesser for same set of conditions, free stream Reynolds number, and angle of attack. Super critical airfoils possess flat suction surfaces which are intended to minimize the wave drag and have high curvature in the aft sections and large leading edge radius compared to NACA 6 series airfoils.^8,9 Although 3D airfoil design remains complex in nature,⁴ the results cannot be directly compared with experimental wind tunnel data. Optimization of wing geometry involves the displacement of slat and flap deflections at the trailing edge to produce high or desired lift. Therefore, they rely on the 2D airfoil design process which is simple to implement in numerical and experimental conditions. The flap orientation angle is important to understand the flow separation caused due to the boundary layer thickness augmentation aft of the airfoil. The flap element and main element overlap for different deflection angles whose position in chord wise direction causes the flow to be energized and followed by subsequent pressure rise aft of the main airfoil element. Figure 2 shows the geometry of airfoils analyzed in this article. Figure 3 depict show flap gap, overhang distance and deflection angle is measured with respect to leading edge radius and chord line of main element of airfoil. The super critical airfoils create increased lift in the aft section, and hence more efficient. It effectively reduces the shock wave drag which is observed in the transonic and supersonic flow conditions. The shock wave occurs on the upper surface of the airfoil and propagates downstream along the chord. When it reaches the aft, the pressure bumps can be noticed in the flap section of the airfoil.¹⁰ The gap between the trailing edge and main element is measured in two directions, normal direction to surface of main element and also in vertical direction. The deflection angle is intended for boundary layer flow over the main element to reenergize and create delayed stall separation in the aft of airfoil. It must be noted that advantage of improved lift on airfoil configuration is created due to the combined effect of high lift on slat section, pressure recovery on the flap sections. The chord length of slat is given in terms % of main chord, maximum thickness, %, of main element airfoil, chord length of flap element expressed in terms of % of main chord which are essential to construct the airfoil geometry.

Figure 2 Geometry of multi-element airfoils.

Figure 3 Schematic of flap gap, overlap distance and deflection angle.¹

Pressure distribution

Computations were performed for multi element airfoil configuration as shown in previous section for angles of attacks 4, 9, 10, 11, 13 and 14 degrees respectively. Since the computations are based upon the numerical panel method, no grid validation studies were involved in the analyses which affect the pressure variable at different angle of attack and free stream Reynolds number. The pressure distribution and contours of tangential velocity for the NLR 7301 supercritical airfoil, GA (W) -1, 30P30N, RAE 16 are discussed.. It utilized Intel 4GHz processor, 2GB RAM, PC at GITAM University high performance computing lab. The pressure contours near the gap region at the trailing edge flap show dense regions with high pressure values distribution and with discontinuous bumps in flow seen near the leading edge of the flap in deflected position. The disturbances are produced near the leading edge slat and translate along the chord wise direction of the airfoil, resulting in nonlinear pressure distribution. The presence of leading edge slats are also intended to reduce the vibration and improve the performance. Table 1 shows the thickness %, chord lengths of main element and the flap gap between the trailing edge of main element and flap. In (Figure 4) (Figure 5), the pressure and velocity contours are shown with the pressure probe on the surface of airfoil are represented on the Y-axis, and angle of attack on the X-axis. The pressure probes are numbered serially starting from upper and lower surfaces and leading edge of slat element to trailing edge of flap element in continuous manner. From Figure 4 it must be noted that the pressure peaks on the suction side of NLR 7301 airfoil are highly over predicted for the flow conditions. The pressure coefficient reaches a maximum of 5 at 40AoA, which represents the transition to turbulent region on the airfoil surface. Further, beyond the 20 % c of main element, the pressure gradient is lowered towards the trailing edge. At 140AoA, the pressure coefficient is over predicted formed on the suction side at leading edge of airfoil. It may further entail in the formation of bubbles due to adverse pressure gradient created on the airfoil surface at higher sonic speeds (Figure 6). Pressure distributions of LS, MS, Natural LaminarFlow (NLF) airfoils at 1 deg Angle of attack. The pressure distribution of GAW-1, GAW-2 airfoils is shown in Figure 6 for 6deg angle of attack. The suction side pressure peak is observed higher in case of GAW-1 due to high thickness and for same set of flow stream and boundary conditions. The pressure coefficient, Cp is ~6 for GAW-2 and ~7 for GAW-1. The trailing edge pressure peak on the pressure side is constrained due to the overlap between the main element and slotted flap with a maximum value of 1.9. It must be noted that the main element is composed of several slots due to which pressure bumps can be observed 20% and 60% chord of airfoil. From Figure 6 a comparison has been made for the series of low speed airfoil with two different thicknesses, ~14 % and 17 % commonly used in cargo or military transport aircrafts where the payload capacity is important criterion than top speed of the aircraft. The flow reversal regime is seen to occur for low to medium speed, MS 313, 317at 1 deg AOA where the suction and pressure surfaces are interchanged. For the 30P30N configuration, a single slotted flap element is chosen at the trailing edge and a retractable leading edge slat. At higher air speeds the slat element ensures the oncoming flow to remain in laminar state which will cause the formation of the boundary layer transitions to occur at large angle of attack. The gap between the slat and main element ensures the pressure recovery at the leading edge of main element where the predominant lift is produced. The flow stream conditions at the trailing edge are altered by the position of the flap and its deflection angle. Further, the gap measured in normal direction between the flap and trailing edge of main element also affects the overall lift produced from the configuration. However, in the present study, no attempt was made to observe the effect of flap deflection and overlap distances between the slotted elements of the airfoil. From Figure 7 it can be noted that the suction side pressure coefficient for flap and slat elements is excessively over predicted at 9 and 11, 13deg due to the turbulent flow nature of aft of airfoil and presence of extended flap, slat gap between the main and flap element. The maximum pressure coefficient obtained is~8 on the flap leading edge. The flow stream past the surface of airfoil increases due to the presence of cove region near the trailing edge on the main element of airfoil. Therefore, the pressure recovery can be observed at high AOA however, the pressure continues to drop at the leading edge of flap section as result of the flap gap in the deflected position (Figure 8).

Figure 4 Pressure distribution of NLR 7301 supercritical airfoil at 14 deg and 4 deg AOA.

Figure 5 Pressure distribution of GAW-1, GAW-2 airfoils at 6deg angle of attack.

Figure 6 Pressure distribution of LS, MS, Natural LaminarFlow (NLF) airfoils at 1 deg Angle of attack.

Figure 7 Pressure distribution of 30P30N airfoil at 9, 11 & 13deg Angle of attack.

Figure 8 Pressure distribution of RAF 15 airfoil at 4 deg and 10deg angle of attack.

The Figure 9 shows the tangential velocity around the RAF 16 airfoil surface for Reynolds number ~ 6.67 x 105. For larger flap deflections, the flow characteristics on the RAF continue with higher values for relative or tangential velocity increments near the trailing edge. Velocity contours indicate that magnitudes are higher towards the leading of flap and main section and show increasing values of relative velocity on the suction side, while continuous decrements on the pressure side of airfoil. Figure 10 shows the pressure contour where the maximum pressure coefficient is observed at the probe number 100-160 and also near 240. The term vortices here refer to the local rotation of fluid and obtained by taking the curl of velocity field vector. Circulation on the other hand can be viewed as the local spinning of the fluid about a fixed center. The vortices contour of the GAW-1 (LS417) airfoil is shown in Figure 11. The probe number is shown on the x-scale beginning from the upper surface of airfoil along the profile boundary and used to measure the pressure on the surface of profile in experimental analysis. For incompressible fluids it can be seen that size of vortices becomes large as the angle of attack (AOA) is increasing. Further, the vortices appear near the trailing edge on upper and lower surfaces with varying size due to augmentation of pressure aft of the airfoil. Figure 12 shows the streamline contour of the GAW-1 airfoil. The relative size of the vortex structure is large near the leading edge of airfoil at higher angle of attack. It can also be seen that the frequency of vortices shed from the trailing edge is higher compared to leading edge for given chord wise position of airfoil. Mathematically vortices can be written as

$\text{ω}=\frac{\partial \text{v}}{\partial \text{x}}-\frac{\partial \text{u}}{\partial \text{y}}=\nabla \times \overline{V}$ (10)
$\text{φ}={{\displaystyle \int }}^{\text{}}\left[\text{udy}-\text{vdx}\right]=\nabla \text{.}\overline{V}$  (11)

Lift and drag characteristics

Models designed based upon the wind tunnel testing and data obtained from measurements are useful to validate the numerical results. Empirical corrections are needed on such wind tunnel data in order to compare airfoil response for efficient characterization at different angle of attack, Mach number and airfoils with variable camber. This data reduction enables to understand the dynamic stall behavior conditions and factors that can delay the stall and improve the lift on the wing span. The relative position of slat and slat gap present between the trailing edge of slat and main element will influence the lift characteristic of airfoil. From Figure 13 the maximum lift coefficient for the airfoil with flap included reached value of 3.67 while the drag for slat increases for large angles of attack, the flap drag reduces due to reattachment of boundary layer occurs with laminar to turbulent flow separation at the trailing edge of main airfoil.

Figure 9 Tangential velocity and Velocity contour of RAF 15 airfoil at Re – 6.6x10⁵.

Figure 10 Pressure contours of RAF 15 airfoil slotted airfoil for -10 to 20 deg AOA.

Figure 11 Vorticity contour of GAW -1 [LS 417] Airfoil.

Figure 12 Stream function contour of GAW -1 Airfoil.

Figure 13 Lift and drag characteristics of 30P30N slotted airfoil for 0 to 20 deg AOA.

Slat and flap elements produce high lift by increasing the effective area of the wing. The pressure distribution for NLR7301, 30P30N, GAW-1, GAW-2, RAF16 airfoils, and lift, drag characteristics are evaluated using the numerical panel method. The flap deflection angle, and the overlapping gap between the main and flap or slat element is important for increasing or decreasing the overall lift of the wing. The flow conditions exhibit the delayed stall characteristics due to the presence of slat element at the leading edge enabling the high lift. The continuity of flow over the main element in near laminar conditions is responsible for producing pressure peaks on the suction side of airfoil and pressure recovery obtained at the leading edge of flap element. At higher angle of attack, the lift coefficient on the flap decreases while it increases on slat. The drag coefficient of flap remains steady until 10 deg and reduces to zero at higher angle of attack. The maximum lift coefficient obtained is 3.67 including the flap and slat element and drag coefficient is 0.36. The lift coefficients for flap and slat elements are under predicted due to the turbulent flow nature aft of the slat element in its deflected position

1. UIUC – University of Illinois Urbana Champaign
2. AOA – Angle of Attack
3. C_L – Lift coefficient
4. C_D – Drag coefficient
5. C_P – Pressure coefficient
6. NLF – Natural Laminar Flow
7. GA – General Aviation,
8. GA (W) – General Aviation (Whitcomb)

None.

Author declares that there is no conflict of interest.

1. Zhang Z, Dong Li. Numerical investigation of flow over multi element airfoils with lift enhancing Tabs. 28th International congress of the Aeronautical Sciences. China; 2012. p. 1‒10.
2. Zahle F, Guanaa M, Sorensen NN, et al. Design and Wind Tunnel Testing of Thick Multi element high lift airfoil. Wind Energy Department, Denmark; 2012.      
3. Narsipur S, Pomeroy BW, Selig M. CFD Analysis of Multi element Airfoils, for Wind Turbines. 30th AIAA Applied Aerodynamics Conference. USA; 2012.
4. Ahaus L, Liggett N, Peters DA, et al. Unsteady aerodynamics of single and multi-element airfoils. ERF. 2010. p. 1‒12.
5. Abbot H, Doenhoff V. A Theory of Wing Sections, Including summary of airfoil data. USA: Dover publications Inc; 2000. p. 1‒704.
6. Haughton, Carpenter. Aerodynamics for Engineering students. 6th ed. USA: Elsevier publications; 2001. p. 1‒740.
7. Harris T, Lowry J. Pressure distribution over NACA 23012 airfoil with fixed slot and slotted flap. National Advisory Committee for Aeronautics. USA: Langley Aeronautical Lab; 1942. p. 85‒98.
8. Weber S, Platzer MF. Computational simulation of dynamic stall on the NLR 7301 airfoil. USA: Naval postgraduate school; 2000. p. 1‒21.
9. B Van den Berg, Gooden JHM. Gooden Low speed surface pressure and boundary layer measurement data for the NLR 7301 airfoil section with trailing edge flap. J Global. 1994;303(2):A9.1‒A9.12.
10. Wenz WH, Fisckom KA. Pressure distributions for the GA(W)-2 airfoil with 20% aileron, 25% slotted flap and 30% Fowler flap. NASA Contractor report 2948. Wichita State Univ, USA; 1978. p. 1‒86.