Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 872 COMPUTATIONAL MODELING OF TURBULENT FLOW AROUND AIRFOIL USING DIFFERENT TURBULENCE MODELS Asst. Lecturer Hakim Tarteeb Kadhim Mechanical Department AL_Dewaniyah Technical Institute Email: Eng.hakim84@yahoo.com Received 4 Jun 2013 Accepted 16 July 2014 ABSTRUCT The analysis of the two dimensional subsonic flow over a National Advisory Committee for Aeronautics (NACA) 0015 airfoil at two angles of attack, (α= 0°, α= 16°) and operating at a Reynolds number of 4.4 x 10 5 is presented. The flow was obtained by solving the steady-state governing equations of continuity and momentum conservation combined with one of three turbulence models [Spalart-Allmaras, Realizable and shear stress transport (SST)] aiming to the validation of these models through the comparison of the predictions and the free field experimental measurements for the selected airfoil. The aim of the work was to show the behavior of the airfoil at these conditions and to establish a verified solution method. Attention is focused on determination of the pressure distribution around the airfoil. Simulations were performed on an map quadratic structured grid with the Fluent (V6.3) software package which solves the Navier–Stokes equations by using finite volume methods. Calculations were done for constant air velocity altering only the angle of attack for every turbulence model tested. Calculations showed that the turbulence models used in commercial CFD codes do not give yet accurate results at high angles of attack and show that the Realizable k- model gave the most closness prediction of pressure distribution when compared with the experimental data [1]. نمذجة حسابية لجريان مضطرب حول مقطع جناح بأستخدام نماذج اضطراب مختلفة م.م حاكم ترتيب كاظم /الديوانية قسم الميكانيك المعهد التقني– الخالصة : لزوايا هجوم (NACA 0015)قدم هذا البحث تحليل لجريان ثنائي األبعاد دون سرعة الصوت لمقطع جناح من نوع (α= 0°, α=16°), (4.4)لعدد رينولدز x 10 5 . الجريان تضمن حل المعادالت الحاكمة لالستمرارية وحفظ الزخم للحالة Spalart-Allmaras, Realizable and shear stress)) المستقره والمتوفره في كل من النماذج المستخدمه في البحث transport SST تم في هذه النتائج المستحصلة مع البيانات التجريبية.ل مقارنة لغرض اختبار صالحية هذه النماذج من خال الدراسة بيان سلوك مقطع الجناح وتثبيت طريقة حل عند هذه الظروف. وكذلك اهتمام رئيسي لتحديد توزيع الضغط حول (Quadratic)ألضالع مقطع الجناح. تم اجراء المحاكاة على شبكة مبنية بشكل منتظم , نوع العنصر المستخدم فيها رباعي ا ستوك بواسطة استخدام طريقة الحجوم المحدودة -( الذي يحل معادالت نافيرFluent V6.3مع استخدام برنامج الحاسوب ) (Finite Volumes) تم تثبيت سرعة الجريان وتغير زوايا الهجوم لكل نموذج اضطراب.النتائج بينت بأن النماذج . mailto:Eng.hakim84@yahoo.com Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 872 بأن المقارنة بين النتائج النظرية والعملية ال تعطي نتائج دقيقه عند زوايا الهجوم العاليه. وكذلك بينت المستخدمه في البرامج التجاريه .[1]يتنبأ توزيع معامالت ضغط مقاربة للنتائج العملية (Realizable k-) نموذج االضطراب NOMENCLATUR  : Density (kg/ m 3 ) u : Velocity component in x direction (m/s) v : Velocity component in y direction (m/s) Re : Reynolds Number X Axial coordinate in the physical domain C : Cord, (m) Cp: Coefficient of pressure x & y; Coordinate direction ABBREVIATIONS CFD: Computational Fluid Dynamics RANS: Reynolds Average Navier-Stokes Equation SST: Shear Stress Transport S-A: Spalart-Allmaras NACA: National Advisory Committee for Aeronautics INTRODUCTION Simulation of the turbulent flow field around airfoil sections is an up to date problem. Many phenomena occur in the flow passing a real airfoil section, such as the transition from the laminar to turbulent region, laminar or turbulent separation, and reattachment of separated flow, etc. This makes this problem complex and difficult to solve numerically. It is important to see how these airfoils can be computationally modelled and what sorts of models produce the most accurate results, in comparison to experimental data [2]. In industrial CFD applications, RANS modelling remains one of the main approaches when dealing with turbulent flows. Over decades, this has consequently facilitated the development of a great variety of RANS turbulence models. Modelling approaches in the context of RANS have shown different degrees of success in various engineering applications, spanning from mixing- length models, linear and nonlinear eddy viscosity models to algebraic and differential Reynolds stress models with a hierarchy of increasing complexity in the modelling formulation and related CFD implementations [3]. Two-equation turbulence models or competitive one-equation models, are considered to be the minimum level of closure in Navier–Stokes solution methods that are capable to cope with complex flows. Two-equation eddy-viscosity models have been credited certain advantages, mostly related to their simplicity and their superior numerical properties, compared to other more sophisticated models, like those based on Reynolds stress closures. In their conventional forms, two-equation models often yield poor predictions for flows that encounter adverse pressure gradients, curvature, rotation or complex strain fields. In order to overcome these deficiencies, without resorting to more sophisticated and more expensive turbulence models, various modifications have been proposed in the literature [4]. Among them, those employed in the present method. Working even with two-equation models, the low-Reynolds terms which are activated close to solid walls are a source of convergence difficulties. For structured grids, a series of papers Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 addressed, in the past, the practical implementation of these models and delineated achievements or suggestions to improve their numerical behavior. In addition, practical problems about how to approximate distances from solid walls or how to compute derivatives normal to the wall are readily solvable on structured grids, with both precision and simplicity. In contrast, the relevant approaches for unstructured grids with triangular elements are, by no means, trivial; they often suffer from ambiguities or lack of precision during the treatment of the corresponding terms [5, 6]. The aim of the present work is to investigate the effect of several turbulent models, available in the commercial code FLUENT (V6.3), on the pressure coefficient distribution around the NACA 0015 airfoil at two angles of attack, (α= 0°, α= 16°) and existing experimental data from [1], are performed to validate the computational results. The velocity distribution with the Realizable k-ɛ turbulence model at (α= 0°, α= 16°) will also be found. GOVERNING EQUATIONS For all flows, the solver solves conservation equations for mass and momentum. Additional transport equations are also solved when the flow is turbulent. The equation for conservation of mass or continuity equation, can be written as follows[7],: m Su t p    )(.   (1) Equation (1) is the general form of the mass conservation equation which valid for incompressible as well as compressible flows. The source m S is the mass added to the continuous phase from the dispersed second phase (for example, due to vaporization of liquid droplets) and any user-defined sources. Conservation of momentum in an inertial reference frame is described by Equation (2) Fgpuuu t       )(.)(.)( (2) Where p is the static pressure,   is the stress tensor (described below) and g   and F  are the gravitational body force and external body forces (for example, that arise from interaction with the dispersed phase), respectively. F  also contains other model-dependent source terms such as porous-media and user-defined sources. The stress tensor is given by:   Iuuu T  . 3 2          (3) Where  is the molecular viscosity, I is the unit tensor, and the second term on the right hand side is the effect of volume dilation. For the 2-D, steady and incompressible flow the continuity equation is: Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 0      y v x u (4) Momentum equations for viscous flow in x and y directions are, respectively: x yxxx f yxx p Dt Du              (5) y yyxy f yxy p Dt Dv              (6) Where due to characteristics of the 2-D flow in continuity equation the term z y   and in momentum equation, z zx   and z zy   drop out. The flow was obtained by solving the governing equations of continuity and momentum conservation combined with one of three turbulence models. COMPUTATIONAL METHOD In this paper, the NACA 0015, the well documented airfoil from the 4-digit series of NACA airfoils, was utilized. The NACA 0015 airfoil is symmetrical; the 00 indicates that it has no camber. The 15 indicates that the airfoil has a 15% thickness to chord length ratio; it is 15% as thick as it is long. Reynolds number for the simulations was 4.4 x 10 5 , same with the reliable experimental data from[1], in order to validate the present simulation. The free stream temperature is 300 K, which is the same as the environmental temperature. The density of the air at the given temperature is ρ=1.225 kg/ m 3 and the viscosity is μ=1.7894×10 -5 kg/ms. For this Reynolds number, the flow can be described as incompressible. This is an assumption close to reality and it is not necessary to resolve the energy equation. A segregated, implicit solver was utilized (ANSYS Fluent 6.3.26., 2006) Calculations were done at two angles of attack (α= 0°, α= 16°). The airfoil was modelled in the software package GAMBIT (2.4.6). The geometry of the airfoil was imported, from an Airfoils [8], and a c-shaped mesh created. The large distance of the domain boundary from the airfoil, in comparison to the airfoil chord length, was so that the boundary conditions had minimal influence on the flow field in the vicinity of the airfoil and allowed the boundary conditions to be simplified [9]. This domain, shown in Figure (1). A fine mesh was needed in the vicinity of the airfoil in order to model the flow field accurately. In the vicinity of the airfoil it was necessary to refine the mesh at the trailing edge in order to ensure accurate modelling, as this would be where events of interest would occur. It was important to maintain good aspect ratios of cells across the domain as well as maintaining smooth changes in cell size, thus quadrilateral elements were used. Figure (2) shows the computational domain and the mesh that was used to solve the problem. Figure (3) shows the mesh around the airfoil. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 828 TURBULANCE MODELS The inlet boundary velocity U, was set to 35 m/sec for all turbulence models for direct comparison with the experimental result [1]. The corresponding Reynolds number is 4.4 x 10 5 based on the chord c of the airfoil (180 mm). A computational grid of 15100 cells was fixed for all models. Three different turbulence models were used, two equation models such as Realizable and SST k- model and Spalart-Allmaras model. These models selected because they are most widely used in aerodynamic industry, and they have well documented strength. Also these models proved to have a superior performance for flows involving strong streamline curvature [10]. All computations have been performed on the same grid to ensure that the presented solution for each model will be compared with each other. Flow conditions around the airfoil were built up by finite volume analysis using FLUENT 6.3 software. RESULTS AND DISCUSSION Simulations for two angles of attack were done in order to be able to compare the results from the different turbulence models and then validate them with existing experimental data from [1]. To do so, the model was solved at zero angle of attack and high angle of attack, α= 16°. Turbulence models are often judged for accuracy based on the comparison of their predictions with the experimentally observed values. Figure (4) shows the variation in the pressure coefficient around the airfoil, for each model and experimental data, at zero degrees angle of attack. The experimental results showed scattered values upon the numerical pressure coefficient. Turbulence models predictions are very close to each other and all it predict high pressure coefficients values at the nose and slightly high pressure coefficients values for the residual part of airfoil when compared with experimental values. Figure (5) shows the variation of the pressure coefficient along the lower surface of the airfoil at angle of attack, α= 16°. All turbulence models predict high pressure coefficients at the nose and low pressure coefficients values for the residual part of airfoil when compared with experimental values and all models show a unsimilar profile behavior between them. Figure (6) shows the variation of the pressure coefficient along the upper surface of the airfoil at angle of attack, α= 16°. Also the experimental results showed scattered values upon the numerical pressure coefficient. The S-A and SST k-ω predictions are very close to each other and their prediction is not acceptable with the experimental data R k-ε model predict close values to experimental data when compared with the other turbulence models. For two angles of attack (α= 0°, α= 16°) the Realizable k-ε model gave the most closeness prediction of pressure distribution when compared with the experimental data for both suction and pressure sides. Figures (7) and (8) show the simulation outcomes of static pressure at angles of attack 0° and 16°with the Realizable k-ε turbulence model. The pressure on the lower surface of the airfoil was greater than that of the incoming flow stream and as a result it effectively “pushed” the airfoil upward, normal to the incoming flow stream. On the other hand, the components of the pressure distribution parallel to the incoming flow stream tended to slow the velocity of the incoming flow relative to the airfoil, as do the viscous stresses. Contours of velocity magnitude at angles of attack (α= 0°, α= 16°) are also shown in figure (9) and (10). At zero angle of attack the velocity distribution is symmetric where the NACA 0015 airfoil is symmetrical while at sixteen angle of attack, the upper surface of the airfoil experienced a higher velocity compared to the lower surface until the separation point where the velocity Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 828 equal zero. Velocity vectors of the flow field around airfoil for different turbulence models are presented in figure (11). The SST k-ω model and S-A model predict very early separation with regard to R k-ε model. That was expected from the pressure distribution where the wall shear stresses occurred on airfoil gave an idea about the flow separation. CONCLOUSION This paper showed the behavior of the 4-digit symmetric airfoil NACA 0015 at two angles of attack. The most appropriate turbulence model for these simulations was the Realizable k-ɛ two equation model, which had the most closeness prediction of pressure distribution when compared with the experimental data for both suction and pressure sides. Turbulence models predictions are very close to each other at zero degrees angle of attack while at high angle of attack, α= 16° the numerical pressure coefficient distributions showed scattered values upon the experimental results. In summary it can be said that the major error in the three turbulence models used for two angles of attack, is at the leading edge. This is apparent from the analysis of the pressure coefficient distribution around the airfoil. REFERENCE [1] Ahmed, Obaid “Study of Aerodynamic Characteristics of An Airfoil with Bumps”, M.Sc. Thesis, Department of Machines and Equipments, University of Technology, 2009. [2] Hills, J. L., “Numerical Modeling of Turbulent Flow Past an Airfoil”, pp.1-10, April, 2005. [3]Shia-Hui, Peng. and Peter, Eliasson. “Examination of the Shear Stress Transport Assumption with a Low-Reynolds Number k − ω Model for Aerodynamic Flows”, Chalmers University of Technology, SE-412 96 Gothenburg, Sweden, 2009. [4] D.G. Koubogiannis. A.N. Athanasiadis. and K.C. Giannakoglou “One- and two-equation turbulence models for the prediction of complex cascade flows using unstructured grids”, Laboratory of Thermal Turbomachines, National Technical University of Athens, P.O. Box 64069, 15710 Athens, Greece,2001. [5] Barth TJ. “Aspects of unstructured grids and finite-volume solvers for the Euler and Navier– Stokes equations”, AGARD Report 787. Special Course on Unstructured Grid Methods for Advection Dominated Flows, 1992. [6] Farhat C. and Lanteri S. “Simulation of compressible viscous flows on a variety of MPPs: computational algorithms for unstructured dynamic meshes and performance results” Comp Meth Appl Mech Engng 94;119:35–60. [7] Douvi C. E., Tsavalos I. A., and Margaris P. D., “Evaluation of the turbulence models for the simulation of the flow over a National Advisory Committee for Aeronautics (NACA) 0012 airfoil”, Journal of Mechanical Engineering Research Vol. 4(3), pp. 100-111, March 2012. [8] Abbtt, I.H. and Von Doenhoff, A.E. “Theory of Wing Sections”, Dover publications, INC. 1959. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 [9] Baxevanou, C.A., and Fidaros, D.K., “Validation of Numerical Schemes and Turbulence Models Combinations for Transient Flow around Airfoil”, Engineering Applications of Computational Fluid Mechanics Vol. 2, No. 2, pp. 208–221 , 2008. [10] Omar Badran, Regis Quadros and Fettah Aldudak, “Two-Equation Turbulence Models for Turbulent Flow over a NACA 4412 Airfoil at Angle of Attack 15 Degree”, p.p 1-8, 2003. Figure (1): Geometry of the Model, Where c = chord length of the airfoil Figure 3 – Computational Domain Figure (2): Computational Domain Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 Figure (3) Mesh around the Airfoil Figure (4) Variation of Pressure Coefficient around Airfoil at Zero Degrees Angle of Attack Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 Figure (5): Variation of the Pressure Coefficient along the Lower Surface of the Airfoil at 16° Angle of Attack. Figure (6): Variation of the Pressure Coefficient along the Upper Surface of the Airfoil at 16° Angle of Attack. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 827 Figure (7): Contours of Static Pressure at Zero Angle of Attack with the Realizable k-ɛ Turbulence Model. Figure (8): Contours of Static Pressure at 16° Angle of Attack with the Realizable k-ɛ Turbulence Model. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 Figure (9): Contours of Velocity Magnitude at Zero Angle of Attack with the Realizable k-ɛ Turbulence Model. Figure (10): Contours of Velocity Magnitude at 16° Angle of Attack with the Realizable k-ɛ Turbulence Model. Al-Qadisiya Journal For Engineering Sciences, Vol. 7……No. 3 ….2014 822 Figure (11): Distribution of Velocity Vector around Airfoil with Different Turbulence Models