

small (250x250 max)
medium (500x500 max)
Large
Extra Large
large ( > 500x500)
Full Resolution


A PARAMETRIC ANALYSIS OF GUSTINDUCED AIRFOIL SURFACEPRESSURE PHASE BEHAVIOR By HARIKISHIN PRAKASH BAKHTIANI Bachelor of Science Oklahoma State University Stillwater, Oklahoma 2002 Submitted to the Faculty of the Graduate College of the Oklahoma State University in partial fulfillment of the requirements for the Degree of MASTER OF SCIENCE July, 2004 ii A PARAMETRIC ANALYSIS OF GUSTINDUCED AIRFOIL SURFACEPRESSURE PHASE BEHAVIOR Thesis Approved: ________________________________________________ Thesis Advisor ________________________________________________ ________________________________________________ _______________________________________________ Dean of the Graduate College iii ACKNOWLEGMENTS A work of this magnitude was made possible only by the assistance of many individuals. First and foremost, I would like to thank my graduate advisor and committee chair Dr. Eric A. Falk. Thank you, sir, for your continuous support, guidance and hard work. You are truly a great friend and mentor. Words cannot express my appreciation of your encouragement, inspiration and patience. Gratitude and appreciation are also extended to my other committee members, Dr. Andrew S. Arena and Dr. Afshin J.Ghajar. To my friends and colleagues at the Turbolab, it’s been both an honor and privilege to have met and worked with you all. Aaron thanks for all your help. Robert thanks for answering all my questions and guiding me when direction was needed. Finally, a special ‘thank you’ to my parents, Prakash and Jaya Bakhtiani for their continuous love, support and words of encouragement. The successes I have achieved did not come without certain sacrifices, which they endured in some form. Finally to my brother, Nicky, who never fails to amaze me. TABLE OF CONTENTS 1 INTRODUCTION..................................................................................................... 1 1.1 OVERVIEW OF AIRFOIL GUST INTERACTIONS....................................... 1 1.2 COMPRESSOR DESIGN INTENT AND OPERATING ENVIRONMENT ... 1 1.3 NORMAL COMPRESSOR OPERATION LEADING TO AIRFOILGUST INTERACTIONS ............................................................................................... 2 1.4 IMPORTANCE OF HCF TO ENGINE COMMUNITY ................................... 4 1.5 IMPORTANCE OF PHASE DISTRIBUTIONS IN PREDICTING AIRFOIL MODAL FORCING ........................................................................................... 5 1.6 NEED FOR CONTINUED SURFACEPRESSURE PHASE ANALYSIS ...... 6 1.7 SCOPE OF CURRENT INVESTIGATION ...................................................... 7 2 PREVIOUS WORK.................................................................................................. 8 2.1 LITERARY REVIEW ........................................................................................ 8 2.1.1 SingleAirfoil Investigations....................................................................... 8 2.1.2 Turbomachine and Cascade Investigations................................................. 9 2.2 UNRESOLVED ISSUES.................................................................................. 14 2.3 CURRENT RESEARCH OBJECTIVE............................................................ 20 3 COMPUTATIONAL METHODOLOGY AND SETUP .................................... 22 3.1 AIRFOIL GEOMETRY AND BOUNDARY CONDITIONS......................... 22 3.2 STATORVANE GEOMETRY AND BOUNDARY CONDITIONS............. 23 3.3 GENERAL FLUENT SOLVER DESCRIPTION............................................ 25 iv 3.3.1 Coupled Solution Method ......................................................................... 27 3.3.2 ReynoldsAveraged NavierStokes (RANS) Equations ........................... 28 3.3.3 k – ε Turbulence Model ............................................................................ 30 3.3.4 FiniteVolume Discretization Methodology ............................................. 31 3.3.5 SecondOrder Upwind Scheme................................................................. 32 3.3.6 SecondOrder Time Discretization ........................................................... 32 3.3.7 Linearization Methodology....................................................................... 34 3.3.8 Periodic Boundary Conditions.................................................................. 35 3.3.9 Standard FLUENT Inlet/Outlet/Wall Boundary Conditions .................... 35 3.3.10 Operating Pressure .................................................................................... 36 3.4 UDF DESCRIPTION........................................................................................ 37 3.4.1 Development Logic/Procedure ................................................................. 37 3.5 COMPUTATIONAL GRID DESCRIPTION .................................................. 39 3.5.1 Gambit Grid Generation Software ............................................................ 39 3.5.2 Grid Methodology..................................................................................... 40 3.6 COMPUTATIONAL SETUP........................................................................... 46 3.6.1 Reference Values ...................................................................................... 48 3.6.2 Boundary Conditions ................................................................................ 48 3.7 GRID INDEPENDENCE ................................................................................. 49 4 TIMEAVERAGED RESULTS ............................................................................ 52 4.1 TIMEAVERAGED METHODOLOGY ......................................................... 52 4.2 BASELINE AIRFOIL TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS ............................................................................................ 53 v 4.3 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (THICKNESS INFLUENCE)................................................................................................... 57 4.4 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (CAMBER INFLUENCE)................................................................................................... 58 4.5 COMPARISON OF TIMEACCURATE BASELINE LIFT DEPENDENCY TO SEAR’S RESULTS.................................................................................... 60 4.6 TIMEAVERAGED RESULTS SUMMARY ................................................. 61 5 RESULTS FOR NACA 0012 BASELINE CONFIGURATION ........................ 63 5.1 DATA REDUCTION METHODOLOGY ....................................................... 63 5.2 AIRFOIL FORCINGFUNCTION TIME DEPENDENCY ............................ 64 5.3 AIRFOIL FORCINGFUNCTION SPECTRAL CONTENT.......................... 67 5.4 AIRFOIL FORCINGFUNCTION PHASE DEPENDENCY ......................... 67 5.5 AIRFOIL SURFACEPRESSURE TIME DEPENDENCY............................ 69 5.6 AIRFOIL SURFACEPRESSURE SPECTRAL CONTENT.......................... 73 5.7 AIRFOIL SURFACEPRESSURE FIRST HARMONIC AMPLITUDE CHORDWISE DEPENDENCY....................................................................... 75 5.8 AIRFOIL SURFACEPRESSURE FIRST HARMONIC CHORDWISE ... PHASE DEPENDENCY.................................................................................. 77 5.9 AIRFOIL SURFACEPRESSURE ANALYTICAL MODEL......................... 79 5.9.1 Interaction Model...................................................................................... 80 5.9.2 Interaction Model Results......................................................................... 86 5.10 SUMMARY...................................................................................................... 88 6 PARAMETRIC ANALYSIS.................................................................................. 90 6.1 DATA REDUCTION METHODOLOGY ....................................................... 90 6.2 AIRFOIL THICKNESS INFLUENCE............................................................. 90 6.3 AIRFOIL CAMBER INFLUENCE.................................................................. 93 vi 6.4 ANGLE OF ATTACK INFLUENCE .............................................................. 96 6.5 SUMMARY...................................................................................................... 98 7 STATORVANE CASCADE RESULTS............................................................ 100 7.1 DATA REDUCTION METHODOLOGY ..................................................... 100 7.1.1 TimeAveraged Results .......................................................................... 100 7.1.2 Unsteady Results..................................................................................... 101 7.2 STATORVANE CASCADE CONFIGURATION....................................... 101 7.3 TIMEAVERAGED RESULTS..................................................................... 101 7.4 UNSTEADY PRESSURE RESULTS............................................................ 103 7.5 SUMMARY.................................................................................................... 107 8 SUMMARY AND CONCLUSIONS ................................................................... 108 8.1 RESULTS ....................................................................................................... 108 8.1.1 TimeAveraged Results .......................................................................... 108 8.1.2 Unsteady Pressure Data .......................................................................... 110 8.2 CORRELATIONS WITH PREVIOUS INVESTIGATIONS........................ 113 8.3 CURRENT CONTRIBUTIONS..................................................................... 114 8.4 RECOMMENDATIONS FOR FUTURE WORK ......................................... 115 A. APPENDIX A........................................................................................................ 120 B. APPENDIX B ........................................................................................................ 122 C. APPENDIX C........................................................................................................ 124 D. APPENDIX D........................................................................................................ 128 E. APPENDIX E ........................................................................................................ 131 F. APPENDIX F ........................................................................................................ 135 vii G. APPENDIX G........................................................................................................ 136 viii LIST OF TABLES Table 3.1 Modeled StatorVane Cascade Geometry......................................................... 24 Table 3.2 Grid Distribution (AirfoilCascade Mesh)........................................................ 41 Table 3.3 BoundaryLayer Mesh Characteristics (AirfoilCascade Mesh). ..................... 43 Table 3.4 StatorVane Cascade Mesh Node Distribution................................................. 45 Table 3.5 FLUENT Configuration for Numerical Simulations. ....................................... 47 Table 3.6 FLUENT Reference Values.............................................................................. 48 Table 3.7 FLUENT Boundary Conditions........................................................................ 49 Table 4.1 Coefficient of Lift: 5° and 10° Angle of Attack. .............................................. 56 Table 4.2 Coefficient of Lift: Various LiftingSurface Cambers. .................................... 60 ix LIST OF FIGURES Figure 21 Relative Phase Variation with Chord: RearwardForced Cascade of Fabian et al. [1996]. .................................................................................................................. 12 Figure 22 Relative Phase Variation with Chord: ForwardForced Cascade of Fabian et al. [1996]. .................................................................................................................. 13 Figure 23 Unsteady Differential SurfacePressure Time Series: (φ (x / c) = 0). .............. 16 Figure 24 Unsteady Differential SurfacePressure Time Series: (φ (x / c) ≠ 0). .............. 16 Figure 25 Unsteady Differential SurfacePressure Variation with Chord: (φ (x / c) = 0). 17 Figure 26 Unsteady Differential SurfacePressure Variation with Chord:(φ (x / c) ≠ 0). 17 Figure 27 RigidBody and FlexibleBody Mode Shapes for a Simply Supported, Infinite Span, TwoDimensional Lifting Surface. ................................................................. 18 Figure 28 Maximum Generalized forces on Structural Modes for Various Phase Distributions.............................................................................................................. 19 Figure 31 AirfoilCascade Computational Boundaries. .................................................. 23 Figure 32 StatorVane Cascade Computational Boundaries. .......................................... 25 Figure 33 Overview of the Coupled Solution Method [Fluent, 2001]. ........................... 28 Figure 34 Rotor Wake Characteristics............................................................................. 38 Figure 35 AirfoilCascade Mesh (NACA 0012 airfoil)................................................... 42 Figure 36 Structured BoundaryLayer Mesh Surrounding NACA 0012 Airfoil Geometry. .................................................................................................................................. 43 Figure 37 StatorVane Cascade. ...................................................................................... 45 Figure 38 Modeled StatorVane Cascade Mesh.............................................................. 46 x Figure 39 Lift Coefficient Time History Showing Convergence: NACA 0012 Profile.. 47 Figure 310 SteadyState Pressure Coefficient Data: NACA 0012. ................................. 50 Figure 311 Pressure Coefficient Data: StatorVane. ....................................................... 51 Figure 41 TimeAveraged TotalPressure Contours (Pa): NACA 0012 Lifting Surface 53 Figure 42 NACA 0012 Lifting Surface Wake Profile..................................................... 54 Figure 43 TimeAveraged Static Pressure : NACA 0012. .............................................. 55 Figure 44 TimeAveraged Static Pressure : 0°, 5° and 10° Angle of Attack .................. 56 Figure 45 TimeAveraged StaticPressure for Various LiftingSurface Thicknesses..... 57 Figure 46 Time averaged Static Pressure Comparison with Experimental data. ............ 58 Figure 47 TimeAveraged StaticPressure Distribution for Various LiftingSurface Cambers .................................................................................................................... 59 Figure 48 LiftTime series Comparison with Sears results. ............................................ 61 Figure 51 Forcing Function TotalPressure Contours, t = 0............................................ 64 Figure 52 Forcing Function TotalPressure Contours, t = T/4........................................ 64 Figure 53 Forcing Function TotalPressure Contours, t = T/2. ....................................... 65 Figure 54 Forcing Function TotalPressure Contours, t = 3T/4. ..................................... 65 Figure 55 TotalPressure Time Series Forward of Airfoil. ............................................. 66 Figure 56 StaticPressure Time series. ............................................................................ 66 Figure 57 Forcing Function TotalPressure Spectral Content. ........................................ 67 Figure 58 Forcing Function StaticPressure Spectral Content. ....................................... 67 Figure 59 Airfoil Forcing Function 1st Harmonic Phase. ................................................ 68 Figure 510 Airfoil Forcing Function Higher Harmonic Phase........................................ 69 Figure 511 Airfoil StaticPressure Contours, t = 0......................................................... 69 xi Figure 512 Airfoil StaticPressure Contours, t = T/4. .................................................... 69 Figure 513 Airfoil StaticPressure Contours, t = T/2. .................................................... 70 Figure 514 Airfoil StaticPressure Contours, t = 3T/4. .................................................. 70 Figure 515 NACA 0012 Upper Surface Unsteady StaticPressure Time series.............. 71 Figure 516 NACA 0012 Lower Surface Unsteady StaticPressure Time series. ............ 72 Figure 517 NACA 0012 UpperSurface Static Pressure Spectral Content. .................... 74 Figure 518 NACA 0012 LowerSurface StaticPressure Spectral Content..................... 75 Figure 519 NACA 0012 Upper Surface 1st Harmonic Pressure Time series. ................. 76 Figure 520 NACA 0012 Lower Surface 1st Harmonic Pressure Time series.................. 76 Figure 521 NACA 0012 1st Harmonic Pressure Amplitude............................................ 77 Figure 522 NACA 0012 1st Harmonic Phase. ................................................................. 78 Figure 523 Wake Induced Pressure Time Series Collected at the Periodic Boundary.... 82 Figure 524 Optimized LiftInduced Pressure Time Series.............................................. 83 Figure 525 Interaction Model 1st Harmonic Time series................................................. 84 Figure 526 Propagating Disturbance Model.................................................................... 85 Figure 527 1st Harmonic Amplitude Comparison. .......................................................... 86 Figure 528 Interaction Model 1st Harmonic Phase.......................................................... 87 Figure 61 1st Harmonic Pressure Amplitude: Various LiftingSurface Thicknesses. .... 92 Figure 62 1st Harmonic Phase : Various LiftingSurface Thicknesses........................... 93 Figure 63 1st Harmonic Amplitude: Various LiftingSurface Cambers. ........................ 94 Figure 64 1st Harmonic Phase: Various LiftingSurface Cambers. ................................. 95 Figure 65 1st Harmonic Amplitude: Various MeanFlow Angles of Attack. .................. 97 Figure 66 1st Harmonic Phase: Various Mean Flow Angles of Attack. .......................... 98 xii Figure 71 TimeAveraged Total Pressure Contours (Pa)............................................. 102 Figure 72 TimeAveraged Static Pressure Distribution. ............................................... 103 Figure 73 TimeAverage Pressure Distribution [Falk, 2000]………………………… 103 Figure 74 RMS Unsteady Pressure Distribution……………………………………... 109 Figure 75 P'RMS Distribution [Falk, 2000]. .................................................................... 104 Figure 76.1st Harmonic Amplitude: StatorVane. ......................................................... 105 Figure 77 1st Harmonic Phase: StatorVane. ................................................................. 106 Figure D1 Unsteady Differential Static Pressure TimeSeries...................................... 128 Figure D2 Unsteady DifferentialPressure Spectral Content. ....................................... 129 Figure D3 1st Harmonic DifferentialPressure TimeSeries. ........................................ 129 Figure D4 Higher Harmonic Amplitudes: NACA 0012................................................ 130 Figure D5 Higher Harmonic Phase: NACA 0012. ........................................................ 130 Figure E6 1st Harmonic TimeSeries: NACA 0010 Upper Surface. ............................. 131 Figure E7 1st Harmonic TimeSeries: NACA 0010 Lower Surface.............................. 131 Figure E8 1st Harmonic TimeSeries: NACA 0015 Upper Surface. ............................. 132 Figure E9 1st Harmonic TimeSeries: NACA 0015 Lower Surface.............................. 132 Figure E10 1st Harmonic TimeSeries: NACA 0020 Upper Surface. ........................... 132 Figure E11 1st Harmonic TimeSeries: NACA 0020 Lower Surface............................ 132 Figure E12 SurfacePressure Spectral Content: NACA 0010 Upper Surface............... 133 Figure E13 SurfacePressure Spectral Content: NACA 0015 Upper Surface............... 133 Figure E14 SurfacePressure Spectral Content: NACA 0020 Upper Surface............... 134 Figure F15 1st Harmonic TimeSeries: 2% Camber Airfoil Upper Surface. ................ 135 Figure F16 1st Harmonic TimeSeries: 2% Camber Airfoil Lower Surface.................. 135 xiii Figure F17 1st Harmonic TimeSeries: 6% Camber Airfoil Upper Surface. ................. 135 Figure F18 1st Harmonic TimeSeries: 6% Camber Airfoil Lower Surface.................. 135 Figure G19 1st Harmonic TimeSeries: 5Degree AOA Upper Surface........................ 136 Figure G20 1st Harmonic TimeSeries: 5Degree AOA Lower Surface. ...................... 136 Figure G21 1st Harmonic TimeSeries: 10Degree AOA Upper Surface……………...136 Figure G22 1st Harmonic TimeSeries: 10Degree AOA Lower Surface. .................... 136 xiv NOMENCLATURE Symbols ΔCp′ = unsteady differentialpressure coefficient Δp′ = differential unsteady pressure n) = normal vector to liftingsurface chord ω = disturbance angular frequency φ = surfacepressure phase distribution ∞ ρ = freestream mean density m ψ r = mth mode shape vector a = forcingdisturbance transverse velocity c = liftingsurface chord m f = mth mode generalized force (2) i H = ith order Hankel function, 2nd kind k = reduced frequency S = Sears function ∞ U = freestream mean velocity x/c = nondimensional distance along chord i u = mean velocity component i u′ = unsteady velocity component xv k G = turbulent kinetic energy due to mean velocity gradients b G = turbulent kinetic energy due to buoyancy ε = dissipation rate vr = velocity vector (=ui vj) ) ) + in 2D A v = surface are vector φ Γ = diffusion coefficient for φ φ S = source of φ per unit volume faces N = number of faces enclosing cell f f f A r ρ υ = mass flux through the face n (∇φ ) = magnitude of ∇φ normal to face f Pt,i = total pressure inlet Pt = total pressure P′ = unsteady pressure P = instantaneous pressure P = timeaveraged pressure TE = trailing edge LE = leading edge Acronyms CFD = computational fluid dynamics AOA = angle of attack RANS = Reynoldsaveraged NavierStokes equations xvi 1 CHAPTER 1 INTRODUCTION 1.1 OVERVIEW OF AIRFOIL GUST INTERACTIONS Liftingsurface response to unsteady aerodynamic forcing is of particular interest in aircraft propulsion applications, primarily due to timeresolved aerodynamic interactions in turbomachinery. In these applications, lifting surfaces (i.e., airfoils) often operate in both randomly turbulent and periodically oscillating fluid environments, including temporally and spatially nonuniform propagating disturbances. Relative unsteady motion between a lifting surface and the fluid results in complex fluidstructure interactions and may produce aerodynamic/structural response. In order to optimize liftingsurface designs, a detailed understanding of inherent fluidstructure interactions is required, where these interactions can be described, in part, by liftingsurface pressure response. Unsteady surfacepressure phase distributions compromise one aspect of liftingsurface fluidstructure interactions, as phase directly affects timeresolved unsteady force/moment behavior, particularly in terms of forcing structural modes. 1.2 COMPRESSOR DESIGN INTENT AND OPERATING ENVIRONMENT For the purpose of emphasizing the continued need for investigating unsteady fluidstructure interactions, consider an axialflow compressor. Axialflow compressors are typically composed of a number of rotating blades for the purpose of turning and 1 adding work to the passing airflow. This turning process accomplishes a desired total enthalpy rise through the device. A single row of rotating blades is often referred as a blade row. A blade row can lead or follow a separate single row of stationary vanes often referred to as a vane row. Vane rows direct rotor inlet/outlet airflow corresponding to the compressor design. In modern compressors, achieving the design rise requires several blade/vane rows, each providing a portion of the enthalpy difference. Each pair of blade/vane rows represents a stage in the compressor [Falk, 2000]. With the overall goal of engine design being often to decrease engine size and weight, compressor size and weight must also be decreased. A lower number of stages and reduced axial spacing between stages helps to accomplish this goal. However, operating at lower number of stages requires a corresponding increase in stage aerodynamic loading to achieve the desired compressor enthalpy rise. In addition, reduced stagetostage spacing leads to greater aerodynamic interactions between vane/blade rows. Such interactions come in the form of disturbances caused by the relative motion between the rotor and stator rows. 1.3 NORMAL COMPRESSOR OPERATION LEADING TO AIRFOILGUST INTERACTIONS The rotational motions between rotor/stator rows, or stages, in a turbineengine compressor generate a designed enthalpy rise across the component. In the process, however, each stage induces propagating aerodynamic disturbances that act as periodic excitations, or forcing functions, for neighboring blade/vane rows. These propagating disturbances are generally grouped as: 2 • Convective downstreampropagating viscous wakes: produced by the frictional interactions between the fluid and lifting surface. • Convective downstreampropagating vortical wakes: produced by vortex shedding in response to bound circulation fluctuations on the lifting surface. • Acoustically propagating potential disturbances: elicited by variations in the velocity potential, or pressure fields, associated with the blades of a given row [Hall, 1991]. Induced potential disturbances may be temporary in nature, decaying exponentially in the near field, or propagate without attenuation into the far field, depending on bladetip Mach numbers. Interactions between a lifting surface and propagating disturbance field induces timedependent angleofattack changes on the body, causing spatially and temporally dependent surfacepressure distributions. Integration of these surfacepressure distributions forms unsteady forces and moments on the body. Moments and forces generate temporally and spatially dependent mechanical stresses, or alternating stresses. If the induced alternating stresses are strong enough, structural fatigue may plague the liftingsurface with the possibility of catastrophic failure. Fatigue is a process of cumulative structural damage caused by repeated load fluctuations, or stresses [Barsom, 1987]. Fatigue occurs in regions deforming plastically under applied loads; thus, under purely elastic stress conditions, localized areas of raised stress must be present to induce fatigue, where these raised stresses exceed the material yield stress. Prolonged exposure to fatigueinducing unsteady loads may cause initiation and subsequent propagation of a crack, or cracks, in plastically deformed structural regions. Eventually, if crack propagation continues, catastrophic fracture and failure of 3 the lifting surface may occur. Typically the number of load fluctuations necessary to initiate a crack within working liftingsurfaces( such as a compressor blade/vane row) is quite large; thus, the cumulative structural damage process is often referred to as high cycle fatigue, or HCF. Aerodynamically induced load fluctuations in a jet engine compressor composed of highstrength structural components are designed to yield only elastic stress fluctuations. However, the occurrence of random material defects, foreign object damage (FOD), or blade rubbing can provide the proper conditions for crack initiation and propagation. Nonetheless, the number of load fluctuations typically required to form a crack within a jet engine blade/vane component is quite large, leading to HCF. 1.4 IMPORTANCE OF HCF TO ENGINE COMMUNITY Highcycle fatigue is of utmost importance in current jet engine design, were small structural failures can greatly affect entire engine operations to the extent of engine failure. Due to its importance, a considerable number of experimental and computational investigations have been conducted with the overall goal of predicting liftingsurface HCF failures in jet engines. Nonetheless, jet engine HCF failures continue to occur and are largely unanticipated. Recent advances in engine technology may only complicate this problem, as current trends toward higher blade loadings, increased operating temperatures, smaller stagetostage spacing, unconventional geometries, and advanced materials reach beyond traditional design domains, complicating HCFresistant technology [Fleeter, 1992]. In practice, HCFrelated engine failure has been identified as a major contributor to enginesafety mishaps in U. S. military fighter aircraft [Thompson, 1998], where as many as 50% of all engine failures have been attributed to HCF. As 4 such, HCF failure presents a major readiness and monetary concern for both the U.S. Air Force and U.S. Navy [Fecke, 1998]. In an attempt to overcome reoccurring HCF problems in military engines, the U.S. Department of Defense established the National Turbine Engine High Cycle Fatigue Program in 1994. The goal of this ongoing program is to develop, implement, and validate damage tolerant design methodologies to avoid HCFrelated engine failures. This goal is to be accomplished by increasing the level of understanding regarding HCF physics, as well as through the development of better HCF predictive capabilities. The specific goal of improving predictions of engine component response to unsteady aerodynamic forcing belongs to the Science and Technology branch of the HCF program. 1.5 IMPORTANCE OF PHASE DISTRIBUTIONS IN PREDICTING AIRFOIL MODAL FORCING In order to avoid structural vibrations and HCFrelated failures, it is important to accurately predict timeresolved generalized forces for each liftingsurface structural mode. Unsteady surfacepressure phase distributions represent one component of such predictions. Chordwisevarying phase distributions influence timeresolved surfacepressure amplitude distributions along the chord. Chordwise varying phase may also affect surfacepressure node locations, where the node locations change positions with different phase distributions. Mispredictions of chordwisevarying phase may result in under predicted modal forces, generating greaterthanexpected mechanical stresses at multiple spatial and temporal frequencies, and therefore HCF. In all, unsteady surfacepressure phase distributions play a very important role in liftingsurface unsteady forcing and thus predictions of modal forcing. A brief example illustrating the importance of 5 considering unsteady surfacepressure phase distributions in liftingsurface forcing is provided in Chapter 2. 1.6 NEED FOR CONTINUED SURFACEPRESSURE PHASE ANALYSIS Given the potential for catastrophic structural failure caused by liftingsurface fatiguelife degradation, accurate fluidstructure interaction predictions have been aggressively sought for aerodynamic unsteady forcing problems. In fact, extensive amounts of information are available regarding the influence of forced response on turbomachine lifting surfaces (due to the propensity of highcycle fatigue failures in modern highperformance gas turbine engines). Representative investigations have predicted, or measured, unsteady surfacepressure distributions on various lifting bodies, and characterized these distributions in terms of amplitude, frequency, and phase. These investigations typically focus on surfacepressure amplitude and frequency, with little attention to the influence of surfacepressure phase. This is not to say that phase has been completely ignored. In fact, the contrary is true. Researchers have reported surfacepressure phase data over lifting surfaces under a variety of forcing conditions, as will be discussed in Chapter #2. Despite the inclusion of phase results in many research investigations, however the dependence of liftingsurface response to variations in chordwise surfacepressure phase distribution remains relatively unexamined. Moreover, no known investigation has developed general “rules of thumb” to act as guidelines in predicting phase distributions for the most common forcing configurations. Lastly, no consistent explanation exists for observed surfacepressure phase variations between different forcing configurations. The 6 role of surfacepressure phase in the production of structural vibrations and HCF failures therefore remains largely unresolved. 1.7 SCOPE OF CURRENT INVESTIGATION The objective of the current research is to conduct a series of numerical simulations to examine the influence of various aerodynamicforcing and liftingsurface configurations on chordwise surfacepressure phase distributions. In particular, the influence of liftingsurface thickness, camber and angle of attack is discussed. The influences of different solidity values as well as the fundamental physics underlying chordwise surfacepressure phase distributions are also explored. Results will assist the interpretation of experimentally and computationally generated chordwise surfacepressure phase data, as well as provide fundamental results to assist future liftingsurface design efforts in resisting aeroelastic modal forcing. Finally, a comparison between the computed phase results and experimental cascade data reported by Fabian et al. and Falk et al. is presented, providing an explanation of the observed experimental phase trends. 7 2 CHAPTER 2 PREVIOUS WORK This chapter reviews previous work in the area of unsteady liftingsurface aerodynamic forcing. The purpose of this review is to place the current research in proper perspective, establishing its motivation, importance and its potential contribution to the aerodynamic forcedresponse community. 2.1 LITERARY REVIEW Many researchers have reported surfacepressure phase distributions over lifting surfaces under a variety of forcing conditions. Numerous experimental investigations have been performed on isolated airfoils, cascades and jet engine blade/vane rows in an attempt to understand the detailed fluidstructure interactions related to liftingsurface forced response. 2.1.1 SingleAirfoil Investigations One of the early works in the field of liftingsurface forced response was conducted by Sears [1938, 1941], who examined unsteady aerodynamic forcing of rigid infinitespan flat plates by convecting chordwise sinusoidal gusts in an incompressible inviscid flow. Sears derived a relationship for the chordwise, unsteady, nondimensional, differentialpressure amplitude distribution on such lifting surfaces as 8 ( ) i t p S k e x c C π x c ω 1 / 2 1 / + − Δ ′ = (2.1) where ( ) [ ] 2 (2) 1 (2) 0 k H iH S k − = π (2.2) Unsteady surfacepressure time series predicted by Eq. (2.1) are found to be synchronous along the chord, indicating no chordwise time delay between pressureamplitude peaks/troughs. Thus, Eq. 2.1 predicts surfacepressure chordwise phase to be independent of gust propagation speed. The results of Sears suggest instantaneous surfacepressure response along the entire liftingsurface chord to the convective sinusoidal gusts. By extension, Sears predicts corresponding surfacepressure phase distributions versus chord would show zero phase slope; i.e., require the entire surface to respond instantaneously to all forcing disturbances. Although not noted by Sears, his results suggest a chordwisevarying phase distribution should correlate with finite surfacepressure propagation speeds over a lifting surface. Thus, the slope of a surfacepressure phase distribution along the chord relates to forcingdisturbance propagation speed. Faster disturbance propagation corresponds to less phase change, or lower slope, with chord and viceversa. 2.1.2 Turbomachine and Cascade Investigations Turbomachine unsteady forcing phenomena were also experimentally investigated to verify existing analytical results, as well as identify new flow physics. For example, Fleeter et. al. conducted an experimental investigation to determine rotorinduced unsteady pressure distributions on downstream stator vanes [Fleeter et. al., 1978]. This was accomplished in the Detroit Diesel Allison (DDA) largescale, low 9 speed, singlestage research compressor. This investigation studied the effects of reduced frequency and incidence angle on statorvane surfacepressure distribution. Measurements were collected with embedded pressure transducers mounted axially along both suction and pressure surfaces of the stator vanes [Fleeter et al., 1978]. Resultant unsteady surfacepressure amplitudes were shown to compare reasonably well with existing analytical results for all reduced frequency values at small incidence values. However, at large negative incidence angles, experimental data correlations with predictions were very poor. This was attributed to convectedwake phenomena not modeled in the analysis. Corresponding surfacepressure phase results were found to be ambiguous (i.e., no clear chordwise trend) leading to the conclusion that rotor wakes travel differently over airfoil suction and pressure surfaces, depending upon their harmonic frequency. Fleeter et al. also observed similar wakepropagation behavior in a later study [Fleeter et al., 1980]. In this investigation, rotorinduced surfacepressure data acquired on cambered statorvanes were compared to flatplate, vortical gust code results to determine the effect of airfoil camber on airfoil unsteady lift. Unsteady surfacepressure amplitudes on cambered stator vanes exhibited amplification at the leading edge decaying in the chordwise direction. As such, amplitude data correlated very well with theoretical predictions, at both zero and negative incidence angles. However, phase data for the cambered statorvanes again exhibited ambiguous characteristics, showing very poor correlation with the theoretical predictions. In particular, phase data were found to correlate with the theoretical predictions only in the leading edge region, varying linearly from the predicted results downstream. Here again, Fleeter et al. attributed this poor 10 phase behavior to unknown convectedwave phenomenon appearing along the camberedairfoil vane row. The convectedwave phenomenon first appeared near the rear of the vane, moving forward as incidence angle decreased. This phenomenon also exhibited different behavior on the vane pressure and suction surfaces. Unfortunately, the observed convectedwave phenomenon was not modeled or investigated fully, in the presented analysis. Liftingsurface chordwise relativephase information was experimentally examined by Fabian et al in a linear transonic cascade [Fabian et. al., 1996]. The cascade consisted of six productionhardware stator vanes collected from the fan stage in a F109 turbofan engine. Stator vanes were placed in a 4 in. × 4 in. crosssection cascade wind tunnel, creating five twodimensional passages; flow turning through the passages induced vane mean aerodynamic loading. Vane unsteady forcing was accomplished via a row of five circular cylinders placed 0.8 vane chords upstream or downstream of the vane row; allowing forward or rearward aerodynamic forcing, respectively. In the rearwardforcing configuration, upstreampropagating potential disturbances, created by shed bound circulation on the downstream cylinders, forced the vane row. Unsteady, phaselocked, surfacepressure measurements were collected on the vanes at various freestream Mach numbers not exceeding 0.59. Surfacepressure results indicated rearward forcing to elicit nearly linear phase behavior with chord, as illustrated in Figure 21, for both first and second harmonic surfacepressure frequencies. Note the slope of the phase data is positive with respect to chord, indicating an upstreampropagating forcing disturbance, as predicted by the superimposed line showing analytical model results. The linear nature of data also 11 suggest phase to be independent of aerodynamic loading, and have a constant propagation speed [Fabian et al., 1996]. Thus, forcing disturbance and surfacepressure propagation correlate in the rearwardforced cascade configuration. Figure 21 Relative Phase Variation with Chord: RearwardForced Cascade of Fabian et al. [1996]. In addition to rearward forcing, forwardforcing investigations conducted by Fabian et al. also examined chordwise phase distributions [Fabian et al., 2001]. By placing forcing cylinders upstream of the vane row, the cascade configuration allowed convective cylinder wakes to propagate across and force the cascade row. Phaselocked surfacepressure measurements, analogous to those collected during rearward forcing, produced chordwise phase distributions such as Figure 22. Note that unlike the rearwardforcing phase data, the forwardforcing data are not linear, have no clear slope or pattern, and do not agree with the analytical phase model results (solid lines). As such, Fabian et al. termed this data behavior as “ambiguous”, attributing the chaotic nature of the data to cylinderwake interaction with downstreampropagating potential disturbances also emanating from the forcing cylinders. The ambiguous phase results of Figure 22 12 raised the question as to whether such phase ambiguity might also be expected in a rotating machine, as opposed to the linear cascade setup. Figure 22 Relative Phase Variation with Chord: ForwardForced Cascade of Fabian et al. [1996]. In order to answer the question raised by the cascade study [Fabian et al., 1996], phaselocked, unsteady, surfacepressure measurements were performed across swept stator vanes in a running F109 turbofan engine [Falk et al., 1997]. Unlike most turbofan engines, the F109 has only a single stage of axial compression. Therefore, no obstructions exist downstream of the stator vanes that might produce upstreampropagating disturbances. The only disturbances propagating across the stator vanes develop from upstream. This forcing configuration is analogous to the forwardforced cascade of Fabian et al. Results from the investigation by Falk again showed liftingsurface phase information to not display a definite propagation direction at either the convected or acoustic disturbance speed along the vane. In fact, the phase ambiguity measured in the F109 supported the previous arguments of Fabian et al. [1996] and provided presumptive evidence of a strong interaction between downstreampropagating potential and 13 convected disturbances. Of these disturbances, one was argued to be the vortical wakes created by the fan blades at the bladepassing frequency, propagating at the local convection velocity. The other disturbance was argued to be a potential disturbance also created at the bladepassing frequency, but propagating downstream at acoustic speeds. Frey and Fleeter [Frey, 1998], performed experiments to investigate and quantify gustgenerated unsteady aerodynamic response of stator blades. In their experiment, 2/revolution unsteady aerodynamic forcing functions were introduced to a first stage rotorblade row, these forcing disturbances having significant vortical and potential components. Obtained results showed unsteady pressure amplitudes to reach a high value near the leading edge, decaying by 75% at midchord and then increasing slightly in the aft chord; such amplitude behavior was also observed by Fabian et al. Results again suggested a strong interaction between vortical and potential disturbances, comparing well with proprietary codes named LINFLO and LINSUB. Unfortunately, no explanation of chordwise phase distributions was given. 2.2 UNRESOLVED ISSUES The aforementioned investigations provide essential improvements toward understanding liftingsurface surfacepressure distributions under a variety of forcing conditions; however, the unexplained behavior of reported phase data remains an open topic. Upon review, a consistent explanation for observed surfacepressure phase variations along examined lifting surfaces under various aerodynamic forcing configurations does not exist. As such, researchers and designers working in the forcedresponse area are largely uneducated about the role of surfacepressure phase distributions in the production of liftingsurface structural vibrations and HCF failures. 14 Therefore, what follows is a brief analytical example illustrating the possible influence of chordwise phase distributions on timeresolved surfacepressure response and modal forcing. It is intended that this example underscores, or even introduces, the importance of accurately considering phase distributions in liftingsurface response models, while also reinforcing the need for continued examination of liftingsurface phase distributions under various, even simplified, forcing conditions. Consider an infinitespan flatplate lifting surface having a chord extending from –1.0 < x/c < 1.0, where the surface is placed in an incompressible inviscid fluid. If the liftingsurface is subjected to a propagating sinusoidal disturbance, and it is assumed that the ingested disturbance is not distorted by interaction with the flat plate, the unsteady differentialpressure distribution along the lifting surface can be described in terms of a periodic function having some amplitude, frequency and phase. This function is given by p' (x / c,t) a(k) U C ' (x / c)sin[ t (x / c)] p Δ = ρ Δ ω +φ ∞ ∞ (2.3) where ΔCṕ is predicted by Sears in Eq. (2.1). Using Eqs. (2.1) and (2.3), nondimensional unsteady surfacepressure time series at various x/c locations are computed, as shown in Figures 23 and 24, for two separate phase distributions. In the nonvarying phase case (i.e.,φ (x / c) = 0) of Figure 23, each x/c timeseries is found to be in phase, reflecting instantaneous chordwise response to each propagating disturbance. This corresponds to the liftingsurface phase distribution of Sears. Conversely, for the chordwisevarying phase case (i.e.,φ (x / c) ≠ 0) of Figure 24, each x/c location responds sequentially to the propagating disturbances. Note that the phase distributions of Figures 23 and 24 assume a constant phase change, or constant disturbance propagation speed, between each chordwise liftingsurface location. 15 Time, t (sec) Differential SurfacePressure, P' 0 50 100 150 200 250 180.0 150.0 120.0 90.0 60.0 30.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 x/c = 0.80 x/c = 0.60 x/c = 0.40 x/c = 0.20 x/c =0.0 x/c =0.20 x/c =0.40 x/c =0.60 x/c =0.80 Figure 23 Unsteady Differential Surface Pressure Time Series: (φ (x / c) = 0). Time, t (sec) Differential SurfacePressure, P' 0 50 100 150 200 250 180.0 150.0 120.0 90.0 60.0 30.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 x/c = 0.80 x/c = 0.60 x/c = 0.40 x/c = 0.20 x/c =0.0 x/c =0.20 x/c =0.40 x/c =0.60 x/c =0.80 Figure 24 Unsteady Differential Surface Pressure Time Series: (φ (x / c) ≠ 0). Figures 25 and 26 exhibit corresponding surfacepressure distributions (computed from Figures 23 and 24) at selected times of t = 1, 45, 89, and 130 μs, for the two examined phase distributions. For the nonvarying phase case of Figure 25, unsteady surfacepressure alternates continuously from positive to negative pressure due to the assumed sinusoidal nature of the forcing disturbance. At no time during the oscillation cycle, however, does the differential surfacepressure have both negative and positive chordwise components. In contrast, unsteady surface pressures corresponding to the chordwisevarying phase case, as shown in Figure 26, have multiple pressurenode locations, where these node locations change chordwise position with time. Furthermore, the chordwisevarying phase case alters the shape of the surfacepressure distribution, particularly along the forward half of the lifting surface. 16 Chord, x/c Differential SurfacePressure, p' 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 150.0 125.0 100.0 75.0 50.0 25.0 0.0 25.0 50.0 75.0 100.0 125.0 150.0 t = 1 t = 45 t = 90 t = 135 Figure 25 Unsteady Differential Surface Pressure Variation with Chord: (φ (x / c) = 0). Chord, x/c Differential SurfacePressure, p' 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 175.0 150.0 125.0 100.0 75.0 50.0 25.0 0.0 25.0 50.0 75.0 100.0 125.0 150.0 175.0 t = 1 t = 45 t = 90 t = 135 Node Locations Figure 26 Unsteady Differential Surface Pressure Variation with Chord:(φ (x / c) ≠ 0). Chordwise integration of the surfacepressures in Figures 25 and 26 provides a maximum unsteady force, or lift, of 140 and 80.0, for the φ (x / c) = 0 and (φ (x / c) ≠ 0) cases, respectively. Given these lift differences, it may be inferred that chordwisevarying phase may be desirable in terms of reducing unsteady aerodynamic loading. However, the integration process ignores modal forcing of the lifting surface. The unsteady generalized force on a particular structural mode, m, can be computed through the integral of the dotproduct between the examined mode shape and surfacepressure distribution over the liftingsurface chord. Thus, for the current example, the generalized force on a particular structural mode can be written as ( ) [ ( / , ) ].[ ( / )] ( / ) 1 1 f t p' x c t n x c d x c m m ψ r ) ∫ + − = Δ (2.4) A simply supported, infinitespan, twodimensional lifting surface has an infinite number of mode shapes grouped into two families. Rigidbody mode shapes correspond to plunging and pitching oscillations of the lifting surface, while flexiblebody mode 17 shapes correspond to elastic structural deflections. Several flexiblebody mode shapes are illustrated in Figure 27, for the first three modes, along with the first rigidbody mode. Chord, x/c Normalized Maximum Deflection 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 2.0 1.5 1.0 0.5 0.0 0.5 1.0 1.5 2.0 Rigid Mode 1 (Lift) Flex Mode 1 Flex Mode 2 Flex Mode 3 Figure 27 RigidBody and FlexibleBody Mode Shapes for a Simply Supported, InfiniteSpan, Two Dimensional Lifting Surface. Accurate predictions of timeresolved generalized forces on each structural mode are important to avoid structural vibrations and HCF failure. To emphasize this fact, the surfacepressure distributions of Figures 25 and 26 are input into Eq. (2.4) along with the mode shapes of Figure 27, producing a generalized force on each examined structural mode. The maximum generalized forces obtained through this exercise are presented in Figure 28 and plotted versus mode number for the two separate chordwise phase distributions. Note that “mode 0” in Figure 28 corresponds to the first rigidbody structural mode, or liftmode, while the higher modes correspond to the first, second and third flexiblebody modes, respectively. Examining Figure 28, the nonvarying phase distribution shows a canonical decay in force with increasing mode number. In contrast, the chordwisevarying phase distribution exhibits decreased force in the lowerorder 18 modes and significantly amplified force in higherorder modes. Therefore, modal force is found to be a function of chordwise phase distribution. Mode # Maximum Generalized Force, fmax 0.0 1.0 2.0 3.0 4.0 0.0 20.0 40.0 60.0 80.0 100.0 120.0 140.0 160.0 φ = 0 φ = 3x/c Figure 28 Maximum Generalized forces on Structural Modes for Various Phase Distributions. The above unsteadyforcing example emphasizes several facts. First, chordwisevarying phase may alter timeresolved surfacepressure amplitude distributions along a liftingsurface chord. Second, chordwisevarying phase distributions also produce surfacepressure node locations, where these node locations may change position with time. Third, generalized modal forces are altered by chordwise phase distribution. In particular, higherorder modal forces may be amplified by chordwisevarying phase distributions, generating greater mechanical stresses at spatial and temporal frequencies not predicted for the nonvarying phase case. Such possible variations in liftingsurface modal forcing may, if inaccurately predicted, lead to decreased fatigue life for the lifting structure. In all, unsteady surfacepressure phase distributions are clearly important to the unsteady forcing problem. This is particularly evident when one considers that some 19 current unsteadyforcing predictive tools employ the results of Sears as a basis of their predictions. Fortunately, an increasing number of forcedresponse predictive tools are not based on the results of Sears, opting rather for direct numerical simulation of the governing fluid dynamic equations. Validation of these computational tools has proven to be laborious and heavily dependent on the availability of properly posed benchmark data. As such, few comparisons between computed an experimentally determined surfacepressure phase distributions have been made. This is not to say that benchmark phase data are not available, but rather that the significance of the data is not well understood or properly examined. In fact, what is intriguing about the previously reported phase data in this chapter is not their lack of inclusion in the open literature, but the almost complete disregard as to their importance and correlation with computed/measured trends. Much of the available computational/experimental phase data do not correspond to the constant disturbancespeed assumptions made in the above example; in fact, certain data sets show almost no discernable trend with chord. Therefore, while the above example emphasizes the importance of considering chordwisevarying phase for a constant propagation speed, the validity of assuming a constant disturbance propagation speed is unclear. Moreover, the effects of phase deviations from the assumed constantspeed disturbance phase on liftingsurface response are unknown. 2.3 CURRENT RESEARCH OBJECTIVE Given the continued need for examining the influence of liftingsurface chordwise phase distribution on surfacepressure response, the current research presents a fundamental study of surfacepressure phase. In particular, twodimensional, timeaccurate, RANS simulations are performed to examine the fundamental physics leading 20 to surfacepressure phase. Attempts are made to reveal the essential dependencies of that phase on forcing configuration. Simulations are performed for a variety of lifting surface geometries and forcing conditions utilizing the commercially available CFD algorithm, Fluent (v. 6.0). 21 3 CHAPTER 3 COMPUTATIONAL METHODOLOGY AND SETUP This chapter discusses the computational methodology and setup employed in the present research. A description of the airfoil and statorvane cascade geometries, as well as associated boundary conditions, is given. Indepth discussions regarding FLUENT and its companion mesh generation software, GAMBIT, are also provided along with a development of the UDF (userdefined function) generating the airfoil forcing function. 3.1 AIRFOIL GEOMETRY AND BOUNDARY CONDITIONS Surfacepressure phase data on forwardforced lifting surfaces are examined using the commercially available CFD algorithm, FLUENT (v. 6.0). Simulations were performed with four symmetric NACA airfoil profiles of 10, 12, 15 and 20% thickness (relative to chord), two cambered airfoils of 2 and 6% camber (relative to chord), two meanflow attack angles of 5 and 10 degrees and two forcingdisturbance frequencies of 150 Hz and 300 Hz. In all simulations, periodic boundary conditions were enforced on the upper and lower computational boundaries located about an otherwise isolated lifting surface, as illustrated in Figure 31. These periodic boundaries simulate the influence of neighboring surfaces, or a series of airfoils in cascade; an airfoil cascade configuration was selected for comparison with previous experimental configurations. A cascade solidity of 4.0, representing weak surfacetosurface pitchwise aerodynamic coupling as 22 compared to modern cascaded blade rows, was selected for the majority of simulations; however, other solidities equaling 2.0 and 6.0 were also briefly investigated. Pressureinlet and pressureoutlet boundary conditions were set for the computational inlet and outlet boundaries, respectively (see Figure 31). Figure 31 AirfoilCascade Computational Boundaries. 3.2 STATORVANE GEOMETRY AND BOUNDARY CONDITIONS In addition to the simplified NACA airfoilcascade configuration, a more complicated cascade configuration was also examined. This configuration employed aerodynamically loaded vanes that mimicked the twodimensional geometry of the statorvane row in the fan compression stage of a F109 turbofan engine (at 87.8% span). The high cascade solidity and double circulararc profile of the statorvane row required definition of additional geometric variables beyond the simplified NACA airfoil cascade. A complete discussion of the cascade geometry and nomenclature is provided in Appendix A. In the present investigation, a vanecentered computational mesh was 23 selected to model the statorvane cascade geometry, with the simulated vane centered in the computational domain, as shown in Figure 32. Periodic boundary conditions were enforced at mesh boundaries above and below the vane, simulating the influence of other vanes in cascade. The periodic boundaries were set at midpitch between vanes, using the vane camberline arc to define the boundary geometry. The statorvane inlet flow angle was set to be 21.9o, while the exit flow was set to be 20.6o, based on previous experimental data [Fabian, 1995]. The stator vanes possess a maximum camber and thickness of 12% and 8% relative to chord, respectively. Vane profile coordinates are listed in Appendix B, reproduced from [Fabian, 1995]. Table 3.1 provides characteristics of the modeled statorvane cascade geometry. Like the NACAairfoil cascade, pressureinlet and pressureoutlet boundary conditions were set for the statorvane computational inlet and outlet boundaries, respectively. Table 3.1 Modeled StatorVane Cascade Geometry. Parameter Value Vane Spacing, S 0.84 m Solidity, σ 1.524 Inlet Flow Angle, α1 21.9o Exit Flow Angle, α2 20.6o 24 Figure 32 StatorVane Cascade Computational Boundaries. 3.3 GENERAL FLUENT SOLVER DESCRIPTION FLUENT is a stateoftheart commercially available flowsolver with capability for modeling unsteady, compressible, viscous flows via numerical solution of the governing fluid dynamic equations. Numerical simulation of the governing fluid dynamic equations in FLUENT is accomplished via a controlvolume based (finitevolume) discretization technique. This technique integrates governing integral equations established within discrete elements (i.e. finite volumes, or cells) of the mesh, resulting in a system of algebraic equations for dependent variables such as velocity and pressure. The discretized algebraic system is then linearized and solved numerically to yield updated variable values at each iteration/time step, using (in the present investigation) implicit linearization schemes. Solution interpolation between adjacent elementface regions is accomplished via one of several userdefined methods, including firstorder upwind, secondorder upwind and powerlaw interpolation. Firstorder or secondorder accurate spatial/temporal discretization is available in FLUENT, where secondorder accuracy is default for coupled solutions. To aid convergenence in highly nonlinear 25 problems, FLUENT allows userdefined controls over underrelaxation and courant (CFL) numbers. Numerical solutions are achieved through one of three userselected solvers, including: segregated, coupledimplicit, or coupledexplicit solvers. The segregated solver linearizes the governing equations implicitly with respect to the dependent variables, solving the resulting set of equations sequentially. Linearized momentum equations are solved individually for fluid velocity, followed by corrective step in which velocity is adjusted based on userselected velocitypressure correlations to satisfy continuity. Conversely, coupled solvers simultaneously solve the set of governing equation defining the dependent variables, where the equation set can be linearized either explicitly or implicitly. For implicit linearization, GaussSeidel solvers are employed in conjunction with an algebraic multigrid (AMG) method to solve the system(s) of equations. Conversely, with explicit linearization, dependentvariable solutions are updated at each time step using a multistep RungeKutta solver, with the additional option of employing a full approximation storage (FAS) multigrid scheme to accelerate convergence. FLUENT allows users to specify several boundary condition types. Supported inlet and outlet boundary conditions include: pressureinlet, velocityinlet, massflowinlet, inletvent, intakefan, pressureoutlet, pressurefarfield, outflow, outletvent, and exhaustfanboundaries. Similarly, wall, repeating, and pole boundary types include: wall, symmetry, periodic and axis boundaries. For the present research, twodimensional numerical simulations of the Reynoldsaveraged NavierStokes equations were accomplished in FLUENT via a finitevolume technique. A coupled solution methodology (Section 3.3.1) was selected, in which the 26 fully coupled systems of equations defining the dependent variables in each cell were discretized, linearized, and solved simultaneously at each iteration/time step. Linearized equation systems were solved using a GaussSeidel solver in conjunction with an algebraic multigrid (AMG) method. Secondorder accurate spatial and temporal discretization was employed for all simulations, with implicit linearization. Simulations were fully viscous, utilizing a standard kε turbulence model (see Section 3.3.3). 3.3.1 Coupled Solution Method The coupled solver used for the current simulations solves the governing equations of continuity, momentum, and (where appropriate) energy and species transport simultaneously (i.e., coupled together). Governing equations for additional scalars (i.e., turbulence, etc.) are solved sequentially using a segregated approach. Since the set of governing equations is nonlinear (and therefore coupled), several subiterations of the solution procedure are performed at each time step before a converged solution is obtained at that time step. Each subiteration consists of the steps illustrated in Figure 33 and outlined below: 1. Fluid properties are updated, based on the current solution. (If the calculation has just started, fluid properties are updated based on an initial solution.) 2. Continuity, momentum, and (where appropriate) energy and species equations are solved simultaneously. 3. Where appropriate, equations for scalars such as turbulence and radiation are solved using the previously updated values of the other variables. 4. A check for convergence of the equation set is made. 27 These steps are continued until the convergence criteria are met at each time step. [FLUENT, 2001] Update properties Solve continuity, momentum, energy, and species equations simultaneously. Solve turbulence and other scalar equations. Converged? Stop Figure 33 Overview of the Coupled Solution Method [Fluent, 2001]. 3.3.2 ReynoldsAveraged NavierStokes (RANS) Equations The Reynoldsaveraged NavierStokes (RANS) equations were selected to represent transport equations for ensembleaveraged, or mean, flow quantities, with all turbulence scales modeled. The approach of permitting a solution for just the meanflow variables greatly reduces the computational effort. If the mean flow is steady, the governing equations will not contain time derivatives and a steadystate solution can be obtained economically. A computational advantage is also provided in required timeaccurate simulations, as time step may be determined by global unsteadiness in the mean flow rather than turbulent unsteadiness. The RANS approach models turbulent flow 28 quantities, through such wellknown models as the SpalartAllmaras, kω, kε, and RSM models [FLUENT, 2001]. 3.3.2.1 Reynolds (Ensemble) Averaging In Reynolds averaging (i.e., ensemble averaging solution variables in the instantaneous (exact) NavierStokes equations are decomposed into mean (ensembleaveraged or timeaveraged) and fluctuating components (about the mean). For velocity components, this decomposition equals i i i u = u + u′ (3.1) Likewise, for pressure and other scalar quantities: φ =φ +φ ′ (3.2) where φ denotes a scalar such as pressure energy or species concentration. Substituting these forms of the flow variables into the instantaneous continuity and momentum equations, and taking a time (or ensemble) average (and dropping the overbar on the mean velocity, u ), yields the ensembleaveraged continuity and momentum equations. These equations can be written in Cartesian tensor form as: ( ) = 0 ∂ ∂ + ∂ ∂ i i u t x ρ ρ (3.3) ) ( ) 3 ( ) ( ) ( 2 1 1 i j j ij i j j i i j i j j i u u x x u x u x u x x u u p x u t − ′ ′ ∂ ∂ + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ − ∂ ∂ + ∂ ∂ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ ρ ρ μ δ ρ (3.4) Equations (3.3) and (3.4) are also known collectively as the RANS equations. 29 3.3.3 k – ε Turbulence Model In the present application, a standard kε turbulence model was used to simulate the effects of turbulence. The standard kε model is a semiempirical model based on model transport equations for turbulence kinetic energy, k, and it dissipation rate. The model transport equation for k was derived from an exact equation, while the model transport equation for ε was obtained using physical reasoning. Equations for k and ε are given as follows k b M k k j t j i i G G Y S x k x ku x k t + + − − + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ ρε σ μ (ρ ) (ρ ) (μ ) (3.5) and ε ε ε ε ε ε ρ ε ε σ μ ρε ρε μ S k G C G C k C x x u t x k b j t j i i + + + − + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ 2 1 3 2 ( ) ( ) ( ) ( ) (3.6) where, 1ε C , , , 2ε C 3ε C k σ and ε σ are model constants nominally having the following values: 1ε C = 1.44 2ε C = 1.92 3ε C = 0.09 k σ = 1.0 ε σ = 1.3 The above values are default in FLUENT, and are employed in the current investigation. 30 3.3.4 FiniteVolume Discretization Methodology FLUENT uses a controlvolumebased technique to convert governing fluid dynamic equations to algebraic equations that can be solved numerically. This controlvolume technique consists of integrating the governing equations about each control volume (or cell) in the computational mesh, yielding discrete equations that conserve each quantity on a controlvolume basis. Discretization of the governing equations can be illustrated most easily by considering the steadystate conservation equation for transport of a scalar quantity,φ . This is demonstrated by the following equation written in integral form for an arbitrary control volume, V, as follows: ∫ =∫Γ ∇ + ∫ V dA dA S dV φ φ φ ρφυ r r r (3.7) The above equation is applied to each control volume, or cell, in the computational domain. Thus, discretization of Eq. (3.7) for a given control volume yields: A A S V faces N faces f n f N f f f f f φ φ Σ ρ υ φ =Σ Γ ∇φ + r r r ( ) (3.8) All governing equations solved by FLUENT take the same general form as Eq. (3.8) and therefore readily apply to unstructured meshes composed of tetrahedra, as in the current investigation. FLUENT stores discrete values of the scalar φ at the cell centers. However, face values, of f φ , are required for convection terms in Eq. 3.8 and therefore must be interpolated from adjacent cellcenter values. This is accomplished using an upwind scheme. 31 3.3.5 SecondOrder Upwind Scheme For the present investigation, a secondorder upwind scheme was selected. When secondorder solution accuracy is desired in FLUENT, flow quantities at cell faces are computed using a multidimensional linear reconstruction approach. In this approach, higherorder accuracy is achieved at cell faces through a Taylor series expansion of the cellcentered solution about the cell centroid. Thus when secondorder upwinding is selected, a face value f φ is computed using the following expression: s f φ =φ + ∇φΔr (3.9) where φ and ∇φ are the cellcentered value and its gradient, and sr Δ is the displacement vector from the cell centroid to the face centroid. This formulation requires the determination of the gradient ∇φ in each cell, where this gradient is computed using the divergence theorem, ∇ = Σ N faces f f A V r φ 1 φ~ (3.10) In Eq. (3.9) face values f φ ~ are computed by averaging φ between cells adjacent to the face in question. Finally, the gradient ∇φ is valuelimited so that no new maxima or minima are introduced in the examined cell region. 3.3.6 SecondOrder Time Discretization For transient simulations, the governing equations must be discretized in both space and time. Spatial discretization for timedependent equations is identical to the steadystate case; however, temporal discretization involves integration of every term in 32 the governing differential equations over one time step,Δt . The integration of transient terms is straightforward, as described below. A generic expression for the time evolution of a variableφ is given by: (φ ) φ F dt d = (3.11) where the function F incorporates both spatial and temporal discretization. If the time derivative is discretized to firstorder approximation, an expression for the discretized derivative may be written as ( ) 1 φ φ φ F t n n = Δ + − (3.12) While a secondorder accurate temporal discretization is given by ( ) 2 3 1 4 1 φ φ φ φ F t n n n = Δ + − + − (3.13) Once the time derivative has been discretized, a choice remains for evaluating F(φ ) . The method employed here to evaluate F(φ ) at a future time level, such as ( 1 ) 1 + + = Δ − n n n F t φ φ φ (3.14) This is referred to as “implicit” integration as in a given cell depends on both sides of Eq. (3.14), giving φ n+1 φ n+1 =φ n + ΔtF(φ n+1 ) (3.15) The implicit relation of Eq. (3.15) can be solved iteratively by initializing to and iterating the equation φ n+1 φ t ( ) 3 2 3 1 3 φ t = 4φ n − φ n−1 + ΔtF φ t (3.16) 33 for secondorder formulation, until stops changing (i.e., converges). At that point, is set equal to . The advantage of a fully implicit scheme is its unconditional stability with respect to timestep size. φ t φ n+1 φ t 3.3.7 Linearization Methodology In the coupledsolution method the discrete, nonlinear governing equations are linearized to produce a system of equations for the dependent variables in every computational cell. The resultant linear system is then solved to yield an updated flowfield solution. The manner in which the governing equations are linearized may take an implicit or explicit form with respect to the dependent variable (or set of variables) of interest. For the present analysis a coupledsolution methodology with implicit linearization was used. This results in a system of linear equations with N equations for each cell in the domain, where N is the number of coupled equations in the set. Because there are N equations per cell, this is sometimes called a block system of equations. A pointimplicit (i.e., block GaussSeidel) linear equation solver is used in conjunction with an algebraic multigrid (AMG) method to solve the resultant block system of equations for all N dependent variables in each cell. For example, linearization of the coupled continuity, x, y, zmomentum, and energy equation set will produce a system of equations in which p, u, v, w, and T are unknowns. Simultaneous solution of this equation system (using the block AMG solver) yields at once updated pressure, velocity, and temperature fields. 34 3.3.8 Periodic Boundary Conditions Periodic boundary conditions are advantageous when the physical geometry of interest and the expected flow solution has a spatially periodic nature. The assumption of spatial flow periodicity implies the velocity components repeat themselves in space as follows: u(r ) = u(r + L) = u(r + 2L) = ........ r r r r r v(r ) = v(r + L) = v(r + 2L) = ........ r r r r r (3.17) where rr is the position vector and L r is the periodic length vector of the domain considered. For viscous flows, the spatially distributed pressure field may not be periodic in the sense of the above Eq. (3.18). Instead, pressure drop between modules maybe periodic in space, giving Δp = p(r ) − p(r + L) = p(r + L) − p(r + 2L) = ........ r r r r r r r (3.18) If the coupled solver is selected in FLUENT, Δp maybe specified as a constant value at any periodic boundary. 3.3.9 Standard FLUENT Inlet/Outlet/Wall Boundary Conditions 3.3.9.1 PressureInlet Boundary Conditions Pressureinlet boundary conditions define fluid pressure at flow inlets, along with all other scalar properties of the flow. They are suitable for both incompressible and compressible flow calculations [Fluent, 2001]. Pressureinlet boundary conditions are typically specified in FLUENT when inlet pressure is known but flow rate and/or velocity is not known. When flow enters through a pressureinlet boundary, FLUENT enforces the input boundary pressure as the fluid totalpressure at the inlet plane, . In 0 p 35 steadystate incompressible flow, the inlet totalpressure and static pressure are related at the inlet velocity via Bernoulli's equation: s p 2 0 2 = + 1 ρν s p p (3.19) In compressible flow, the totalpressure and static pressure are related by 2 1 0 ) 2 (1 1 − − = + γ γ γ p p M s (3.20) 3.3.9.2 PressureOutlet Boundary Condition: As specified earlier, a pressureoutlet boundary condition was selected at the computational outlet plane for all simulations. Pressureoutlet boundary conditions in FLUENT require the specification of a static (gauge) pressure at the outlet boundary. For a subsonic flow, FLUENT enforces the input boundary pressure as the fluid static pressure at the outlet plane, and extrapolates all other flow conditions to the outlet from the interior of the domain. 3.3.10 Operating Pressure In the present investigation, the operating pressure condition in FLUENT was set to 0 Pa, allowing all pressure calculations to be treated as gauge pressures. Setting the computational operating pressure equal to 0 Pa is significant as it reduces roundoff error problems. In lowMach number flows, overall pressure changes are typically small compared to the atmospheric static pressure, and therefore can be significantly affected by numerical roundoff error. Moreover, FLUENT always uses gauge pressure for all its calculations, as such setting the operating pressure equal to zero makes gauge and absolute pressures equivalent. 36 3.4 UDF DESCRIPTION Aerodynamic forcing of the examined lifting surfaces was achieved by modeling waveform characteristics similar to a passing wake, as might be found in a rotorstator compression stage. In doing so, a userdefined function (UDF) was written in FLUENT. A UDF can be dynamically loaded into the FLUENT solver to enhance standard features of the available code. Such functions are written in the C programming language, allowing standard C library functions to be used, as well as predefined macros (provided by Fluent Inc.). A UDF can be implemented as an interpreted or compiled function. An interpreted UDF is read in and interpreted at run time. Alternatively, a compiled UDF is grouped into shared libraries when it is compiled and linked with the standard FLUENT executable. An interpreted UDF is simple to use but imparts coding and speed limitations. A compiled UDF executes faster and has no coding limitations, but requires more effort to develop and implement. The current investigation employs a compiled UDF. 3.4.1 Development Logic/Procedure The developed inlet UDF established an inlet totalpressure profile, simulating the movement of a rotorblade wake passing across the computational inlet. The following steps were employed to code and implement the UDF in FLUENT for the present investigation. 1. Define the wake model. 2. Create a C source code file 3. Start FLUENT and setup the case file. 37 4. Compile the C source code 5. Activate the UDF in FLUENT. The first step to generate and implement the inlet UDF was to develop an analytical model of a typical rotorblade wake; this wake representing the aerodynamic forcing function for the examined lifting surfaces. Experimental lowspeed compressor investigations by Reynolds et. al (1979) showed that rotorblade wake profiles approximately follow a Gaussian profile. As such, a simple aerodynamic forcingfunction model was developed using a Gauss function. This model is characterized by the wake centerline defect, , and the semiwake width, dc W δ (see Figure 34). An equation for totalpressure deficit produced by a rotorblade wake function was defined as [ ( t )2 ] c y ti dc P W e − − = δ 3.21 Pitch, y/c Total Pressure, Pt (Pa) 4.0 2.0 0.0 2.0 4.0 101265 101270 101275 101280 101285 101290 101295 101300 101305 101310 101315 101320 101325 101330 Wake Centerline Wdc δ Figure 34 Rotor Wake Characteristics. 38 Values for the wake centerline defect, , and the semiwake width, dc W δ , were optimized to obtain the waveform shown in Figure 34. Having developed the aerodynamic forcingfunction model, a C source code was written. This source code incorporated standard C library functions as well as predefined macros (provided by FLUENT). These include the DEFINE_PROFILE macro, the F_CENTROID macro, begin_f_loop macro and F_PROFILE macro. The UDF (source code) written was compiled in FLUENT as an interpreted UDF and applied to the pressureinlet boundary. A complete description of the developed UDF with a detailed discussion of each macro is provided in Appendix C. 3.5 COMPUTATIONAL GRID DESCRIPTION The FLUENT software package contains an available preprocessor for geometry modeling and mesh generation known as GAMBIT. As such, the computational grids for the present investigation were created with GAMBIT grid generation software [FLUENT, 2001]. 3.5.1 Gambit Grid Generation Software Mesh generation for FLUENT simulations may be provided through its companion software, GAMBIT; however, select thirdparty software is also supported. GAMBIT is a graphically oriented software package mimicking many computer aided drawing (CAD) programs. GAMBIT allows users to graphically reproduce complex physical geometries of interest, and mesh these geometries using several userselected meshgeneration algorithms. GAMBIT geometries are created by manipulating components such as edges, faces, and volumes that represent the physical system being 39 considered. Vertex coordinates for these components can be directly generated in GAMBIT or imported from existing data files. Created components can be rotated, translated, copied, split, united and subtracted from one another, to name a few operations, as is typical of most CAD software. GAMBIT geometries can also be created in different coordinate systems, or multiple local coordinate systems, depending on the particular application. The software provides complete meshing flexibility in both two and three dimensions, featuring algorithms for structured quadrilateral and hexahedral meshes, unstructured triangular, tetrahedral, pyramid and wedge meshes, and mixed structuredunstructured meshes. GAMBIT also supports boundarylayer subroutines to construct structured meshes in close proximity to geometric walls, providing increased grid resolution in viscous boundary layers. Socalled “boundarylayer meshes” can be coupled with unstructured meshes in the adjacent potential flow, forming a hybrid mesh. Mesh density and distribution are easily controlled via direct input of nodal positions, or through appropriate sizing functions that concentrate mesh points around userdefined regions. User input during mesh construction allows highdensity elements to be concentrated in areas of highflow gradients. 3.5.2 Grid Methodology 3.5.2.1 AirfoilCascade Mesh The methodology followed to create meshes for the parametric NACA airfoilcascade analysis was as follows. Airfoil geometry (obtained from XFOIL) was read into GAMBIT as vertex data, and interpolation was performed between vertex data to create the basic airfoil shape. Once airfoil shape was defined, inlet and outlet boundary planes 40 were created. The inlet boundary was located one airfoil chord forward of the airfoil leading edge, while the outlet boundary was placed 16 chords downstream of the airfoil trailing edge. In order to reduce computational expense, only one airfoil was considered with periodic boundary conditions simulating the influence of neighboring cascaded airfoils. A twodimensional hybrid mesh, consisting of a structured Hgrid immediately surrounding the airfoil, and an unstructured triangular mesh representing the potential region outside of the airfoil boundary layer, was created. The resulting grid system, for the NACA 0012 airfoil profile, is represented in Figure 35. A detailed description of the boundarylayer mesh and its characteristics is provided in Section 3.5.2.1.1. Grid sizes employed for the airfoilcascade meshes are provided in Section 3.5.2.1.2. 3.5.2.1.1 Grid Density and Grid Distribution Grid node distributions employed for the airfoilcascade parametric analysis are outlined in Table 3.2. The total number of mesh elements for the airfoilcascade mesh equaled 26446, where 24246 were triangular mesh elements and 2200 were quadrahedral mesh elements. The TriPave meshing scheme provided by GAMBIT was used to create the triangular mesh elements. Table 3.2 Grid Distribution (AirfoilCascade Mesh). Surface # of Nodes Inlet 100 Outlet 25 Periodic Boundary 150 Airfoil Suction Side 55 Airfoil Pressure Side 55 41 Figure 35 AirfoilCascade Mesh (NACA 0012 airfoil) 3.5.2.1.2 Boundary Layer (ViscousFlow Region) Grid As mentioned in Section 3.5.2.1, in the immediate region surrounding the airfoil a structured boundarylayer grid was generated. Boundarylayer grids in FLUENT define mesh node spacing in regions immediately adjacent to an edge and/or face. These grids control mesh density and location in a specific region of interest, here the upper (suction) and lower (pressure) surfaces of the airfoil. In the case of the NACA airfoilcascade simulations, the first boundarylayer nodes were positioned 0.12% (relative to chord) from the airfoil surface with a growth rate of 1.1, with the total number of boundarylayer mesh levels equaling 20. Table 3.3 outlines the boundarylayer mesh characteristics employed for the NACA airfoilcascade simulations. Figure 36 depicts the structured boundarylayer mesh surrounding the NACA 0012 airfoil. 42 Table 3.3 BoundaryLayer Mesh Characteristics (AirfoilCascade Mesh). Parameter Value Algorithm Uniform First Row Height 0.12% (relative to chord) Growth Factor 1.1 No. of Rows 20 Total Depth 0.06873 Internal Continuity On Wedge Corner Shape On Transition Pattern 1:1 x(m) y(m) 0 0.2 0.4 0.6 0.8 1 0.4 0.2 0 0.2 0.4 Figure 36 Structured BoundaryLayer Mesh Surrounding NACA 0012 Airfoil Geometry. 43 3.5.2.2 StatorVane Cascade Mesh The steps followed to create the statorvane cascade mesh are similar to those followed for the airfoilcascade mesh. First, statorvane coordinates (reproduced from Fabian [1995]) were read into GAMBIT as vertex data and interpolation was performed between vertices to create the statorvane shape. Second, inlet and outlet boundary planes were created. Similar to the airfoilcascade mesh, the inlet boundary was located one statorvane chord forward of the vane leading edge, while the outlet boundary was placed 16 chords downstream of the vane trailing edge. The inlet and outlet boundaries were set exactly perpendicular to the xaxis with an inlet flow angle of 21.9o and an outlet flow angle of 20.6o. In order to reduce computational expense, only one statorvane was considered with periodic boundary conditions simulating the influence of neighboring cascaded statorvanes. A structured Hgrid mesh was created for the entire statorvane computational domain. 3.5.2.2.1 Grid Density and Grid Distribution Grid node distributions employed for the statorvane cascade mesh are outlined in Table 3.4. A successive ratio of 0.96 was used for periodic boundary edgenode grading, while biexponent type grading, with a ratio of 0.75, was selected for the statorvane upper and lower surfaces. The total number of mesh elements for the statorvane cascade mesh equaled 28000, all of which were quadrahedral mesh elements. The Map meshing scheme provided by GAMBIT was used to create the quadrahedral mesh elements. Figures 37 and 38 display the statorvane cascade geometry and the modeled statorvane mesh, respectively. 44 Table 3.4 StatorVane Cascade Mesh Node Distribution. Surface # of Nodes Inlet 80 Outlet 80 Periodic Boundary 200 Airfoil Suction Side 75 Airfoil Pressure Side 75 x (m) y (m) 0 2 4 1 0 1 2 3 4 Figure 37 StatorVane Cascade. 45 x (m) y (m) 0 1 2 3 1.5 1 0.5 0 0.5 1 1.5 Figure 38 Modeled StatorVane Cascade Mesh. 3.6 Computational Setup Table 3.5 summarizes the common inputs for all FLUENT simulations. All unsteady simulations for the airfoil parametric analysis were performed using 100 time steps per disturbance forcing period, corresponding to a time step of 5x103 s. Local convergence was achieved at each time step, with as many as 35 subiterations per time step. Global convergence for all unsteady simulations was based on lift periodicity; convergence was achieved when lift periodicity reached less than 0.1% variability between forcing periods. Figure 39 illustrates liftconvergence history for the NACA 0012 airfoil; this history is representative of all simulations conducted in the present research. 46 Table 3.5 FLUENT Configuration for Numerical Simulations. Parameter Value Fluent Solver CoupledImplicit Discretization Scheme SecondOrder (Momentum and Viscosity) Material Properties Standard Day Air as an IdealGas Operating Pressure 0 Pa Viscosity Model Standard kε Turbulence Model Inlet Boundary Condition Pressure Inlet Outlet Boundary Condition Pressure Outlet Periodic Boundary Condition Periodic Time, t (sec) Coefficient of Lift, Cl 0.0 1.0 2.0 3.0 4.0 5.0 6.0 0.04 0.02 0 0.02 0.04 0.06 0.08 0.1 0.12 Figure 39 Lift Coefficient Time History Showing Convergence: NACA 0012 Profile. 47 3.6.1 Reference Values During calculations, FLUENT employs nondimensionalized variables. All variables are nondimensionalized using the reference values listed in Table 3.6. Table 3.6 FLUENT Reference Values Parameter Value Reference Density 1.225 kg/m3 Reference Length 1 m Reference Pressure 101325 Pa Reference Temperature 300 K Reference Velocity 10.0001 m/s Reference Viscosity 1.7894e5 kg/ms Ratio of Specific Heats 1.40 Operating Pressure 0 Pa 3.6.2 Boundary Conditions Pressureinlet and pressureoutlet boundary conditions were employed for the computational inlet and outlet, respectively, with inlet pressure and temperature boundary conditions being set to standardday atmospheric values for the steadystate simulation. However, for the unsteady simulations a UDF was employed at the pressureinlet boundary, mimicking the convecting wake from a rotating rotor blade. Standard wall boundary conditions were used for airfoil/statorvane suction and pressure surfaces. A detailed description of the UDF is presented in Section 3.4. Timeaveraged freestream velocity was set to 10 m/s for all NACA airfoilcascade simulations, while timeaveraged freestream velocity for the statorvane cascade simulations was set to 200 m/s. Table 3.7 48 summarizes the inlet and outlet boundary conditions employed for all FLUENT simulations performed in this research. Table 3.7 FLUENT Boundary Conditions Parameter Value Inlet Boundary Conditions Total Temperature 300 K Direction Specification Normal to Boundary Turbulence Intensity (%) 1 Turbulence Viscosity Ratio 1 Outlet Boundary Conditions Backflow Total Temperature 300 K Turbulence Intensity (%) 1 Turbulence Viscosity Ratio 1 3.7 Grid Independence Gird independence was examined by obtaining NACA 0012 steadystate solutions for higher density (fine) mesh, as compared to the nominal mesh of Figure 35. The fine mesh contained 42018 elements, where 37618 were triangular mesh elements and 4400 were quadrahedral mesh elements. Pressure coefficient results from the NACA 0012 airfoil upper and lower surfaces were compared for each mesh. Figure 310 shows steadystate pressure coefficient data along the NACA 0012 airfoil upper and lower surfaces for the nominal and fine meshes. As seen in the figure, the obtained pressure data compare favorably well, showing minimal discrepancies. As such, the nominal grid is deemed to be of sufficient density for the present investigation. 49 Chord, x/c Coefficient of Pressure, Cp 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.4 1.2 1 0.8 0.6 0.4 0.2 0 Nominal Mesh Fine Mesh Figure 310 SteadyState Pressure Coefficient Data: NACA 0012. Similar grid independence studies were performed for the statorvane cascade analysis. In this case, the high density (fine) mesh contained 38000 quadrahedral elements as compared to the 28000 quadrahedral elements contained by the nominal mesh shown in Figure 38. Figure 311 shows steadystate pressure coefficient data along the statorvane upper and lower surfaces for the nominal and fine meshes. Here again, the pressure data obtained compare favorably showing negligible discrepancies. In particular, along the statorvane upper surface small differences in pressure coefficient values are seen near the midchord region. However, near the leading and trailing edges the pressure data for the two meshes match perfectly well. Given the agreement of data and the additional CPU time required per computation by the fine mesh as compared to the nominal mesh, the nominal mesh is deemed to be of sufficient density for the present investigation. 50 Chord, x/c Coefficient of Pressure, Cp 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 2.0 1.8 1.6 1.4 1.2 1.0 0.8 0.6 0.4 0.2 0.0 Fine Mesh  Upper Nominal Mesh  Upper Fine Mesh  Lower Nominal Mesh  Lower Figure 311 Pressure Coefficient Data: StatorVane. 51 4 Chapter 4 SOLUTION VALIDATION This chapter describes the timeaveraged staticpressure results for the baseline NACA 0012 airfoil case. The timeaveraged flowfield represents an important aspect of the overall unsteadyforcing simulations for two reasons. First, unsteady results develop by subtracting timeaveraged parameters from the instantaneous solution at each time step; thus, giving parameter perturbations about the timeaveraged field. Second, the timeaveraged flowfield facilitates examination to determine the accuracy and applicability of the FLUENT simulation results. Note that the timeaveraged computational results described herein are primarily compared with experimental data from Abbott and Von Doenhoff [1959]. 4.1 TIMEAVERAGED METHODOLOGY Timeaveraged parameter distributions for each simulation are computed by summing flow parameters (i.e., pressure, velocity, etc.) at each lifting surface grid location over 50 time steps (one aerodynamic forcing period), respectively. Resulting summations are then divided by 50, giving a timeaveraged value for each solution parameter at each grid point. 52 4.2 BASELINE AIRFOIL TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS A timeaveraged totalpressure contour plot for the NACA 0012 airfoil is presented in Figure 41, with attached numerical values indicating respective pressure contour levels. 101102192075 101290 101320 101305 101305 101315 101320 101320 101320 Chord, x/c Pitch, y/c 0.4 0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 Figure 41 TimeAveraged TotalPressure Contours (Pa): NACA 0012 Lifting Surface As expected for this liftingsurface profile, the contours illustrate good qualitative characteristics of the simulated flowfield, showing a wellbehaved symmetric total pressure distribution about the profile and in the liftingsurface wake. The flowfield on the airfoil surface shown in Figure 41 shows no indication of largescale separation, exhibiting smooth attached flow. Thus, the timeaveraged flow about the baseline liftingsurface is argued to be attached and producing a similar wake to that which would be expected. 53 Figure 42 exhibits the totalpressure wake profile of the NACA 0012 lifting surface. Pitchwise total pressure wake profiles at x/c = 0.005 and 0.01 downstream of the airfoil are exhibited. As discussed in Section 3.4 wake depth is defined as the maximum total pressure deficit, while wake width represents the maximum pitchwise effect of the total pressure deficit. Symmetric wake profiles are observed at both x/c locations downstream of the airfoil. Thus, confirming the inferences made earlier from Figure 41. Pitch, y/c Total Pressure, Pt (Pa) 0.2 0.1 0.0 0.1 0.2 101280 101284 101288 101292 101296 101300 101304 101308 101312 101316 101320 101324 x/c = 0.005 aft x/c = 0.01 aft Figure 42 NACA 0012 Lifting Surface Wake Profile. Figure 43 displays timeaveraged static pressures computed along the NACA 0012 lifting surface chord. Pressure distributions along the suction and pressure surfaces are individually displayed. As expected, the symmetric NACA 0012 profile produces equivalent timeaveraged staticpressure distributions on the suction (upper) and pressure (lower) surfaces. Thus, no timeaveraged aerodynamic loading exists on the airfoil. Staticpressure distributions predicted by FLUENT along the NACA 0012 profile are also compared with available experimental data [Abbott and Von Doenhoff, 1959], as 54 shown in Figure 43. The computational data compared favorably with the experimental data, except near the leading and trailing edge regions where slight discrepancies were found. In particular, the minimum pressure near the leading edge was underestimated, while a slightly lower pressure was predicted in the trailingedge region. Given the reasonable agreement of the timeaveraged data, the simulations presented herein correctly predict the timeaveraged flowfield about the examined lifting surface. Chord, x/c Time Averaged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101235 101250 101265 101280 101295 101310 101325 NACA0012 Suction NACA0012 Pressure NACA0012 Suction (Abbot & Von Doenhof) NACA0012 Pressure (Abbott & Von Doenhoff) Figure 43 TimeAveraged Static Pressure : NACA 0012. Figure 44 illustrates timeaveraged staticpressure distributions computed for flow over the NACA 0012 airfoil at 0°, 5° and 10° angles of attack relative to the freestream. The nonequivalent surfacepressure distributions on the suction and pressure surfaces indicate aerodynamic loading caused by the nonzero angle of attack of the flow. The static pressure at the suctionsurface leading edge decreases as angle of attack increases, while the static pressure at the pressuresurface leading edge increases. As 55 angle of attack increases, the difference in the pressure between the suction and pressure surfaces increases, producing higher aerodynamic loading. Chord, x/c Time Averaged Static Pressure 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101100 101125 101150 101175 101200 101225 101250 101275 101300 101325 5AOA Suction 10AOA Suction NACA0012 Suction 5AOA Pressure 10AOA Pressure NACA0012 Pressure Figure 44 TimeAveraged Static Pressure : 0°, 5° and 10° Angle of Attack Table 4.1 compares the computed lift coefficient values for the NACA 0012 airfoil at 0°, 5° and 10° angles of attack (relative to the freestream) with available experimental values from Abbott and Von Doenhoff [1959]. Note that the experimental values were measured at Re = 6x106, while Reynolds number for the current simulations was calculated to be Re = 5.88x106. As such, the discrepancy in the lift coefficient value is attributed to the lower Re used for the current investigation. Table 4.1 Coefficient of Lift: 5° and 10° Angle of Attack. Angle of Attack Computed Value Experimental Value 5° 0.54 0.59 10° 0.82 0.87 56 4.3 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (THICKNESS INFLUENCE) As described later in this thesis, the influence of liftingsurface thickness on timeaccurate surfacepressure phase is examined via simulation of three NACA profiles of 10, 12, 15 and 20% thickness (based on chord). Figure 45 illustrates the corresponding timeaveraged staticpressure distributions from these simulations, collected over one aerodynamic forcing period. The data in Figure 45 show good agreement between suction and pressure surfaces for each thickness, as anticipated for symmetric profiles. The equivalent surfacepressure distributions on both suction and pressure surfaces for each thickness case are indicative of no timeaveraged aerodynamic loading. Each profile has distinct timeaveraged pressure gradients along the chord, with higher thickness values resulting in more severe chordwise gradients. Chord, x/c Time Average Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101220 101240 101260 101280 101300 101320 NACA0010 Suction NACA0012 Suction NACA0015 Suction NACA0020 Suction NACA0010 Pressure NACA0012 Pressure NACA0015 Pressure NACA0020 Pressure Figure 45 TimeAveraged StaticPressure for Various LiftingSurface Thicknesses. 57 The staticpressure distributions predicted by FLUENT for the 10 and 15% thickness profiles were compared with available experimental data from Abbott and Von Doenhoff [1959] as shown in Figure 46. Here again, as in the NACA 0012 case, the computational data compared favorably with the experimental data, except in the leading and trailing edge regions where slight discrepancies are observed. In particular, the minimum staticpressure near the leading edge is underestimated and a slightly lower staticpressure is predicted in the trailing edge region. Chord, x/c TimeAveraged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101235 101250 101265 101280 101295 101310 101325 NACA0010 (Computational) NACA0010 (Abbott & Von Doenhoff) NACA0015 (Computational) NACA0015 (Abbott & Von Doenhoff) Figure 46 Time averaged Static Pressure Comparison with Experimental data. 4.4 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (CAMBER INFLUENCE) In addition to the liftingsurface cases, three cambered liftingsurface geometries were also considered, including 0, 2, and 6% camber (relative to chord). Figure 47 shows timeaveraged staticpressure distributions along the suction and pressure surfaces 58 of each cambered profile. Note from Figure 47 that while each profile exhibits a unique chordwise surfacepressure gradient, as for various profile thicknesses, the cambered profiles also exhibit timeaveraged aerodynamic loading. This is exhibited by the nonequivalent surface pressure distributions on the suction and pressure surfaces in Figure 47. Differential pressure across the liftingsurface increases as the percentage of camber relative to chord increases. Chord, x/c Time Averaged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101220. 101240. 101260. 101280. 101300. 101320. NACA0012 Suction 2%Camber Suction 6% Camber suction NACA0012 Pressure 2%Camber Pressure 6% Camber Pressure Figure 47 TimeAveraged StaticPressure Distribution for Various LiftingSurface Cambers Table 4.2 compares the computed lift coefficient values for 2, and 6% cambered (relative to chord) airfoils with available experimental values from Abbott and Von Doenhoff [1959]. Here again, as was seen earlier for the AOA case, the computed lift coefficient values were slightly lower than the experimental values. 59 Table 4.2 Coefficient of Lift: Various LiftingSurface Cambers. Camber) Computed Value Experimental Value 2 % 0.21 0.25 6 % 0.74 0.78 4.5 COMPARISON OF TIMEACCURATE BASELINE LIFT DEPENDENCY TO SEAR’S RESULTS Figure 48 presents the firstharmonic (i.e., fundamental frequency) lift time series obtained from FLUENT for the baseline NACA 0012 airfoil, as well as known lift prediction for airfoils experiencing a sinusoidal transverse gust, as described by Sears [1938]. As can be seen in Figure 47, the computed lift compares perfectly well with the Sears predicted lift series. Note that, Sears unsteady aerodynamic forcing was modeled as a perfect sinewave transverse gust as opposed to the rotorwake model (transverse and chordwise varying velocity components) used for unsteady aerodynamic forcing in the present analysis. Given this agreement of the lift time series data, the unsteady results obtained from FLUENT are further validated. 60 Time (sec) Lift 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.4 2 1.6 1.2 0.8 0.4 0 0.4 0.8 1.2 1.6 2 2.4 Sears Computational Figure 48 LiftTime series Comparison with Sears results. 4.6 TIMEAVERAGED RESULTS SUMMARY Timeaveraged results obtained from FLUENT for the baseline NACA 0012 airfoil showed wellbehaved characteristics. Total pressure contours around the airfoil exhibited smooth attached flow with no signs of largescale separation and a symmetric wake profile. Pitchwise wake profiles developed downstream of the NACA 0012 lifting surface at x/c = 0.005 and x/c = 0.01 exhibited symmetry as well, further ascertaining this inference. Timeaveraged static pressure distributions along the upper and lower surface of the NACA 0012 airfoil were perfectly symmetric, as expected for this profile. Thus, indicating no aerodynamic loading exists on the airfoil. Comparison of the timeaveraged pressures with experimentally obtained values from Abbott and Von Doenhoff [1959] exhibited reasonable agreement with very slight discrepancies observed near the leading edge and trailing edge regions. Given this agreement of timeaveraged data, the simulations presented correctly predict the timeaveraged flowfield about the examined lifting surface. 61 Timeaveraged pressure distributions obtained for various mean flow angles of attack (AOA), various thickness profiles and various camber profiles exhibited characteristics similar to that, which would be expected for the respective profiles. Nonequivalent timeaveraged pressure distributions were exhibited at various mean flow angles of attack and by the various camber profiles, indicating aerodynamic loading on these profiles. Differential pressure (aerodynamic loading) was observed to increase as both AOA and percentage of camber relative to chord increased. Conversely, the various thickness profiles exhibited equivalent timeaveraged pressure distributions along the upper and lower surfaces indicative of no timeaveraged aerodynamic loading. Each thickness profile exhibits distinct timeaveraged pressure gradients along the chord, with higher thickness values resulting in more severe chordwise gradients. Computationally obtained timeaveraged results for the parametric analyses were compared with those obtained experimentally by Abbott and Von Doenhoff [1959]. The computationally obtained timeaveraged results compared favorably with the experimental data. Thereby, validating the steady state solution obtained from FLUENT for the current simulations. To validate the unsteady solution, firstharmonic lift timeseries obtained from FLUENT for the NACA 0012 airfoil was compared with the Sears predicted lift series. The timeseries matched perfectly well, even though Sears unsteady aerodynamic forcing was modeled as a perfect sinewave transverse gust as opposed to the rotorwake model (transverse and chordwise varying velocity components) used for unsteady aerodynamic forcing in the present analysis. Thus, the unsteady solution is validated. 62 5 Chapter 5 RESULTS FOR NACA 0012 BASELINE CONFIGURATION This chapter describes the aerodynamic forcing function, as well as resulting airfoil forced response, for the baseline NACA 0012 airfoil simulations. Airfoil results are reported in terms of chordwise unsteady surfacepressure distributions, spectral content, and phase. Examination of the surfacepressure phase data reveals characteristics indicative of multiple disturbance interactions, similar to that discussed by Fabian and Jumper [1996]. An analytical model is developed to reproduce these observed interaction characteristics, and compared to simulation results. 5.1 DATA REDUCTION METHODOLOGY Unsteady pressure results were obtained by removing the timeaveraged pressure from the instantaneous pressure via P′ = P − P Unsteady pressure data were further reduced into elements of amplitude, frequency, and phase for first, second and third harmonics (i.e. one, two and three times the fundamental frequency) components. This was accomplished via FastFourier Transform (FFT) techniques. 63 5.2 AIRFOIL FORCINGFUNCTION TIME DEPENDENCY Figure 51 through Figure 54 exhibit totalpressure contours of the airfoil aerodynamic forcing function, or wake, moving downstream over the NACA 0012 airfoil. The plots illustrate totalpressure contours at four time instances during a single aerodynamic forcing period, T; illustrated times correlate to t = 0, T/4, T/2, and 3T/4. The totalpressure disturbance is periodic and repeats with every forcing period. As the totalpressure disturbance translates, it directly impacts the airfoil leading edge. At the airfoil leading edge, the disturbance splits and propagates downstream along both airfoil lower and upper surfaces. This disturbance splitting process is most evident from t = T/4 to t = T/2. At impact on the airfoil upper surface, the disturbance accelerates (as compared to lower surface) before propagating downstream along the airfoil, decaying in strength during the process. On the airfoil lower surface, disturbance impact is much less prominent; disturbance chordwise propagation is delayed (i.e. it remains further upstream) relative to the uppersurface disturbance. Figure 51 Forcing Function TotalPressure Contours, t = 0. Figure 52 Forcing Function TotalPressure Contours, t = T/4. 64 Figure 53 Forcing Function TotalPressure Contours, t = T/2. Figure 54 Forcing Function TotalPressure Contours, t = 3T/4. Figure 55 displays totalpressure time series relating the airfoil forcingfunction at five equally spaced x/c locations (x/c = 0.2, 0.4, 0.6, 0.8, and –1.0) forward of the airfoil leading edge. Note that the time series have been arbitrarily shifted vertically by four units at each sequential x/c location, to provide better viewing. As expected, Figure 55 indicates periodic totalpressure variations corresponding to the aerodynamicforcing frequency. Disturbance amplitude changes as the forcing wake travels downstream, indicating expected wake decay with convective distance. A monotonic phase shift is also observed between the time series, indicating constantspeed disturbance propagation downstream, or convection. Note that while the pressure fluctuations in Figure 55 are periodic, they are not purely sinusoidal, indicating the existence of forcingfunction harmonic content. 65 Time, t (sec) "Shifted"Total Pressure, Pt (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 35.0 30.0 25.0 20.0 15.0 10.0 5.0 0.0 5.0 10.0 15.0 20.0 25.0 x/c = 0.2 x/c = 0.4 x/c = 0.6 x/c = 0.8 x/c = 1.0 Figure 55 TotalPressure Time Series Forward of Airfoil. Figure 56 displays forcingfunction staticpressure time series at five equally spaced x/c locations (20, 40, 60, 80 and 100%) forward of the airfoil leading edge. Similar amplitude trends (decreasing downstream) and periodicity are seen as in the totalpressure time series; however, no phase shift is observed between the series. Here again, the static pressure fluctuations are periodic but not purely sinusoidal, indicating the existence of harmonic content. Time, t (sec) "Shifted" Static Pressure, Ps (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 10.0 5.0 0.0 5.0 10.0 15.0 20.0 25.0 30.0 x/c = 0.2 x/c = 0.4 x/c = 0.6 x/c = 0.8 x/c = 1.0 Figure 56 StaticPressure Time series. 66 5.3 AIRFOIL FORCINGFUNCTION SPECTRAL CONTENT Figure 57 and Figure 58 depict spectral content corresponding to the forcingfunction totalpressure and staticpressure time series in Figure 55 and Figure 56, respectively. Note the xaxis in the spectral plots is reversed in order to indicate flow direction from the simulation inlet to airfoil leading edge. Significant harmonic frequencies reaching three times the primary aerodynamicforcing frequency are observed in both figures. Both total and static pressures exhibit increased harmonic content near the inlet, with totalpressure showing greater amplitude. In each case, harmonic content decays as the wake travels downstream towards the airfoil leading edge. Such harmonic amplitude decay is expected, as wake mixing with convective distance inherently causes rapid spectral content loss. Chord, x/c Unsteady Total Pressure, Pt' (Pa) 1.0 0.8 0.6 0.4 0.2 0.0 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 1st Harmonic 2ndHarmonic 3rdHarmonic Inlet Airfoil LE Figure 57 Forcing Function TotalPressure Spectral Content. Chord, x/c Unsteady Static Pressure, Ps' (Pa) 1.0 0.8 0.6 0.4 0.2 0.0 0.0 1.0 2.0 3.0 4.0 5.0 1st Harmonic 2ndHarmonic 3rdHarmonic Inlet Airfoil LE Figure 58 Forcing Function StaticPressure Spectral Content. 5.4 AIRFOIL FORCINGFUNCTION PHASE DEPENDENCY Figure 59 shows relative phase data for the firstharmonic static pressure at five equally spaced x/c locations (x/c = 0.2, 0.4, 0.6, 0.8, and –1.0) forward of the airfoil leading edge. Note the xaxis direction on the phase plot is reversed corresponding to 67 disturbance propagation direction. The firstharmonic phase data in Figure 59 shows consistent downstream disturbance propagation (i.e., convection) towards the airfoil leading edge, as indicated by the constant negative slope of the phase line. However, unsteady disturbance propagation changes dramatically as the airfoil leading edge is approached, as indicated by the positive phase slope beginning x/c = 0.2 forward of the leading edge. This change in phase slope correlates to rapid disturbance deformation near the leading edge, as can be observed in Figure 52. Chord, x/c Relative Phase (Radian) 1.0 0.8 0.6 0.4 0.2 0.0 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 1st Harmonic Figure 59 Airfoil Forcing Function 1st Harmonic Phase. Relativephase data at the second and thirdharmonic frequencies are shown in Figure 510, for the forcingfunction staticpressure data at five equally spaced x/c locations forward of the airfoil leading edge. It is observed that the second harmonic has a positively increasing phase slope; the increase in slope is very small until x/c = 0.4, after which it increases significantly. In comparison, the third harmonic phase slope is approximately zero, indicating very little propagation of the thirdharmonic static pressure component. It is interesting to note the relative phase differences between the 68 static pressure first, second and third harmonics. Each harmonic component exhibits a different phase behavior, but collectively they combine to produce a waveform showing almost no propagation characteristics, as noted in Figure 56. Chord, x/c Relative Phase (Radian) 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2nd Harmonic 3rd Harmonic Figure 510 Airfoil Forcing Function Higher Harmonic Phase. 5.5 AIRFOIL SURFACEPRESSURE TIME DEPENDENCY Figure 511 through Figure 514 show contours of static pressure for the NACA 0012 airfoil at four time instances during a single aerodynamic forcing period, T; illustrated times correlate to t = 0, T/4, T/2, and 3T/4. Figure 511 Airfoil StaticPressure Contours, t = 0. Figure 512 Airfoil StaticPressure Contours, t = T/4. 69 Figure 513 Airfoil StaticPressure Contours, t = T/2. Figure 514 Airfoil StaticPressure Contours, t = 3T/4. As desired in construction of the forcing function, staticpressure disturbances repeat in time, mimicking the effect of a wake convecting from an upstream rotor in a turbomachine compressor. At the airfoil leading edge, the wakeinduced pressure disturbances split and propagate downstream along both airfoil lower and upper surfaces. This disturbance splitting process is most evident from t = T/4 to t = T/2. At impact on the airfoil upper surface, the pressure waves also reflect back into the oncoming disturbance field. The impacted wave and its reflection travel together downstream along the airfoil, decaying in strength during the process. On the airfoil lower surface, pressure wave impact and reflection is much less prominent. Lowersurface chordwise disturbance propagation is delayed (i.e. further downstream) relative to the upper surface, as illustrated from t = T/4 to t = T/2. Figure 515 shows airfoil surfacepressure time series at various chordwise locations along the airfoil upper surface. Note, the staticpressure series have been arbitrarily shifted vertically by two units at each sequential x/c location, to provide better viewing. As expected, these figures indicate periodic staticpressure variations corresponding to the aerodynamicforcing frequency. While the pressure fluctuations are 70 periodic, they are not purely sinusoidal, indicating the existence of harmonic content. Higher amplitude unsteady pressure fluctuations exist near the airfoil leading edge, decaying rapidly downstream. Figure 5.15 also indicates almost no phase delay between time series at each x/c locations, a curious finding given the convective nature of the forcing function. This finding corresponds to the freestream staticpressure time series of Figure 56, and will be examined further in following sections. Time, t (sec) "Shifted"Unsteady Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 515 NACA 0012 Upper Surface Unsteady StaticPressure Time series. Figure 516 shows airfoil surfacepressure time series at various chordwise locations along the airfoil lower surface. Similar to the upper surface, the lower surface exhibits periodic timeseries behavior corresponding to the aerodynamicforcing frequency, unsteady leadingedge pressure amplification, and very little phase shift between chordwise locations. The lowersurface pressure fluctuations also display an almost reversed behavior as compared to the uppersurface, showing pressure increases at 71 the same instances as the upper surface experiences a pressure decrease. However, overall timeseries amplitudes are greater for the lower surface as compared to the upper surface. Time, t (sec) "Shifted"Unsteady Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 24.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 516 NACA 0012 Lower Surface Unsteady StaticPressure Time series. Differential surfacepressure time series data at various chordwise locations along the airfoil are provided in Appendix D. Differential surface pressure represents uppersurface minus lowersurface unsteady pressure, as presented in Eq. (5.1). u l Δp′ = p′ − p′ 5.1 Here again, periodic timeseries behavior corresponding to the aerodynamicforcing frequency is obtained. The differential pressure is highest in amplitude at the leading edge, decays in strength downstream along the airfoil, and shows almost no phase shift between chordwise locations. Differentialpressure fluctuations are periodic in nature, but not purely sinusoidal, indicating the existence of harmonic content. 72 5.6 AIRFOIL SURFACEPRESSURE SPECTRAL CONTENT Figure 517 depicts spectral content related to the airfoil uppersurface pressure time series shown in Figure 515. Relevant surfacepressure harmonic frequencies reaching four times the aerodynamicforcing frequency are observed. The airfoil upper surface exhibits increased higherorder pressure harmonics near the airfoil leading edge, generally decaying downstream along the upper surface. Specifically, the firstharmonic amplitude decays until x/c = 0.4, then increases toward the airfoil trailing edge. The second harmonic amplitude also shows a decaying trend downstream, but between x/c = 0.4 and x/c = 0.6 it exhibits higher amplitude than the first harmonic. The third and fourth harmonic amplitudes generally show a decaying trend with chordwise distance, almost vanishing near the airfoil trailing edge. Note that harmonic content can be related to the “sinelike” behavior of a time series. For example, it would be expected that a wake timeseries profile, such as in Figure 55, would exhibit significant harmonic content, while a pure sine wave would not. Thus, the relative change in harmonic content with chord position, as seen in Figure 517, can be directly attributed to disturbance deformation and interaction over the airfoil. Near the leading edge the disturbance is “wakelike” while at the trailing edge it is much more “sinelike”. This statement is supported by Figure 515. 73 Chord, x/c Surface Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.5 1.0 1.5 2.0 2.5 1st Harmonic 2nd Harmonic 3rd Harmonic 4th Harmonic Figure 517 NACA 0012 UpperSurface Static Pressure Spectral Content. Spectral content related to the airfoil lowersurface pressure time series is depicted in Figure 518, corresponding to the timeseries data of Figure 516. Similar to the airfoil upper surface, the lower surface shows relevant surfacepressure frequencies reaching four times the aerodynamicforcing frequency, as well as increased higherorder pressure harmonics near the airfoil leading edge. The firstharmonic amplitude shows a similar trend to the upper surface, decaying downstream until x/c = 0.6 and then increasing. Unlike the upper surface, the second harmonic amplitude decays constantly downstream, always having lower amplitude than the first harmonic. Higher order harmonics decay downstream along the airfoil lower surface, approximately reaching a zero value near the trailing edge. Again, the change in harmonic content between the upper and lower surfaces indicates a difference in disturbance deformation and interaction over the airfoil lower surface. This statement is supported by Figure 516. 74 Chord, x/c Surface Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.5 1.0 1.5 2.0 2.5 1st Harmonic 2nd Harmonic 3rd Harmonic 4th Harmonic Figure 518 NACA 0012 LowerSurface StaticPressure Spectral Content. 5.7 AIRFOIL SURFACEPRESSURE FIRST HARMONIC AMPLITUDE CHORDWISE DEPENDENCY Figure 519 and Figure 520 display firstharmonic staticpressure time series at various chordwise locations along the airfoil upper and lower surfaces, respectively. Note, the staticpressure series have been arbitrarily shifted vertically by two units at each successive x/c location, to provide better viewing. As expected, these figures indicate periodic pressure variations corresponding to the aerodynamic forcing frequency. A phase shift is also observed between chord locations, (as indicated by the arrows) indicating disturbance propagation direction along the chord. On both upper and lower surfaces, large unsteady pressure fluctuations exist near the airfoil leading edge, decaying rapidly downstream. The overall timeseries amplitudes are greater for the lower surface as compared to the upper surface. 75 Time, t (sec) 1st Harmonic Unsteady Pressure, (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 x/c = 0.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 519 NACA 0012 Upper Surface 1st Harmonic Pressure Time series. Time, t (sec) 1st Harmonic Unsteady Pressure 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 24.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 520 NACA 0012 Lower Surface 1st Harmonic Pressure Time series. 76 Figure 521 exhibits firstharmonic pressure amplitude dependence on chord; upper, lower and differential unsteady pressure amplitudes are shown. The airfoil lower surface exhibits higheramplitude unsteady pressures as compared to the upper surface; thus, confirming inferences made earlier related to Figures 519 and 520. Upper, lower and differential unsteady pressure amplitudes follow similar chordwise trends; i.e., amplified at leading edge and decaying downstream along the chord. Chord, x/c 1st Harmonic Amplitude, (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 Upper Lower Differential Figure 521 NACA 0012 1st Harmonic Pressure Amplitude. 5.8 AIRFOIL SURFACEPRESSURE FIRST HARMONIC CHORDWISE PHASE DEPENDENCY Figure 522 shows relative phase data for the firstharmonic unsteady surface pressures along the airfoil. The figure displays airfoil uppersurface, lowersurface and differential pressure phase along the airfoil chord. 77 Chord, x/c 1st Harmonic Relative Phase (Radians) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 3.0 2.0 1.0 0.0 1.0 2.0 3.0 4.0 Upper Lower Differential Figure 522 NACA 0012 1st Harmonic Phase. Several important results can be immediately obtained from Figure 522. First, the uppersurface and lowersurface phase data exhibit different chordwise trends. Second, the uppersurface data show a significantly varying slope, particularly near the leading and trailing edges. Third, the lowersurface phase exhibits almost no slope near the midchord, increasing in slope towards the trailing edge. Finally, phase changes appear rapidly near the leading edges on both upper and lower surfaces. Fabian and Jumper [1996] established phasemap analysis as a viable means of determining airfoil disturbancepropagation characteristics in an unsteady flow. Fabian and Jumper argued that a negative slope in a phase map, such as near the leading edge in Figure 522, indicates downstreampropagating surfacepressure waves, while a positive slope predicts upstreampropagating disturbances. Clearly, if these arguments are true, the data in Figure 522 exhibit nonphysical disturbancepropagation characteristics. 78 Similar nonphysical phase data were also discussed by Fabian and Jumper, from which they indicated “ambiguous” phase maps, like the one in Figure 522, most likely result from multiple pressure disturbances interacting over the airfoil chord. 5.9 AIRFOIL SURFACEPRESSURE ANALYTICAL MODEL The phase data of Figure 522 display ambiguous disturbance propagation, similar to that discussed by Fabian and Jumper [1996]. The phase data do not follow the convected speed along the airfoil, from which the downstream propagation of an airfoil forcing disturbance would exhibit a negative phase slope. In fact, based on the results of Fabian and Jumper, the uppersurface phase data of Figure 522 suggest an initially downstreampropagating disturbance from x/c = 0.0 to x/c = 0.1, with the disturbance gaining speed from x/c = 0.1 to x/c = 0.3, as evidenced by the flat phase slope. At x/c = 0.3, the disturbance presumably changes direction, propagating upstream from x/c = 0.55 to x/c = 0.8. Finally, the disturbance propagates upstream at a lower speed near the airfoil trailing edge. Conversely, the lowersurface phase map of Figure 522 suggests an initially downstreampropagating disturbance gaining speed from x/c = 0.1 to x/c = 0.65, finally propagating upstream near the airfoil trailing edge. These inferences are reinforced by examining the fundamental (firstharmonic) frequency time series in Figure 519 and Figure 520, for the airfoil upper and lower surfaces, respectively. In Figures 519 and 520, the probable propagation direction of the waveforms from the leading edge to each successive x/c location has been indicated, where it has been assumed that the minimum time between each successive pressure peak provides the correct wave propagation direction. When compared to the uppersurface 79 and lowersurface phase maps of Figure 522, the wave propagation directions in Figure 519 and 520 are consistent, showing the same overall chordwise characteristics. 5.9.1 INTERACTION MODEL Although the phase information in Figures 519, 520 and 522 correlate, the disturbance propagation directions implied by the data are clearly nonphysical and incorrect. In the investigation of Fabian and Jumper, ambiguous phase maps similar to Figure 522 were argued to result from convectivelypropagating vortical and acousticallyradiating potential disturbance int
Click tabs to swap between content that is broken into logical sections.
Rating  
Title  Parametric Analysis of Gustinduced Airfoil Surfacepressure Phase Behavior 
Date  20040701 
Author  Bakhtiani, Harikishin Prakash 
Keywords  Cfd, Gustinduced, Phase, Fluent, Airfoil, Cascade 
Department  Mechanical Engineering 
Document Type  
Full Text Type  Open Access 
Abstract  The current study investigates chordwise phase distributions on twodimensional lifting surfaces using computational fluid dynamic (CFD) analysis, particularly through the CFD software FLUENT 6.0. A series of CFD simulations are performed in which convecting twodimensional fluid dynamic disturbances are allowed to propagate over a variety of lifting surfaces (or airfoils). These unsteady forcing scenarios examine the most basic chordwise phase behavior, focusing on the influence of airfoil thickness, camber, mean loading, and neighboring surfaces (e.g., a cascade of airfoils). Finally, a comparison between the computed phase results and experimental data is conducted, in an attempt to provide a definitive explanation of the observed experimental phase trends. Unsteady surfacepressure results correlate with previous investigations indicating accurately modeled flow physics. Unsteady pressure amplitudes on the lifting surfaces exhibit leading edge amplification decaying rapidly downstream. Spectral analyses show surfacepressure harmonic frequencies reaching four times the aerodynamic forcing frequency. Firstharmonic phase data exhibit ambiguous disturbance positionvs.phase maps, from which disturbance propagation direction can not be inferred. Based on the results of a previous investigation, the ambiguous phase results suggest the existence of multiple disturbance interactions across the airfoil, these disturbances having the same primary frequency, but different amplitudes and propagation speeds. A firstorder disturbanceinteraction model developed in an attempt to reproduce the ambiguous phase behavior, provided evidence suggesting phasemap ambiguity results from interaction between local convectively propagating and global acoustically propagating pressure waves. Parametric results indicate increase in midchord phase slope as liftingsurface thickness and camber were increased. However, as angleofattack increased a mere shift in the midchord phase slope variation towards the leading edge was observed. Finally, the computed results compared favorably to results collected in a previous experiment. Phase data exhibited similar ambiguity to that seen in the experimental results, further suggesting the existence of convectively propagating and acoustically propagating wave interaction. 
Note  Thesis 
Rights  © Oklahoma Agricultural and Mechanical Board of Regents 
Transcript  A PARAMETRIC ANALYSIS OF GUSTINDUCED AIRFOIL SURFACEPRESSURE PHASE BEHAVIOR By HARIKISHIN PRAKASH BAKHTIANI Bachelor of Science Oklahoma State University Stillwater, Oklahoma 2002 Submitted to the Faculty of the Graduate College of the Oklahoma State University in partial fulfillment of the requirements for the Degree of MASTER OF SCIENCE July, 2004 ii A PARAMETRIC ANALYSIS OF GUSTINDUCED AIRFOIL SURFACEPRESSURE PHASE BEHAVIOR Thesis Approved: ________________________________________________ Thesis Advisor ________________________________________________ ________________________________________________ _______________________________________________ Dean of the Graduate College iii ACKNOWLEGMENTS A work of this magnitude was made possible only by the assistance of many individuals. First and foremost, I would like to thank my graduate advisor and committee chair Dr. Eric A. Falk. Thank you, sir, for your continuous support, guidance and hard work. You are truly a great friend and mentor. Words cannot express my appreciation of your encouragement, inspiration and patience. Gratitude and appreciation are also extended to my other committee members, Dr. Andrew S. Arena and Dr. Afshin J.Ghajar. To my friends and colleagues at the Turbolab, it’s been both an honor and privilege to have met and worked with you all. Aaron thanks for all your help. Robert thanks for answering all my questions and guiding me when direction was needed. Finally, a special ‘thank you’ to my parents, Prakash and Jaya Bakhtiani for their continuous love, support and words of encouragement. The successes I have achieved did not come without certain sacrifices, which they endured in some form. Finally to my brother, Nicky, who never fails to amaze me. TABLE OF CONTENTS 1 INTRODUCTION..................................................................................................... 1 1.1 OVERVIEW OF AIRFOIL GUST INTERACTIONS....................................... 1 1.2 COMPRESSOR DESIGN INTENT AND OPERATING ENVIRONMENT ... 1 1.3 NORMAL COMPRESSOR OPERATION LEADING TO AIRFOILGUST INTERACTIONS ............................................................................................... 2 1.4 IMPORTANCE OF HCF TO ENGINE COMMUNITY ................................... 4 1.5 IMPORTANCE OF PHASE DISTRIBUTIONS IN PREDICTING AIRFOIL MODAL FORCING ........................................................................................... 5 1.6 NEED FOR CONTINUED SURFACEPRESSURE PHASE ANALYSIS ...... 6 1.7 SCOPE OF CURRENT INVESTIGATION ...................................................... 7 2 PREVIOUS WORK.................................................................................................. 8 2.1 LITERARY REVIEW ........................................................................................ 8 2.1.1 SingleAirfoil Investigations....................................................................... 8 2.1.2 Turbomachine and Cascade Investigations................................................. 9 2.2 UNRESOLVED ISSUES.................................................................................. 14 2.3 CURRENT RESEARCH OBJECTIVE............................................................ 20 3 COMPUTATIONAL METHODOLOGY AND SETUP .................................... 22 3.1 AIRFOIL GEOMETRY AND BOUNDARY CONDITIONS......................... 22 3.2 STATORVANE GEOMETRY AND BOUNDARY CONDITIONS............. 23 3.3 GENERAL FLUENT SOLVER DESCRIPTION............................................ 25 iv 3.3.1 Coupled Solution Method ......................................................................... 27 3.3.2 ReynoldsAveraged NavierStokes (RANS) Equations ........................... 28 3.3.3 k – ε Turbulence Model ............................................................................ 30 3.3.4 FiniteVolume Discretization Methodology ............................................. 31 3.3.5 SecondOrder Upwind Scheme................................................................. 32 3.3.6 SecondOrder Time Discretization ........................................................... 32 3.3.7 Linearization Methodology....................................................................... 34 3.3.8 Periodic Boundary Conditions.................................................................. 35 3.3.9 Standard FLUENT Inlet/Outlet/Wall Boundary Conditions .................... 35 3.3.10 Operating Pressure .................................................................................... 36 3.4 UDF DESCRIPTION........................................................................................ 37 3.4.1 Development Logic/Procedure ................................................................. 37 3.5 COMPUTATIONAL GRID DESCRIPTION .................................................. 39 3.5.1 Gambit Grid Generation Software ............................................................ 39 3.5.2 Grid Methodology..................................................................................... 40 3.6 COMPUTATIONAL SETUP........................................................................... 46 3.6.1 Reference Values ...................................................................................... 48 3.6.2 Boundary Conditions ................................................................................ 48 3.7 GRID INDEPENDENCE ................................................................................. 49 4 TIMEAVERAGED RESULTS ............................................................................ 52 4.1 TIMEAVERAGED METHODOLOGY ......................................................... 52 4.2 BASELINE AIRFOIL TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS ............................................................................................ 53 v 4.3 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (THICKNESS INFLUENCE)................................................................................................... 57 4.4 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (CAMBER INFLUENCE)................................................................................................... 58 4.5 COMPARISON OF TIMEACCURATE BASELINE LIFT DEPENDENCY TO SEAR’S RESULTS.................................................................................... 60 4.6 TIMEAVERAGED RESULTS SUMMARY ................................................. 61 5 RESULTS FOR NACA 0012 BASELINE CONFIGURATION ........................ 63 5.1 DATA REDUCTION METHODOLOGY ....................................................... 63 5.2 AIRFOIL FORCINGFUNCTION TIME DEPENDENCY ............................ 64 5.3 AIRFOIL FORCINGFUNCTION SPECTRAL CONTENT.......................... 67 5.4 AIRFOIL FORCINGFUNCTION PHASE DEPENDENCY ......................... 67 5.5 AIRFOIL SURFACEPRESSURE TIME DEPENDENCY............................ 69 5.6 AIRFOIL SURFACEPRESSURE SPECTRAL CONTENT.......................... 73 5.7 AIRFOIL SURFACEPRESSURE FIRST HARMONIC AMPLITUDE CHORDWISE DEPENDENCY....................................................................... 75 5.8 AIRFOIL SURFACEPRESSURE FIRST HARMONIC CHORDWISE ... PHASE DEPENDENCY.................................................................................. 77 5.9 AIRFOIL SURFACEPRESSURE ANALYTICAL MODEL......................... 79 5.9.1 Interaction Model...................................................................................... 80 5.9.2 Interaction Model Results......................................................................... 86 5.10 SUMMARY...................................................................................................... 88 6 PARAMETRIC ANALYSIS.................................................................................. 90 6.1 DATA REDUCTION METHODOLOGY ....................................................... 90 6.2 AIRFOIL THICKNESS INFLUENCE............................................................. 90 6.3 AIRFOIL CAMBER INFLUENCE.................................................................. 93 vi 6.4 ANGLE OF ATTACK INFLUENCE .............................................................. 96 6.5 SUMMARY...................................................................................................... 98 7 STATORVANE CASCADE RESULTS............................................................ 100 7.1 DATA REDUCTION METHODOLOGY ..................................................... 100 7.1.1 TimeAveraged Results .......................................................................... 100 7.1.2 Unsteady Results..................................................................................... 101 7.2 STATORVANE CASCADE CONFIGURATION....................................... 101 7.3 TIMEAVERAGED RESULTS..................................................................... 101 7.4 UNSTEADY PRESSURE RESULTS............................................................ 103 7.5 SUMMARY.................................................................................................... 107 8 SUMMARY AND CONCLUSIONS ................................................................... 108 8.1 RESULTS ....................................................................................................... 108 8.1.1 TimeAveraged Results .......................................................................... 108 8.1.2 Unsteady Pressure Data .......................................................................... 110 8.2 CORRELATIONS WITH PREVIOUS INVESTIGATIONS........................ 113 8.3 CURRENT CONTRIBUTIONS..................................................................... 114 8.4 RECOMMENDATIONS FOR FUTURE WORK ......................................... 115 A. APPENDIX A........................................................................................................ 120 B. APPENDIX B ........................................................................................................ 122 C. APPENDIX C........................................................................................................ 124 D. APPENDIX D........................................................................................................ 128 E. APPENDIX E ........................................................................................................ 131 F. APPENDIX F ........................................................................................................ 135 vii G. APPENDIX G........................................................................................................ 136 viii LIST OF TABLES Table 3.1 Modeled StatorVane Cascade Geometry......................................................... 24 Table 3.2 Grid Distribution (AirfoilCascade Mesh)........................................................ 41 Table 3.3 BoundaryLayer Mesh Characteristics (AirfoilCascade Mesh). ..................... 43 Table 3.4 StatorVane Cascade Mesh Node Distribution................................................. 45 Table 3.5 FLUENT Configuration for Numerical Simulations. ....................................... 47 Table 3.6 FLUENT Reference Values.............................................................................. 48 Table 3.7 FLUENT Boundary Conditions........................................................................ 49 Table 4.1 Coefficient of Lift: 5° and 10° Angle of Attack. .............................................. 56 Table 4.2 Coefficient of Lift: Various LiftingSurface Cambers. .................................... 60 ix LIST OF FIGURES Figure 21 Relative Phase Variation with Chord: RearwardForced Cascade of Fabian et al. [1996]. .................................................................................................................. 12 Figure 22 Relative Phase Variation with Chord: ForwardForced Cascade of Fabian et al. [1996]. .................................................................................................................. 13 Figure 23 Unsteady Differential SurfacePressure Time Series: (φ (x / c) = 0). .............. 16 Figure 24 Unsteady Differential SurfacePressure Time Series: (φ (x / c) ≠ 0). .............. 16 Figure 25 Unsteady Differential SurfacePressure Variation with Chord: (φ (x / c) = 0). 17 Figure 26 Unsteady Differential SurfacePressure Variation with Chord:(φ (x / c) ≠ 0). 17 Figure 27 RigidBody and FlexibleBody Mode Shapes for a Simply Supported, Infinite Span, TwoDimensional Lifting Surface. ................................................................. 18 Figure 28 Maximum Generalized forces on Structural Modes for Various Phase Distributions.............................................................................................................. 19 Figure 31 AirfoilCascade Computational Boundaries. .................................................. 23 Figure 32 StatorVane Cascade Computational Boundaries. .......................................... 25 Figure 33 Overview of the Coupled Solution Method [Fluent, 2001]. ........................... 28 Figure 34 Rotor Wake Characteristics............................................................................. 38 Figure 35 AirfoilCascade Mesh (NACA 0012 airfoil)................................................... 42 Figure 36 Structured BoundaryLayer Mesh Surrounding NACA 0012 Airfoil Geometry. .................................................................................................................................. 43 Figure 37 StatorVane Cascade. ...................................................................................... 45 Figure 38 Modeled StatorVane Cascade Mesh.............................................................. 46 x Figure 39 Lift Coefficient Time History Showing Convergence: NACA 0012 Profile.. 47 Figure 310 SteadyState Pressure Coefficient Data: NACA 0012. ................................. 50 Figure 311 Pressure Coefficient Data: StatorVane. ....................................................... 51 Figure 41 TimeAveraged TotalPressure Contours (Pa): NACA 0012 Lifting Surface 53 Figure 42 NACA 0012 Lifting Surface Wake Profile..................................................... 54 Figure 43 TimeAveraged Static Pressure : NACA 0012. .............................................. 55 Figure 44 TimeAveraged Static Pressure : 0°, 5° and 10° Angle of Attack .................. 56 Figure 45 TimeAveraged StaticPressure for Various LiftingSurface Thicknesses..... 57 Figure 46 Time averaged Static Pressure Comparison with Experimental data. ............ 58 Figure 47 TimeAveraged StaticPressure Distribution for Various LiftingSurface Cambers .................................................................................................................... 59 Figure 48 LiftTime series Comparison with Sears results. ............................................ 61 Figure 51 Forcing Function TotalPressure Contours, t = 0............................................ 64 Figure 52 Forcing Function TotalPressure Contours, t = T/4........................................ 64 Figure 53 Forcing Function TotalPressure Contours, t = T/2. ....................................... 65 Figure 54 Forcing Function TotalPressure Contours, t = 3T/4. ..................................... 65 Figure 55 TotalPressure Time Series Forward of Airfoil. ............................................. 66 Figure 56 StaticPressure Time series. ............................................................................ 66 Figure 57 Forcing Function TotalPressure Spectral Content. ........................................ 67 Figure 58 Forcing Function StaticPressure Spectral Content. ....................................... 67 Figure 59 Airfoil Forcing Function 1st Harmonic Phase. ................................................ 68 Figure 510 Airfoil Forcing Function Higher Harmonic Phase........................................ 69 Figure 511 Airfoil StaticPressure Contours, t = 0......................................................... 69 xi Figure 512 Airfoil StaticPressure Contours, t = T/4. .................................................... 69 Figure 513 Airfoil StaticPressure Contours, t = T/2. .................................................... 70 Figure 514 Airfoil StaticPressure Contours, t = 3T/4. .................................................. 70 Figure 515 NACA 0012 Upper Surface Unsteady StaticPressure Time series.............. 71 Figure 516 NACA 0012 Lower Surface Unsteady StaticPressure Time series. ............ 72 Figure 517 NACA 0012 UpperSurface Static Pressure Spectral Content. .................... 74 Figure 518 NACA 0012 LowerSurface StaticPressure Spectral Content..................... 75 Figure 519 NACA 0012 Upper Surface 1st Harmonic Pressure Time series. ................. 76 Figure 520 NACA 0012 Lower Surface 1st Harmonic Pressure Time series.................. 76 Figure 521 NACA 0012 1st Harmonic Pressure Amplitude............................................ 77 Figure 522 NACA 0012 1st Harmonic Phase. ................................................................. 78 Figure 523 Wake Induced Pressure Time Series Collected at the Periodic Boundary.... 82 Figure 524 Optimized LiftInduced Pressure Time Series.............................................. 83 Figure 525 Interaction Model 1st Harmonic Time series................................................. 84 Figure 526 Propagating Disturbance Model.................................................................... 85 Figure 527 1st Harmonic Amplitude Comparison. .......................................................... 86 Figure 528 Interaction Model 1st Harmonic Phase.......................................................... 87 Figure 61 1st Harmonic Pressure Amplitude: Various LiftingSurface Thicknesses. .... 92 Figure 62 1st Harmonic Phase : Various LiftingSurface Thicknesses........................... 93 Figure 63 1st Harmonic Amplitude: Various LiftingSurface Cambers. ........................ 94 Figure 64 1st Harmonic Phase: Various LiftingSurface Cambers. ................................. 95 Figure 65 1st Harmonic Amplitude: Various MeanFlow Angles of Attack. .................. 97 Figure 66 1st Harmonic Phase: Various Mean Flow Angles of Attack. .......................... 98 xii Figure 71 TimeAveraged Total Pressure Contours (Pa)............................................. 102 Figure 72 TimeAveraged Static Pressure Distribution. ............................................... 103 Figure 73 TimeAverage Pressure Distribution [Falk, 2000]………………………… 103 Figure 74 RMS Unsteady Pressure Distribution……………………………………... 109 Figure 75 P'RMS Distribution [Falk, 2000]. .................................................................... 104 Figure 76.1st Harmonic Amplitude: StatorVane. ......................................................... 105 Figure 77 1st Harmonic Phase: StatorVane. ................................................................. 106 Figure D1 Unsteady Differential Static Pressure TimeSeries...................................... 128 Figure D2 Unsteady DifferentialPressure Spectral Content. ....................................... 129 Figure D3 1st Harmonic DifferentialPressure TimeSeries. ........................................ 129 Figure D4 Higher Harmonic Amplitudes: NACA 0012................................................ 130 Figure D5 Higher Harmonic Phase: NACA 0012. ........................................................ 130 Figure E6 1st Harmonic TimeSeries: NACA 0010 Upper Surface. ............................. 131 Figure E7 1st Harmonic TimeSeries: NACA 0010 Lower Surface.............................. 131 Figure E8 1st Harmonic TimeSeries: NACA 0015 Upper Surface. ............................. 132 Figure E9 1st Harmonic TimeSeries: NACA 0015 Lower Surface.............................. 132 Figure E10 1st Harmonic TimeSeries: NACA 0020 Upper Surface. ........................... 132 Figure E11 1st Harmonic TimeSeries: NACA 0020 Lower Surface............................ 132 Figure E12 SurfacePressure Spectral Content: NACA 0010 Upper Surface............... 133 Figure E13 SurfacePressure Spectral Content: NACA 0015 Upper Surface............... 133 Figure E14 SurfacePressure Spectral Content: NACA 0020 Upper Surface............... 134 Figure F15 1st Harmonic TimeSeries: 2% Camber Airfoil Upper Surface. ................ 135 Figure F16 1st Harmonic TimeSeries: 2% Camber Airfoil Lower Surface.................. 135 xiii Figure F17 1st Harmonic TimeSeries: 6% Camber Airfoil Upper Surface. ................. 135 Figure F18 1st Harmonic TimeSeries: 6% Camber Airfoil Lower Surface.................. 135 Figure G19 1st Harmonic TimeSeries: 5Degree AOA Upper Surface........................ 136 Figure G20 1st Harmonic TimeSeries: 5Degree AOA Lower Surface. ...................... 136 Figure G21 1st Harmonic TimeSeries: 10Degree AOA Upper Surface……………...136 Figure G22 1st Harmonic TimeSeries: 10Degree AOA Lower Surface. .................... 136 xiv NOMENCLATURE Symbols ΔCp′ = unsteady differentialpressure coefficient Δp′ = differential unsteady pressure n) = normal vector to liftingsurface chord ω = disturbance angular frequency φ = surfacepressure phase distribution ∞ ρ = freestream mean density m ψ r = mth mode shape vector a = forcingdisturbance transverse velocity c = liftingsurface chord m f = mth mode generalized force (2) i H = ith order Hankel function, 2nd kind k = reduced frequency S = Sears function ∞ U = freestream mean velocity x/c = nondimensional distance along chord i u = mean velocity component i u′ = unsteady velocity component xv k G = turbulent kinetic energy due to mean velocity gradients b G = turbulent kinetic energy due to buoyancy ε = dissipation rate vr = velocity vector (=ui vj) ) ) + in 2D A v = surface are vector φ Γ = diffusion coefficient for φ φ S = source of φ per unit volume faces N = number of faces enclosing cell f f f A r ρ υ = mass flux through the face n (∇φ ) = magnitude of ∇φ normal to face f Pt,i = total pressure inlet Pt = total pressure P′ = unsteady pressure P = instantaneous pressure P = timeaveraged pressure TE = trailing edge LE = leading edge Acronyms CFD = computational fluid dynamics AOA = angle of attack RANS = Reynoldsaveraged NavierStokes equations xvi 1 CHAPTER 1 INTRODUCTION 1.1 OVERVIEW OF AIRFOIL GUST INTERACTIONS Liftingsurface response to unsteady aerodynamic forcing is of particular interest in aircraft propulsion applications, primarily due to timeresolved aerodynamic interactions in turbomachinery. In these applications, lifting surfaces (i.e., airfoils) often operate in both randomly turbulent and periodically oscillating fluid environments, including temporally and spatially nonuniform propagating disturbances. Relative unsteady motion between a lifting surface and the fluid results in complex fluidstructure interactions and may produce aerodynamic/structural response. In order to optimize liftingsurface designs, a detailed understanding of inherent fluidstructure interactions is required, where these interactions can be described, in part, by liftingsurface pressure response. Unsteady surfacepressure phase distributions compromise one aspect of liftingsurface fluidstructure interactions, as phase directly affects timeresolved unsteady force/moment behavior, particularly in terms of forcing structural modes. 1.2 COMPRESSOR DESIGN INTENT AND OPERATING ENVIRONMENT For the purpose of emphasizing the continued need for investigating unsteady fluidstructure interactions, consider an axialflow compressor. Axialflow compressors are typically composed of a number of rotating blades for the purpose of turning and 1 adding work to the passing airflow. This turning process accomplishes a desired total enthalpy rise through the device. A single row of rotating blades is often referred as a blade row. A blade row can lead or follow a separate single row of stationary vanes often referred to as a vane row. Vane rows direct rotor inlet/outlet airflow corresponding to the compressor design. In modern compressors, achieving the design rise requires several blade/vane rows, each providing a portion of the enthalpy difference. Each pair of blade/vane rows represents a stage in the compressor [Falk, 2000]. With the overall goal of engine design being often to decrease engine size and weight, compressor size and weight must also be decreased. A lower number of stages and reduced axial spacing between stages helps to accomplish this goal. However, operating at lower number of stages requires a corresponding increase in stage aerodynamic loading to achieve the desired compressor enthalpy rise. In addition, reduced stagetostage spacing leads to greater aerodynamic interactions between vane/blade rows. Such interactions come in the form of disturbances caused by the relative motion between the rotor and stator rows. 1.3 NORMAL COMPRESSOR OPERATION LEADING TO AIRFOILGUST INTERACTIONS The rotational motions between rotor/stator rows, or stages, in a turbineengine compressor generate a designed enthalpy rise across the component. In the process, however, each stage induces propagating aerodynamic disturbances that act as periodic excitations, or forcing functions, for neighboring blade/vane rows. These propagating disturbances are generally grouped as: 2 • Convective downstreampropagating viscous wakes: produced by the frictional interactions between the fluid and lifting surface. • Convective downstreampropagating vortical wakes: produced by vortex shedding in response to bound circulation fluctuations on the lifting surface. • Acoustically propagating potential disturbances: elicited by variations in the velocity potential, or pressure fields, associated with the blades of a given row [Hall, 1991]. Induced potential disturbances may be temporary in nature, decaying exponentially in the near field, or propagate without attenuation into the far field, depending on bladetip Mach numbers. Interactions between a lifting surface and propagating disturbance field induces timedependent angleofattack changes on the body, causing spatially and temporally dependent surfacepressure distributions. Integration of these surfacepressure distributions forms unsteady forces and moments on the body. Moments and forces generate temporally and spatially dependent mechanical stresses, or alternating stresses. If the induced alternating stresses are strong enough, structural fatigue may plague the liftingsurface with the possibility of catastrophic failure. Fatigue is a process of cumulative structural damage caused by repeated load fluctuations, or stresses [Barsom, 1987]. Fatigue occurs in regions deforming plastically under applied loads; thus, under purely elastic stress conditions, localized areas of raised stress must be present to induce fatigue, where these raised stresses exceed the material yield stress. Prolonged exposure to fatigueinducing unsteady loads may cause initiation and subsequent propagation of a crack, or cracks, in plastically deformed structural regions. Eventually, if crack propagation continues, catastrophic fracture and failure of 3 the lifting surface may occur. Typically the number of load fluctuations necessary to initiate a crack within working liftingsurfaces( such as a compressor blade/vane row) is quite large; thus, the cumulative structural damage process is often referred to as high cycle fatigue, or HCF. Aerodynamically induced load fluctuations in a jet engine compressor composed of highstrength structural components are designed to yield only elastic stress fluctuations. However, the occurrence of random material defects, foreign object damage (FOD), or blade rubbing can provide the proper conditions for crack initiation and propagation. Nonetheless, the number of load fluctuations typically required to form a crack within a jet engine blade/vane component is quite large, leading to HCF. 1.4 IMPORTANCE OF HCF TO ENGINE COMMUNITY Highcycle fatigue is of utmost importance in current jet engine design, were small structural failures can greatly affect entire engine operations to the extent of engine failure. Due to its importance, a considerable number of experimental and computational investigations have been conducted with the overall goal of predicting liftingsurface HCF failures in jet engines. Nonetheless, jet engine HCF failures continue to occur and are largely unanticipated. Recent advances in engine technology may only complicate this problem, as current trends toward higher blade loadings, increased operating temperatures, smaller stagetostage spacing, unconventional geometries, and advanced materials reach beyond traditional design domains, complicating HCFresistant technology [Fleeter, 1992]. In practice, HCFrelated engine failure has been identified as a major contributor to enginesafety mishaps in U. S. military fighter aircraft [Thompson, 1998], where as many as 50% of all engine failures have been attributed to HCF. As 4 such, HCF failure presents a major readiness and monetary concern for both the U.S. Air Force and U.S. Navy [Fecke, 1998]. In an attempt to overcome reoccurring HCF problems in military engines, the U.S. Department of Defense established the National Turbine Engine High Cycle Fatigue Program in 1994. The goal of this ongoing program is to develop, implement, and validate damage tolerant design methodologies to avoid HCFrelated engine failures. This goal is to be accomplished by increasing the level of understanding regarding HCF physics, as well as through the development of better HCF predictive capabilities. The specific goal of improving predictions of engine component response to unsteady aerodynamic forcing belongs to the Science and Technology branch of the HCF program. 1.5 IMPORTANCE OF PHASE DISTRIBUTIONS IN PREDICTING AIRFOIL MODAL FORCING In order to avoid structural vibrations and HCFrelated failures, it is important to accurately predict timeresolved generalized forces for each liftingsurface structural mode. Unsteady surfacepressure phase distributions represent one component of such predictions. Chordwisevarying phase distributions influence timeresolved surfacepressure amplitude distributions along the chord. Chordwise varying phase may also affect surfacepressure node locations, where the node locations change positions with different phase distributions. Mispredictions of chordwisevarying phase may result in under predicted modal forces, generating greaterthanexpected mechanical stresses at multiple spatial and temporal frequencies, and therefore HCF. In all, unsteady surfacepressure phase distributions play a very important role in liftingsurface unsteady forcing and thus predictions of modal forcing. A brief example illustrating the importance of 5 considering unsteady surfacepressure phase distributions in liftingsurface forcing is provided in Chapter 2. 1.6 NEED FOR CONTINUED SURFACEPRESSURE PHASE ANALYSIS Given the potential for catastrophic structural failure caused by liftingsurface fatiguelife degradation, accurate fluidstructure interaction predictions have been aggressively sought for aerodynamic unsteady forcing problems. In fact, extensive amounts of information are available regarding the influence of forced response on turbomachine lifting surfaces (due to the propensity of highcycle fatigue failures in modern highperformance gas turbine engines). Representative investigations have predicted, or measured, unsteady surfacepressure distributions on various lifting bodies, and characterized these distributions in terms of amplitude, frequency, and phase. These investigations typically focus on surfacepressure amplitude and frequency, with little attention to the influence of surfacepressure phase. This is not to say that phase has been completely ignored. In fact, the contrary is true. Researchers have reported surfacepressure phase data over lifting surfaces under a variety of forcing conditions, as will be discussed in Chapter #2. Despite the inclusion of phase results in many research investigations, however the dependence of liftingsurface response to variations in chordwise surfacepressure phase distribution remains relatively unexamined. Moreover, no known investigation has developed general “rules of thumb” to act as guidelines in predicting phase distributions for the most common forcing configurations. Lastly, no consistent explanation exists for observed surfacepressure phase variations between different forcing configurations. The 6 role of surfacepressure phase in the production of structural vibrations and HCF failures therefore remains largely unresolved. 1.7 SCOPE OF CURRENT INVESTIGATION The objective of the current research is to conduct a series of numerical simulations to examine the influence of various aerodynamicforcing and liftingsurface configurations on chordwise surfacepressure phase distributions. In particular, the influence of liftingsurface thickness, camber and angle of attack is discussed. The influences of different solidity values as well as the fundamental physics underlying chordwise surfacepressure phase distributions are also explored. Results will assist the interpretation of experimentally and computationally generated chordwise surfacepressure phase data, as well as provide fundamental results to assist future liftingsurface design efforts in resisting aeroelastic modal forcing. Finally, a comparison between the computed phase results and experimental cascade data reported by Fabian et al. and Falk et al. is presented, providing an explanation of the observed experimental phase trends. 7 2 CHAPTER 2 PREVIOUS WORK This chapter reviews previous work in the area of unsteady liftingsurface aerodynamic forcing. The purpose of this review is to place the current research in proper perspective, establishing its motivation, importance and its potential contribution to the aerodynamic forcedresponse community. 2.1 LITERARY REVIEW Many researchers have reported surfacepressure phase distributions over lifting surfaces under a variety of forcing conditions. Numerous experimental investigations have been performed on isolated airfoils, cascades and jet engine blade/vane rows in an attempt to understand the detailed fluidstructure interactions related to liftingsurface forced response. 2.1.1 SingleAirfoil Investigations One of the early works in the field of liftingsurface forced response was conducted by Sears [1938, 1941], who examined unsteady aerodynamic forcing of rigid infinitespan flat plates by convecting chordwise sinusoidal gusts in an incompressible inviscid flow. Sears derived a relationship for the chordwise, unsteady, nondimensional, differentialpressure amplitude distribution on such lifting surfaces as 8 ( ) i t p S k e x c C π x c ω 1 / 2 1 / + − Δ ′ = (2.1) where ( ) [ ] 2 (2) 1 (2) 0 k H iH S k − = π (2.2) Unsteady surfacepressure time series predicted by Eq. (2.1) are found to be synchronous along the chord, indicating no chordwise time delay between pressureamplitude peaks/troughs. Thus, Eq. 2.1 predicts surfacepressure chordwise phase to be independent of gust propagation speed. The results of Sears suggest instantaneous surfacepressure response along the entire liftingsurface chord to the convective sinusoidal gusts. By extension, Sears predicts corresponding surfacepressure phase distributions versus chord would show zero phase slope; i.e., require the entire surface to respond instantaneously to all forcing disturbances. Although not noted by Sears, his results suggest a chordwisevarying phase distribution should correlate with finite surfacepressure propagation speeds over a lifting surface. Thus, the slope of a surfacepressure phase distribution along the chord relates to forcingdisturbance propagation speed. Faster disturbance propagation corresponds to less phase change, or lower slope, with chord and viceversa. 2.1.2 Turbomachine and Cascade Investigations Turbomachine unsteady forcing phenomena were also experimentally investigated to verify existing analytical results, as well as identify new flow physics. For example, Fleeter et. al. conducted an experimental investigation to determine rotorinduced unsteady pressure distributions on downstream stator vanes [Fleeter et. al., 1978]. This was accomplished in the Detroit Diesel Allison (DDA) largescale, low 9 speed, singlestage research compressor. This investigation studied the effects of reduced frequency and incidence angle on statorvane surfacepressure distribution. Measurements were collected with embedded pressure transducers mounted axially along both suction and pressure surfaces of the stator vanes [Fleeter et al., 1978]. Resultant unsteady surfacepressure amplitudes were shown to compare reasonably well with existing analytical results for all reduced frequency values at small incidence values. However, at large negative incidence angles, experimental data correlations with predictions were very poor. This was attributed to convectedwake phenomena not modeled in the analysis. Corresponding surfacepressure phase results were found to be ambiguous (i.e., no clear chordwise trend) leading to the conclusion that rotor wakes travel differently over airfoil suction and pressure surfaces, depending upon their harmonic frequency. Fleeter et al. also observed similar wakepropagation behavior in a later study [Fleeter et al., 1980]. In this investigation, rotorinduced surfacepressure data acquired on cambered statorvanes were compared to flatplate, vortical gust code results to determine the effect of airfoil camber on airfoil unsteady lift. Unsteady surfacepressure amplitudes on cambered stator vanes exhibited amplification at the leading edge decaying in the chordwise direction. As such, amplitude data correlated very well with theoretical predictions, at both zero and negative incidence angles. However, phase data for the cambered statorvanes again exhibited ambiguous characteristics, showing very poor correlation with the theoretical predictions. In particular, phase data were found to correlate with the theoretical predictions only in the leading edge region, varying linearly from the predicted results downstream. Here again, Fleeter et al. attributed this poor 10 phase behavior to unknown convectedwave phenomenon appearing along the camberedairfoil vane row. The convectedwave phenomenon first appeared near the rear of the vane, moving forward as incidence angle decreased. This phenomenon also exhibited different behavior on the vane pressure and suction surfaces. Unfortunately, the observed convectedwave phenomenon was not modeled or investigated fully, in the presented analysis. Liftingsurface chordwise relativephase information was experimentally examined by Fabian et al in a linear transonic cascade [Fabian et. al., 1996]. The cascade consisted of six productionhardware stator vanes collected from the fan stage in a F109 turbofan engine. Stator vanes were placed in a 4 in. × 4 in. crosssection cascade wind tunnel, creating five twodimensional passages; flow turning through the passages induced vane mean aerodynamic loading. Vane unsteady forcing was accomplished via a row of five circular cylinders placed 0.8 vane chords upstream or downstream of the vane row; allowing forward or rearward aerodynamic forcing, respectively. In the rearwardforcing configuration, upstreampropagating potential disturbances, created by shed bound circulation on the downstream cylinders, forced the vane row. Unsteady, phaselocked, surfacepressure measurements were collected on the vanes at various freestream Mach numbers not exceeding 0.59. Surfacepressure results indicated rearward forcing to elicit nearly linear phase behavior with chord, as illustrated in Figure 21, for both first and second harmonic surfacepressure frequencies. Note the slope of the phase data is positive with respect to chord, indicating an upstreampropagating forcing disturbance, as predicted by the superimposed line showing analytical model results. The linear nature of data also 11 suggest phase to be independent of aerodynamic loading, and have a constant propagation speed [Fabian et al., 1996]. Thus, forcing disturbance and surfacepressure propagation correlate in the rearwardforced cascade configuration. Figure 21 Relative Phase Variation with Chord: RearwardForced Cascade of Fabian et al. [1996]. In addition to rearward forcing, forwardforcing investigations conducted by Fabian et al. also examined chordwise phase distributions [Fabian et al., 2001]. By placing forcing cylinders upstream of the vane row, the cascade configuration allowed convective cylinder wakes to propagate across and force the cascade row. Phaselocked surfacepressure measurements, analogous to those collected during rearward forcing, produced chordwise phase distributions such as Figure 22. Note that unlike the rearwardforcing phase data, the forwardforcing data are not linear, have no clear slope or pattern, and do not agree with the analytical phase model results (solid lines). As such, Fabian et al. termed this data behavior as “ambiguous”, attributing the chaotic nature of the data to cylinderwake interaction with downstreampropagating potential disturbances also emanating from the forcing cylinders. The ambiguous phase results of Figure 22 12 raised the question as to whether such phase ambiguity might also be expected in a rotating machine, as opposed to the linear cascade setup. Figure 22 Relative Phase Variation with Chord: ForwardForced Cascade of Fabian et al. [1996]. In order to answer the question raised by the cascade study [Fabian et al., 1996], phaselocked, unsteady, surfacepressure measurements were performed across swept stator vanes in a running F109 turbofan engine [Falk et al., 1997]. Unlike most turbofan engines, the F109 has only a single stage of axial compression. Therefore, no obstructions exist downstream of the stator vanes that might produce upstreampropagating disturbances. The only disturbances propagating across the stator vanes develop from upstream. This forcing configuration is analogous to the forwardforced cascade of Fabian et al. Results from the investigation by Falk again showed liftingsurface phase information to not display a definite propagation direction at either the convected or acoustic disturbance speed along the vane. In fact, the phase ambiguity measured in the F109 supported the previous arguments of Fabian et al. [1996] and provided presumptive evidence of a strong interaction between downstreampropagating potential and 13 convected disturbances. Of these disturbances, one was argued to be the vortical wakes created by the fan blades at the bladepassing frequency, propagating at the local convection velocity. The other disturbance was argued to be a potential disturbance also created at the bladepassing frequency, but propagating downstream at acoustic speeds. Frey and Fleeter [Frey, 1998], performed experiments to investigate and quantify gustgenerated unsteady aerodynamic response of stator blades. In their experiment, 2/revolution unsteady aerodynamic forcing functions were introduced to a first stage rotorblade row, these forcing disturbances having significant vortical and potential components. Obtained results showed unsteady pressure amplitudes to reach a high value near the leading edge, decaying by 75% at midchord and then increasing slightly in the aft chord; such amplitude behavior was also observed by Fabian et al. Results again suggested a strong interaction between vortical and potential disturbances, comparing well with proprietary codes named LINFLO and LINSUB. Unfortunately, no explanation of chordwise phase distributions was given. 2.2 UNRESOLVED ISSUES The aforementioned investigations provide essential improvements toward understanding liftingsurface surfacepressure distributions under a variety of forcing conditions; however, the unexplained behavior of reported phase data remains an open topic. Upon review, a consistent explanation for observed surfacepressure phase variations along examined lifting surfaces under various aerodynamic forcing configurations does not exist. As such, researchers and designers working in the forcedresponse area are largely uneducated about the role of surfacepressure phase distributions in the production of liftingsurface structural vibrations and HCF failures. 14 Therefore, what follows is a brief analytical example illustrating the possible influence of chordwise phase distributions on timeresolved surfacepressure response and modal forcing. It is intended that this example underscores, or even introduces, the importance of accurately considering phase distributions in liftingsurface response models, while also reinforcing the need for continued examination of liftingsurface phase distributions under various, even simplified, forcing conditions. Consider an infinitespan flatplate lifting surface having a chord extending from –1.0 < x/c < 1.0, where the surface is placed in an incompressible inviscid fluid. If the liftingsurface is subjected to a propagating sinusoidal disturbance, and it is assumed that the ingested disturbance is not distorted by interaction with the flat plate, the unsteady differentialpressure distribution along the lifting surface can be described in terms of a periodic function having some amplitude, frequency and phase. This function is given by p' (x / c,t) a(k) U C ' (x / c)sin[ t (x / c)] p Δ = ρ Δ ω +φ ∞ ∞ (2.3) where ΔCṕ is predicted by Sears in Eq. (2.1). Using Eqs. (2.1) and (2.3), nondimensional unsteady surfacepressure time series at various x/c locations are computed, as shown in Figures 23 and 24, for two separate phase distributions. In the nonvarying phase case (i.e.,φ (x / c) = 0) of Figure 23, each x/c timeseries is found to be in phase, reflecting instantaneous chordwise response to each propagating disturbance. This corresponds to the liftingsurface phase distribution of Sears. Conversely, for the chordwisevarying phase case (i.e.,φ (x / c) ≠ 0) of Figure 24, each x/c location responds sequentially to the propagating disturbances. Note that the phase distributions of Figures 23 and 24 assume a constant phase change, or constant disturbance propagation speed, between each chordwise liftingsurface location. 15 Time, t (sec) Differential SurfacePressure, P' 0 50 100 150 200 250 180.0 150.0 120.0 90.0 60.0 30.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 x/c = 0.80 x/c = 0.60 x/c = 0.40 x/c = 0.20 x/c =0.0 x/c =0.20 x/c =0.40 x/c =0.60 x/c =0.80 Figure 23 Unsteady Differential Surface Pressure Time Series: (φ (x / c) = 0). Time, t (sec) Differential SurfacePressure, P' 0 50 100 150 200 250 180.0 150.0 120.0 90.0 60.0 30.0 0.0 30.0 60.0 90.0 120.0 150.0 180.0 x/c = 0.80 x/c = 0.60 x/c = 0.40 x/c = 0.20 x/c =0.0 x/c =0.20 x/c =0.40 x/c =0.60 x/c =0.80 Figure 24 Unsteady Differential Surface Pressure Time Series: (φ (x / c) ≠ 0). Figures 25 and 26 exhibit corresponding surfacepressure distributions (computed from Figures 23 and 24) at selected times of t = 1, 45, 89, and 130 μs, for the two examined phase distributions. For the nonvarying phase case of Figure 25, unsteady surfacepressure alternates continuously from positive to negative pressure due to the assumed sinusoidal nature of the forcing disturbance. At no time during the oscillation cycle, however, does the differential surfacepressure have both negative and positive chordwise components. In contrast, unsteady surface pressures corresponding to the chordwisevarying phase case, as shown in Figure 26, have multiple pressurenode locations, where these node locations change chordwise position with time. Furthermore, the chordwisevarying phase case alters the shape of the surfacepressure distribution, particularly along the forward half of the lifting surface. 16 Chord, x/c Differential SurfacePressure, p' 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 150.0 125.0 100.0 75.0 50.0 25.0 0.0 25.0 50.0 75.0 100.0 125.0 150.0 t = 1 t = 45 t = 90 t = 135 Figure 25 Unsteady Differential Surface Pressure Variation with Chord: (φ (x / c) = 0). Chord, x/c Differential SurfacePressure, p' 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 175.0 150.0 125.0 100.0 75.0 50.0 25.0 0.0 25.0 50.0 75.0 100.0 125.0 150.0 175.0 t = 1 t = 45 t = 90 t = 135 Node Locations Figure 26 Unsteady Differential Surface Pressure Variation with Chord:(φ (x / c) ≠ 0). Chordwise integration of the surfacepressures in Figures 25 and 26 provides a maximum unsteady force, or lift, of 140 and 80.0, for the φ (x / c) = 0 and (φ (x / c) ≠ 0) cases, respectively. Given these lift differences, it may be inferred that chordwisevarying phase may be desirable in terms of reducing unsteady aerodynamic loading. However, the integration process ignores modal forcing of the lifting surface. The unsteady generalized force on a particular structural mode, m, can be computed through the integral of the dotproduct between the examined mode shape and surfacepressure distribution over the liftingsurface chord. Thus, for the current example, the generalized force on a particular structural mode can be written as ( ) [ ( / , ) ].[ ( / )] ( / ) 1 1 f t p' x c t n x c d x c m m ψ r ) ∫ + − = Δ (2.4) A simply supported, infinitespan, twodimensional lifting surface has an infinite number of mode shapes grouped into two families. Rigidbody mode shapes correspond to plunging and pitching oscillations of the lifting surface, while flexiblebody mode 17 shapes correspond to elastic structural deflections. Several flexiblebody mode shapes are illustrated in Figure 27, for the first three modes, along with the first rigidbody mode. Chord, x/c Normalized Maximum Deflection 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.4 0.6 0.8 1.0 2.0 1.5 1.0 0.5 0.0 0.5 1.0 1.5 2.0 Rigid Mode 1 (Lift) Flex Mode 1 Flex Mode 2 Flex Mode 3 Figure 27 RigidBody and FlexibleBody Mode Shapes for a Simply Supported, InfiniteSpan, Two Dimensional Lifting Surface. Accurate predictions of timeresolved generalized forces on each structural mode are important to avoid structural vibrations and HCF failure. To emphasize this fact, the surfacepressure distributions of Figures 25 and 26 are input into Eq. (2.4) along with the mode shapes of Figure 27, producing a generalized force on each examined structural mode. The maximum generalized forces obtained through this exercise are presented in Figure 28 and plotted versus mode number for the two separate chordwise phase distributions. Note that “mode 0” in Figure 28 corresponds to the first rigidbody structural mode, or liftmode, while the higher modes correspond to the first, second and third flexiblebody modes, respectively. Examining Figure 28, the nonvarying phase distribution shows a canonical decay in force with increasing mode number. In contrast, the chordwisevarying phase distribution exhibits decreased force in the lowerorder 18 modes and significantly amplified force in higherorder modes. Therefore, modal force is found to be a function of chordwise phase distribution. Mode # Maximum Generalized Force, fmax 0.0 1.0 2.0 3.0 4.0 0.0 20.0 40.0 60.0 80.0 100.0 120.0 140.0 160.0 φ = 0 φ = 3x/c Figure 28 Maximum Generalized forces on Structural Modes for Various Phase Distributions. The above unsteadyforcing example emphasizes several facts. First, chordwisevarying phase may alter timeresolved surfacepressure amplitude distributions along a liftingsurface chord. Second, chordwisevarying phase distributions also produce surfacepressure node locations, where these node locations may change position with time. Third, generalized modal forces are altered by chordwise phase distribution. In particular, higherorder modal forces may be amplified by chordwisevarying phase distributions, generating greater mechanical stresses at spatial and temporal frequencies not predicted for the nonvarying phase case. Such possible variations in liftingsurface modal forcing may, if inaccurately predicted, lead to decreased fatigue life for the lifting structure. In all, unsteady surfacepressure phase distributions are clearly important to the unsteady forcing problem. This is particularly evident when one considers that some 19 current unsteadyforcing predictive tools employ the results of Sears as a basis of their predictions. Fortunately, an increasing number of forcedresponse predictive tools are not based on the results of Sears, opting rather for direct numerical simulation of the governing fluid dynamic equations. Validation of these computational tools has proven to be laborious and heavily dependent on the availability of properly posed benchmark data. As such, few comparisons between computed an experimentally determined surfacepressure phase distributions have been made. This is not to say that benchmark phase data are not available, but rather that the significance of the data is not well understood or properly examined. In fact, what is intriguing about the previously reported phase data in this chapter is not their lack of inclusion in the open literature, but the almost complete disregard as to their importance and correlation with computed/measured trends. Much of the available computational/experimental phase data do not correspond to the constant disturbancespeed assumptions made in the above example; in fact, certain data sets show almost no discernable trend with chord. Therefore, while the above example emphasizes the importance of considering chordwisevarying phase for a constant propagation speed, the validity of assuming a constant disturbance propagation speed is unclear. Moreover, the effects of phase deviations from the assumed constantspeed disturbance phase on liftingsurface response are unknown. 2.3 CURRENT RESEARCH OBJECTIVE Given the continued need for examining the influence of liftingsurface chordwise phase distribution on surfacepressure response, the current research presents a fundamental study of surfacepressure phase. In particular, twodimensional, timeaccurate, RANS simulations are performed to examine the fundamental physics leading 20 to surfacepressure phase. Attempts are made to reveal the essential dependencies of that phase on forcing configuration. Simulations are performed for a variety of lifting surface geometries and forcing conditions utilizing the commercially available CFD algorithm, Fluent (v. 6.0). 21 3 CHAPTER 3 COMPUTATIONAL METHODOLOGY AND SETUP This chapter discusses the computational methodology and setup employed in the present research. A description of the airfoil and statorvane cascade geometries, as well as associated boundary conditions, is given. Indepth discussions regarding FLUENT and its companion mesh generation software, GAMBIT, are also provided along with a development of the UDF (userdefined function) generating the airfoil forcing function. 3.1 AIRFOIL GEOMETRY AND BOUNDARY CONDITIONS Surfacepressure phase data on forwardforced lifting surfaces are examined using the commercially available CFD algorithm, FLUENT (v. 6.0). Simulations were performed with four symmetric NACA airfoil profiles of 10, 12, 15 and 20% thickness (relative to chord), two cambered airfoils of 2 and 6% camber (relative to chord), two meanflow attack angles of 5 and 10 degrees and two forcingdisturbance frequencies of 150 Hz and 300 Hz. In all simulations, periodic boundary conditions were enforced on the upper and lower computational boundaries located about an otherwise isolated lifting surface, as illustrated in Figure 31. These periodic boundaries simulate the influence of neighboring surfaces, or a series of airfoils in cascade; an airfoil cascade configuration was selected for comparison with previous experimental configurations. A cascade solidity of 4.0, representing weak surfacetosurface pitchwise aerodynamic coupling as 22 compared to modern cascaded blade rows, was selected for the majority of simulations; however, other solidities equaling 2.0 and 6.0 were also briefly investigated. Pressureinlet and pressureoutlet boundary conditions were set for the computational inlet and outlet boundaries, respectively (see Figure 31). Figure 31 AirfoilCascade Computational Boundaries. 3.2 STATORVANE GEOMETRY AND BOUNDARY CONDITIONS In addition to the simplified NACA airfoilcascade configuration, a more complicated cascade configuration was also examined. This configuration employed aerodynamically loaded vanes that mimicked the twodimensional geometry of the statorvane row in the fan compression stage of a F109 turbofan engine (at 87.8% span). The high cascade solidity and double circulararc profile of the statorvane row required definition of additional geometric variables beyond the simplified NACA airfoil cascade. A complete discussion of the cascade geometry and nomenclature is provided in Appendix A. In the present investigation, a vanecentered computational mesh was 23 selected to model the statorvane cascade geometry, with the simulated vane centered in the computational domain, as shown in Figure 32. Periodic boundary conditions were enforced at mesh boundaries above and below the vane, simulating the influence of other vanes in cascade. The periodic boundaries were set at midpitch between vanes, using the vane camberline arc to define the boundary geometry. The statorvane inlet flow angle was set to be 21.9o, while the exit flow was set to be 20.6o, based on previous experimental data [Fabian, 1995]. The stator vanes possess a maximum camber and thickness of 12% and 8% relative to chord, respectively. Vane profile coordinates are listed in Appendix B, reproduced from [Fabian, 1995]. Table 3.1 provides characteristics of the modeled statorvane cascade geometry. Like the NACAairfoil cascade, pressureinlet and pressureoutlet boundary conditions were set for the statorvane computational inlet and outlet boundaries, respectively. Table 3.1 Modeled StatorVane Cascade Geometry. Parameter Value Vane Spacing, S 0.84 m Solidity, σ 1.524 Inlet Flow Angle, α1 21.9o Exit Flow Angle, α2 20.6o 24 Figure 32 StatorVane Cascade Computational Boundaries. 3.3 GENERAL FLUENT SOLVER DESCRIPTION FLUENT is a stateoftheart commercially available flowsolver with capability for modeling unsteady, compressible, viscous flows via numerical solution of the governing fluid dynamic equations. Numerical simulation of the governing fluid dynamic equations in FLUENT is accomplished via a controlvolume based (finitevolume) discretization technique. This technique integrates governing integral equations established within discrete elements (i.e. finite volumes, or cells) of the mesh, resulting in a system of algebraic equations for dependent variables such as velocity and pressure. The discretized algebraic system is then linearized and solved numerically to yield updated variable values at each iteration/time step, using (in the present investigation) implicit linearization schemes. Solution interpolation between adjacent elementface regions is accomplished via one of several userdefined methods, including firstorder upwind, secondorder upwind and powerlaw interpolation. Firstorder or secondorder accurate spatial/temporal discretization is available in FLUENT, where secondorder accuracy is default for coupled solutions. To aid convergenence in highly nonlinear 25 problems, FLUENT allows userdefined controls over underrelaxation and courant (CFL) numbers. Numerical solutions are achieved through one of three userselected solvers, including: segregated, coupledimplicit, or coupledexplicit solvers. The segregated solver linearizes the governing equations implicitly with respect to the dependent variables, solving the resulting set of equations sequentially. Linearized momentum equations are solved individually for fluid velocity, followed by corrective step in which velocity is adjusted based on userselected velocitypressure correlations to satisfy continuity. Conversely, coupled solvers simultaneously solve the set of governing equation defining the dependent variables, where the equation set can be linearized either explicitly or implicitly. For implicit linearization, GaussSeidel solvers are employed in conjunction with an algebraic multigrid (AMG) method to solve the system(s) of equations. Conversely, with explicit linearization, dependentvariable solutions are updated at each time step using a multistep RungeKutta solver, with the additional option of employing a full approximation storage (FAS) multigrid scheme to accelerate convergence. FLUENT allows users to specify several boundary condition types. Supported inlet and outlet boundary conditions include: pressureinlet, velocityinlet, massflowinlet, inletvent, intakefan, pressureoutlet, pressurefarfield, outflow, outletvent, and exhaustfanboundaries. Similarly, wall, repeating, and pole boundary types include: wall, symmetry, periodic and axis boundaries. For the present research, twodimensional numerical simulations of the Reynoldsaveraged NavierStokes equations were accomplished in FLUENT via a finitevolume technique. A coupled solution methodology (Section 3.3.1) was selected, in which the 26 fully coupled systems of equations defining the dependent variables in each cell were discretized, linearized, and solved simultaneously at each iteration/time step. Linearized equation systems were solved using a GaussSeidel solver in conjunction with an algebraic multigrid (AMG) method. Secondorder accurate spatial and temporal discretization was employed for all simulations, with implicit linearization. Simulations were fully viscous, utilizing a standard kε turbulence model (see Section 3.3.3). 3.3.1 Coupled Solution Method The coupled solver used for the current simulations solves the governing equations of continuity, momentum, and (where appropriate) energy and species transport simultaneously (i.e., coupled together). Governing equations for additional scalars (i.e., turbulence, etc.) are solved sequentially using a segregated approach. Since the set of governing equations is nonlinear (and therefore coupled), several subiterations of the solution procedure are performed at each time step before a converged solution is obtained at that time step. Each subiteration consists of the steps illustrated in Figure 33 and outlined below: 1. Fluid properties are updated, based on the current solution. (If the calculation has just started, fluid properties are updated based on an initial solution.) 2. Continuity, momentum, and (where appropriate) energy and species equations are solved simultaneously. 3. Where appropriate, equations for scalars such as turbulence and radiation are solved using the previously updated values of the other variables. 4. A check for convergence of the equation set is made. 27 These steps are continued until the convergence criteria are met at each time step. [FLUENT, 2001] Update properties Solve continuity, momentum, energy, and species equations simultaneously. Solve turbulence and other scalar equations. Converged? Stop Figure 33 Overview of the Coupled Solution Method [Fluent, 2001]. 3.3.2 ReynoldsAveraged NavierStokes (RANS) Equations The Reynoldsaveraged NavierStokes (RANS) equations were selected to represent transport equations for ensembleaveraged, or mean, flow quantities, with all turbulence scales modeled. The approach of permitting a solution for just the meanflow variables greatly reduces the computational effort. If the mean flow is steady, the governing equations will not contain time derivatives and a steadystate solution can be obtained economically. A computational advantage is also provided in required timeaccurate simulations, as time step may be determined by global unsteadiness in the mean flow rather than turbulent unsteadiness. The RANS approach models turbulent flow 28 quantities, through such wellknown models as the SpalartAllmaras, kω, kε, and RSM models [FLUENT, 2001]. 3.3.2.1 Reynolds (Ensemble) Averaging In Reynolds averaging (i.e., ensemble averaging solution variables in the instantaneous (exact) NavierStokes equations are decomposed into mean (ensembleaveraged or timeaveraged) and fluctuating components (about the mean). For velocity components, this decomposition equals i i i u = u + u′ (3.1) Likewise, for pressure and other scalar quantities: φ =φ +φ ′ (3.2) where φ denotes a scalar such as pressure energy or species concentration. Substituting these forms of the flow variables into the instantaneous continuity and momentum equations, and taking a time (or ensemble) average (and dropping the overbar on the mean velocity, u ), yields the ensembleaveraged continuity and momentum equations. These equations can be written in Cartesian tensor form as: ( ) = 0 ∂ ∂ + ∂ ∂ i i u t x ρ ρ (3.3) ) ( ) 3 ( ) ( ) ( 2 1 1 i j j ij i j j i i j i j j i u u x x u x u x u x x u u p x u t − ′ ′ ∂ ∂ + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ − ∂ ∂ + ∂ ∂ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ ρ ρ μ δ ρ (3.4) Equations (3.3) and (3.4) are also known collectively as the RANS equations. 29 3.3.3 k – ε Turbulence Model In the present application, a standard kε turbulence model was used to simulate the effects of turbulence. The standard kε model is a semiempirical model based on model transport equations for turbulence kinetic energy, k, and it dissipation rate. The model transport equation for k was derived from an exact equation, while the model transport equation for ε was obtained using physical reasoning. Equations for k and ε are given as follows k b M k k j t j i i G G Y S x k x ku x k t + + − − + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ ρε σ μ (ρ ) (ρ ) (μ ) (3.5) and ε ε ε ε ε ε ρ ε ε σ μ ρε ρε μ S k G C G C k C x x u t x k b j t j i i + + + − + ⎥ ⎥⎦ ⎤ ⎢ ⎢⎣ ⎡ ∂ ∂ + ∂ ∂ = ∂ ∂ + ∂ ∂ 2 1 3 2 ( ) ( ) ( ) ( ) (3.6) where, 1ε C , , , 2ε C 3ε C k σ and ε σ are model constants nominally having the following values: 1ε C = 1.44 2ε C = 1.92 3ε C = 0.09 k σ = 1.0 ε σ = 1.3 The above values are default in FLUENT, and are employed in the current investigation. 30 3.3.4 FiniteVolume Discretization Methodology FLUENT uses a controlvolumebased technique to convert governing fluid dynamic equations to algebraic equations that can be solved numerically. This controlvolume technique consists of integrating the governing equations about each control volume (or cell) in the computational mesh, yielding discrete equations that conserve each quantity on a controlvolume basis. Discretization of the governing equations can be illustrated most easily by considering the steadystate conservation equation for transport of a scalar quantity,φ . This is demonstrated by the following equation written in integral form for an arbitrary control volume, V, as follows: ∫ =∫Γ ∇ + ∫ V dA dA S dV φ φ φ ρφυ r r r (3.7) The above equation is applied to each control volume, or cell, in the computational domain. Thus, discretization of Eq. (3.7) for a given control volume yields: A A S V faces N faces f n f N f f f f f φ φ Σ ρ υ φ =Σ Γ ∇φ + r r r ( ) (3.8) All governing equations solved by FLUENT take the same general form as Eq. (3.8) and therefore readily apply to unstructured meshes composed of tetrahedra, as in the current investigation. FLUENT stores discrete values of the scalar φ at the cell centers. However, face values, of f φ , are required for convection terms in Eq. 3.8 and therefore must be interpolated from adjacent cellcenter values. This is accomplished using an upwind scheme. 31 3.3.5 SecondOrder Upwind Scheme For the present investigation, a secondorder upwind scheme was selected. When secondorder solution accuracy is desired in FLUENT, flow quantities at cell faces are computed using a multidimensional linear reconstruction approach. In this approach, higherorder accuracy is achieved at cell faces through a Taylor series expansion of the cellcentered solution about the cell centroid. Thus when secondorder upwinding is selected, a face value f φ is computed using the following expression: s f φ =φ + ∇φΔr (3.9) where φ and ∇φ are the cellcentered value and its gradient, and sr Δ is the displacement vector from the cell centroid to the face centroid. This formulation requires the determination of the gradient ∇φ in each cell, where this gradient is computed using the divergence theorem, ∇ = Σ N faces f f A V r φ 1 φ~ (3.10) In Eq. (3.9) face values f φ ~ are computed by averaging φ between cells adjacent to the face in question. Finally, the gradient ∇φ is valuelimited so that no new maxima or minima are introduced in the examined cell region. 3.3.6 SecondOrder Time Discretization For transient simulations, the governing equations must be discretized in both space and time. Spatial discretization for timedependent equations is identical to the steadystate case; however, temporal discretization involves integration of every term in 32 the governing differential equations over one time step,Δt . The integration of transient terms is straightforward, as described below. A generic expression for the time evolution of a variableφ is given by: (φ ) φ F dt d = (3.11) where the function F incorporates both spatial and temporal discretization. If the time derivative is discretized to firstorder approximation, an expression for the discretized derivative may be written as ( ) 1 φ φ φ F t n n = Δ + − (3.12) While a secondorder accurate temporal discretization is given by ( ) 2 3 1 4 1 φ φ φ φ F t n n n = Δ + − + − (3.13) Once the time derivative has been discretized, a choice remains for evaluating F(φ ) . The method employed here to evaluate F(φ ) at a future time level, such as ( 1 ) 1 + + = Δ − n n n F t φ φ φ (3.14) This is referred to as “implicit” integration as in a given cell depends on both sides of Eq. (3.14), giving φ n+1 φ n+1 =φ n + ΔtF(φ n+1 ) (3.15) The implicit relation of Eq. (3.15) can be solved iteratively by initializing to and iterating the equation φ n+1 φ t ( ) 3 2 3 1 3 φ t = 4φ n − φ n−1 + ΔtF φ t (3.16) 33 for secondorder formulation, until stops changing (i.e., converges). At that point, is set equal to . The advantage of a fully implicit scheme is its unconditional stability with respect to timestep size. φ t φ n+1 φ t 3.3.7 Linearization Methodology In the coupledsolution method the discrete, nonlinear governing equations are linearized to produce a system of equations for the dependent variables in every computational cell. The resultant linear system is then solved to yield an updated flowfield solution. The manner in which the governing equations are linearized may take an implicit or explicit form with respect to the dependent variable (or set of variables) of interest. For the present analysis a coupledsolution methodology with implicit linearization was used. This results in a system of linear equations with N equations for each cell in the domain, where N is the number of coupled equations in the set. Because there are N equations per cell, this is sometimes called a block system of equations. A pointimplicit (i.e., block GaussSeidel) linear equation solver is used in conjunction with an algebraic multigrid (AMG) method to solve the resultant block system of equations for all N dependent variables in each cell. For example, linearization of the coupled continuity, x, y, zmomentum, and energy equation set will produce a system of equations in which p, u, v, w, and T are unknowns. Simultaneous solution of this equation system (using the block AMG solver) yields at once updated pressure, velocity, and temperature fields. 34 3.3.8 Periodic Boundary Conditions Periodic boundary conditions are advantageous when the physical geometry of interest and the expected flow solution has a spatially periodic nature. The assumption of spatial flow periodicity implies the velocity components repeat themselves in space as follows: u(r ) = u(r + L) = u(r + 2L) = ........ r r r r r v(r ) = v(r + L) = v(r + 2L) = ........ r r r r r (3.17) where rr is the position vector and L r is the periodic length vector of the domain considered. For viscous flows, the spatially distributed pressure field may not be periodic in the sense of the above Eq. (3.18). Instead, pressure drop between modules maybe periodic in space, giving Δp = p(r ) − p(r + L) = p(r + L) − p(r + 2L) = ........ r r r r r r r (3.18) If the coupled solver is selected in FLUENT, Δp maybe specified as a constant value at any periodic boundary. 3.3.9 Standard FLUENT Inlet/Outlet/Wall Boundary Conditions 3.3.9.1 PressureInlet Boundary Conditions Pressureinlet boundary conditions define fluid pressure at flow inlets, along with all other scalar properties of the flow. They are suitable for both incompressible and compressible flow calculations [Fluent, 2001]. Pressureinlet boundary conditions are typically specified in FLUENT when inlet pressure is known but flow rate and/or velocity is not known. When flow enters through a pressureinlet boundary, FLUENT enforces the input boundary pressure as the fluid totalpressure at the inlet plane, . In 0 p 35 steadystate incompressible flow, the inlet totalpressure and static pressure are related at the inlet velocity via Bernoulli's equation: s p 2 0 2 = + 1 ρν s p p (3.19) In compressible flow, the totalpressure and static pressure are related by 2 1 0 ) 2 (1 1 − − = + γ γ γ p p M s (3.20) 3.3.9.2 PressureOutlet Boundary Condition: As specified earlier, a pressureoutlet boundary condition was selected at the computational outlet plane for all simulations. Pressureoutlet boundary conditions in FLUENT require the specification of a static (gauge) pressure at the outlet boundary. For a subsonic flow, FLUENT enforces the input boundary pressure as the fluid static pressure at the outlet plane, and extrapolates all other flow conditions to the outlet from the interior of the domain. 3.3.10 Operating Pressure In the present investigation, the operating pressure condition in FLUENT was set to 0 Pa, allowing all pressure calculations to be treated as gauge pressures. Setting the computational operating pressure equal to 0 Pa is significant as it reduces roundoff error problems. In lowMach number flows, overall pressure changes are typically small compared to the atmospheric static pressure, and therefore can be significantly affected by numerical roundoff error. Moreover, FLUENT always uses gauge pressure for all its calculations, as such setting the operating pressure equal to zero makes gauge and absolute pressures equivalent. 36 3.4 UDF DESCRIPTION Aerodynamic forcing of the examined lifting surfaces was achieved by modeling waveform characteristics similar to a passing wake, as might be found in a rotorstator compression stage. In doing so, a userdefined function (UDF) was written in FLUENT. A UDF can be dynamically loaded into the FLUENT solver to enhance standard features of the available code. Such functions are written in the C programming language, allowing standard C library functions to be used, as well as predefined macros (provided by Fluent Inc.). A UDF can be implemented as an interpreted or compiled function. An interpreted UDF is read in and interpreted at run time. Alternatively, a compiled UDF is grouped into shared libraries when it is compiled and linked with the standard FLUENT executable. An interpreted UDF is simple to use but imparts coding and speed limitations. A compiled UDF executes faster and has no coding limitations, but requires more effort to develop and implement. The current investigation employs a compiled UDF. 3.4.1 Development Logic/Procedure The developed inlet UDF established an inlet totalpressure profile, simulating the movement of a rotorblade wake passing across the computational inlet. The following steps were employed to code and implement the UDF in FLUENT for the present investigation. 1. Define the wake model. 2. Create a C source code file 3. Start FLUENT and setup the case file. 37 4. Compile the C source code 5. Activate the UDF in FLUENT. The first step to generate and implement the inlet UDF was to develop an analytical model of a typical rotorblade wake; this wake representing the aerodynamic forcing function for the examined lifting surfaces. Experimental lowspeed compressor investigations by Reynolds et. al (1979) showed that rotorblade wake profiles approximately follow a Gaussian profile. As such, a simple aerodynamic forcingfunction model was developed using a Gauss function. This model is characterized by the wake centerline defect, , and the semiwake width, dc W δ (see Figure 34). An equation for totalpressure deficit produced by a rotorblade wake function was defined as [ ( t )2 ] c y ti dc P W e − − = δ 3.21 Pitch, y/c Total Pressure, Pt (Pa) 4.0 2.0 0.0 2.0 4.0 101265 101270 101275 101280 101285 101290 101295 101300 101305 101310 101315 101320 101325 101330 Wake Centerline Wdc δ Figure 34 Rotor Wake Characteristics. 38 Values for the wake centerline defect, , and the semiwake width, dc W δ , were optimized to obtain the waveform shown in Figure 34. Having developed the aerodynamic forcingfunction model, a C source code was written. This source code incorporated standard C library functions as well as predefined macros (provided by FLUENT). These include the DEFINE_PROFILE macro, the F_CENTROID macro, begin_f_loop macro and F_PROFILE macro. The UDF (source code) written was compiled in FLUENT as an interpreted UDF and applied to the pressureinlet boundary. A complete description of the developed UDF with a detailed discussion of each macro is provided in Appendix C. 3.5 COMPUTATIONAL GRID DESCRIPTION The FLUENT software package contains an available preprocessor for geometry modeling and mesh generation known as GAMBIT. As such, the computational grids for the present investigation were created with GAMBIT grid generation software [FLUENT, 2001]. 3.5.1 Gambit Grid Generation Software Mesh generation for FLUENT simulations may be provided through its companion software, GAMBIT; however, select thirdparty software is also supported. GAMBIT is a graphically oriented software package mimicking many computer aided drawing (CAD) programs. GAMBIT allows users to graphically reproduce complex physical geometries of interest, and mesh these geometries using several userselected meshgeneration algorithms. GAMBIT geometries are created by manipulating components such as edges, faces, and volumes that represent the physical system being 39 considered. Vertex coordinates for these components can be directly generated in GAMBIT or imported from existing data files. Created components can be rotated, translated, copied, split, united and subtracted from one another, to name a few operations, as is typical of most CAD software. GAMBIT geometries can also be created in different coordinate systems, or multiple local coordinate systems, depending on the particular application. The software provides complete meshing flexibility in both two and three dimensions, featuring algorithms for structured quadrilateral and hexahedral meshes, unstructured triangular, tetrahedral, pyramid and wedge meshes, and mixed structuredunstructured meshes. GAMBIT also supports boundarylayer subroutines to construct structured meshes in close proximity to geometric walls, providing increased grid resolution in viscous boundary layers. Socalled “boundarylayer meshes” can be coupled with unstructured meshes in the adjacent potential flow, forming a hybrid mesh. Mesh density and distribution are easily controlled via direct input of nodal positions, or through appropriate sizing functions that concentrate mesh points around userdefined regions. User input during mesh construction allows highdensity elements to be concentrated in areas of highflow gradients. 3.5.2 Grid Methodology 3.5.2.1 AirfoilCascade Mesh The methodology followed to create meshes for the parametric NACA airfoilcascade analysis was as follows. Airfoil geometry (obtained from XFOIL) was read into GAMBIT as vertex data, and interpolation was performed between vertex data to create the basic airfoil shape. Once airfoil shape was defined, inlet and outlet boundary planes 40 were created. The inlet boundary was located one airfoil chord forward of the airfoil leading edge, while the outlet boundary was placed 16 chords downstream of the airfoil trailing edge. In order to reduce computational expense, only one airfoil was considered with periodic boundary conditions simulating the influence of neighboring cascaded airfoils. A twodimensional hybrid mesh, consisting of a structured Hgrid immediately surrounding the airfoil, and an unstructured triangular mesh representing the potential region outside of the airfoil boundary layer, was created. The resulting grid system, for the NACA 0012 airfoil profile, is represented in Figure 35. A detailed description of the boundarylayer mesh and its characteristics is provided in Section 3.5.2.1.1. Grid sizes employed for the airfoilcascade meshes are provided in Section 3.5.2.1.2. 3.5.2.1.1 Grid Density and Grid Distribution Grid node distributions employed for the airfoilcascade parametric analysis are outlined in Table 3.2. The total number of mesh elements for the airfoilcascade mesh equaled 26446, where 24246 were triangular mesh elements and 2200 were quadrahedral mesh elements. The TriPave meshing scheme provided by GAMBIT was used to create the triangular mesh elements. Table 3.2 Grid Distribution (AirfoilCascade Mesh). Surface # of Nodes Inlet 100 Outlet 25 Periodic Boundary 150 Airfoil Suction Side 55 Airfoil Pressure Side 55 41 Figure 35 AirfoilCascade Mesh (NACA 0012 airfoil) 3.5.2.1.2 Boundary Layer (ViscousFlow Region) Grid As mentioned in Section 3.5.2.1, in the immediate region surrounding the airfoil a structured boundarylayer grid was generated. Boundarylayer grids in FLUENT define mesh node spacing in regions immediately adjacent to an edge and/or face. These grids control mesh density and location in a specific region of interest, here the upper (suction) and lower (pressure) surfaces of the airfoil. In the case of the NACA airfoilcascade simulations, the first boundarylayer nodes were positioned 0.12% (relative to chord) from the airfoil surface with a growth rate of 1.1, with the total number of boundarylayer mesh levels equaling 20. Table 3.3 outlines the boundarylayer mesh characteristics employed for the NACA airfoilcascade simulations. Figure 36 depicts the structured boundarylayer mesh surrounding the NACA 0012 airfoil. 42 Table 3.3 BoundaryLayer Mesh Characteristics (AirfoilCascade Mesh). Parameter Value Algorithm Uniform First Row Height 0.12% (relative to chord) Growth Factor 1.1 No. of Rows 20 Total Depth 0.06873 Internal Continuity On Wedge Corner Shape On Transition Pattern 1:1 x(m) y(m) 0 0.2 0.4 0.6 0.8 1 0.4 0.2 0 0.2 0.4 Figure 36 Structured BoundaryLayer Mesh Surrounding NACA 0012 Airfoil Geometry. 43 3.5.2.2 StatorVane Cascade Mesh The steps followed to create the statorvane cascade mesh are similar to those followed for the airfoilcascade mesh. First, statorvane coordinates (reproduced from Fabian [1995]) were read into GAMBIT as vertex data and interpolation was performed between vertices to create the statorvane shape. Second, inlet and outlet boundary planes were created. Similar to the airfoilcascade mesh, the inlet boundary was located one statorvane chord forward of the vane leading edge, while the outlet boundary was placed 16 chords downstream of the vane trailing edge. The inlet and outlet boundaries were set exactly perpendicular to the xaxis with an inlet flow angle of 21.9o and an outlet flow angle of 20.6o. In order to reduce computational expense, only one statorvane was considered with periodic boundary conditions simulating the influence of neighboring cascaded statorvanes. A structured Hgrid mesh was created for the entire statorvane computational domain. 3.5.2.2.1 Grid Density and Grid Distribution Grid node distributions employed for the statorvane cascade mesh are outlined in Table 3.4. A successive ratio of 0.96 was used for periodic boundary edgenode grading, while biexponent type grading, with a ratio of 0.75, was selected for the statorvane upper and lower surfaces. The total number of mesh elements for the statorvane cascade mesh equaled 28000, all of which were quadrahedral mesh elements. The Map meshing scheme provided by GAMBIT was used to create the quadrahedral mesh elements. Figures 37 and 38 display the statorvane cascade geometry and the modeled statorvane mesh, respectively. 44 Table 3.4 StatorVane Cascade Mesh Node Distribution. Surface # of Nodes Inlet 80 Outlet 80 Periodic Boundary 200 Airfoil Suction Side 75 Airfoil Pressure Side 75 x (m) y (m) 0 2 4 1 0 1 2 3 4 Figure 37 StatorVane Cascade. 45 x (m) y (m) 0 1 2 3 1.5 1 0.5 0 0.5 1 1.5 Figure 38 Modeled StatorVane Cascade Mesh. 3.6 Computational Setup Table 3.5 summarizes the common inputs for all FLUENT simulations. All unsteady simulations for the airfoil parametric analysis were performed using 100 time steps per disturbance forcing period, corresponding to a time step of 5x103 s. Local convergence was achieved at each time step, with as many as 35 subiterations per time step. Global convergence for all unsteady simulations was based on lift periodicity; convergence was achieved when lift periodicity reached less than 0.1% variability between forcing periods. Figure 39 illustrates liftconvergence history for the NACA 0012 airfoil; this history is representative of all simulations conducted in the present research. 46 Table 3.5 FLUENT Configuration for Numerical Simulations. Parameter Value Fluent Solver CoupledImplicit Discretization Scheme SecondOrder (Momentum and Viscosity) Material Properties Standard Day Air as an IdealGas Operating Pressure 0 Pa Viscosity Model Standard kε Turbulence Model Inlet Boundary Condition Pressure Inlet Outlet Boundary Condition Pressure Outlet Periodic Boundary Condition Periodic Time, t (sec) Coefficient of Lift, Cl 0.0 1.0 2.0 3.0 4.0 5.0 6.0 0.04 0.02 0 0.02 0.04 0.06 0.08 0.1 0.12 Figure 39 Lift Coefficient Time History Showing Convergence: NACA 0012 Profile. 47 3.6.1 Reference Values During calculations, FLUENT employs nondimensionalized variables. All variables are nondimensionalized using the reference values listed in Table 3.6. Table 3.6 FLUENT Reference Values Parameter Value Reference Density 1.225 kg/m3 Reference Length 1 m Reference Pressure 101325 Pa Reference Temperature 300 K Reference Velocity 10.0001 m/s Reference Viscosity 1.7894e5 kg/ms Ratio of Specific Heats 1.40 Operating Pressure 0 Pa 3.6.2 Boundary Conditions Pressureinlet and pressureoutlet boundary conditions were employed for the computational inlet and outlet, respectively, with inlet pressure and temperature boundary conditions being set to standardday atmospheric values for the steadystate simulation. However, for the unsteady simulations a UDF was employed at the pressureinlet boundary, mimicking the convecting wake from a rotating rotor blade. Standard wall boundary conditions were used for airfoil/statorvane suction and pressure surfaces. A detailed description of the UDF is presented in Section 3.4. Timeaveraged freestream velocity was set to 10 m/s for all NACA airfoilcascade simulations, while timeaveraged freestream velocity for the statorvane cascade simulations was set to 200 m/s. Table 3.7 48 summarizes the inlet and outlet boundary conditions employed for all FLUENT simulations performed in this research. Table 3.7 FLUENT Boundary Conditions Parameter Value Inlet Boundary Conditions Total Temperature 300 K Direction Specification Normal to Boundary Turbulence Intensity (%) 1 Turbulence Viscosity Ratio 1 Outlet Boundary Conditions Backflow Total Temperature 300 K Turbulence Intensity (%) 1 Turbulence Viscosity Ratio 1 3.7 Grid Independence Gird independence was examined by obtaining NACA 0012 steadystate solutions for higher density (fine) mesh, as compared to the nominal mesh of Figure 35. The fine mesh contained 42018 elements, where 37618 were triangular mesh elements and 4400 were quadrahedral mesh elements. Pressure coefficient results from the NACA 0012 airfoil upper and lower surfaces were compared for each mesh. Figure 310 shows steadystate pressure coefficient data along the NACA 0012 airfoil upper and lower surfaces for the nominal and fine meshes. As seen in the figure, the obtained pressure data compare favorably well, showing minimal discrepancies. As such, the nominal grid is deemed to be of sufficient density for the present investigation. 49 Chord, x/c Coefficient of Pressure, Cp 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.4 1.2 1 0.8 0.6 0.4 0.2 0 Nominal Mesh Fine Mesh Figure 310 SteadyState Pressure Coefficient Data: NACA 0012. Similar grid independence studies were performed for the statorvane cascade analysis. In this case, the high density (fine) mesh contained 38000 quadrahedral elements as compared to the 28000 quadrahedral elements contained by the nominal mesh shown in Figure 38. Figure 311 shows steadystate pressure coefficient data along the statorvane upper and lower surfaces for the nominal and fine meshes. Here again, the pressure data obtained compare favorably showing negligible discrepancies. In particular, along the statorvane upper surface small differences in pressure coefficient values are seen near the midchord region. However, near the leading and trailing edges the pressure data for the two meshes match perfectly well. Given the agreement of data and the additional CPU time required per computation by the fine mesh as compared to the nominal mesh, the nominal mesh is deemed to be of sufficient density for the present investigation. 50 Chord, x/c Coefficient of Pressure, Cp 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 2.0 1.8 1.6 1.4 1.2 1.0 0.8 0.6 0.4 0.2 0.0 Fine Mesh  Upper Nominal Mesh  Upper Fine Mesh  Lower Nominal Mesh  Lower Figure 311 Pressure Coefficient Data: StatorVane. 51 4 Chapter 4 SOLUTION VALIDATION This chapter describes the timeaveraged staticpressure results for the baseline NACA 0012 airfoil case. The timeaveraged flowfield represents an important aspect of the overall unsteadyforcing simulations for two reasons. First, unsteady results develop by subtracting timeaveraged parameters from the instantaneous solution at each time step; thus, giving parameter perturbations about the timeaveraged field. Second, the timeaveraged flowfield facilitates examination to determine the accuracy and applicability of the FLUENT simulation results. Note that the timeaveraged computational results described herein are primarily compared with experimental data from Abbott and Von Doenhoff [1959]. 4.1 TIMEAVERAGED METHODOLOGY Timeaveraged parameter distributions for each simulation are computed by summing flow parameters (i.e., pressure, velocity, etc.) at each lifting surface grid location over 50 time steps (one aerodynamic forcing period), respectively. Resulting summations are then divided by 50, giving a timeaveraged value for each solution parameter at each grid point. 52 4.2 BASELINE AIRFOIL TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS A timeaveraged totalpressure contour plot for the NACA 0012 airfoil is presented in Figure 41, with attached numerical values indicating respective pressure contour levels. 101102192075 101290 101320 101305 101305 101315 101320 101320 101320 Chord, x/c Pitch, y/c 0.4 0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 Figure 41 TimeAveraged TotalPressure Contours (Pa): NACA 0012 Lifting Surface As expected for this liftingsurface profile, the contours illustrate good qualitative characteristics of the simulated flowfield, showing a wellbehaved symmetric total pressure distribution about the profile and in the liftingsurface wake. The flowfield on the airfoil surface shown in Figure 41 shows no indication of largescale separation, exhibiting smooth attached flow. Thus, the timeaveraged flow about the baseline liftingsurface is argued to be attached and producing a similar wake to that which would be expected. 53 Figure 42 exhibits the totalpressure wake profile of the NACA 0012 lifting surface. Pitchwise total pressure wake profiles at x/c = 0.005 and 0.01 downstream of the airfoil are exhibited. As discussed in Section 3.4 wake depth is defined as the maximum total pressure deficit, while wake width represents the maximum pitchwise effect of the total pressure deficit. Symmetric wake profiles are observed at both x/c locations downstream of the airfoil. Thus, confirming the inferences made earlier from Figure 41. Pitch, y/c Total Pressure, Pt (Pa) 0.2 0.1 0.0 0.1 0.2 101280 101284 101288 101292 101296 101300 101304 101308 101312 101316 101320 101324 x/c = 0.005 aft x/c = 0.01 aft Figure 42 NACA 0012 Lifting Surface Wake Profile. Figure 43 displays timeaveraged static pressures computed along the NACA 0012 lifting surface chord. Pressure distributions along the suction and pressure surfaces are individually displayed. As expected, the symmetric NACA 0012 profile produces equivalent timeaveraged staticpressure distributions on the suction (upper) and pressure (lower) surfaces. Thus, no timeaveraged aerodynamic loading exists on the airfoil. Staticpressure distributions predicted by FLUENT along the NACA 0012 profile are also compared with available experimental data [Abbott and Von Doenhoff, 1959], as 54 shown in Figure 43. The computational data compared favorably with the experimental data, except near the leading and trailing edge regions where slight discrepancies were found. In particular, the minimum pressure near the leading edge was underestimated, while a slightly lower pressure was predicted in the trailingedge region. Given the reasonable agreement of the timeaveraged data, the simulations presented herein correctly predict the timeaveraged flowfield about the examined lifting surface. Chord, x/c Time Averaged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101235 101250 101265 101280 101295 101310 101325 NACA0012 Suction NACA0012 Pressure NACA0012 Suction (Abbot & Von Doenhof) NACA0012 Pressure (Abbott & Von Doenhoff) Figure 43 TimeAveraged Static Pressure : NACA 0012. Figure 44 illustrates timeaveraged staticpressure distributions computed for flow over the NACA 0012 airfoil at 0°, 5° and 10° angles of attack relative to the freestream. The nonequivalent surfacepressure distributions on the suction and pressure surfaces indicate aerodynamic loading caused by the nonzero angle of attack of the flow. The static pressure at the suctionsurface leading edge decreases as angle of attack increases, while the static pressure at the pressuresurface leading edge increases. As 55 angle of attack increases, the difference in the pressure between the suction and pressure surfaces increases, producing higher aerodynamic loading. Chord, x/c Time Averaged Static Pressure 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101100 101125 101150 101175 101200 101225 101250 101275 101300 101325 5AOA Suction 10AOA Suction NACA0012 Suction 5AOA Pressure 10AOA Pressure NACA0012 Pressure Figure 44 TimeAveraged Static Pressure : 0°, 5° and 10° Angle of Attack Table 4.1 compares the computed lift coefficient values for the NACA 0012 airfoil at 0°, 5° and 10° angles of attack (relative to the freestream) with available experimental values from Abbott and Von Doenhoff [1959]. Note that the experimental values were measured at Re = 6x106, while Reynolds number for the current simulations was calculated to be Re = 5.88x106. As such, the discrepancy in the lift coefficient value is attributed to the lower Re used for the current investigation. Table 4.1 Coefficient of Lift: 5° and 10° Angle of Attack. Angle of Attack Computed Value Experimental Value 5° 0.54 0.59 10° 0.82 0.87 56 4.3 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (THICKNESS INFLUENCE) As described later in this thesis, the influence of liftingsurface thickness on timeaccurate surfacepressure phase is examined via simulation of three NACA profiles of 10, 12, 15 and 20% thickness (based on chord). Figure 45 illustrates the corresponding timeaveraged staticpressure distributions from these simulations, collected over one aerodynamic forcing period. The data in Figure 45 show good agreement between suction and pressure surfaces for each thickness, as anticipated for symmetric profiles. The equivalent surfacepressure distributions on both suction and pressure surfaces for each thickness case are indicative of no timeaveraged aerodynamic loading. Each profile has distinct timeaveraged pressure gradients along the chord, with higher thickness values resulting in more severe chordwise gradients. Chord, x/c Time Average Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101220 101240 101260 101280 101300 101320 NACA0010 Suction NACA0012 Suction NACA0015 Suction NACA0020 Suction NACA0010 Pressure NACA0012 Pressure NACA0015 Pressure NACA0020 Pressure Figure 45 TimeAveraged StaticPressure for Various LiftingSurface Thicknesses. 57 The staticpressure distributions predicted by FLUENT for the 10 and 15% thickness profiles were compared with available experimental data from Abbott and Von Doenhoff [1959] as shown in Figure 46. Here again, as in the NACA 0012 case, the computational data compared favorably with the experimental data, except in the leading and trailing edge regions where slight discrepancies are observed. In particular, the minimum staticpressure near the leading edge is underestimated and a slightly lower staticpressure is predicted in the trailing edge region. Chord, x/c TimeAveraged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101235 101250 101265 101280 101295 101310 101325 NACA0010 (Computational) NACA0010 (Abbott & Von Doenhoff) NACA0015 (Computational) NACA0015 (Abbott & Von Doenhoff) Figure 46 Time averaged Static Pressure Comparison with Experimental data. 4.4 TIMEAVERAGED STATICPRESSURE DISTRIBUTIONS (CAMBER INFLUENCE) In addition to the liftingsurface cases, three cambered liftingsurface geometries were also considered, including 0, 2, and 6% camber (relative to chord). Figure 47 shows timeaveraged staticpressure distributions along the suction and pressure surfaces 58 of each cambered profile. Note from Figure 47 that while each profile exhibits a unique chordwise surfacepressure gradient, as for various profile thicknesses, the cambered profiles also exhibit timeaveraged aerodynamic loading. This is exhibited by the nonequivalent surface pressure distributions on the suction and pressure surfaces in Figure 47. Differential pressure across the liftingsurface increases as the percentage of camber relative to chord increases. Chord, x/c Time Averaged Static Pressure (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 101220. 101240. 101260. 101280. 101300. 101320. NACA0012 Suction 2%Camber Suction 6% Camber suction NACA0012 Pressure 2%Camber Pressure 6% Camber Pressure Figure 47 TimeAveraged StaticPressure Distribution for Various LiftingSurface Cambers Table 4.2 compares the computed lift coefficient values for 2, and 6% cambered (relative to chord) airfoils with available experimental values from Abbott and Von Doenhoff [1959]. Here again, as was seen earlier for the AOA case, the computed lift coefficient values were slightly lower than the experimental values. 59 Table 4.2 Coefficient of Lift: Various LiftingSurface Cambers. Camber) Computed Value Experimental Value 2 % 0.21 0.25 6 % 0.74 0.78 4.5 COMPARISON OF TIMEACCURATE BASELINE LIFT DEPENDENCY TO SEAR’S RESULTS Figure 48 presents the firstharmonic (i.e., fundamental frequency) lift time series obtained from FLUENT for the baseline NACA 0012 airfoil, as well as known lift prediction for airfoils experiencing a sinusoidal transverse gust, as described by Sears [1938]. As can be seen in Figure 47, the computed lift compares perfectly well with the Sears predicted lift series. Note that, Sears unsteady aerodynamic forcing was modeled as a perfect sinewave transverse gust as opposed to the rotorwake model (transverse and chordwise varying velocity components) used for unsteady aerodynamic forcing in the present analysis. Given this agreement of the lift time series data, the unsteady results obtained from FLUENT are further validated. 60 Time (sec) Lift 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 2.4 2 1.6 1.2 0.8 0.4 0 0.4 0.8 1.2 1.6 2 2.4 Sears Computational Figure 48 LiftTime series Comparison with Sears results. 4.6 TIMEAVERAGED RESULTS SUMMARY Timeaveraged results obtained from FLUENT for the baseline NACA 0012 airfoil showed wellbehaved characteristics. Total pressure contours around the airfoil exhibited smooth attached flow with no signs of largescale separation and a symmetric wake profile. Pitchwise wake profiles developed downstream of the NACA 0012 lifting surface at x/c = 0.005 and x/c = 0.01 exhibited symmetry as well, further ascertaining this inference. Timeaveraged static pressure distributions along the upper and lower surface of the NACA 0012 airfoil were perfectly symmetric, as expected for this profile. Thus, indicating no aerodynamic loading exists on the airfoil. Comparison of the timeaveraged pressures with experimentally obtained values from Abbott and Von Doenhoff [1959] exhibited reasonable agreement with very slight discrepancies observed near the leading edge and trailing edge regions. Given this agreement of timeaveraged data, the simulations presented correctly predict the timeaveraged flowfield about the examined lifting surface. 61 Timeaveraged pressure distributions obtained for various mean flow angles of attack (AOA), various thickness profiles and various camber profiles exhibited characteristics similar to that, which would be expected for the respective profiles. Nonequivalent timeaveraged pressure distributions were exhibited at various mean flow angles of attack and by the various camber profiles, indicating aerodynamic loading on these profiles. Differential pressure (aerodynamic loading) was observed to increase as both AOA and percentage of camber relative to chord increased. Conversely, the various thickness profiles exhibited equivalent timeaveraged pressure distributions along the upper and lower surfaces indicative of no timeaveraged aerodynamic loading. Each thickness profile exhibits distinct timeaveraged pressure gradients along the chord, with higher thickness values resulting in more severe chordwise gradients. Computationally obtained timeaveraged results for the parametric analyses were compared with those obtained experimentally by Abbott and Von Doenhoff [1959]. The computationally obtained timeaveraged results compared favorably with the experimental data. Thereby, validating the steady state solution obtained from FLUENT for the current simulations. To validate the unsteady solution, firstharmonic lift timeseries obtained from FLUENT for the NACA 0012 airfoil was compared with the Sears predicted lift series. The timeseries matched perfectly well, even though Sears unsteady aerodynamic forcing was modeled as a perfect sinewave transverse gust as opposed to the rotorwake model (transverse and chordwise varying velocity components) used for unsteady aerodynamic forcing in the present analysis. Thus, the unsteady solution is validated. 62 5 Chapter 5 RESULTS FOR NACA 0012 BASELINE CONFIGURATION This chapter describes the aerodynamic forcing function, as well as resulting airfoil forced response, for the baseline NACA 0012 airfoil simulations. Airfoil results are reported in terms of chordwise unsteady surfacepressure distributions, spectral content, and phase. Examination of the surfacepressure phase data reveals characteristics indicative of multiple disturbance interactions, similar to that discussed by Fabian and Jumper [1996]. An analytical model is developed to reproduce these observed interaction characteristics, and compared to simulation results. 5.1 DATA REDUCTION METHODOLOGY Unsteady pressure results were obtained by removing the timeaveraged pressure from the instantaneous pressure via P′ = P − P Unsteady pressure data were further reduced into elements of amplitude, frequency, and phase for first, second and third harmonics (i.e. one, two and three times the fundamental frequency) components. This was accomplished via FastFourier Transform (FFT) techniques. 63 5.2 AIRFOIL FORCINGFUNCTION TIME DEPENDENCY Figure 51 through Figure 54 exhibit totalpressure contours of the airfoil aerodynamic forcing function, or wake, moving downstream over the NACA 0012 airfoil. The plots illustrate totalpressure contours at four time instances during a single aerodynamic forcing period, T; illustrated times correlate to t = 0, T/4, T/2, and 3T/4. The totalpressure disturbance is periodic and repeats with every forcing period. As the totalpressure disturbance translates, it directly impacts the airfoil leading edge. At the airfoil leading edge, the disturbance splits and propagates downstream along both airfoil lower and upper surfaces. This disturbance splitting process is most evident from t = T/4 to t = T/2. At impact on the airfoil upper surface, the disturbance accelerates (as compared to lower surface) before propagating downstream along the airfoil, decaying in strength during the process. On the airfoil lower surface, disturbance impact is much less prominent; disturbance chordwise propagation is delayed (i.e. it remains further upstream) relative to the uppersurface disturbance. Figure 51 Forcing Function TotalPressure Contours, t = 0. Figure 52 Forcing Function TotalPressure Contours, t = T/4. 64 Figure 53 Forcing Function TotalPressure Contours, t = T/2. Figure 54 Forcing Function TotalPressure Contours, t = 3T/4. Figure 55 displays totalpressure time series relating the airfoil forcingfunction at five equally spaced x/c locations (x/c = 0.2, 0.4, 0.6, 0.8, and –1.0) forward of the airfoil leading edge. Note that the time series have been arbitrarily shifted vertically by four units at each sequential x/c location, to provide better viewing. As expected, Figure 55 indicates periodic totalpressure variations corresponding to the aerodynamicforcing frequency. Disturbance amplitude changes as the forcing wake travels downstream, indicating expected wake decay with convective distance. A monotonic phase shift is also observed between the time series, indicating constantspeed disturbance propagation downstream, or convection. Note that while the pressure fluctuations in Figure 55 are periodic, they are not purely sinusoidal, indicating the existence of forcingfunction harmonic content. 65 Time, t (sec) "Shifted"Total Pressure, Pt (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 35.0 30.0 25.0 20.0 15.0 10.0 5.0 0.0 5.0 10.0 15.0 20.0 25.0 x/c = 0.2 x/c = 0.4 x/c = 0.6 x/c = 0.8 x/c = 1.0 Figure 55 TotalPressure Time Series Forward of Airfoil. Figure 56 displays forcingfunction staticpressure time series at five equally spaced x/c locations (20, 40, 60, 80 and 100%) forward of the airfoil leading edge. Similar amplitude trends (decreasing downstream) and periodicity are seen as in the totalpressure time series; however, no phase shift is observed between the series. Here again, the static pressure fluctuations are periodic but not purely sinusoidal, indicating the existence of harmonic content. Time, t (sec) "Shifted" Static Pressure, Ps (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 10.0 5.0 0.0 5.0 10.0 15.0 20.0 25.0 30.0 x/c = 0.2 x/c = 0.4 x/c = 0.6 x/c = 0.8 x/c = 1.0 Figure 56 StaticPressure Time series. 66 5.3 AIRFOIL FORCINGFUNCTION SPECTRAL CONTENT Figure 57 and Figure 58 depict spectral content corresponding to the forcingfunction totalpressure and staticpressure time series in Figure 55 and Figure 56, respectively. Note the xaxis in the spectral plots is reversed in order to indicate flow direction from the simulation inlet to airfoil leading edge. Significant harmonic frequencies reaching three times the primary aerodynamicforcing frequency are observed in both figures. Both total and static pressures exhibit increased harmonic content near the inlet, with totalpressure showing greater amplitude. In each case, harmonic content decays as the wake travels downstream towards the airfoil leading edge. Such harmonic amplitude decay is expected, as wake mixing with convective distance inherently causes rapid spectral content loss. Chord, x/c Unsteady Total Pressure, Pt' (Pa) 1.0 0.8 0.6 0.4 0.2 0.0 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 1st Harmonic 2ndHarmonic 3rdHarmonic Inlet Airfoil LE Figure 57 Forcing Function TotalPressure Spectral Content. Chord, x/c Unsteady Static Pressure, Ps' (Pa) 1.0 0.8 0.6 0.4 0.2 0.0 0.0 1.0 2.0 3.0 4.0 5.0 1st Harmonic 2ndHarmonic 3rdHarmonic Inlet Airfoil LE Figure 58 Forcing Function StaticPressure Spectral Content. 5.4 AIRFOIL FORCINGFUNCTION PHASE DEPENDENCY Figure 59 shows relative phase data for the firstharmonic static pressure at five equally spaced x/c locations (x/c = 0.2, 0.4, 0.6, 0.8, and –1.0) forward of the airfoil leading edge. Note the xaxis direction on the phase plot is reversed corresponding to 67 disturbance propagation direction. The firstharmonic phase data in Figure 59 shows consistent downstream disturbance propagation (i.e., convection) towards the airfoil leading edge, as indicated by the constant negative slope of the phase line. However, unsteady disturbance propagation changes dramatically as the airfoil leading edge is approached, as indicated by the positive phase slope beginning x/c = 0.2 forward of the leading edge. This change in phase slope correlates to rapid disturbance deformation near the leading edge, as can be observed in Figure 52. Chord, x/c Relative Phase (Radian) 1.0 0.8 0.6 0.4 0.2 0.0 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 1st Harmonic Figure 59 Airfoil Forcing Function 1st Harmonic Phase. Relativephase data at the second and thirdharmonic frequencies are shown in Figure 510, for the forcingfunction staticpressure data at five equally spaced x/c locations forward of the airfoil leading edge. It is observed that the second harmonic has a positively increasing phase slope; the increase in slope is very small until x/c = 0.4, after which it increases significantly. In comparison, the third harmonic phase slope is approximately zero, indicating very little propagation of the thirdharmonic static pressure component. It is interesting to note the relative phase differences between the 68 static pressure first, second and third harmonics. Each harmonic component exhibits a different phase behavior, but collectively they combine to produce a waveform showing almost no propagation characteristics, as noted in Figure 56. Chord, x/c Relative Phase (Radian) 1.0 0.8 0.6 0.4 0.2 0.0 0.2 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 2.2 2nd Harmonic 3rd Harmonic Figure 510 Airfoil Forcing Function Higher Harmonic Phase. 5.5 AIRFOIL SURFACEPRESSURE TIME DEPENDENCY Figure 511 through Figure 514 show contours of static pressure for the NACA 0012 airfoil at four time instances during a single aerodynamic forcing period, T; illustrated times correlate to t = 0, T/4, T/2, and 3T/4. Figure 511 Airfoil StaticPressure Contours, t = 0. Figure 512 Airfoil StaticPressure Contours, t = T/4. 69 Figure 513 Airfoil StaticPressure Contours, t = T/2. Figure 514 Airfoil StaticPressure Contours, t = 3T/4. As desired in construction of the forcing function, staticpressure disturbances repeat in time, mimicking the effect of a wake convecting from an upstream rotor in a turbomachine compressor. At the airfoil leading edge, the wakeinduced pressure disturbances split and propagate downstream along both airfoil lower and upper surfaces. This disturbance splitting process is most evident from t = T/4 to t = T/2. At impact on the airfoil upper surface, the pressure waves also reflect back into the oncoming disturbance field. The impacted wave and its reflection travel together downstream along the airfoil, decaying in strength during the process. On the airfoil lower surface, pressure wave impact and reflection is much less prominent. Lowersurface chordwise disturbance propagation is delayed (i.e. further downstream) relative to the upper surface, as illustrated from t = T/4 to t = T/2. Figure 515 shows airfoil surfacepressure time series at various chordwise locations along the airfoil upper surface. Note, the staticpressure series have been arbitrarily shifted vertically by two units at each sequential x/c location, to provide better viewing. As expected, these figures indicate periodic staticpressure variations corresponding to the aerodynamicforcing frequency. While the pressure fluctuations are 70 periodic, they are not purely sinusoidal, indicating the existence of harmonic content. Higher amplitude unsteady pressure fluctuations exist near the airfoil leading edge, decaying rapidly downstream. Figure 5.15 also indicates almost no phase delay between time series at each x/c locations, a curious finding given the convective nature of the forcing function. This finding corresponds to the freestream staticpressure time series of Figure 56, and will be examined further in following sections. Time, t (sec) "Shifted"Unsteady Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 515 NACA 0012 Upper Surface Unsteady StaticPressure Time series. Figure 516 shows airfoil surfacepressure time series at various chordwise locations along the airfoil lower surface. Similar to the upper surface, the lower surface exhibits periodic timeseries behavior corresponding to the aerodynamicforcing frequency, unsteady leadingedge pressure amplification, and very little phase shift between chordwise locations. The lowersurface pressure fluctuations also display an almost reversed behavior as compared to the uppersurface, showing pressure increases at 71 the same instances as the upper surface experiences a pressure decrease. However, overall timeseries amplitudes are greater for the lower surface as compared to the upper surface. Time, t (sec) "Shifted"Unsteady Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 24.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 516 NACA 0012 Lower Surface Unsteady StaticPressure Time series. Differential surfacepressure time series data at various chordwise locations along the airfoil are provided in Appendix D. Differential surface pressure represents uppersurface minus lowersurface unsteady pressure, as presented in Eq. (5.1). u l Δp′ = p′ − p′ 5.1 Here again, periodic timeseries behavior corresponding to the aerodynamicforcing frequency is obtained. The differential pressure is highest in amplitude at the leading edge, decays in strength downstream along the airfoil, and shows almost no phase shift between chordwise locations. Differentialpressure fluctuations are periodic in nature, but not purely sinusoidal, indicating the existence of harmonic content. 72 5.6 AIRFOIL SURFACEPRESSURE SPECTRAL CONTENT Figure 517 depicts spectral content related to the airfoil uppersurface pressure time series shown in Figure 515. Relevant surfacepressure harmonic frequencies reaching four times the aerodynamicforcing frequency are observed. The airfoil upper surface exhibits increased higherorder pressure harmonics near the airfoil leading edge, generally decaying downstream along the upper surface. Specifically, the firstharmonic amplitude decays until x/c = 0.4, then increases toward the airfoil trailing edge. The second harmonic amplitude also shows a decaying trend downstream, but between x/c = 0.4 and x/c = 0.6 it exhibits higher amplitude than the first harmonic. The third and fourth harmonic amplitudes generally show a decaying trend with chordwise distance, almost vanishing near the airfoil trailing edge. Note that harmonic content can be related to the “sinelike” behavior of a time series. For example, it would be expected that a wake timeseries profile, such as in Figure 55, would exhibit significant harmonic content, while a pure sine wave would not. Thus, the relative change in harmonic content with chord position, as seen in Figure 517, can be directly attributed to disturbance deformation and interaction over the airfoil. Near the leading edge the disturbance is “wakelike” while at the trailing edge it is much more “sinelike”. This statement is supported by Figure 515. 73 Chord, x/c Surface Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.5 1.0 1.5 2.0 2.5 1st Harmonic 2nd Harmonic 3rd Harmonic 4th Harmonic Figure 517 NACA 0012 UpperSurface Static Pressure Spectral Content. Spectral content related to the airfoil lowersurface pressure time series is depicted in Figure 518, corresponding to the timeseries data of Figure 516. Similar to the airfoil upper surface, the lower surface shows relevant surfacepressure frequencies reaching four times the aerodynamicforcing frequency, as well as increased higherorder pressure harmonics near the airfoil leading edge. The firstharmonic amplitude shows a similar trend to the upper surface, decaying downstream until x/c = 0.6 and then increasing. Unlike the upper surface, the second harmonic amplitude decays constantly downstream, always having lower amplitude than the first harmonic. Higher order harmonics decay downstream along the airfoil lower surface, approximately reaching a zero value near the trailing edge. Again, the change in harmonic content between the upper and lower surfaces indicates a difference in disturbance deformation and interaction over the airfoil lower surface. This statement is supported by Figure 516. 74 Chord, x/c Surface Static Pressure, Ps' (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.5 1.0 1.5 2.0 2.5 1st Harmonic 2nd Harmonic 3rd Harmonic 4th Harmonic Figure 518 NACA 0012 LowerSurface StaticPressure Spectral Content. 5.7 AIRFOIL SURFACEPRESSURE FIRST HARMONIC AMPLITUDE CHORDWISE DEPENDENCY Figure 519 and Figure 520 display firstharmonic staticpressure time series at various chordwise locations along the airfoil upper and lower surfaces, respectively. Note, the staticpressure series have been arbitrarily shifted vertically by two units at each successive x/c location, to provide better viewing. As expected, these figures indicate periodic pressure variations corresponding to the aerodynamic forcing frequency. A phase shift is also observed between chord locations, (as indicated by the arrows) indicating disturbance propagation direction along the chord. On both upper and lower surfaces, large unsteady pressure fluctuations exist near the airfoil leading edge, decaying rapidly downstream. The overall timeseries amplitudes are greater for the lower surface as compared to the upper surface. 75 Time, t (sec) 1st Harmonic Unsteady Pressure, (Pa) 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 x/c = 0.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 519 NACA 0012 Upper Surface 1st Harmonic Pressure Time series. Time, t (sec) 1st Harmonic Unsteady Pressure 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 1.6 1.8 2.0 0.0 2.0 4.0 6.0 8.0 10.0 12.0 14.0 16.0 18.0 20.0 22.0 24.0 x/c = 0.1 x/c = 0.2 x/c = 0.3 x/c = 0.4 x/c = 0.5 x/c = 0.6 x/c = 0.7 x/c = 0.8 x/c = 0.9 x/c = 1.0 Figure 520 NACA 0012 Lower Surface 1st Harmonic Pressure Time series. 76 Figure 521 exhibits firstharmonic pressure amplitude dependence on chord; upper, lower and differential unsteady pressure amplitudes are shown. The airfoil lower surface exhibits higheramplitude unsteady pressures as compared to the upper surface; thus, confirming inferences made earlier related to Figures 519 and 520. Upper, lower and differential unsteady pressure amplitudes follow similar chordwise trends; i.e., amplified at leading edge and decaying downstream along the chord. Chord, x/c 1st Harmonic Amplitude, (Pa) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 0.0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 Upper Lower Differential Figure 521 NACA 0012 1st Harmonic Pressure Amplitude. 5.8 AIRFOIL SURFACEPRESSURE FIRST HARMONIC CHORDWISE PHASE DEPENDENCY Figure 522 shows relative phase data for the firstharmonic unsteady surface pressures along the airfoil. The figure displays airfoil uppersurface, lowersurface and differential pressure phase along the airfoil chord. 77 Chord, x/c 1st Harmonic Relative Phase (Radians) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 3.0 2.0 1.0 0.0 1.0 2.0 3.0 4.0 Upper Lower Differential Figure 522 NACA 0012 1st Harmonic Phase. Several important results can be immediately obtained from Figure 522. First, the uppersurface and lowersurface phase data exhibit different chordwise trends. Second, the uppersurface data show a significantly varying slope, particularly near the leading and trailing edges. Third, the lowersurface phase exhibits almost no slope near the midchord, increasing in slope towards the trailing edge. Finally, phase changes appear rapidly near the leading edges on both upper and lower surfaces. Fabian and Jumper [1996] established phasemap analysis as a viable means of determining airfoil disturbancepropagation characteristics in an unsteady flow. Fabian and Jumper argued that a negative slope in a phase map, such as near the leading edge in Figure 522, indicates downstreampropagating surfacepressure waves, while a positive slope predicts upstreampropagating disturbances. Clearly, if these arguments are true, the data in Figure 522 exhibit nonphysical disturbancepropagation characteristics. 78 Similar nonphysical phase data were also discussed by Fabian and Jumper, from which they indicated “ambiguous” phase maps, like the one in Figure 522, most likely result from multiple pressure disturbances interacting over the airfoil chord. 5.9 AIRFOIL SURFACEPRESSURE ANALYTICAL MODEL The phase data of Figure 522 display ambiguous disturbance propagation, similar to that discussed by Fabian and Jumper [1996]. The phase data do not follow the convected speed along the airfoil, from which the downstream propagation of an airfoil forcing disturbance would exhibit a negative phase slope. In fact, based on the results of Fabian and Jumper, the uppersurface phase data of Figure 522 suggest an initially downstreampropagating disturbance from x/c = 0.0 to x/c = 0.1, with the disturbance gaining speed from x/c = 0.1 to x/c = 0.3, as evidenced by the flat phase slope. At x/c = 0.3, the disturbance presumably changes direction, propagating upstream from x/c = 0.55 to x/c = 0.8. Finally, the disturbance propagates upstream at a lower speed near the airfoil trailing edge. Conversely, the lowersurface phase map of Figure 522 suggests an initially downstreampropagating disturbance gaining speed from x/c = 0.1 to x/c = 0.65, finally propagating upstream near the airfoil trailing edge. These inferences are reinforced by examining the fundamental (firstharmonic) frequency time series in Figure 519 and Figure 520, for the airfoil upper and lower surfaces, respectively. In Figures 519 and 520, the probable propagation direction of the waveforms from the leading edge to each successive x/c location has been indicated, where it has been assumed that the minimum time between each successive pressure peak provides the correct wave propagation direction. When compared to the uppersurface 79 and lowersurface phase maps of Figure 522, the wave propagation directions in Figure 519 and 520 are consistent, showing the same overall chordwise characteristics. 5.9.1 INTERACTION MODEL Although the phase information in Figures 519, 520 and 522 correlate, the disturbance propagation directions implied by the data are clearly nonphysical and incorrect. In the investigation of Fabian and Jumper, ambiguous phase maps similar to Figure 522 were argued to result from convectivelypropagating vortical and acousticallyradiating potential disturbance int 



A 

B 

C 

D 

E 

F 

I 

J 

K 

L 

O 

P 

R 

S 

T 

U 

V 

W 


