Week 7 - Mini project - Flow over an Airfoil
AIM : To setup a flow over NACA-2412 airfoil for different angle of attack and calculate the values of lift and drag coefficient (Cl and Cd), hence the compare the result of turbulence.
SIMULATION-SETUP :
- MODEL DESIGN : CATIA V5
- Case setup : CONVERGE CFD STUDIO
- Commanad editor : Cygwin
- Post-processor : Paraview
THEORY (AIRFOIL AND WING) :
As most of them might wonder, how an aircraft fly???. Well there are many important concepts that are responsible for an aircraft to fly. One of those is the pressure difference created on an aircraft wing which is responsible for generating an upward force called as “ LIIFT ”.
Well, it should be understood that all the forces should be balanced in order for the body to be in equilibrium. Lift generated is balanced by the ‘ Weight ‘ of the aircraft. In order for the aircraft to move forward , the engine generates a force called ‘ thrust ’ which is countered by the aircraft’s ‘ drag ’.
Now let us try to understand, what is an airfoil.
An airfoil is just a 2D of a wing or a cross-section of a wing. Let us see the technical aspects and terms related to airfoil.
The various terms related to airfoils are defined below.
- The suction surface( upper surface) is generally associated with higher velocity and lower static pressure.
- The pressure surface(lower surface) has a comparatively higher static pressure than the suction surface. The pressure gradient between these two surfaces contributes to the lift force generated for a given airfoil.
The geometry of the airfoil is described with a variety of terms.
- The leading edgeis the point at the front of the airfoil that has maximum curvature (minimum radius).
- The trailing edgeis defined similarly as the point of maximum curvature at the rear of the airfoil.
- The chord lineis the straight line connecting leading and trailing edges. The chord length, or simply chord, {\displaystyle c}is the length of the chord line. That is the reference dimension of the airfoil section.
The shape of the airfoil is defined using the following geometrical parameters:
- The mean camber lineor mean line is the locus of points midway between the upper and lower surfaces. Its shape depends on the thickness distribution along the chord;
- The thicknessof an airfoil varies along the chord. It may be measured in either of two ways:
- Thickness measured perpendicular to the camber line.
- Thickness measured perpendicular to the chord line.
Now comes the technical aspects such as ;
- Centre of pressure (COP) : which is the chord-wise location about which the pitching moment is zero. In other words the resultant of all the forces vanishes at this point
- Aerodynamic centre : The point in a chord-wise location where the aerodynamic moment vanishes.
- Angle of Attack (AOA) : The angle b’n the free stream or relative stream velocity and the chord line.
EFFECT OF DIFFERENT AOA ON COP.
The AOA changes as the aircraft’s attitude changes. Since the AOA has a great deal to do with determining lift, it is given primary consideration when designing airfoils. In a properly designed airfoil, the lift increases as the AOA is increased.
When the AOA is increased gradually toward a positive AOA, the lift component increases rapidly up to a certain point and then suddenly begins to drop off. During this action the drag component increases slowly at first, then rapidly as lift begins to drop off.
When the AOA increases to the angle of maximum lift, the burble point is reached. This is known as the critical angle. When the critical angle is reached, the air ceases to flow smoothly over the top surface of the airfoil and begins to burble or eddy. This means that air breaks away from the upper camber line of the wing. What was formerly the area of decreased pressure is now filled by this burbling air. When this occurs, the amount of lift drops and drag becomes excessive. The force of gravity exerts itself, and the nose of the aircraft drops. This is a stall. Thus, the burble point is the stalling angle.
As previously seen, the distribution of the pressure forces over the airfoil varies with the AOA. The application of the resultant force, or COP, varies correspondingly. As this angle increases, the COP moves forward; as the angle decreases, the COP moves back. The unstable travel of the COP is characteristic of almost all airfoils.
HOW DOES AIRFOIL GENERATES LIFT ??
When the flow encounters the leading edge of an airfoil as shown in the above figure, the flow velocity initially decreases and as it climbs the upper curvature (upper or suction surface) the velocity starts increasing due to the negative slope or geometry of the airfoil. With increase in the velocity obviously the pressure decreases on the upper surface, creating a negative pressure. On the other hand ie,,, on the lower side of an airfoil the velocity decreases , because the air molecules consume more of their internal energy to travel from the lower side, when compare to the upper side of an airfoil, again this is due to the geometry of the airfoil. Hence pressure on the lower side increases, thus creating the positive pressure on the lower side.
Thus it is this pressure difference that is responsible for generating an upward force called “ LIFT ”.
MODEL-DESIGN :
The modeling of the 2D airfoil into 3D wing was done with the knowledge of CAD package ‘ CATIA ’. The airfoil NACA-2412 coordinates were imported from the excel sheet to the catia softwere.
Later using the tools available in the cad package the airfoil (2D) was extruded to a 3D wing. The chord length was captured to be 0.099871 m and the length of the wing is 1 m.
CASE-SETUP :
The designed model in Catia V5 was then imported Converge studio as a .stl file, where the entire surface was divided into triangles.
OBSERVATIONS :
Chord-length : 0.0999981 m.
Length of the wing : 1m.
The imported model to converge studio is as shown.
The diagnosis was done and was checked for any surface errors like non-manifold, intersection and normal-orientation errors. Surface clean-up was done and all the errors was resolved.
CREATING VIRTUAL WIND-TUNNEL :
Since for a 2D simulation we require an inlet, outlet , front and back surface it is necessary for us to design a box that is capable of holding the above created wing in a fixed position, which is nothing but a ‘ virtual wind-tunnel ’. The dimensions of the virtual wind-tunnel for aerodynamic simulation is as shown below.
For an aerodynamic airfoil simulation usually the length of the virtual wind-tunnel should be 50 times the chord length and the height should be 30 times the chord length. ie,,,,, length = 50*C and height = 30*C
Above figure shows the schematic representation of a Virtual-Wind Tunnel.
From CFD point of view for a good wind-tunnel the inflow length should be smaller (L1) and the outflow length should be larger (L2) and also the height should of the wind-tunnel should be smaller.
VIRTUAL WIND-TUNNEL DESIGNED IN CONVERGE
Keeping the above concept in mind, the wind-tunnel was designed in converge itself. Once the design was complete surface clean-up was done and all the surface errors were resolved.
The length of the wind tunnel is 1 m long (in a X-Y plane) and height is about 0.35 m. As we are running the 2D simulation the dimension in Z-plane is irrelevant for us, as it involves only the front 2D surface and back 2d surface.
BOUNDARY FLAGGING :
In order to assign the boundary conditions the model which is divided into triangles must be grouped into boundaries. The given model is divided into 6 different boundaries which is shown below.
BOUNDARY CONDITIONS AND REGION INTIALIZATION
The inlet is an in-flow boundary with an caluclated inlet velocity of 30 m/s and pressure being a Nueman boundary condition (NE).
The outlet is an outflow boundary with the velocity being Nueman boundary condition and pressure is a Dirchlet boundary condition with a atmospheric pressure of 101.325 KPa.
Front and back surfaces are just a 2D walls, where the simulation can be observed.
The top and bottom surfaces are symmetric walls. The meaning is that, whatever are the flow properties above the chord line (for reference), the property will remain same below the chord line too.
Airfoil is just a wall-type boundary condition assigned such that it is stationary and surface movement to be fixed. Temperature condition being ‘ law of wall ’ @ 300K.
REGIONS AND INTIALIZATION.
In the manner how we group the boundaries from triangles, we have to group these boundaries into volumetric regions. Hence all the above boundaries were grouped into a volumetric region with velocity in x-direction to be 30m/s.
SIMULATION TIME PARAMETERS :
The entire length of the wind tunnel is 1 m and we know the velocity (calculated) is 30m/s. Therefore if the flow starts from inlet and travels till the outlet the time taken by the flow to reach exit point once is,
t= dx/v
t=1/30 = 0.033s.
Say that we want to simulate the flow for 20-30 times, since it’s an aerodynamic simulation the flow can change any instant of time and hence the starting time for the simulation is 0 s and end time is 0.033*30 = 1s.
Now let us take an initial and minimum time step of 1e-5 and maximum time step of 1s to reach the simulation end, although the time step of the simulation does not reaches 1sec.
SOLVER :
We are running the simulation for transient based full-hydrodynamic mode, compressible flow simulation where the flow density varies.
Density-based ‘ Navier-Stokes ’ solver is used.
TURBULENCE-MODELLING :
The turbulence model taken here is RNG K-Ɛ.
GRID CONTROL AND MESHING :
For any aerodynamic simulation the mesh size should be very accurate and finer mesh should be refined. Hence around 0.1 million cells were created in the entire model with a base grid of dx = dy = dz = 0.0036m. which gave an estimation of 1,02,113 or 0.102 million cells.
Since it’s an aerodynamic simulation we need to capture the data very accurately and precisely. Hence fixed embedding was added around the airfoil, with scale being 3 and number of embedding layer being 2. Hence new size of the cells created can be calculated as…
Grid size = (base grid)/2^(embedding level)
= 0.0036/2^2
= 9e-4 m. (size of the new cells or elements in the fixed embedding).
OUTPUT FILES :
Around 100 3d output files were created and was dumped by the command editor ‘ cygwin ’ to the output folder of the specified location which then can be converted into a format that the post-processor ‘ paraview ’ can understand.
The text output was written at every time step given ie,,,, at every 1e-5 sec the text output for plotting the graphs was being printed and the graphs were observed.
RESULTS AND CALUCLATIONS :
VELOCITY :
The given Reynolds number is 2e5. Therefore wkt the Reynolds number Re is given by,
R e = ρ V L /µ
For air at 300K and 1bar pressure, ρ=1.225 kg/m^3, µ=1.85e-5
Therefore, V = (Re*µ)/(ρ*L)
V = 30.08 m/s
Lift-coefficient (Cl) and Drag-coefficient (Cd) :
For any airfoil the characteristic parameter Cl and Cd plays an important role in accounting for effect of lift and drag generation.
Cl and Cd is defined as lift and drag co-efficients respectively, which is dimensionless number, and is given as
Where L is the lift force in Newton (N)
D is the drag force in Newton (N)
is the density of the fluid or continuum
V is the free-stream or relative stream velocity
S is the surface area
The value of ‘L’ and ‘D’ is obtained from the simulation for different angle of attack.
We know density of air at 300K, which is 1.225 kg/m^3
We also have velocity, V = 30m/s.
The only unknown here is area of the wing, which can be calculated, we know the chord length C = 0.099981m and we have the length of the wing as 1m, therefore the area,
A=0.099981*1=0.099981 m^2=0.1m^2
As calculate only the projected surface area of the wing. Hence now the value of Cl and Cd can be calculated
For AOA of 1º :
Above graph shows the pressure force in Y-direction ie,,, nothing but the lift generated for a free-stream velocity of 30m/s.
Above graph shows the pressure force in X-direction ie,,, nothing but the drag generated for a free-stream velocity of 30m/s.
Hence from the above simulated results the value of Lift (L) and Drag (D) for AOA of 1º is,
L = 20.0676 N.
D = 1.00305 N.
Hence,
For AOA of 5º :
Above graph shows the pressure force in Y-direction ie,,, nothing but the lift generated for a free-stream velocity of 30m/s.
Above graph shows the pressure force in X-direction ie,,, nothing but the drag generated for a free-stream velocity of 30m/s.
Hence from the above simulated results the value of Lift (L) and Drag (D) for AOA of 5º is,
L = 43.0298 N.
D = 1.0751 N.
For AOA of 10º :
Above graph shows the pressure force in Y-direction ie,,, nothing but the lift generated for a free-stream velocity of 30m/s.
Above graph shows the pressure force in X-direction ie,,, nothing but the drag generated for a free-stream velocity of 30m/s.
Hence from the above simulated results the value of Lift (L) and Drag (D) for AOA of 10º is,
L = 69.2883 N.
D = 1.8172N.
For AOA of 15º :
Above graph shows the pressure force in Y-direction ie,,, nothing but the lift generated for a free-stream velocity of 30m/s.
Above graph shows the pressure force in X-direction ie,,, nothing but the drag generated for a free-stream velocity of 30m/s.
NOTE : The above plot results for 15º is almost converged, since we ran the simulation for 1sec which almost took 3hrs. If the simulation had ran for end time of 2s the solution would have been converged
Hence from the above simulated results the value of Lift (L) and Drag (D) for AOA of 15º is,
L = 51.8851 N.
D = 3.5092 N.
Now let us try to plot these as graphs ie,,, Cl for different AOA and Cd for different AOA.
Plot of Cl v/s AOA :
Above figure shows the plot of lift coefficient (Cl) for different angle of attack (AOA).
Initially when the angle of attack is 1º , the flow is very smooth and steady over the airfoil body and also the flow is attached to the body . When the flow encounters the leading edge of the airfoil it’s velocity decreases initially and the flow starts accelerating on the suction side (upper) and deceleration in the pressure side (lower) , this leads to a higher pressure on the lower side and lesser pressure on the upper side, which creates a pressure difference and is responsible for generation of ‘ lift ’. With the increase in AOA the lift starts increasing, because more pressure difference is created by turning the flow in the downward direction at the trailing edge and hence more lift force is generated.
The increase in lift co-efficient is observed only upto a certain point called “ stall limit” where the phenomenon called ‘ stalling’ happens. Stall is the maximum AOA that any airfoil can have where the lift starts decreasing from maximum lift generated. Once stalling happens the lift coefficient starts decreasing, because the lift starts decreasing with increasing angle of attack as the flow separation takes place and leads to the formation of eddies and vortices in the trailing edge. The formation of this is shown below
In the above figure one can see, at the trailing edge the flow is getting separated and results in the formation of vortices, which leads to decrease in lift generated and hence Cl decreases.
Plot of Cd v/s AOA :
Above graph shows the plot of drag coefficient for different angle of attack. With the increase in AOA the value of drag coefficient increases initially slowly upto a certain point (around 10-12º) and once the ‘stalling’ happens the drag coefficient increases rapidly, as the drag produced will be in greater values during stalling when compared to the initial AOA’s.
MESHING
AOA 1º
Base grid of dx = dy = dz = 0.036m was provided and converge estimated cells around 102013 cells or elements without any embedding. Since for us it is necessary to capture data around airfoil in accurate way, fixed embedding was added.
Above figure shows fixed embedding enabled around the airfoil with 2 levels of embedding layer and scale of 3.
The size of the cells or elements in the fixed embedding can be calculated as shown below.
Grid size = (base grid)/2^(embedding level)
= 0.0036/2^2
= 9e-4 m. (size of the new cells or elements in the fixed embedding).
AOA 5º
AOA 10º
AOA 15º
EFFECT OF TURBULENCE (TURBULENT VISCOSITY) ON DIFFERENT AOA
AOA 1º
AOA 5º
AOA 10º
AOA 15º
Above plots shows the variation of turbulent viscosity for different AOA.
Turbulent viscosity is the term used to describe how well the fluid is mixing or how well the flow is around the body. With the increase in AOA the turbulent viscosity decreases, this means that the flow over an airfoil is slowly separating due to the growth of boundary layer and formation of vortices.
PRESSURE AND VELOCITY PLOTS AND CONTOUR
AOA 1º
PRESSURE
As discussed above, in order to generate lift there should be a pressure diffrence created. In the above plot the pressure on the upper surface of the airfoil decreases, at the same time the lower surface pressure is
increasing which leads to diffrence in pressure that generates lift.
PRESSURE CONTOUR
The theory of the pressure profile can be observed here. In the above contour plot the pressure on the upper side is lesser than the pressure on the lower side, which leads to
a pressure diffrence and generation of lift.
VELOCITY
VELOCITY CONTOUR
In the above contour plot the velocity is clearly high on the upper surface of the airfoil and the velocity is lesser on the ower side of it. This is beacuse of the airfoil geometry where the upper surface
has a negative slope (curved downwards) where the air-molecules can move down easily and also they have their internal kinetic energy which increasess the velocity. On the other hand the lower side of airfoil
is somewhat curved upwards, where the air-molecules consumes their internal energy to climb it and hence their velocity decreases.
GLYPH (VECTOR CONTOUR TO UNDERSTAND FLOW BEHAVIOUR)
Above contour plot is for AOA of 1º , since there is no much AOA the flow is smooth and attatched to the body, there is no disturbance
in the flow over an airfoil.
AOA 5º
PRESSURE
PRESSURE CONTOUR
VELOCITY
VELOCITY CONTOUR
GLYPH (VECTOR CONTOUR TO UNDERSTAND FLOW BEHAVIOUR
There is slight disturbance in the flow, where the flow vectors at the leading edge is pileing up and moving upward and downard of anairfoil, but
still the flow is attatched to the body.
AOA 10º
PRESSURE
PRESSURE CONTOUR
VELOCITY
VELOCITY CONTOUR
GLYPH (VECTOR CONTOUR TO UNDERSTAND FLOW BEHAVIOUR)
A slight seperation in the flow can be observed here, where the flow is seperating from the leading edge of an airfoil. The maximum velocity over an airfoil for 10º
AOA is 54 m/s.
AOA 15º
PRESSURE
PRESSURE CONTOUR
VELOCITY
VELOCITY CONTOUR
GLYPH (VECTOR CONTOUR TO UNDERSTAND FLOW BEHAVIOUR)
Above contour plot shows the huge disturbances in the flow and how the flow is seperating in the leading and trailing edges. The velocity vectors or the flow
at the leading edge is very chaotic and is deviating too much. The flow is not fully attatched to the body, as stalling is taking place, where lift is decreasing. At the trailing edge
the flow is recirculating back , this means that the flow velocity is increasing on the pressure or lower side of the airfoil, where the pressure zone is decreasing and as the presssure diffrence
decreases the lift generation decreases and hence the lift coefficient decreases, also at the trailing edge there is a formation of vortices, which increases the drag on the airfoil, and hence thrust generated
in an airplane decreases.
CONCLUSION
The behaviour of flow over an airfoil is a very important concept in order to understand how the lift is generated and how an aircraft fly.
The flow is attatched to the body and the lift generation increases till the effect of stallling comes into picture. Once stalling happens the lift generation
decreases and flow seperation takes place leading to the formation of vortices and eddies at the trailing edge.
Hence the results that are caluclated from the simulation is as follows.
For AOA of 1º :
Cl = 0.364038 Cd = 0.01819For AOA of 5º :
Cl = 0.78058 Cd = 0.01958For AOA of 10º :
Cl = 1.25693 Cd = 0.03055For AOA of 15º :
Cl = 0.94122 Cd = 0.06365ANIMATION LINKS
for 15º :
Velocity contour
Pressure contour
GLYPH PLOT