Document
International Mechanical Engineering Congress & Exposition IMECE15 November 13-19, 2015, Houston, Texas, USA
IMECE2015 - 53052
NUMERICAL CFD SIMULATION AND TEST CORRELATION IN A FLI GHT PROJECT ENVIRONMENT K. K. Gupta S. F. Lung NASA Armstrong Flight Research Center Jacobs Technology Edwards, CA, USA Edwards, CA, USA A. H. Ibrahim Norfolk State University Norfolk, VA, USA ABSTRACT flow represented by the Navier - Stokes formulation. An unstructured gr id is used for domain decomposition.
This paper presents detailed description of a novel CFD The one equation model (Ref. 16) has been adapted for procedure and comparison of its solution results to that turbulence modelling. In this process both the viscous stresses obtained by other available CFD codes as well as actual flight pertaining to the linear viscous flow and the flux in the energy and wind tunnel test data pertaining to the GIII aircraft, equation, are duly modified.
currently undergoing flight testing at AFRC . It is t hen followed by detailed results of analyses which are next compared with actual flight test and wind tunnel simulation results. These results indicate that most CFD INTRODUCTION solutions compare reasonably well with the test data. The FE solutions in particular prove to be efficient and accurate and the Two in - house software as well as a number of related software are available for public use.
2 - 6 Finally, some summarizations and discussions of the commercially available CFD codes were used to analyze the problem, for comparison purposes. In this process both finite current effor t is given in the ‘Concluding Remarks’ section.
volume and finite element discretization were used for Euler and Navier - Stokes simulations. Both unstructured and structured grids were employed, as appropriate and solutions NOMENCLATURE were derived for Mach 0.701 and angle of attack = 3.92 AFRC = Armstrong Flight Research Center degree.
CFD = Computational Fluid Dynamics Extensive flight tests were performed for validation purposes.
FE = Finite Element Also these tests were complimented with detailed wind tunnel simulations. All such test results are c ompared with the t = time step numerical solution data obtained by the various CFD codes. = Density 7 8,9 Associated finite difference and finite volume techniques are = Dynamic viscosity 10,11 well described in the literature . The finite element = Viscous stress tensor technique for the discretization of fluid flow employs u = free stream velocity unstructured mesh and is based on a Taylor - Galerkin E = Total energy 13 - 15 procedure .
a = Shape function A description of the finite element fluids solver is = Conservation variable presented in some detail. It pertains to the solution of viscous f = Convection 𝒈 = Diffusion k = Thermal conductivity 𝜕𝒂 𝑇 𝑇 𝑴 = ∫ 𝒂 𝒂𝑑𝑉 ; 𝑲 = ∫ 𝒂 𝑢̅ 𝑑𝑉 ; 𝒊 p = Pressure 𝜕𝑥 𝑖 𝑉 𝑉 M = Mass matrix 𝜕𝑒 𝜕𝑒 𝑖 𝑖 𝑇 𝑇 K = Convection matrix 𝒇̂ = ∫ 𝒂 𝑝̅ 𝑑𝑉; 𝒇̂ = ∫ 𝒂 𝑒̅ 𝑑𝑉; 1 𝑖 2 𝑖 𝜕𝑥 𝜕𝑥 Re = Reynolds number 𝑉 𝑖 𝑉 𝑖 𝑇 𝑇 Pr = Prandtl number 𝜕𝒂 𝜕𝒂 𝑲 = − ∫ 𝑒 𝜎 𝑑𝑉 − ∫ 𝒎 𝑞 𝑑𝑉 ; 𝜎 𝑗 𝑖𝑗 𝑗 𝑗 𝜕𝑥 𝜕𝑥 𝑗 𝑗 𝑉 𝑉 PROCEDURE 𝑇 𝑇 (11) 𝒇 = ∫ 𝒂 𝑒 𝜎 𝑛̂𝑑Γ + ∫ 𝒂 𝒎 𝑞 𝑛̂𝑑Γ 𝜎 𝑗 𝑖𝑗 𝑗 𝑗 Γ Γ In these equations, 𝑝̅ , 𝑢̅ , 𝒆̅ are the average values; 𝒆 = The Navier - Sto kes equation can be written as 𝑖 𝑖 𝑖 1 𝑇 𝑇 𝑇 𝜕𝒗 𝜕𝒇 𝜕𝒈 [0 1 0 0 𝑢 ] ] ] 𝒊 𝒊 , 𝒆 = [0 0 1 0 𝑢 , 𝒆 = [0 0 0 1 𝑢 , 𝑹̂ is the 1 2 2 3 3 + + = 0 𝑖 = 1, 2, 3 (1) 𝑇 𝜕𝑡 𝜕𝑥 𝜕𝑥 𝑖 𝑖 artificial dissipation, and 𝒎 = 𝒎 = 𝒎 = [0 0 0 0 1] .
1 2 3 i n which the conservation variables, flux, and body force Turbulence terms are included by modifying the viscous column vectors, as well as the viscous stress are defined as effects.
𝑇 𝜌𝑢 𝜌𝐸] 17 𝒗 = [𝜌 𝑗 , 𝑗 = 1,2,3 (2) A novel two - step solution procedure is adopted for the 𝑇 (𝜌𝑢 𝑢 + 𝑝𝛿 ) 𝑢 (𝑝 + 𝜌𝐸)] 𝑓 = [𝜌𝑢 , 𝑗 = 1,2,3 (3) 𝑗 𝑖 𝑗 𝑖𝑗 𝑗 𝑗 flow equation, the inviscid solution being a ugmented with the 𝑇 𝜕𝑇 viscous term and stabilized with artificial dissipation terms.
𝜎 (𝑢 𝜎 + 𝑘 )] 𝐸 = [0 (4) 𝑖𝑗 𝑖 𝑖𝑗 𝑗 𝜕𝑥 Assuming, 𝑗 𝑇 𝑓 𝑢 𝑓 ] ∆𝒗̃ = 𝒗̃ − 𝒗̃ (12) 𝑓 = [0 (5) 𝑏 𝑖 𝑏 𝒏+𝟏 𝒏 𝑏 𝑖 𝑖 𝜕𝑢 𝜕𝑢 2 𝜕𝑢 then, 𝑖 𝑗 𝑙 𝜎 = 𝜇 [ + − 𝛿 ] 𝑙 = 1,2,3 (6) 𝑖𝑗 𝑖𝑗 −∆𝑡 𝜕𝑥 𝜕𝑥 3 𝜕𝑥 𝑗 𝑖 𝑙 𝑴(𝒗̃ − 𝒗̃ ) = [𝑐𝑴 + 𝑲](𝒗̃ + 𝒗̃ ) − ∆𝑡(𝒇̂ + 𝑛+1 𝑛 𝑛+1 𝑛 1 w h ere 𝑢 are velocity components in the 𝑥 coordinate system; 𝑖 𝑖 𝒇̂ ) (13) , p , and E are the density, pressure, and total energy which becomes respectively; is the dynamic viscosity; k is the thermal ∆𝑡 ∆𝑡 ∆𝑡 ∆𝑡 [(1 + 𝑐) 𝑴 + 𝑲] 𝒗̃ = [(1 − 𝑐) 𝑴 − 𝑲] 𝒗̃ + conductivity, the heat flux 𝑞 being −𝑘𝜕𝑇/𝜕𝑥 ; T is the 𝑛+1 𝑛 𝑗 𝑗 2 2 2 2 (14) temperature; 𝒇 represents the body forces. ∆𝑡𝑹 𝑏 or The preceding equations are nondimensionalised for [𝑴 ]𝒗̃ ]𝒗̃ numerical calculations. In this process the governing equations = [𝑴 + ∆𝑡𝑹 (15) + 𝑛+1 − 𝑛 where remain in the same form excepting 𝑔 , which becomes 𝑗 𝑇 𝑹 = −(𝒇̂ + 𝒇̂ ) (16) 1 2 𝑔 = [0 𝜎 (𝑢 𝜎 − 𝑞 )] (7) 𝑗 𝑖𝑗 𝑖 𝑖𝑗 𝑗 Let and also the viscous stress tensor and he at flux take the 𝑴 = 𝑫 + 𝑴′ (17) + + + following form : the matrix 𝑫 having diagonal elements. Equation (8) may then + 𝜇 𝜕𝑢 𝜕𝑢 2 𝜕𝑢 𝑖 𝑗 𝑙 𝜎 = [ + − 𝛿 ] be solved as follows.
𝑖𝑗 𝑖𝑗 𝑅𝑒 𝜕𝑥 𝜕𝑥 3 𝜕𝑥 𝑗 𝑖 𝑙 Step 1: Form 1 𝜕𝑇 [𝑫 ]𝒗̃ ]𝒗̃ ]𝒗̃ (18) = [𝑴 − [𝑴′ + ∆𝑡𝑹 𝑞 = (8) + 𝑛+1 − 𝑛 + 𝑛+1 𝑗 𝑅𝑒𝑃𝑟 𝜕𝑥 𝑗 Step 2: Solve 𝒗̃ iteratively 𝑛+1 (𝒊+𝟏) (𝒊) i n which the Reynolds number is defined as 𝑅𝑒 = 𝑢 𝐿/𝜐 ; (−𝟏) ′ ∞ ∞ 𝝊̃ = [𝑫 ] {[𝑴 ]𝝊̃ − [𝑴 ]𝝊̃ + ∆𝑡(𝑹 + 𝑹̂ + + − 𝒏 + 𝒏+𝟏 𝒏+𝟏 𝜐 = 𝜇 /𝜌 is termed the kinematic viscosity; P r is the ∞ ∞ ∞ 𝑲 + 𝒇 )} (19) 𝝈 𝝈 Prandtl number, 𝑃𝑟 = 𝜐 /𝛼 , with 𝛼 = 𝑘/(𝜌 𝑐 ) is the (𝒊+𝟏) (𝒊) ∞ ∞ ∞ ∞ 𝑝 Step 3: If ‖𝒗̃ ‖ ≠ EPS1‖𝒗̃ ‖ go to Step 2.
𝑛+1 𝑛+1 thermal diffusivity.
(𝒊+𝟏) (𝒊) Step 4: If ‖𝒗̃ ‖ ≠ EPS2‖𝒗̃ ‖ go to Step 1.
𝑛+1 𝑛+1 The Taylor’s expansion of the solution 𝒗(𝑥, 𝑡) in the t ime Step 5: Repeat Steps 1 to 4 NITER times until desired domain, neglecting second order term and body forces, yields 𝜕𝒇 𝜕𝒈 convergence is achieved, that is until 𝒗̃ ≈ 𝒗̃ ; EPS1 and 𝒊 𝒊 𝑛+1 𝑛 ∆𝒗 = −∆𝑡 [ + ] (9) 𝜕𝑥 𝜕𝑥 EPS2 are suitable convergence criteria factors, specified by the 𝑖 𝑖 (𝑡) users.
in which ∆𝒗 = 𝒗(𝑡 + ∆𝑡) − 𝒗(𝑡) . Applying Galerkin’s spatial The iterative process in Step 2 requires a small number of idealization 𝒗 = 𝒂𝒗̃ , 𝒗̃ being the nodal values and 𝒂 the shape steps, usually 1, and achieves a stabl e, convergent solution.
functions vector, the flow equation can be expressed as 𝜕𝑢 In regions of high pressure gradients, artificial dissipation term 𝑖 𝑴∆𝒗̃ = −∆𝑡 [ 𝑴 + 𝑲] 𝒗̃ − ∆𝑡(𝒇̂ + 𝒇̂ ) + ∆𝑡𝑹̂ + 1 2 𝜕𝑥 𝑖 is applied to prevent oscillations near discontinuities. This is ] (10) ∆𝑡 [𝑲 + 𝒇 implemented by incorporating pressure - switched diffusion 𝜎 𝜎 in which 𝑴 is the consistent mass matrix, 𝑲 the convection coefficients as appropriate. Thus, 𝐶 𝑆 matrix, 𝒇̂ , 𝒇̂ the pressure matrices, 𝑲 the second - order 𝑆 𝑒 1 2 𝜎 −1 [𝑀 ]𝜐̃ 𝑅̂ = 𝑀 − 𝑀 (20) 𝐿 𝑐 𝐿 𝑛 matrix that includes viscous and heat flux effects, and 𝒇 the 𝜎 ∆𝑡 boundary integral matrix from second - order terms. Then, i n which 𝐶 is a shock capturing constant, 𝑆 is the averaged 𝑆 𝑒 element value of the nodal pressure switch defined as right wing section was used for CFD analysis. Total number of |∑ 𝑝 − 𝑝 | 𝑖 𝑗 𝑆 = (21) 𝑖 CFD mesh using triangular element on wing surface is 31k for ∑(|𝑝 − 𝑝 |) 𝑖 𝑗 coarse mesh and 59k for finer mesh. Total number of 3 - D CFD and 𝑀 and 𝑀 are the consistent and lumped mass matrices 𝑐 𝐿 mesh using tetrahedron element in aerodynamic domain is 1.2m respectively; l is the node under consideration and j are the for co arse mesh and 2.8m for finer mesh.
nodes connected to i .
. The flight condition was for Mach 0.7 and angle of attack = To obtain the viscous components, 𝜎 in Eq. (4) is written 𝑖𝑗 3.92 degrees. Table 2 provides extensive description of relevant as analyses hardware employed for each of the participative code 𝜕𝑢 2 𝜇 𝜕𝑢 𝜇 𝜕𝑢 𝑗 𝑙 𝑖 𝜎 = − 𝛿 + ( + ) (22) and solution CPU time for a co nverged solution. The STARS 𝑖𝑗 𝑖𝑗 3 𝑅𝑒 𝜕𝑥 𝑅𝑒 𝜕𝑥 𝜕𝑥 𝑙 𝑗 𝑖 has two solution option modules, namely CFDSOL and MG a nd the diffusion flux of the Navier - Stokes equation being and both appear to be competitive in terms of solution time, 𝑇 1 𝜕𝑇 accuracy, grid size and CPU numbers.
𝑔 = (0 𝜎 𝜎 𝜎 𝑢 𝜎 + ) 𝑖 𝑖1 𝑖2 𝑖3 𝑗 𝑖𝑗 𝑅𝑒𝑃𝑟 𝜕𝑥 Figures 4 to 6 depict pressure (C ) distribution around the 𝑖 p 𝑖 = 1,2,3; 𝑗 = 1,2,3 (23) wing ai rfoil cr oss section at the wing 368.3 cm, 584.2 cm and is the nondimensional viscosity term, whereas R e and P r are 1003.5 cm span wise locations. Further, the wind tunnel and the Reynolds and Pra n dtl numbers, respectively. Next actual flight test results are also shown for comparison and components of 𝜕𝑔 /𝜕𝑥 are evaluated term by term and then validation purpose. Due to the proprietary nature of the wind 𝑖 𝑖 discretized by Galerkin approximation. tunnel and fight test data , actual sca les on the figures cannot be This procedure is adopted in the STARS - CFDSOL code shown. Each of the codes shows reasonable correlation; that enables effective solution of the Naviar - Stokes equation in solution of the CFDSOL and MG codes appear to be rather most flight regimes. close to the two test results.
Figure 7 depicts the C distribution along the airfoil at p different span locations.
NUMERICAL AND TEST R ESULTS CONCLUDING REMARKS Accuracy of the STARS CFD code was verified pertainin g to the Hyper - X flight vehicle, carrying the X - 43 ve hicle for The paper presents detailed comparison of solutions of the subsequent hypersonic flight at Mach 5.0 and 7.0. Table 1 GIII aircraft wing obtained by a number of commercially provides such a comparison of computational results and actual available CFD codes as well as two AFRC in - house codes that flight test d ata at various sensor locations; these data pertain to use a finite element fluids discretization employing the ascent state of Hyper - X at Mach 0.9 and an altitude o f unstructured grids; related formulations of the novel CFDSOL 22,500 ft. Figure 1 provides a graphical depiction of code are also presented in detail. Importantly these solutions comparison of the two sets of results, signifying accuracy of the are compared with actual flight and also wind tunnel test data.
relevant procedures. Also Table 1 shows the numerical values Each of the codes shows reasonable correlation; solution o f the of fli ght test and computed aerodynamic pressures ; excellent CFDSOL and MG codes appear to be rather close to the two correlation is observed for primary data values; the last three test results, particularly around the leading edge; further, use of values in the Table are comparatively small and hence pron e to a single CPU to derive solutions testifies to their cost measurement inaccuracy. This code was next used, along with a effectiveness.
variety of existing commercially available programs, to solve a practical proje ct problem. The results of which were also compared to that obtained by actual flight and wind tunnel REFERENCES tests.
The Gulfstream GIII airplane (Gulfstream Aerospace 1. Gupta, K. K. and Meek, J. L., A Primer on Corporation, Savannah, Georgia) , currently undergoing flight Multidisciplinary Finite Element Engineering tests at NASA AFRC, was chosen as the example problem for Analysis , AIAA Education Series, Sept. 2013.
ver ification purposes. The GIII business jet as shown in Figure 2. Prinkey, M., Shahnam, M., and Rogers, W. A., “SOFC 2 is being modified and instrumented by NASA's Armstrong FLUENT Model Theory Guide and User Manual,” Flight Research Center to serve as a test bed for a variety of Release Version 1.0, FLUENT, Inc., 2004.
flight research exp eriments, in support of the Environmentally 3. FUN3D Manual: 12.4 - 69883, B iedron, R. T., Derlaga, Responsible Aviation (ERA) project. The twin - turbofan aircraft J. M., Gnoffo, P. A., Hammond D. P., Jones, W. T., provides long - term capability for efficient testing of subsonic Kleb B., Lee - Rausch E. M., Nielse E. J., Park M.l A., flight experiments for NASA, the U.S. Air Force, other Rumsey C. L., Thomas J. L., and Wood W. A., NASA - government agencies, academia, and private industry.
TM - 2014 - 218179, 2014.
The wing span of the GIII aircraft is 23.7226 meter with 4. Johnson, F. T., Tinoco, E. N., and Yu, N. J., “Thirty sweep angle 27.66 degree. The airfoil section is a NACA 0012 Ye ars of Development and Application of CFD at modification. The aerodynamic model of the GIII wing used in Boeing Commercial Airplanes, Seattle,” Computers & the CFDSOL and MG solutions is shown in Figure 3; onl y the Fluids, Vol. 34, Issue 10, December 2005, pp. 1115 - 1151.
5. Pandya, M. J., Frink, N. T., and Noack, R. W, “ Progress Toward Overset - Grid Moving Body Capability for USM3D Unstructured Flow Solver,” th 17 AIAA CFD Conference, June 8 - 9 2005, Toronto, Canada.
6. CD - adapco, STAR - CCM + 7. Smith, G. D., Numerical Solution of Partial rd Differential Equations: Finite Difference Methods, 3 edn, Clarendon Press, Oxford, 1985.
8. Versteeg , H. K., and Malalasekera, W., An introduction to Computational Fluid Dynamics - The Finite Volume nd Method, 2 edn., Pearson/Prentice Hall, 2007.
9. Jameson, A., and Caughey, D., “ A Finite Volume Method for T ransonic Potential Flow,” AIAA p aper rd 77 - 635, 3 AIAA Computational Fluid Dynamics Conference, Albuquerque, New Mwxico, 1977.
10. Ferziger, Joel H., and Peric, Milovan, Computational rd Methods for Fluid Dynamics, 3 , rev. edition, Springer, 2002.
11. Fletcher, C. A. J., Computational Techniques for Fluid Dynamic s, Vols I and II, springer - Verlag, Berlin, 1991.
12. Zienkiewicz, O. C. and Taylor, R. L., The Finite Element Method – Vol. 2: Solid and Fluid Mechanics, McGraw - Hill, New York, 1991.
13. Donea, J., “A Taylor - Galerkin Method for Convective Transport Problems”, Inte rnational Journal of Numerical Methods in engineering, Vol. 20, No. 1, pp.
101 - 119, 1984.
14. Peraire, J., Peiro, J., Formaggia, L., Morgan, K., and Zienkiewicz, O. C., “Finite Element Euler Computations in three dimensions,” International Journal of Numerica l Methods in engineering, Vol.
26, No. 10, pp. 2135 - 2159, 1988.
15. Morgan, K., Peraire, J., and Peiro, J., “The Computation of Three Dimensional Flows using Unstructured Grids”, Computational Methods in Applied mechanics and engineering, Vol. 87, No. 3, pp. 335 - 332, 1991.
16. Spalart, P. R., and Allwaras, S. R.,”A One - Equation Turbulance Model for Aerodynamic Flows,” AIAA Paper 92 - 0439, 1992.
17. Gupta, K. K. and Meek, J. L., Finite Element Multidisciplinary Analysis , AIAA Education Series, Sept, 2003.
18. Gupta, K. K. and Bach, C., “Computational Fluid Dynamics – Based Aeroservoelastic Analysis with Hyper - X Applications”, AIAA Journal, Vol. 45, No. 7, pp 1459 - 1471, 2007.
19. Baumann, E., Hernandez, J., and Ruhf, J., “An Overview of NASA’s SubsoniC Research Aircraft Testbed (SCRAT),” AIAA - 2013 - 5083, 2013.
Table 1 Comparison of computed and flight test measured pressure data for the Hyper - X /X - 43 vehicle Pressure, Mpa Sensor Percent point Flight test CFD computed difference 001 0 . 01165 0.01193 2.34 003 0 . 01227 0.01164 6.12 007 - 0.00167 - 0. 00 096 42.12 085 - 0.00108 - 0.00268 147.99 090 0.00048 - 0.00055 2.56 Table 2 CFD Solvers Comparison No. of CFD Solver Flow Equation Platform Total CPU time Grid Size Note CPU number of processers is 7.2M an estimate, RANS, finite 6hr, 40min (533 volume, K - polyhedra/prismatic and the time STARCCM+ Cluster ~80 cpu hours) - omega SST for half model is an estimate 3000 iterations turbulence without T - tail for that number of processors 1 Intel 2.8 hr, (100 Euler, finite Dell M620 8GB Core i7 1.2 M Tetrahedrons STARS (MG) step s , 25 inner element Ram, 64 bit @2.67 for wing only cycles) GHz 1 Intel STARS Full N - S, finite Dell M620 8GB Core i7 13.8 hr (10000 2.8 M Tetrahedrons (CFDSOL) element Ram, 64 bit @2.67 step s ) for wing only GHz Full N - S, finite 1.9 M cells for half USM3D Mac 64 bit 2 CPUs 16 hr volume model without T - tail Full potential + Linux 1 CPU 2h, 28min 1.7M cells for full TRANAIR viscosity Workstation model (boundary layer) Fig. 1 Comparison of flight measured and calculated (CFD) pressure on Hyper - X/X - 43 vehicle Fig. 2 Grumman Gulfstream III ( G III) business jet .
(a) Domain Discretization (b) Surface Mesh Fig. 3 Aerodynamic model of the GIII aircraft wing.
STARS (MG) Wind Tunnel Flight Test USM3D STARS (CFDSOL) STARCCM TRANAIR Cp x/c Fig. 4 C plot at span station 145 p STARS (MG) Wing Tunnel Flight Test USM3D STARS (CFDSOL) STARCCM TRANAIR Cp x/c Fig. 5 C plot at span station 230 p STARS (MG) Wing Tunnel Flight Test USM3D STARS (CFDSOL) STARCCM TRANAIR Cp x/c Fig. 6 C plot at span station 385.
p (a) C distribution on wing surface (b) C at station 145 p p (c) C at station 230 (d) C at station 385 p p Fig.7 Typical C plots at various locations p