4. 

NASA Contractor Report 185217 

Investigation of Advanced Counterrotation 

Blade Configuration Concepts for High Speed 

Turboprop Systems 

Task 1 - Ducted Propfan Analysis 

Final Report 

Edward J. Hall, Robert A. Delaney, and James L. Bettner 

Allison Gas Turbine Division of General Motors 

Indianapolis, Indiana 

April 1990 

Prepared for 

Lewis Research Center 

Under Contract NAS3-25270

ORIGINAL CONTAINS 

#rIr !'. NASA 

National Aeronautics and 

Space Administration 

NAA -L'.- l e' I 7) iNVT (CAT Ifl'1 ',j r- ALVA!Ct 

COUNTRPLTATION aLAnE CONFIGURATION CflCLPTS 

FOR HIGH PFEO TURBOPROP SYSTEMS, TASK I: 

DUCTtn PROPFAN ANALYSIS Final Report 

(General Motors Corp.) 18 p CSCI 2IE 3/O7

unci as 

0261598
Errata for NASA Contract NAS3-25270

Investigation of Advanced Counterrotation Blade Configuration

Concepts for High Speed Turboprop Systems

Taski-Ducted Propfan Analysis

Final Report 

Page 29: Second line from the bottom should read - Pref is determined from the freestream 

total density 

Page 78: Figure title should read - Axisymmetric plane projection of 160x41x15 H-grid for 

Kruger

ORIGINAL CONTAII4S 

COLOR ILLDJSjj10
PRECEDING PAGE BLANK NOT FILMED 

111 

Preface 

This manual was prepared by E. J. Hall and R. A. Delaney of the Advanced 

Turbomachinery Department, Allison Gas Turbine Division, General Motors Corporation, Indianapolis, IN. The work was performed under NASA Contract NAS3-25270 

from September, 1988 to January, 1990. The grid generation and flow code theory 

and programming modifications necessary for the analysis of ducted propfans were 

performed by E. J. Hall based on existing NASA software. The Allison program 

manager for this contract was J. L. Bettner. 

Acknowledgements 

The authors would like to express their appreciation to the following NASA 

personnel who contributed to this program: 

M. Celestina, K. Curtley, and T. Beach for their helpful suggestions during the 

early phases of this program; 

Dr. J. J. Adamczyk for the many helpful technical discussions concerning the 

development of the computer codes, and for providing the initial 3D Euler solver. 

Dr. C. H. Miller and Dr. L. Bober for their suggestions and critical review of 

the program. 

The services of the NASA Numerical Aerodynamic Simulation (NAS) facilities 

and personnel are gratefully acknowledged.
V 

TABLE OF CONTENTS 

NOTATION . 

1. SUMMARY .................................1 

2. INTRODUCTION .............................3 

3. GRID GENERATION ALGORITHM ................9 

3.1 Computational Domain ...........................9 

3.2 H-Grid Algorithm .............................12 

3.3 C-Grid Algorithm .............................20 

4. 3D EULER NUMERICAL ALGORITHM .............. 29 

4.1 Nondimensionalizatjon .......................... 29 

4.2 Governing Equations ............................ 30 

4.3 Runge-Kutta Time Integration ...................... 33 

4.4 Artificial Dissipation ........................... 34 

4.5 Implicit Residual Smoothing ....................... 35 

4.6 Boundary Conditions ........................... 37 

4.7 Multiple-Block Coupling ......................... 39 

4.8 Solution Procedure ............................ 40

, 4L.JHILNT%0IISAO PRECEDING PAGE BLANK NOT FILMED 
ril 

5. RESULTS • 41 

5.1 SR7 8-Bladed Propfan .......................... 42 

5.2 SR7 2-Bladed Propfan Modane Tests .................. 43 

5.3 Symmetric Circular Arc Cowl Test Case ................ 50 

5.4 Ducted SR7 Demonstration ....................... 58 

5.5 Low Speed Ducted Propeller ....................... 68 

5.6 NASA 1.15 Pressure Ratio Fan ..................... 79

6. CONCLUSIONS ................................97 

REFERENCES ...................................99 

APPENDIX A. GRID GENERATOR OPERATING INSTRUCTIONS .....................................103 

A.1 Introduction ................................103 

A.2 Input/Output File Description ......................105 

A.3 Subroutine Description ..........................117 

APPENDIX B. 3D EULER SOLVER OPERATING INSTRUCTIONS .....................................125 

B.1 Introduction ................................125 

B.2 Input/Output Files ............................128 

B.3 H-grid Euler Solver Subroutine Description ...............133 

B.4 Multiple-Block C-Grid Euler Solver Subroutine Description ......138 
vu 

LIST OF FIGURES 

Figure 2.1: Propfan aerodynamic characteristics . 5

Figure 2.2: Ducted propfan (very high bypass fan) geometry .......6 

Figure 3.1: Ducted propfan analysis computational domain ........11 

Figure 3.2: Axisymmetric plane projection grid generation subregions for 

an unducted propeller ......................13 

Figure 3.3: Axisymmetric plane projection grid generation subregions for 

a ducted propeller ........................14 

Figure 3.4: Sample H-grid for a ducted propfan ..............21 

Figure 3.5: Axisymmetric plane projection subregions for multiple-block 

C-grid generation ........................23 

Figure 3.6: Sample multiple-block C-grid for a ducted propfan ......28 

Figure 4.1: Three-dimensional finite volume cell ..............32 

Figure 5.1: SR7 propfan design characteristics.................44

Figure 5.2: Axisymmetric plane projection of 70x35x15 SR7 grid .....45 

Figure 5.3: Comparison of predicted, experimental, and design power coefficient distributions for 8-bladed SR7 propfan (M=0.8) . . . 46
yin 

Figure 5.4: Predicted propfan surface static/total pressure ratio contours 

for 8-bladed SR7 propfan (M=0.8) ...............47 

Figure 5.5: Comparison of predicted and experimental blade spanwise elemental power coefficient loading distributions for the 8-bladed 

SR7 propfan (M=0.7) ......................48 

Figure 5.6: 2-bladed SR7 propfan geometry and 120x45x25 grid for Modane 

test comparison .........................51 

Figure 5.7: Comparison of predicted and experimental airfoil surface static 

pressure coefficient distributions for 2-bladed SR7 propfan Modane 

test (9.9% span) ..........................52 

Figure 5.8: Comparison of predicted and experimental airfoil surface static 

pressure coefficient distributions for 2-bladed SR7 propfan Modane 

test (55.6% span) ........................53 

Figure 5.9: Comparison of predicted and experimental airfoil surface static 

pressure coefficient distributions for 2-bladed SR7 propfan Modane 

test (94.4% span) ........................54 

Figure 5.10: Comparison of predicted and experimental airfoil pressure surface static/total pressure ratio contours for 2-bladed SR7 propfan................................55 

Figure 5.11: Comparison of predicted and experimental airfoil suction surface static/total pressure ratio contours for 2-bladed SR7 propfan................................ 56
Ix 

Figure 5.12: Symmetric double circular arc cowl ducted propeller test case 

geometry.............................59 

Figure 5.13: Comparison of 91x21 2D grid and axisymmetric plane projection of 91x42x5 3D H-grid for the symmetric double circular 

arc cowl test case ........................60 

Figure 5.14: Comparison of 2D and 3D H-grid predicted cowl surface grid 

line Mach number distributions (M=0.5) ...........61 

Figure 5.15: Comparison of 2D and 3D H-grid predicted cowl surface grid 

line Mach number distributions (M=0.7) ...........62 

Figure 5.16: Comparison of 91x21 2D grid and axisymmetric plane projection of multiple-block C-grid for the symmetric double circular 

arc cowl test case ........................63 

Figure 5.17: Comparison of 2D and 3D multiple-block C-grid predicted 

cowl surface grid line Mach number distributions (M=0.5) . . 64 

Figure 5.18: Comparison of 2D and 3D multiple-block C-grid predicted 

cowl surface grid line Mach number distributions (M=0.7) . . 65 

Figure 5.19: Comparison of the predicted static pressure contours for the 

axisymmetric plane projection of the 3D H-grid (top) and 

multiple-block C-grid (bottom) Euler analyses of the symmetric double circular arc cowl test case (M=0.7) .........66 

Figure 5.20: Ducted 8-bladed SR7 geometry and 69x45x15 H-grid .....69 

Figure 5.21: Axisymmetric plane projection of 69x45x15 H-grid for ducted 

SR7propfan ...........................70
x 

Figure 5.22: Comparison of H-grid predicted elemental blade power coefficient distribution with frequency domain panel method results 

for ducted 8-bladed SR7 ....................71 

Figure 5.23: Comparison of H-grid predicted blade elemental thrust coefficient distributions with frequency domain panel method results for ducted 8-bladed SR7 ..................72 

Figure 5.24: Predicted H-grid propfan surface static/total pressure ratio 

contours for ducted SR7 .....................73 

Figure 5.25: Axisymmetric plane projection of multiple-block C-grid for 

ducted SR7 propfan ......................74 

Figure 5.26: Comparison of multiple-block C-grid predicted elemental blade 

power coefficient distribution with frequency domain panel 

method results for ducted 8-bladed SR7 ............75 

Figure 5.27: Comparison of multiple-block C-grid predicted blade elemental thrust coefficient distributions with frequency domain panel 

method results for ducted 8-bladed SR7 ............76 

Figure 5.28: Geometry and 127x48x15 grid for Kruger low speed ducted 

propeller.............................77 

Figure 5.29: Aid symmetric plane projection of 127x48x15 H-grid for Kruger 

low speed ducted propeller ...................78 

Figure 5.30: Comparison of H-grid predicted and experimental cowl surface 

static pressure coefficient distributions for Kruger low speed 

propeller.............................80
ni 

Figure 5.31: Axisymmetric plane projection of multiple-block C-grid for 

Kruger low speed ducted propeller ...............81 

Figure 5.32: Comparison of multiple-block C-grid predicted and experimental cowl surface static pressure coefficient distributions for 

Kruger low speed propeller ...................82 

Figure 5.33: NASA 1.15 pressure ratio fan stage geometry (dimensions in 

cm) ................................83 

Figure 5.34: NASA 1.15 pressure ratio fan stage and 3D grid overlay . . . 84 

Figure 5.35: Axisymmetric plane projection of NASA 1.15 pressure ratio 

fan stage 225x52x15 H-grid ...................86 

Figure 5.36: Comparison of predicted and experimental cowl leading edge 

static pressure coefficient distributions for NASA 1.15 pressure 

ratio fan (H-grid, rotor only, M=0.75) .............87 

Figure 5.37: Comparison of predicted and experimental cowl leading edge 

static pressure coefficient distributions for NASA 1.15 pressure 

ratio fan (H-grid, rotor only, M=0.85) .............88 

Figure 5.38: Comparison of predicted and experimental pressure data for 

NASA 1.15 pressure ratio fan stage (H-grid, full stage, M=0.75) 90 

Figure 5.39: Predicted surface static pressure ratio contours for NASA 1.15 

pressure ratio fan stage (H-grid, full stage, M=0.75) .....91 

Figure 5.40: Multiple-block C-grid for NASA 1.15 pressure ratio fan stage 

(rotor only) ...........................93
XII 

Figure 5.41: Comparison of predicted cowl leading edge static pressure coefficient distributions for NASA 1.15 pressure ratio fan stage 

(C-grid, rotor only, M=0.75) ..................94 

Figure 5.42: Comparison of predicted cowl leading edge static pressure coefficient distributions for NASA 1.15 pressure ratio fan stage 

(C-grid, rotor only, M=0.85) ...................95 

Figure A.1: Sample input data file for grid generation ...........106 

Figure A.2: Blade setting angle reference plane description ........116 

Figure A.3: Sample grid based on H-grid generation .............118 

Figure A.4: Sample grid based on multiple-block C-grid generation . . . . 119 

Figure A.5: Program flowchart for ducted propfan analysis grid generation 123 

Figure B.1: Sample input data file for ducted propfan 3D Euler solver . . 129 

Figure B.2: Program flowchart for ducted propfan analysis H-grid Euler 

solver...............................139 

Figure B.3: Program flowchart for ducted propfan analysis multiple-block 

C-grid Euler solver .......................142
)Q11 

LIST OF TABLES 

Table A.1: Description of input/output files and UNIX-based filenames 

for H-grid generation ......................117 

Table A.2: Description of input/output files and UNIX-based filenames 

for C-grid generation ......................120 

Table B.1: Description of input/output files and UNIX-based filenames 

for H-grid Euler solver ......................134 

Table B.2: Description of input/output files and UNIX-based filenames 

for multiple-block C-grid Euler solver .............135
xv

NOTATION 

A list of the symbols used throughout this document and their definitions is 

provided below for convenience.

Roman Symbols 

a... speed of sound 

Cp... specific heat at constant pressure 

CV'... specific heat at constant volume 

e... internal energy 

i... z index of numerical solution 

j... r index of numerical solution 

k... 9 index of numerical solution 

n... time step index of numerical solution or rotational speed (rpm) 

il... outward unit normal vector 

p... pressure 

r... radius or radial coordinate 

t... time 

V... velocity 

z... axial coordinate 

ML4LL INTINT10149tV.

 SAW.

PRECEDING PAGE BLANK NOT FILMED
xv' 

A... surface area 

Cp... pressure coefficient (C = p - poo/2poouCIO 

Cp... power coefficient (C = P/pn3D5) 

Ct... thrust coefficient (Ct = T/pn2D4) 

CFL... Courant-Frei drichs- Levy 

D... dissipation flux vector or diameter 

DELT... adjacent cell grid spacing 

DELTLE... adjacent cell grid spacing at leading edge 

DELTTE... adjacent cell grid spacing at trailing edge 

F... flux vector in . z direction 

G... flux vector in r direction 

H... flux vector in 9 direction 

Ht... total enthalpy 

J... advance ratio (J = U/nD) 

K... source term flux vector 

L ... length 

N... number of propeller blades 

o... orthogonality 

P... point distribution function or power 

Q... vector of dependent variables 

R... gas constant or residual or maximum radius 

RAT... maximum ratio of adjacent cell grid spacings 

S... arc length or pertaining to surface area normal
xvii 

T... temperature or torque 

U... flight velocity 

V... volume

Greek Symbols 

a... time-stepping factor 

3... 3/4 radius propfan blade setting angle 

modified second-order damping coefficient 

64 ... modified fourth-order damping coefficient 

p... density 

ic2 ... second order damping coefficient 

K4 ... fourth order damping coefficient 

7... specific heat ratio 

S... spatial second-order central difference operator 

\... blockage factor 

ii... radial transformed variable for C-grid generation 

axial transformed variable for C-grid generation 

ii... damping factor 

T... relative pressure gradient coefficient 

increment of change

Special Symbols 

V... spatial vector gradient operator 

spatial forward difference operator
xviii 

v... spatial backward difference operator 

Superscripts 

[]... averaged variable 

[}... dimensional variable 

[]... implicitly smoothed variable 

[]... vector variable 

1*... intermediate variable 

time step index of variable 

Subscripts 

grid point index of variable 

[ ]max ... maximum value 

[imin... minimum value 

]p ... related to pressure 

pressure (high pressure) surface 

Jss ... suction (low pressure) surface 

[it... total quantity 

Liz... derivative or value with respect to z 

[}r... derivative or value with respect to r 

1 9 ... derivative or value with respect to 8 

[ioo ... freestream value 

[iref... reference value
1. SUMMARY 

The time-dependent three-dimensional Euler equations of gas dynamics have 

beensolved numerically to study the steady compressible transonic flow about ducted 

propfan propulsion systems. Aerodynamic calculations were based on a four-stage 

Runge-Kutta time-marching finite volume solution technique with added numerical 

dissipation. An implicit residual smoothing operator was used to aid convergence. 

Two calculation grids were employed in this study. The first grid utilized an H-type 

mesh network with a branch cut opening to represent the axisymmetric cowl. The 

second grid utilized a multiple-block mesh system with a C-type grid about the cowl. 

The individual blocks were numerically coupled in the Euler solver. Grid systems were 

generated by a combined algebraic/elliptic algorithm developed specifically for ducted 

propfans. Numerical calculations were initially performed for unducted propfans to 

verify the accuracy of the three-dimensional Euler formulation. The Euler analyses 

were then applied for the calculation of ducted propfan flows, and predicted results 

were compared with experimental data for two cases. The three-dimensional Euler 

analyses displayed exceptional accuracy, although certain parameters were observed 

to be very sensitive to geometric deflections. Both solution schemes were found to 

be very robust and demonstrated nearly equal efficiency and accuracy, although it 

was observed that the multiple-block C-grid formulation provided somewhat better

/
resolution of the cowl leading edge region.
3

2. INTRODUCTION 

Modern turbopropeller (propfan) aircraft propulsion systems have demonstrated 

the potential for significant fuel savings over existing high bypass turbofan engines. 

By incorporating innovative airfoil design concepts utilizing a high degree of blade 

lean, sweep, and twist, propfans can perform effectively at relative Mach numbers 

extending well into the transonic regime. A number of aerodynamic phenomena common in propfan installations are illustrated on Fig. 2.1. The process of designing 

advanced propfans requires a complicated iterative interaction between the disciplines of aerodynamics, structural analysis, and aeroacoustics. Reliable attainment 

of performance goals necessitates the use of proven aerodynamic analyses for the 

highly three-dimensional blades. In this respect, the realization of the performance 

improvements available in advanced propfan designs is directly related to advances 

in computational fluid dynamics. 

Several researchers have recently published advanced numerical analyses for propfan flows. Barton et al. [1] and Usab et al. [2] presented results from several threedimensional Euler codes for single-rotation propellers which were verified through 

comparison with experimental data. Celestina et al. [3] performed a time-averaged 

solution of three-dimensional counter-rotating propeller flowfields using an averagepassage system of equations. Whitfield et al. [4] performed calculations for unsteady
4 

three-dimensional single-rotation and counter-rotation propfans using a multipleblock grid formulation. Saito et al. [5], Matsuo et al. [6], and Kobayawa and Hatano 

[7] have all demonstrated the application of the three-dimensional Navier-Stokes equations for the analysis of high speed propfan designs. Many of the aerodynamic features 

relevant to propfan flowfields have been successfully captured by these and other analyses, including blade leading edge and tip vortices, transonic/shock flow phenomena, 

and counter-rotating blade row unsteady interaction. 

In light of the desire to optimize propfan designs, the concept of a ducted propfan configuration (very high bypass fan, see Fig. 2.2) has recently been explored. By 

ducting the propulsor, a significant increase in thrust and efficiency can be obtained, 

particularly at static, or low forward flight velocities. This thrust improvement could 

significantly impact future designs by permitting reduced diameters and underwing 

installation. Ducted propellers may also have advantages related to community noise 

levels and enhanced public acceptance. The concept of a ducted propeller has existed 

for more than fifty years, although research was initially directed at low-speed applications due to the drag penalty associated with the cowl at higher flight velocities. 

Early tests of ducted propellers are given in a review by Sacks and Burnell [8]. A 

more recent discussion on the technical issues relating to ducted propfans was given 

by Groeneweg and Bober [9]. In order to assess the relative aerodynamic merits 

of advanced ducted propfan design configurations, detailed computational analyses 

similar to those described above for unducted propulsors are required. 

This document contains the Final Report and User's Manual for the 3D Euler Ducted Propfan Analysis codes developed by the Allison Gas Turbine Division
5 

— 

U,

3

CC CC 

U) 

U.'

C))

<z -w 

LU

c 

to 0 Z 

Z cc-0 

LU L)N>

U) CD 

cn cc 

U) 

z 

LU —j U) 

U) U) 0 Z LU Z 

CD 3 

0 LU 

U) U- —;

—J 

U)

cr 

co 

C.) 

<U) 

z >-

= >-

U.' 

z 

z = 

0 C 

U)

Cl) 

C.) 

Cl) 

U 

U 

U 

E 

cd 

C 

V 

ce 

0 

UJUJ 

cc

V 

C-C 

bO 

H 

C) 

Li 

-J 

U 

C/) 

Cl) 

uJ 

C) 

= 

C-) 

C-) 

>-

LU 

LL 

CL–

Cl-
ORIGINAL PAGE

COLOR PHOTOGRAPH 

6 

Figure 2.2: Ducted propfan (very high bypass fan) geometry
7 

of the General Motors Corporation under Task 1 of NASA Contract NAS3-25270. 

This study was directed at the development of a three-dimensional Euler analysis for 

high-speed ducted propfans. A finite volume multiple-block four-stage Run ge-Kutta

numerical algorithm was utilized to predict the aerodynamics of ducted fan flows. Of 

particular interest was the ability to accurately predict the aerodynamics of the cowl. 

Two different grid arrangements were tested. The first system involved a sheared 

H-type grid which utilized a branch cut opening for the cowl. The second system was 

based on a numerically coupled multiple-block grid arrangement with a body-centered 

C-type grid about the cowl. Predicted results using both grid systems were compared 

with experimental data for two cases: a single-rotation low speed ducted propeller, 

and a high-speed 1.15 pressure ratio fan (both single-rotation and counter-rotation 

arrangements). This program is considered the initial step in developing a complete 

aerodynamic model of the flowfield through a ducted fan and the external flow over 

the engine nacelle and fan cowl. 

To predict the flow about a ducted propfan using the analyses described in this 

document, the necessary sequence is: 

1. Generate a numerical grid for the domain of interest. 

2. Run the Euler code to predict the steady aerodynamics. 

3. Process the results as needed. 

Separate sections are provided in the chapters which follow to describe the basis 

and operation of the codes used in the first two steps. A theoretical development of 

the grid generation algorithms for both the H-grid and C-grid architectures is given in
8 

Chapter 2. The 3D time-marching Euler analysis is detailed in Chapter 3. In Chapter 

4, a summary of the predicted results and verification studies performed to validate 

the accuracy of the analysis is presented. A summary of the conclusion of this study 

is given in Chapter 5. A user's manual for the grid generation code is provided in 

Appendix A. The operation of both the H-grid and multiple-block C-grid 3D Euler 

computer codes is outlined in Appendix B. 

It is worthwhile mentioning that most of the plotting and graphical postprocessing of the solutions are performed on graphics workstations. Presently, the PLOT3D 

[10] graphics software developed at NASA Ames Research Center is being extensively 

used for this purpose, and plot output has been tailored for this software.
9

3. GRID GENERATION ALGORITHM 

In this chapter, the numerical algorithm forming the basis of the ducted propfan 

analysis grid generation schemes is described. The geometry and computational domain are briefly described in the first section below. During the course of this study, 

two types of grids were tested for the numerical calculations. The grid systems are 

divided into two categories: H-grid, and multiple-block C-grid. Separate sections are 

provided to describe each category. 

3.1 Computational Domain 

The problem of interest is a ducted propfan geometry in steady operation as 

shown in Fig 2.2. Geometric parameters are expressed in a cylindrical coordinate 

system. The axial coordinate lies along the propeller axis of rotation, the radial 

coordinate is perpendicular to the axis of rotation, and the circumferential coordinate 

sweeps in the counter-clockwise direction when viewed down the axis of rotation. This 

coordinate system is illustrated in Fig. 3.1. The hub and duct contours are assumed 

to be axisymmetric surfaces (no circumferential variation). It is further assumed that 

the airfoil tips are fully enclosed within the duct (although this limitation could be 

easily relaxed). For any number of blades, under the assumptions of steady flow and
10 

a periodic geometry, the computational domain may be reduced to a single blade 

passage, circumferentially, as shown in Fig. 3.1. For counter-rotating configurations, 

a single blade passage is used for each blade row, and the multiple passages are then 

numerically coupled in the flow solver. The problem is further simplified by fixing the 

computational domain to the rotating blade elements. The resulting flow predictions 

are therefore based on the steady flow relative to the rotating blade. 

In order to implement a finite volume numerical solution of the steady aerodynamics for this computational domain, the region of interest must be subdivided 

into a finite number of smaller elements within the overall region of interest. To 

enhance the quality of the numerical solution, we further specify that these subdivisions be constructed in a relatively smooth and orderly manner with some arbitrary 

limitation on both the relative spacing between cells and the total number of cells 

in the computational domain. Each of these subdivisions will be used to define an 

elemental computational cell which will form the geometric basis of the finite volume 

aerodynamic solver to be described in a later chapter. 

In this study, two methods for constructing the computational cells were developed. The first method, referred to as the H-grid, is based on a single-block array of 

computational cells ordered in an H-type fashion. The second method, referred to as 

the multiple-block C-grid, utilizes a coupled system of five blocks each containing an 

array of computational cells. The five-block system is constructed by first subdividing 

the overall computational domain into five subregions, then further subdividing the 

individual subregions. The multiple-block arrangement permits the use of a special 

C-type grid wrapped about the leading edge of the duct to enhance the accuracy of
11 

Figure 3.1: Ducted propfan analysis computational domain
12 

the aerodynamic calculations in this critical region. 

Both the H-grid and C-grid algorithms are described in the sections below.

3.2 H-Grid Algorithm 

The H-grid generation scheme is based on a procedure originally developed by 

Mulac [11] for multiple blade row flow calculations using the average-passage [12] 

system of equations. In solving the average-passage system of equations, it is a requirement that the individual blade row grids have identical axisymmetric projections. 

This same technique was also applied in this study. The grid generation technique 

applies equally to single-rotation or multiple-rotation geometries. The approach here 

is to first construct a two-dimensional axisymmetric mesh in the z,r plane, followed 

by separate constructions for each of the individual blade rows in the circumferential 

direction to determine the final three-dimensional meshes. 

The two-dimensional meshes are constructed by dividing the axisymmetric plane 

into several subregions as shown in Fig. 3.2 for an unducted propfan and Fig. 3.3 for a 

ducted propfan. The regions of interest include the inlet, mid-inlet, blade, mid-exit, 

exit, between blade, tip gap, and outer flow regions. The mid-inlet and mid-exit 

regions result from the assumption that the blade leading and trailing edges are fully 

contained within the axial span of the duct. These regions are individually gridded 

in the two-dimensional plane to satisfy the conditions of common grid points along 

region interface boundaries, and a suitable distribution of points within each region.
13 

Unducted Propfan H-Grid Generation Subregions 

Between Blade Region 

Figure 3.2: Axisymmetric plane projection grid generation subregions for an tinducted propeller
14 

Ducted Propfan H-Grid Generation Subregions 

Between Blade Region 

Figure 3.3: Axisymmetric plane projection grid generation subregions for a ducted 

propeller
15 

The various regions differed in their individual construction. The distribution of 

points in any region is given axially and radially by one of three point distribution 

routines. These routines are described in the paragraphs below. All interpolations of 

coordinates were performed using a cubic spline interpolation routine. 

The distribution of the points on the axisymmetric blade plane in both the radial 

and axial direction is determined by a simple geometric progression radiating from a 

symmetric centerline referred to as packing algorithm #1 described below. 

Packing algorithm 1 

Given:

D total projected length of cubic spline curve 

RAT maximum ratio of adjacent cell projected lengths 

M number of points distributed across D (must be odd) 

the initial cell length of the progression is calculated as 

D 

DELT= 7

(3.1) 

(RAT 2 ' 

The point distribution is then given by 

P(i) = L (DELT)(RAT)2 ; i = 1,(M/2) + 1 (3.2) 

j=1

Packing algorithm 2 

The second construction utilizes a slightly different packing algorithm. In this 

case, D, RAT, and DELT are specified along with N, the total number of parallel
16 

spline curves, leaving M to be determined. The value of M is calculated iteratively 

until the the following conditions are satisfied: 

DM 

(DELT)max = E(RATr 

i=1 

(DEL T

D\ M 

T)mjn = (RAT)' 

i=1 

where

/ D \ / D \ 

DELT) max = max (DELT)j=1,N 

i D \ i D \ 

( DELT) min

= mm 

%DELT)=1,N

(3.3)

These conditions guarantee that neither the adjacent cell ratios, nor their inverses, exceed RAT. The resulting point distribution is then given by 

P(i) = E (DELT)(RAT) 1 ,i = 1,M (3.4) 

j1 

This approach is used axially in the inlet and exit regions, and radially in the 

outer boundary region.

Packing algorithm A3 

The third packing algorithm is used to determine the axial point distribution 

between blade rows. 

Given:

D total projected length of cubic spline curve 

RAT maximum ratio of adjacent cell projected lengths 

DELTTE initial cell width (forward blade row trailing edge)
17 

DELTLE initial cell width (alt blade row leading edge) 

N total number of parallel cubic spline curves 

the following values are then calculated: 

DELTTE DTE 

= (DELTTE+DELTLE) (D) 

DELTLE DLE = ( DELTTE + DELTLE) (D) 

The value of M is determined iteratively until the following conditions are satisfied: 

ID' M 

( DELT) max ^ RAT' (3.7)

i=1 

ID' M 

DELT) min E RAT' (3.8)

i=1 

where:

=max 1,N (3.9)

I D' (DTE)(DLE) ,j \DELT/max DELTTE DELTLE = 

ID'

= rain ( DTE

)(DELTLE D)DLE

DELTTE \ELTj rain = 1N (3.10) 

Again, this ensures that neither the adjacent cell ratios, nor their inverses exceed 

RAT. The point distribution is then given by: 

P(i) = L (DELTTE)(RAT)J, i =1, M (3.11)

j=1 

P(i) = L (DELTLE)(RAT)J , i = 1,M (3.12)

j=1 

The blade region is generated initially. This region is defined by the axisymmetric 

projection of the blade hub and tip boundaries, and the blade leading and trailing 

(3.5) 

edges. The number of points defining the blade region axially and radially is specified.
18 

The point distributions are determined by packing algorithm #1. Once the blade 

points are determined, the full inlet region pictured in Fig. 3.2 is computed next using 

packing algorithm #2, regardless of whether a ducted or unducted configuration is 

being generated. Once this region has been discretized, a simple linear interpolation is 

performed to determine the coordinates of the radial grid line which exactly intersects 

the duct leading edge. This grid line then determines the boundary between the midinlet region and the inlet region in Fig. 3.3. Once this boundary is determined, the 

upper boundary of the mid-inlet region is constructed by extrapolating the blade tip 

grid line axially, and maintaining the known tip gap spacing. Packing algorithm #1 

is then used to distribute points axially in exactly the same manner in the mid-inlet 

region as the blade region, with the additional constraint that the grid lines must 

match at the region interface. Once this region is constructed, the inlet region can be 

constructed using packing algorithm #2. The point distributions for additional blade 

rows (if necessary) are then determined using packing algorithm #1. For multiple 

blade rows, the next step is to calculate the region between blades using packing 

algorithm #3. The mid-exit and exit regions are constructed in the same manner as 

the mid-inlet and inlet regions by locating the grid line which exactly intersects the 

cowl trailing edge. The tip gap region is constructed algebraically with a specified 

number of equally spaced points in the radial direction between the blade tip grid 

line and the cowl lower surface. Upstream and downstream of the cowl, the cowl 

surface grid line is extended axially from the leading and trailing edges. It should be 

noted that this may impose some limitation on the maximum diameter of the hub 

upstream and downstream of the cowl. Finally, the radial distribution of points in
19 

the outer region is determined from packing algorithm #2. The grid lines which form 

the upper and lower surface of the cowl overlap upstream and downstream of the 

cowl. A duplicate grid line is placed on top of the cowl lower surface grid line to 

provide a natural storage location for the phantom points used in the application of 

the cowl surface boundary conditions (see Section 4.6). Once the axisymmetric plane 

grid has been calculated, the full three-dimensional grid is constructed individually. 

for each blade row. The blade twist and thickness distributions are used to determine 

the circumferential coordinates for the blade regions. The remaining circumferential 

coordinates are based on utilizing packing algorithm #1 in the blade-to-blade plane. 

The circumferential variation of grid points upstream and downstream of the blade 

regions is adjusted to provide a smooth transition such that the grid lines become 

parallel to the axis of rotation away from the blade region. 

A sample H-grid for a ducted propfan based on this grid generation technique 

is given in Fig. 3.4. The tip gap has been enlarged for this geometry in order to 

illustrate the tip gap grid. The emphasis in this development, as in Mulac [11], was 

to maintain a reasonable grid quality, rather than a specific number of grid points. 

Due to the algebraic nature of the scheme, the construction is performed relatively 

easily and inexpensively. Unfortunately, this approach can have drawbacks under 

certain conditions. The cubic spline interpolation . routine is subject to errors when 

a boundary or grid line has a discontinuous or sudden change in slope (particularly 

troublesome spots are the spinner or cowl leading edges). In addition, blades or cowls 

with a large thickness or blunt leading edge may not be well represented by this 

scheme when a relatively low number of grid points is used. Finally, the construction
20 

may occasionally have limits based on the number of grid points in a particular region 

due to the inability to satisfy one or more of the adjacent cell ratio constraints in the 

packing algorithms.

3.3 C-Grid Algorithm 

The C-grid generation scheme utilizes a number of concepts developed for the 

H-grid scheme. Again, the mesh is generated by first constructing a two-dimensional 

grid in the axisymmetric plane. In this case, a slightly different subdivision approach 

is utilized. The individual subregions and the five-block structure for the C-grid 

approach are listed and displayed in Fig. 3.5. Several of the subdomains are similar 

to those used in the H-grid construction, so the extension of the algorithm to the 

C-grid is obvious. In the case of the H-grid, the individual subregions are eventually 

combined into a single grid for use with the 3D H-grid Euler analysis. In the multipleblock C-grid generation, the subdomains are left divided and the individual flow 

solutions are numerically coupled. The blade, mid-inlet, inlet, mid-exit, exit, and 

outer regions are all constructed in basically the same manner for the C-grid as they 

were for the H-grid. In this case, the axial locations of the upstream boundary of 

the mid-inlet and the downstream boundary of the mid-exit regions are specified as 

a fraction of the cowl axial chord upstream and downstream of the cowl leading and 

trailing edges, respectively. The tip gap is extended upstream and downstream of the 

cowl, so there are no overlapping grid lines (resulting in the need for blocks 2 and 4), 

and the outer boundary extends radially from a grid line which is displaced from the 

cowl outer surface by a distance equal to the inner surface tip gap spacing. The near-
21 

Figure 3.4: Sample H-grid for a ducted propfan
cowl region is then discretized by the C-grid as shown in Fig. 3.5. The construction of 

blocks 1,2,4, and 5 follows directly from the development of the H-grid, so no further 

discussion is presented here. Instead, the underlying algorithm for the C-grid (block 

3) is presented in the paragraphs below. 

Grid points for the C-grid block surrounding the cowl are determined in a twostep procedure. In the first step, the inner and outer boundary points are specified. 

In the second step, the interior points are determined. Each step has several options 

each of which is described below. 

The exterior boundary points of the C-grid are determined by the surrounding 

four blocks. This is based on the requirement of coincident points along all block 

boundaries. The number of points along the inlet plane is directly specified, and 

must be an odd number. The number of points from the cowl surface to the outer 

boundary of the C-grid must also be specified, which fixes the number of points along 

the C-grid block exit plane. The cowl surface points are determined through one 

of two procedures. The first method performs a simple interpolation of the input 

coordinates based on arc length around the airfoil such that the relative arc length 

distributions of the airfoil surface points and the corresponding outer boundary points 

are the same. Unfortunately, this can often lead to highly skewed mesh lines near the 

airfoil surface. To circumvent this problem, an option was provided which allows the 

cowl surface grid points to "float" along the contour of the airfoil in a manner that 

would impose orthogonality at the cowl surface from a line extending from the cowl 

point to the corresponding outer boundary point. The movement of the airfoil grid 

points is controlled by a secant iteration procedure which optimizes the orthogonality
23 

Ducted Propfan Multiple-Block C-Grid Subregions 

Blade Region 

Figure 3.5: Axisymmetric plane projection subregions for multiple-block C-grid generation
24 

of the surface grid point location as a function of arc length along the airfoil while 

maintaining a smooth transition through neighboring surface grid points. The new 

airfoil coordinates are determined from the updated value of arc length through a 

linear interpolation of the arc length and airfoil coordinates originally specified. The 

resulting grid thus possesses more desirable orthogonality characteristics along the 

airfoil contour, which is desirable in terms of solution accuracy. 

The secant iteration procedure described above for the cowl surface grid points 

is expressed as:

-

+ (0.0 -O)(S' - sk_1) 

- (QC - Ok_i) (3.13) 

where Si is arc length measured clockwise around the cowl from a fixed reference 

location (1) to the point (i), k is the secant iteration count, and Oi is the measure of 

nonorthogonality. The orthogonality measure was based on the change of arc length 

Sgrjd between the outer point and the corresponding inner point. The optimized 

point occurs when O = dS9Td/dSk = 0. 

The starting values for the secant iteration are determined by using an initial 

point distribution obtained from the relative arc length cowl surface grid point interpolation routine (method one, above). 

In order to avoid overlapping grid lines and to maintain stability, the new surface 

grid point locations were never allowed to migrate more than one third of the distance 

between the previous surface location and that of the neighboring grid points. 

Interior grid points are determined through one of two methods. The first method 

is a simple algebraic interpolation of coordinates between the inner and outer C-grid 

boundary point distributions. In the second method, the interior points are generated
25 

through the numerical solution of a set of elliptic equations controlling a weighted 

distribution of grid smoothness, orthogonality, and grid point density based on a 

variational formulation originally developed by Brackbill and Saltzman [13]. A brief 

description of this scheme is given below. 

For a two-dimensional grid, the following integral expressions may be derived to 

evaluate critical aspects of the overall grid quality in physical space: 

Grid smoothness:

= + (V7)2Jdzdr (3.14) 

Grid orthogonality:

Io 

= f f [(V^) . (Vr)]dzdr (315) 

Grid point density:

Iw = f f w(z,r)Jdzdr (3.16) 

where:

=(z,r) i7=77(z,r) (3.17) 

O(,77) ( 8(z, r) 3.18)

and where V is the Cartesian gradient vector operator. The term w(z, r) is a userspecified function the magnitude of which is proportional to the desired grid point 

density in physical space. 

Obviously, the smoothest possible grid is obtained when 1s is minimized, the 

most orthogonal grid is obtained when lo is minimized, and the grid with the most 

desirable point density is obtained when Iw is minimized. By minimizing a weighted 

sum of these terms, i.e:

113+C010+CWIW (3.19)
26 

the constants Co, Cw may be used to control the relative importance of orthogonality 

and point density, respectively, in the overall grid point distribution. By exchanging dependent and independent variables, and applying the concepts of variational 

calculus for minimizing functions using the Euler-Lagrange equations, the following 

nonlinear coupled set of equations results: 

b1 z + b2 z + b3z77,q + a1 r + a2r + a3rr = 2w Oz (3.20) 

a1 z + a2 z 77 + a3 z + c1 r + c2 r 71 + C3?1777 = _ j2 ... L OW 

2w Or (3.21)

where the coefficients (a1 , b1 , c1 , i = 1,3) are all functions of the coordinate derivatives 

as: 

a1 = a3 1+C0a01+Cwavi a2 = a32 +Coao2+Cwav2 a3 = a33+Coa03+Cwav3 

= 1i3 1 + C0b01 + Cwbj b2 = b32 + C0602 + Cwb 2 b3 = b33 + C0b03 + Civbv3 

Cl = c31 + C0c01 + Cwc1 C2 = C32 + C0c02 + Cwc2 C3 = C3 3 + C0c03 + Cwc3 

(3.22) 

where:

a31 = —(aa)a a32 = 2(aa)f3 a33 = —(aa)y 

b31 = (bb)a b32 = —2(bb)fi b33 = (bb)7 

c3j = (cc)a c32 = —2(cc)13 c33 = —(cc)'y 

a01 = Z7,?j a02 = zr17 + z77re a03 = zr 

bol = Z

2

b02 = 2(2zz17 + rrij) b03 =

77 ZZ 

C01 = r2 CO2 = 2(zz17 + 2r r17 ) CO3 = 
27 

avi = —z71r77 a.2 =

+ zllr av3 _ze 

bVj = r b,2 = _2rErv b3 = rE2 

Cvi = z77v2 = —2zE s,7 Cv3 =ZZ 

(aa)=zEr.+z,?r,, (bb)=r+r (cc)=z+z 

a = (z + y)/J3 /3 = (zx,7 + y)/J2 y = (z + y)/J3 (3.23) 

When CO 0 and C, 0 the more standard Laplace or Poisson grid generation 

schemes result. This system is more complex than the usual Poisson-based grid 

generation schemes [14], but is still solvable using standard relaxation techniques. 

In this case, an iterative successive overrelaxation Gauss-Seidel solution technique is 

applied to solve the finite-difference equations resulting from a second-order centraldifference approximation of the resulting equations. For example: 

- + Zj...4j

()2 (3.24) 

j+1 i+l,3 _1 + i_i,j_1 - i1,j+1 (3.25) 

() 

A sample C-grid mesh network is given in Fig. 3.6 for a ducted propfan geometry. While this arrangement is not remarkably different from the H-grid network, 

the elimination of the highly skewed grid lines from the cowl leading edge region is 

certainly a desirable quality.
28 

Figure 3.6: Sample multiple-block C-grid for a ducted propfan
29 

4. 3D EULER NUMERICAL ALGORITHM 

In this chapter, the numerical algorithm forming the basis of the ducted propfan aerodynamic analysis is described. Separate sections are presented for noñdimensionalization, governing equations, boundary conditions, and discretization and 

time-integration schemes.

4.1 Nondimensjonaljzatjon 

To simplify the numerical treatment, all variables in the numerical solution are 

nondimensionalized by reference values as follows: 

z=r V9 

Lref Lref Vref Vref Vref 

P p= ____ Lref 

Pref Pref ref 

The reference quantities are defined as follows: 

Lref is determined from the maximum diameter of all the blade rows 

Pref is determined from the freestrearn total pressure 

Pref is determined from the freestrearn total density 

aref is determined from the freestream conditions

(4.1)
30 

= s/7Pref/Pref 

Vref is determined from the freestream acoustic velocity as 

are Vref =

4.2 Governing Equations 

The Euler equations of gas dynamics form the basis for the numerical solution 

procedure. In cylindrical coordinates, the strong conservation form of the inviscid 

equations of gas dynamics may be expressed in a reference frame rotating with the 

propeller as: 

f" (,XQ)dV

Yt

+ fdA [)tPdAz + )dA + )(fl - rwQ)dA9] = / )tSdV + / AKdV 

(4.2) 

where:

p 

pVz 

Q = pvr (4.3) 

rpv9 

pet 

pVz Pvr Pt,9 

PVz+ 

F = PVZVT G = pv +P H = pVrV9 (4.4) 

?pVZV9 ?pVrV9 r(pv + p) 

pvH L pvH pv9H 

F = F(Q) 0 = G(Q) ff = H(Q) (4.5)
31 

101 

101 

K __

(1.0, outside blade row; 

- r = 1+ (085—Ops)N inside blade row

(4.6) 

0 11 2ir 

L• 0 .1 

Here, V, V, v9 represent the velocity components in the axial, radial, and circumferential directions, respectively. The total energy function, e, is defined as: 

et= ()+(1z+Vr+V9) (4.7) 

The total enthalpy, H, is related to the total energy by: 

P

(4.8) 

These equations are based on the average-passage system of equations developed by 

Adamczyk [15] for the prediction of time-averaged flows in counter-rotating and multiple blade row machines. The system of equations is closed through the specification 

of the blockage \) and source terms (S). The closure model developed by Adamczyk 

[16] was utilized in this study to define these terms. Further details are given in the 

references cited. 

The governing equations are applied to a generalized finite volume in physical 

space as shown in Fig. 4.1. The cell surface areas dAt, dAt, and dA9 are calculated 

using the cross product of the diagonals of a cell face, and the cell volume is determined 

by a procedure outlined by Hung and Kordulla [17] for generalized nonorthogonal 

cells.

For convenience, a special operator L(Q) is defined as: 

L() = - UA [)PdAz + )tdAr + k - rwQ)dA9] ) +1 .xsdv+J )KdV (4.9)
I ,J 

i,j,k+1

i+1,j+1,k+ 

32

i,j,k 

Figure 4.1: Three-dimensional finite volume cell
33 

The flux variables P, U, and A are determined at each grid cell interface by averaging

the cell-centered dependent variables from the individual finite volumes adjoining the

interface. The L operator represents the summation of fluxes of mass, momentum,

and energy for the control volume V under the influence of blockage and source terms.

4.3 Runge-Kutta Time Integration 

The time-stepping scheme used to advance the discretized equations is a fourstage Runge-Kutta integration. The solution proceeds as: 

Qi = Qfl - a1 L&t[L(Q') + D(Q')] 

Q2 = Qfl - cx2Lt[L(Q1 ) + D(Q")J 

Q = Qfl - a3Lt[L(Q2 ) + D(Q')J 

Q = Qfl - cz4Et[L(Q3 ) + D(Q'2)] 

Qfl+l = (4.10) 

where:

1 1 1 

a

1 = 1 a2 = a3 = a4 =1 (4.11) 

An artificial dissipation operator (D) has been added to control oscillations in the 

numerical solution. A further discussion on this term is given in the following section. 

Jameson [18] determined that this scheme is stable for all time steps satisfying the 

CFL-related time step limitation

CFL < 2v(2- (4.12)
34 

An acceleration technique known as local time stepping is used to enhance convergence 

to the steady-state solution. Local time stepping utilizes the maximum allowable time 

increment at each point during the course of the solution. While this destroys the 

physical nature of the transient solution, the steady-state solution is unaffected and 

can be obtained more efficiently. 

4.4 Artificial Dissipation 

The Runge-Kutta algorithm described above utilizes an artificial numerical dissipation term D(Q) to suppress the odd-even point decoupling often observed for 

central-difference formulations. This problem is especially prevalent near shock waves, 

and it has been observed that the formulation of the dissipative term can have a significant influence on the final numerical solution. Jameson [18] demonstrated that a 

dissipative system combining second and fourth difference smoothing terms can effectively eliminate undesirable numerical oscillations without destroying the accuracy of 

the solution. This dissipation operator is constructed in the following manner: 

1 . Dz(Q) = dj+ —d. (4.13) jk z—j,k 

- 

V.

[2 

1

tZQ. 1 . - dj+jk - (st) 1i+,j,k s+,k i+,j,ki+,j,k] (4.14)

where:

62 i IC2 max(vj+l,j,k, 'ij,k) (4.15)

i+ ,j,k 

C4 

1 = max (0, r.4- f2 (4.16) 

0
35 

IPi+1,j,k - 2Pi,j,k + Pi— 1,j,kI

"i,j,k (4.17) 

= IPi+1j,k + 2p1 ,j k + Pi-1,j,kI 

Typical values for the second and fourth difference damping constants are: 

21 1 

K X4 = (4.18) 

The complete dissipation operator D1,3k is constructed as the sum of the dissipation 

operators in each of the respective coordinate directions as: 

Dj,j, k = (Dz) j,j,k + (Dr) i,j,k + ( D9)1,j,k (4.19) 

4.5 Implicit Residual Smoothing 

An effective technique for extending the stability limit of many explicit schemes is 

implicit residual smoothing. Residual smoothing attempts to accelerate the propagation of changes in the dependent variables by filtering the residuals of the calculation 

(which may also be interpreted as the local time derivative of the computational solution) at each time step. By enhancing the transfer of information between grid 

points, calculation time steps much larger than the stability-limited values may be 

utilized. Residual smoothing was originally introduced by Lerat (see e.g. Hollanders, 

et al. [19]) for use with the Lax-Wendroff scheme and later applied to Runge-Kutta 

schemes by Jameson and Baker [20] as a technique to acclerate convergence for steadystate calculations. Jorgenson and Chima [21] applied this technique for unsteady flows 

and obtained a speedup of about 5 for their calculations without any change in the 

unsteady solution.
) T)*** 

ILi,j,k = (4.25) 

36 

Since the time rate of change of the dependent variables 8Q/Ot is in essence 

controlled by a residual operator R(Q) = L(Q) - D(Q), it would follow that any 

measure which accelerates the propagation of changes in the residual throughout 

the domain would ultimately enhance convergence. The implicit residual smoothing 

operator used in this study can be written as: 

(1 - Ez6zz)(1 - ErSrr)(l - = R,j,k (4.20) 

where the differencing operator S is expressed as: 

Szz RI,j,k = - 2Ri,j,k + R _1,j,k (4.21) 

A value of EZ = Err = coo = 2 is typically used. 

The reduction is applied along each coordinate direction separately as: 

= (1 - E z8zz) . 'RI,j,k (4.22) 

1,3

R"• - (1 - 6z5rr) 4R ,,k (4.23) zj7k- 

- (1 - 6zE99)'R (4.24 

' ) -

where each of the first three steps above requires the inversion of a scalar tridiagonal 

matrix. The residual smoothing operator is applied at the first and third stage during 

the four-stage Runge-Kutta algorithm. The time-marching scheme then becomes: 

Qi =Q"—a1(Q') 

Q2 = Qfl - a2R(Q1)
37 

Q3=Q—a3k(Q2) 

= a4R(Q3) 

Qfl+l 

= Q (4.26) 

The implicit residual smoothing operator applied in this context can increase the 

effective time step by a factor of 3 or more, and can result in a 50% reduction in 

overall CPU time.

4.6 Boundary Conditions 

Inflow and exit boundary conditions are applied numerically using characteristic 

theory. A one-dimensional isentropic system of equations is utilized to derive the 

following characteristic equations at an axial inflow/outflow boundary: 

ac__ - 

__ __ _ - (VZ - a) 

Oz 

= 0 (4.27) 

oc+ oc+ 

at + (t

9z 

+ a) = 0 (4.28) 

where:

- 2a 2a 

C = Vz - C' = V + (4.29) 7-1 

In order to efficiently process boundary information in the numerical solution, phantom cells are located just outside the computational domain to permit the unmodified 

application of the interior point scheme at near boundary cells. Boundary condition 

information is effectively introduced into the solution by properly controlling the 

dependent variables in the phantom cells while permitting the application of the 

standard interior point scheme at near boundary cells.
Op Pv 

Or r

(4.30) 

38 

For subsonic normal inflow, the upstream running invariant C is extrapolated 

to the inlet, and along with the equation of state, specified total pressure, total 

temperature, and flow angle, the flow variables at the boundary may be determined. 

A second boundary condition was available as an option which specifies the mass flow 

through the inlet, although this option is not recommended. At the exit, a static 

pressure is specified at the hub for internal flows, and at the outer boundary for 

external flows. The remaining pressures along the outflow boundary are calculated 

by integrating the radial momentum equation: 

In this case, the downstream running invariant C+ is used to update the phantom 

cells at the exit boundary. Far-field boundaries also use this characteristic technique 

based on whether the local flow normal to the boundary passes into or out of the 

domain. The solid surfaces (hub, cowl, airfoils) must satisfy flow tangency:

(4.31) 

In this case, we specify no convective flux through the boundary (an impermeable 

surface), and hence, only pressure is needed at the phantom cell. The pressure may 

be extrapolated, or updated using a variant of the normal momentum equation. In 

this case, extrapolation was found to be the most effective technique based on rapid 

convergence and adequate results.
39

4.7 Multiple-Block Coupling 

For the multiple-block C-grid scheme, the Euler solution is performed on a single 

grid block at a time. Special boundary conditions along block boundaries are therefore 

required to provide some transport of information between blocks. This transport is 

provided through a simple procedure which relies on the fact that the grid block 

boundaries have coincident grid points. Since the standard interior point scheme 

is performed at the first grid cell off the block boundary, all that is required is to 

determine the fluxes along the block boundaries themselves. Phantom points are 

provided along each side of a block boundary to accommodate the flux calculation at 

the block boundary cell face. The flow variables in the phantom cells are provided 

by interrogating the corresponding cell in the adjacent block when a particular value 

is required. After each stage of the Runge-Kutta integration for each block, the 

block boundary phantom cells are updated by utilizing the new dependent variable 

values from the adjacent blocks. This direct specification technique was compared 

to a characteristic-based method similar to the inlet and exit boundary condition 

routines, and was found to provide enhanced convergence behavior. 

Artificial damping is applied at the block boundaries by neglecting the fourth 

order derivative term. Implicit residual smoothing. is applied at the block boundary 

by imposing a zero residual gradient (i.e. (dR/dz) = 0.0) condition at the boundary.
40

4.8 Solution Procedure 

The numerical solution is performed in an identical manner for both the U-grid 

and multiple-block C-grid Euler analyses. Assuming that the numerical grid and flow 

parameters are known, the time-marching procedure may begin from some set of initial data. This initial data is initially specified as a uniform flow, or may be implied 

from a previous solution. The time-marching procedure is applied iteratively to update the flow variables as the solution proceeds. The solution is deemed converged 

when the average residual Rhas been reduced by a factor of io. 

In the case of a counter-rotation calculation, the overall process of obtaining a 

converged solution is applied iteratively for each blade row as the source terms for the 

adjacent blade rows are updated. Thus, the time-marching loop is contained in an 

overall cycle which updates these source terms. The counter-rotation solutions were 

deemed converged when the differences in residuals between cycles were less than 

io_2.
41

5. RESULTS 

Several numerical results from the 3D Euler analyses described in Chapter 3 are 

presented in the sections which follow. Initial results are given to verify the accuracy 

of the 3D Euler formulation for steady flow about unducted propfans. Calculations are 

presented for both the 8-bladed SRi propfan geometry, and a 2-bladed SR7 propfan 

for which experimental airfoil surface static pressure distribution data were available. 

Several calculations for ducted propfans are presented. The numerics of the ducted 

propfan Euler analyses are first verified through comparisons with a two-dimensional 

numerical prediction of a circular arc cowl about a zero thickness nonrotating propeller 

aligned with the incoming flow. A second demonstration of the analyses is then 

performed for a fictitious ducted SR7 propfan. Predicted results are compared with 

numerical data obtained from a frequency domain panel method analysis for ducted 

propellers. A third set of calculations is given for a low speed ducted propeller. The 

predicted cowl surface pressure coefficient distribution is compared with experimental 

data for this geometry. Finally, the analyses are extensively compared with high 

speed experimental data for a 1.15 pressure ratio ducted fan configuration. Each case 

is discussed separately in the sections which follow.
42

5.1 SR7 8-Bladed Propfan 

In order to verify the numerics of the 3D Euler solver, several initial calculations 

were performed for steady flows about unducted propfans. A number of experimental 

studies have been performed for the unducted SR7 propfan geometry ranging from 

scaled wind tunnel tests to in-flight measurements. The SR7 propfan utilizes 8 blades 

with 41 degrees of sweep at the tip. A special hub contour is employed to eliminate 

flow choking at the hub. A description of the SR7 design parameters is given in 

Fig. 5.1. A particularly useful set of experimentally determined performance data 

for the SR7 was given by Stefko et al. [22]. In this study, a 2/9 scale model of the 

SR7 propfan was tested in a wind tunnel over a wide range of blade setting angles 

and operating conditions. A 70x35x15 (axial z radial x circumferential) numerical 

H-grid (as shown in Fig. 5.2) was generated to simulate this flow. Calculations were 

performed over a wide range of advance ratios and blade setting angles for flight Mach 

numbers of 0.7 and 0.8. All of the calculations for the unducted propfans described 

below were performed using the 3D H-grid Euler analysis operating in the unducted 

mode (no cowl boundary conditions are applied). 

A comparison of the predicted and experimental power coefficient map for the 

SR7 propfan at a flight Mach number of 0.8 is given in Fig. 5.3. During the investigation of the SR7 geometry, it was immediately observed that the predicted results were 

extremely sensitive to the blade setting angle and required an accurate assessment 

of the deflected blade shape. Since no reasonable means was available to adjust the 

blade shape for each run, the blade angle was adjusted in a manner which permitted 

the predicted power coefficient to match the experimental power coefficient at the
43 

design operating conditions. All other runs were made relative to this blade setting 

angle. The predicted power coefficient map illustrated in Fig. 5.3 is in fairly good 

agreement with the experimental curves for each blade setting angle. The comparison 

appears to worsen for larger blade setting angles as the power coefficient was consistently underpredicted by about 10%. This may be due to the inviscid flow assumption 

in the present analysis. 

A plot of the predicted propfan surface static/total pressure ratio contours for a 

Mach number of 0.8, 3/4 blade setting angle of 60.2 degrees, and an advance ratio of 

3.06 is given in Fig. 5.4. This contour plot clearly illustrates the spinner and blade 

root stagnation regions, and a shock which extends along the suction surface of the 

propfan blade. 

A similar series of calculations was performed at a flight Mach number of 0.7. 

In this case, an examination of the blade spanwise loading profiles was performed. A 

comparison of the predicted spanwise blade elemental power coefficient distributions 

with experimental data is given for several loading factors for M = 0.7 for the 8-bladed 

SR7 propfan in Fig. 5.5. The predicted power loading factors were somewhat different 

from the experimentally measured values; however, it is obvious that the experimental 

data brackets the predicted spanwise distributions for each power loading value. 

5.2 SR7 2-Bladed Propfan Modane Tests 

A second, more extensive comparison of predicted results for steady unducted 

propfan flows was performed based on the test of the SR7 airfoil at the Modane wind 

tunnel test facility reported by Bushnell [23]. In this case, the propfan driver did
44 

Characteristic 

Number of blades 

Tip sweep angle, degrees 

Tip Speed (m/s) 

Power loading (kW/m**2) 

Activiy factor 

Integrated design lift coefficient 

Airfoils 

Ratio of nacelle maximum diameter to propeller diameter 

Cruise design Mach number 

Cruise design advance ratio 

Cruise design power coefficient

SR7 

8

4].

243.8

256.85

227

0.202

NACA 16 and 65/CA

0.35

0.8

3.06

1.45 

Figure 5.1: SR7 propfan design characteristics
45 

Figure 5.2: Axisymmetric plane projection of 70x35x15 SR7 grid
46 

2.50 

2.00 

1.50 

Cp

1.00 

0.500 

SR7 Power Coefficient Performance (M=0.8) 

A Design Point (60.2 degrees) 

0 Experiment (57.7 degrees) 

VVV

Experiment (60.2 degrees) 

V Experiment (63.3 degrees) 

V N 3D Euler (57.7 degrees) 

V A 3D Euler (60.2 degrees) 

3D Euler (63.3 degrees)

0.000+ 

2.50

A V

V A

V 

V 

A

V A

A 0 A 

0 A 

0 A 0

A 0

A 0

A

A 

0

P I 

3.00 3.50 4.00 

Advance Ratio (J)

V

V

V

V 

I V 

4.50 5.00 

Figure 5.3: Comparison of predicted, experimental, and design power coefficient distributions for 8-bladed SR7 propfan (M=0.8) 
V 

bO HI 

C) - r-

'— 0 0 0

It M w— LC) 0 

r ((U) 

CD 0 0 0 d 

0 

4—. 

cz 

cr 

CD 

C') 

U) 

cz 

Cl)

ORIGINAL PAGE

COLOR PHO1OGRAPJ1 

47

PAGE IS 

QUALITY 

0 

0 U 

0 

V 

'-4 

V 

1 

0 

.4-

U 

.4-

U 

'4-

'-I

U 
5 

Radius Ratio (r/R) 

48 

5.00

4.00

3.00

d(Cp)/d(r/F

2.00

1.00 

0.000 

0.

0 0 Experiment: Power Loading (Cp/J**3) 74.1% 

Experiment: Power Loading (Cp/J**3) 957% 

v Experiment: Power Loading (Cp/J**3) 124.0% 

p Experiment: Power Loading (Cp/J**3) 147.0% 

3D Euler: Power Loading (Cp/J**3) 66.0% 

3D Euler: Power Loading (Cp/J**3) 109.0% fl Pi iIr PrlA,r I r d,rri I('r/ I\ I 73 nO!

Figure 5.5: Comparison of predicted and experimental blade spanwise elemental 

power coefficient loading distributions for the 8-bladed SR7 propfan 

(M=O.7)
49 

not have enough power to drive the 8-bladed configuration, so a 2-bladed version was 

tested instead. The airfoil surfaces were instrumented to permit the measurement of 

steady airfoil surface static pressures. The geometry and 120x45x25 H-grid used for 

these calculations are given in Fig. 5.6. This grid system utilized 49 points axially and 

21 points radially on the airfoil surface, which should be considered near minimum 

to accurately predict airfoil surface pressure distributions. 

Predictions were compared for a 3/4 radius blade setting angle of 30.58 degrees. 

The freestream Mach number was 0.2, with an advance ratio of 0.881. Predicted airfoil 

surface static pressure coefficient distributions are compared with experimental data 

at spanwise locations of 9.9%, 55.6%, and 94.4% (r/R) in Figs. 5.7-5.9, respectively. 

The agreement between experiment and prediction is generally very good, with some 

small disagreement near the leading edge of each cross section. This discrepancy may 

be due to several factors. As in the 8-bladed SR7 calculations, the exact deflected 

airfoil shape and blade setting angle are unknown (experimental uncertainty is on the 

order of 1 degree), which introduces a significant source of error in the calculation. 

In addition, this particular flow condition gives rise to a leading edge vortex which 

bends across the suction surface of the blade. This feature is evident in both the 

experimental data and the calculation. 

In order to qualitatively compare the predicted airfoil surface pressure distributions, a comparison of the predicted and experimental airfoil surface static pressure 

contours is given in Figs. 5.10-5.11 for the pressure (face) and suction (camber) surfaces, respectively. The pressure surface shows little variation; however, the suction 

surface displays a noticable region of high pressure (presumably a result of the leading
50 

edge vortex) which bends across the airfoil surface near the leading edge from roughly 

50%-100% span. This feature is also evident in the experimental static pressure contours, and illustrates the remarkable similarity between prediction and experiment. 

5.3 Symmetric Circular Arc Cowl Test Case 

In order to verify the accuracy of the Euler formulation for a ducted propfan 

geometry, a numerical test case was initiated based on a symmetric double circular 

arc airfoil cowl about a zero thickness, nonrotating propeller with a large hub diameter. A single passage of this geometry is given in Fig. 5.12. The propeller blades are 

aligned with the incoming flow and are therefore "invisible" to the inviscid numerical 

solution. These aspects, combined with the large hub diameter and large number of 

propeller airfoils (320) essentially reduce the aerodynamics to two dimensions. Therefore, the 3D Euler results were compared with an extensively verified 2D Euler code 

based on the hopscotch [24] time-marching technique for a similar two-dimensional 

geometry. A comparison with numerical, rather than experimental, data should be 

more enlightening in this instance since discrepancies between experimental data and 

calculations may not necessarily be due to boundary condition errors (i.e., viscous 

effects, grid induced errors, unsteadiness, etc.). 

The 3D H-grid calculations were performed on a 91x41x5 grid, while the 2D 

calculations were performed on a 91x21 grid. A comparison of the 2D grid and 

the 3D grid axisymmetric projection is given in Fig. 5.13. Since the 2D problem is 

perfectly symmetric, only half of the fiowfield was computed. Each grid utilized 31
51 

Figure 5.6: 2-bladed SR7 propfan geometry and 120x45x25 grid for Modane test 

comparison
52 

2-Bladed SR7 Modane Test Comparison (M=0.2, 9.9% span) 

3.00 

2.00 

1.00 

Cp

0.000 

IsI,J 

-2.00 

0.000 0.200 0.400 0.600 0.800 1.00 

Axial Chord (x/Cx) 

Figure 5.7: Comparison of predicted and experimental airfoil surface static pressure 

coefficient distributions for 2-bladed SR7 propfan Modane test (9.9% 

span)
3.00 

2.00 

1.00 

Cp

0.00C 

-1.00 

-2.00 

53 

2-Bladed SR7 Modane Test Comparison (M=0.2, 55.6% Span) 

0.000 0.200 0.400 0.600 0.800 1.00 

Axial Chord (x/Cx) 

Figure 5.8: Comparison of predicted and experimental airfoil surface static pressure 

coefficient distributions for 2-bladed SR7 propfan Modane test (55.6% 

span)
-2 

Cp

54 

2-Bladed SR7 Modane Test Comparison (M=0.2, 94.4% span) 

nn 

u.'uu U.4UU 0.600 0.800 1.00 

Axial Chord (x/Cx) 

Figure 5.9: Comparison of predicted and experimental airfoil surface static pressure 

coefficient distributions for 2-bladed SR7 propfan Modane test (94.4% 

span)
ORIGINAL PAGE

COLOR PHOTOGRAPH 

55 

w 

a cq)

1i4AL PAGE IS 

OOR QUALITY 

E 

a) 

a-

>< 

Ui

D =J = 

= = D = 

N 

I- U)U1U ----

_ C C C C C C C = 

—J 

LU 

LU

C C = = = C C C 

CD C C C C C C CD C C C = = C 

C C N 

1 

C = - = = = 

Ci

C C

1= C C = 

= C 

= 

= C = = 

= C = = 

N C =

V 

Q 

'-4 

V 

cI 'I-I 

V0 

V 

--

- 

'4 

ftu '4 

0 

cd 

- 0 V •-

- 

V '4 

it cs —I 

V 

'4 

bO
ORIGINAL PAGE

COLOR PHOTOGRAPH 56 

ci) 

LU 

0 

C)

ORiGINAL PA' 

OF POOR QU 

a) 

E 

a) 

0 

>< 

LU

= = = -j 

-

C = C = = C = = C = C 

CDD 

cx C C C C = = = C C C C = = = C C C C === ====== 

:a-

....... 

CDCCD CCC 

C.,

V 

11 

I -

V 

— 

V . 

Vo 

-2 

co 0 

U H 

V V 

I-I P-4 

'-I-I 

o.I 

o 

n -+. 

. 0 

d•.. 

qrA 

—I 

V 

to 

Pk 

LW 

MC30 ===MM C 

N = 

L--------

- ii u r'- j = = C•• -C C = =
57 

points to define the cowl upper and lower surfaces. 

Calculations were compared for inflow Mach numbers of 0.5 and 0.7. A comparison of the predicted cowl surface grid line Mach number distributions is given in 

Fig. 5.14. The 3D Euler results for both the inner and outer cowl surfaces display 

excellent agreement with the 2D Euler code results. A similar comparison is given in 

Fig. 5.15 for the 0.7 Mach number case. For this case, a strong symmetric shock wave 

pattern forms on the cowl surface, affording an examination of the shock-capturing 

characteristics of the cowl boundary algorithm. In this case, small discrepancies exist 

between the 2D solution, and the upper and lower surfaces of the 3D solution, although the agreement is still generally very good. The rapid decrease in Mach number 

across the shock wave is evident in each calculation. The small discrepancies between 

the solutions, and between the upper and lower surfaces of the 3D Euler solution, are 

easily explained by the differences in the mesh spacing near the cowl. 

A similar comparison was performed for the multiple-block C-grid Euler solver. 

In this case, the five-block grid arrangement shown in Fig. 5.16 was generated for the 

numerical predictions. Again, the 91x21 2D grid is also displayed for comparison. 

The total number of grid points is roughly similar to the H-grid calculations with 

individual block sizes of (65x13x5), (16x13x5), (81x7x5), (16x13x5), and (65x14x5). 

A comparison of the predicted cowl surface grid line Mach number distribution from 

the 3D multiple-block C-grid Euler analysis and the 2D solution for an inflow Mach 

number of 0.5 is given in Fig. 5.17. Again, good agreement is observed between 

the 3D and 2D results. A second comparison of results is given for an inflow Mach 

number of 0.7 in Fig. 5.18. Again, the predicted Mach number distributions are in
58 

good agreement, and the effect of the shock which results from the higher inflow 

Mach number is clearly displayed. Static pressure contour plots for the axisymmetric 

plane projection of the 3D H-grid and multiple-block C-grid Euler results for the 0.7 

inflow Mach number case are given in Fig. 5.19. The individual block boundaries are 

outlined in this case to demonstrate the smooth transition of contour lines across the 

block interfaces. The H-grid and multiple-block C-grid results are very similar in spite 

of the differences in the solution technique. The asymmetry of the contour lines is a 

result of the differences in boundary conditions at the hub and outer boundaries. This 

figure serves to validate the accuracy of the block coupling procedure and indicates 

that no nonphysical behavior is introduced by this procedure. 

5.4 Ducted SR7 Demonstration 

An interesting comparison of results was based on the calculation of the flowfield 

about a fictitious ducted 8-bladed SR7 propfan. In spite of the fact that there are 

no experimental data for this configuration, these calculations were useful for preliminary debugging of the analysis for ducted configurations, and also provided an 

opportunity to compare the predicted Euler results with an existing frequency domain 

panel method also developed for ducted propfan flows. 

The geometry was based on the original 8-bladed SR7 geometry given in Section 

4.1 with a thin profile duct about the propeller. A relatively large tip gap (2% of the 

duct diameter) was imposed to permit an examination of the tip flow. 

The geometry and grid surfaces are shown in Fig. 5.20. The axisymmetric projection of the 69x45x15 H-grid generated for this case is illustrated in Fig. 5.21. The
59 

Figure 5.12: Symmetric double circular arc cowl ducted propeller test case geometry

2 ;n 

- i .uu -.uu 0.000 0.500 1.00 1.50 

MacI

61 

Circular Arc Cowl Test Case (M=0.5) 

Cowl Axial Chord (x/Cx) 

Figure 5.14: Comparison of 2D and 3D H-grid predicted cowl surface grid line Mach 

number distributions (M=0.5)
2 

iJ 

LII 

Mach I

62 

Circular Arc Cowl Test Case (M=0.7) 

2 cn 

-1.00 -.500 0.000 0.500 1.00 1.50

Cowl Axial Chord (x/Cx) 

Figure 5.15: Comparison of 2D and 3D H-grid predicted cowl surface grid line Mach 

number distributions (M=0.7)

64 

Circular Arc Cowl Test Case (M=0.5) 

2.50 

2.00 

1.50

Mach No.

1.00 

0.500 

-1.00 -.500 0.000 0.500 1.00 1.50

Cowl Axial Chord (x/Cx) 

Figure 5.17: Comparison of 2D and 3D multiple-block C-grid predicted cowl surface 

grid line Mach number distributions (M=0.5)
65 

2.50 

2.00

Circular Arc Cowl Test Case (M=0.7) 

3D Euler Cowl Surface 

2D Hopscotch 

1.50

Mach No.

1.00 

0.500.4 • 

['1.1.1']

I • 

'I 

I. 

I 

-1.00 -.500 0.000 0.500 1.00 1.50

Cowl Axial Chord (x/Cx) 

Figure 5.18: Comparison of 2D and 3D multiple-block C-grid predicted cowl surface 

grid line Mach number distributions (M=0.7)
66 

Figure 5.19: Comparison of the predicted static pressure contours for the axisymmetric plane projection of the 3D H-grid (top) and multiple-block Cgrid (bottom) Euler analyses of the symmetric double circular arc cowl 

test case (M=0.7)
67 

case tested was for a true 3/4 blade setting angle of 60.2 degrees, flight Mach number of 0.7, and an advance ratio of 3.06. Predicted results from the 3D Euler code 

were compared with predicted results from a frequency domain panel code developed 

for ducted propellers by Williams et al. [25]. A comparison of the predicted airfoil 

spanwise elemental power coefficient loading distributions is given in Fig. 5.22 for 

this case. A similar comparison is given for the predicted airfoil spanwise elemental 

thrust coefficient loading distributions in Fig. 5.23. The predicted power and thrust 

coefficient loading distributions for the unducted case are also given in Figs. 5.22-5.23 

for comparison. Unfortunately, the panel code does not accurately model the hub 

and cowl contour; however, the loading distributions are largely similar, although the 

Euler code predicts a somewhat larger aerodynamic loading near the blade tip region. 

The predicted propfan surface static pressure contours for the 3D Euler results are 

pictured in Fig. 5.24. The effect of adding the duct to the normally unducted blade 

is immediately evident in the predicted results. The ducted blade becomes highly 

loaded at the tip, and generates a rather strong shock which extends to the inner 

duct surface. The spanwise loading indicates that the addition of the duct results in 

a 30% increase in thrust and a 25% increase in power. 

A second 3D solution was performed for this geometry using the multiple-block 

C-grid Euler analysis. In this case, the block grid pictured in Fig. 5.25 was utilized. 

The individual grid blocks contained (69x21x15), (11x9x15), (103x5x15), (12x9x15), 

and (69x18x15) points, respectively. Again, the predicted results were compared to 

the frequency-domain panel method of Williams et al. [25]. A comparison of the 

predicted airfoil spanwise elemental power coefficient loading distributions is given in
68 

Fig. 5.26 for this case. A similar comparison is given for the predicted airfoil spanwise 

elemental thrust coefficient loading distributions in Fig. 5.27. The results exhibited 

by the multiple-block C-grid Euler analysis are largely similar to the H-grid Euler 

analysis results.

5.5 Low Speed Ducted Propeller 

The third ducted configuration tested was based on the low speed ducted propeller originally tested by Kruger [26]. Kruger's geometry was based on an axisymmetric elliptical centerbody with an 8-bladed NACA 23009-23012 (hub to tip) section 

propeller. Data were taken for the static pressure distributions on the centerbody and 

cowl surfaces. 

A 160x41x15 grid was generated for this geometry, as shown in Fig. 5.28. The 

axisymmetric projection of this grid is given in Fig. 5.29. The case presented corresponds to Kruger's cowl #1, propeller #1, and nacelle #1. The calculation was based 

on an advance ratio of 0.95 and a flight Mach number of 0.4275. This experimental 

arrangement was particularly convenient for numerical simulation since the duct was 

supported separately from the centerbody, and thus mounting struts which might 

cause aerodynamic interference were not required. 

A comparison of the H-grid predicted and experimental cowl surface static pressure coefficient distributions is given in Fig. 5.30. The somewhat jagged appearance 

of the computational results is due to a lack of smoothness in the coordinates used to 

define Kruger's cowl geometry (the cowl coordinates were obtained by manually digitizing a drawing of the cowl cross section). In spite of the crudeness of the geometric
ORIGINAL PAGE 

COLOR PHOTOGRApJl 69 

- . - - 

. .

' ./..; 

"1 • •'' 1• 

Er" "i, .' .. MJ 

$ 4 

•

-

/ 

Figure 5.20: Ducted 8-bladed SR7 geometry and 69x45x15 H-grid 

ORIGINAL PAGE IS 

OF POOR QUALITY 
70 

Figure 5.21: Axisymmetric plane projection of 69x45x15 H-grid for ducted SR7 

propfan
6.00.

dCp/d(r/R)

4.00. 

2.00. 

71 

Ducted SR7 Spanwise Load Analysis M=0.7 J=3.06 

10.0- a a 3D Euler H-Grid Ducted Propfan Analysis (Unducted) 

D----4J 3D Euler H-Grid Ducted Propfan Analysis (Ducted) 

Frequency Domain Panel Method (Unducted) 

---- Frequency Domain Panel Method (Ducted) 

8.00. 

[11.1.11] 

0.000 0.200 0.400 0.600 0.800 1.00 

Radius r/R 

Figure 5.22: Comparison of 11-grid predicted elemental blade power coefficient distribution with frequency domain panel method results for ducted 8-bladed 

SR7
Ducted SR7 Spanwise Load Analysis M=0.7 J=3.06 

2.50- p 3D Euler H-Grid Ducted Proptan Analysis (Unducted) 

B-----a 3D Euler H-Grid Ducted Propfan Analysis (Ducted) 

A Frequency Domain Panel Method (Unducted) 

----. Frequency Domain Panel Method (Ducted) 

2.00 

1.50

F - 

dCt/d(r/R) 0 

/ #' 

1.00

El 

['IoI.II

I I I I.- 0.000 0.200 0.400 0.600 0.800 1.00 

Radius nA 

Figure 5.23: Comparison of H-grid predicted blade elemental thrust coefficient distributions with frequency domain panel method results for ducted 8-

bladed SR7
0 

0 

U 

0 

4 

Cl) 

U) 

V 

0 

U 

U) 

V 

U 

Cl) 

0 

'-4 

-I

CID 

ORIGt4AL is

1 

OF pOOR

U 

V

1 

0 

ca

k 

V I

P-I 

U) 

CD C)) bO 

ci) 

4.-

IL

c 194, r- - - aD - L() 0 

cz

OCaD 

'- 0 0 0 0 0 d o d
74 

Figure 5.25: Axisymmetric plane projection of multiple-block C-grid for ducted SR7 

propfan
75 

Ducted SR7 Spanwise Load Analysis M=0.7 J=3.06 

10.0 .., o D3D Euler C-Grid Ducted Propfan Analysis (Unducted) 

B-----a 3D Euler C-Grid Ducted Proptan Analysis (Ducted) 

Frequency Domain Panel Method (Unducted) 

----. Frequency Domain Panel Method (Ducted) 

8.00. 

6.00.

dCp/d( r/R)

4.00. 

2.00- 

0.000_ 

0.000 0.200 0.400 0.600 0.800 1.00 

Radius nA 

Figure 5.26: Comparison of multiple-block C-grid predicted elemental blade power 

coefficient distribution with frequency domain panel method results for 

ducted 8-bladed SR7
76 

Ducted SR7 Spanwise Load Analysis M=0.7 J=3.06 

2.50 .., o a 3D Euler C-Grid Ducted Propfan Analysis (Unducted) 

13----D 3D Euler C-Grid Ducted Propfan Analysis (Ducted) 

Frequency Domain Panel Method (Unducted) 

---- Frequency Domain Panel Method (Ducted) 

2.00 

1.50 

dCt/d(r/R) 

1.00 

0.500 

[.1.1.1'] 

0.000 0.200 0.400 0.600 0.800 1.00 

Radius OR 

Figure 5.27: Comparison of multiple-block C-grid predicted blade elemental thrust 

coefficient distributions with frequency domain panel method results 

for ducted 8-bladed SR7
ORIGINAL PAGE

COLOR PHOTOGRAPH 

Figure 5.28: Geometry and 127x48x15 grid for Kruger low speed ducted propeller 

ORIGINAL PAGE Is 

OF POOR QUALITY
78 

Figure 5.29: Axisymmetric plane projection of 127x48x15 H-grid for Kruger low 

speed ducted propeller
79 

model, good agreement was obtained with the experimental data. 

A second calculation was performed for this case using the 3D multiple-block Cgrid Euler analysis. The axisyminetric plane projection of the calculation grid is given 

in Fig. 5.31. The respective grid block sizes for this calculations were (150x21x15), 

(65x13x15), (147x7x15), (19x13x15), and (150x15x15). Again, a comparison of the 

multiple-block C-grid predicted and experimental cowl surface local static pressure 

coefficient distributions is given in Fig. 5.32. A similar level of agreement was obtained 

from the multiple-block C-grid calculation, in spite of the geometric uncertainties. 

5.6 NASA 1.15 Pressure Ratio Fan 

Final calculations were performed for a 1.15 pressure ratio fan stage originally 

tested at NASA for a wide range of flows [27]-[30]. A description of the geometry is 

given in Fig. 5.33. This geometry is representative of a 25:1 bypass ratio turbofan engine fan stage. Both single-rotation and counter-rotation calculations were performed 

for this geometry. For each calculation, it was not possible to exactly match the experimentally measured mass flow through the fan due to the expense of iterating on 

the three-dimensional solution; therefore, predictions are based on the flight Mach 

number and estimated fan rotational speed alone. A three-dimensional view of the 

geometry and grid system overlay are given in Fig. 5.34. 

H-grid predictions were based on the 225x52x15 H-grid pictured in the axisymmetric plane projection given in Fig. 5.35. The blades were represented by 23 points 

in the axial direction and 21 points in the radial direction. Both blade rows are represented in the axisymmetric projection of the grid. Single-rotation calculations were
80 

Kruger Ducted Propeller Test (M=0.4275, J=0.95) 

3.00 - 

2.00- 

1.00_

0 

Cp

0.000_ \ \

3D Euler Cowl Pressure Coefficient Distribution 

Experimental Cowl Pressure Coefficient Distribution 

-1.00—_ 

0.000 0.250 0.500 0.750 1.00 1.25

Cowl Axial Chord (xICx) 

Figure 5.30: Comparison of H-grid predicted and experimental cowl surface static 

pressure coefficient distributions for Kruger low-speed propeller
81 

Figure 5.31: Axisymmetric plane projection of multiple-block C-grid for Kruger low 

speed ducted propeller
82 

Kruger Ducted Propeller Test (M=0.4275, J=0.95) 

3.00 

2.00 

1.00 

Op

0.000 

-1.00

0 

-2.00

3D Euler Cowl Pressure Coefficient Distribution 

0 Experimental Cowl Pressure Coefficient Distribution 

0.000 0.250 0.500 0.750 1.00 1.25

Cowl Axial Chord (xICx) 

Figure 5.32: Comparison of multiple-block C-grid predicted and experimental cowl 

surface static pressure coefficient distributions for Kruger low speed 

propeller
83

DMAX 54.61 (21.5) C 50.80 (20.0) Wind shieldD flID MAX 0.930 XH!DM 0.847 

X OA/DMAX 2.05 DFH!Dfl. 0.4° Total external 

wetted are, 

XFC!DMAX 1.51 DCORE/Dfl 0.425 20387 cm 

Pylon tic 0.1125

Pylon 

Dorsal 

DFI LJ 1 DMAX I 

D Hd L ___ 

Load cell 

XH 

XFC

XOA

Figure 5.33: NASA 1.15 pressure ratio fan stage geometry (dimensions in cm) 

M''IM 
ORIGINAL PAGE

COLOR PHOTOGRAPH 

84

ORIGINAL PAGE Is 

OF POOR QUALITY 

Figure 5.34: NASA 1.15 pressure ratio fan stage and 3D grid overlay
85 

performed by neglecting the thickness of the deswirl vanes. 

For the initial calculations, only the single-rotation rotor was modelled. A comparison of predicted and experimental cowl leading edge surface static pressure coefficient distributions for a flight Mach number of 0.75 and an advance ratio (J) of 

2.86 for the rotor alone is given in Fig. 5.36. The data are presented in an unwrapped 

form to expand the critical region of interest at the cowl leading edge. The agreement between prediction and experiment is outstanding. Each of the features of 

the experimental pressure distribution is well defined by the predicted results. Any 

discrepancies between experiment and prediction are thought to be due to small differences between the predicted and experimentally measured mass flow rate through 

the fan or grid inadequacies. A similar comparison of results is given in Fig. 5.37 for a 

flight Mach number of 0.85 and an advance ratio of 3.23. In this case, good agreement 

is observed over much of the cowl surface, except for a rather large discrepancy on 

the outer cowl surface downstream of the leading edge region. Experimental evidence 

indicates that a flow separation occurs for the outer cowl surface in this case, and the 

disagreement is a result of the inviscid flow assumption in the present analysis. This 

would suggest that a viscous flow analysis will ultimately be required for high speed 

ducted fan design and off-design analysis. 

In order to assess the impact of the deswirl vanes on the aerodynamic solution, the 

0.75 Mach number solution was repeated for a counter-rotation configuration using 

the average-passage approach. This demonstration was performed under a no-cost 

contract amendment to verify the average-passage formulation in the ducted propfan 

Euler solver. The solution procedure involves iteratively updating representative body
86 

Figure 5.35: Axisymmetric plane projection of NASA 1.15 pressure ratio fan stage 

225x52x15 H-grid
: 

87 

1.15 Pressure Ratio Fan Cowl Pressure Coefficient (M=0.75) 

3.00

I

13 Experimental Data 

I

I H-Grid Analysis 

I Throat 

2.00] Highlight 

I S 

1.00 

CP

0.000 

-1.00 

-2.00 

-.400 -.300 -.200 -.100E+00 0.000 0.100E+00

Area Ratio (A/A(MAX)-1 .0) 

Figure 5.36: Comparison of predicted and experimental cowl leading edge static 

pressure coefficient distributions for NASA 1.15 pressure ratio fan (Hgrid, rotor only, M=0.75)
1.00 

CP

-1.00

88 

1.15 Pressure Ratio Fan Cowl Pressure Coefficient (M=0.85) 

3.00

o Exp. Data 

H-Grid Analysis 

Throat 

2.00 Highlight 

I I 

$ S 

I I U 

I I 

I I 

I I 

I S 

-2.00 

-.400 -.300 -.200 -.100E+00 0.000 0.100E+00

Area Ratio (A/A(MAX)-1 .0) 

Figure 5.37: Comparison of predicted and experimental cowl leading edge static 

pressure coefficient distributions for NASA 1.15 pressure ratio fan (Hgrid, rotor only, M=0.85)
89 

forces for the adjacent blade in each individual blade row solution. The introduction 

of the axisymmetric cowl has no bearing on this solution procedure. 

Again, the grid pictured in Fig. 5.35 was used although now the thickness of the 

deswirl vanes is included. A comparison of predicted and experimental data for the 

cowl leading edge surface static pressure coefficient distributions, nozzle boattail surface static pressure coefficient distributions, and nozzle inlet radial total pressure ratio 

distributions is given in Fig. 5.38. Again good agreement is obtained for each of these 

stations, demonstrating the capability of accurately simulating the aerodynamics of 

a complete ducted fan engine configuration. 

An illustration of the predicted fan surface static pressure ratio contours for this 

calculation is given in Fig. 5.39. The extreme flow gradients at the cowl leading edge 

and the aerodynamic interaction between the blade tips and cowl inner surface are 

well displayed in the predicted contours. 

The single-rotation fan geometry flowfield was also calculated using the multipleblock C-grid Euler solver. Predictions were based on the multiple-block grid pictured 

in Fig. 5.40. Grid block sizes were generated to be compatible with the total number 

of points used in the H-grid calculation, with block sizes of (172x25x15), (68x17x15), 

(193x9x15), (17x17x15), and (172x21x15), respectively. A comparison of the predicted and experimental cowl leading edge surface static pressure ratios is given in 

Fig. 5.41 for a Mach number of 0.75 and an advance ratio of 2.86. Again, excellent 

agreement between experiment and calculation is observed. The embedded C-grid 

block permits a large number of grid points in the region of the cowl leading edge 

without incurring an excessive number of points in the remainder of the flowfield. The
cJ — 

dO 

C 

C, 

C.) 

C) 0 

C) 

C, 

U) a 

U) ac x 

C) < W

a 

= 

.! uixZo

CC 

o 

co a 

C 

.c C.) 

a) 

o 

CD 

0 

U) U) 

0 ca 

h 

0.. o

C)

a H 

o o 

dO

•0 

U) CY CU C;

. 

in - o

a 

o a 

° 8

.2 

I.

8 

a 

2 

co

-

E 

I.-

U)

a a 

•

.2 co 

C) 

'-

a 0 

0..

a

0 a 

w 

C) 

C)

•

0 

CY 

U, 

Cl) a I I

a 8 8 8° 

jU!d4Od

90 

L() 0

II 

z 0 

Fa 

w 

rr 

az 

U0 

p 

w 

rr 

C') 

0) 

w 

cr 

0- 

Lo 

1 

C/) z
V 

I-I 

U) 

U) 

V 

p-I 

z 

U) 

s-I 

o8d 

o 

1-Il.': 

V bO 

1.4 cd 

U) 0) 

U) 

cd

V 

p-I 

U 

+bO U) 

U - 

S-I bO 

o 

cd 

U) 4 

-d 

V 

.4-

U 

.4p.4 

GNAL PAGE IS 

POOR QUALITY 

0 

crci 

V 

p-I 

bO 

a) 

U) 

U) 

ci) 

ci

ORIGINAL PAGE

COLOR PHOTOGRAPH 

91

C) C) (0 (0 U) C') - C) 

i.o c'
92 

high shear associated with the H-grid in this region is also eliminated. The solution 

was found to have a smooth continuous transition across block boundaries, and a 

relatively rapid convergence (1500 iterations) was achieved. 

A second calculation was performed for the C-grid algorithm at a Mach number of 

0.85. A comparison of the predicted and experimental cowl leading edge surface static 

pressure coefficient distributions is given in Fig. 5.42. Again, excellent agreement 

is observed over most of the cowl surface, except near the downstream portion of 

the cowl outer surface. This is again thought to be due to an experimental flow 

separation, although the deviation is less pronounced in this case. The C-grid results 

do not differ significantly from the H-grid results for either of the cases examined, 

although a detailed grid resolution study has not been performed. In addition, the 

H-grid calculation was performed with a relatively large number of grid points about 

the cowl leading edge to enhance the accuracy in this critical region. As the number 

of grid points decreases, the C-grid analysis will likely be superior.
Figure 5.40: Multiple-block C-grid for NASA 1.15 pressure ratio fan stage (rotor 

only)
1.00 

OP

['IIIiIi 

.1.00 

94 

1. 15 Pressure Ratio Fan Cowl Pressure Coefficient (M=0.75) 

o Experimental Data 

I C-Grid Analysis 

I Throat 

2.00] Highlight 

-2.00 

-.400 -.300 -.200 -.100E+00 0.000 0.100E+00

Area Ratio (A/A(MAX)-1 .0) 

Figure 5.41: Comparison of predicted cowl leading edge static pressure coefficient 

distributions for NASA 1.15 pressure ratio fan stage (C-grid, rotor only, 

M=O.75)
1.00 

CP

0.000 

-1.00. 

95 

1.15 Pressure Ratio Fan Cowl Pressure Coefficient (M=0.85) 

3.00

D Exp. Data I C-Grid Analysis I Throat 

2.00] Highlight 

-2.00 

-.400 -.300 -.200 -.100E+00 0.000 0.100E+00 

Area Ratio (A/A(MAX)-1 .0) 

Figure 5.42: Comparison of predicted cowl leading edge static pressure coefficient 

distributions for NASA 1.15 pressure ratio fan stage (C-grid, rotor only, 

M=O.85)
PRECEDING PAGE BLANK NOT FILMED 

Nrl 

6. CONCLUSIONS 

A three-dimensional Euler analysis and grid generation scheme have been developed for the numerical analysis of ducted propfan flowfields. Both H-grid and 

multiple-block C-grid solution structures have been examined. The analysis has been 

verified through comparisons with experimental results from both a low speed ducted 

propeller and a high speed 1.15 pressure ratio fan stage for both single-rotation and 

counter-rotation configurations. At least one case was observed which indicated that 

a viscous flow analysis might improve the accuracy of the predictions. A comparison 

of the predicted results for both the H-grid and multiple-block C-grid analysis did 

not display significant differences, although as the number of grid points decreases, 

the C-grid analysis will likely be better. The capability of accurately simulating the 

aerodynamics about a complete ducted fan propulsor system has been demonstrated. 

Clearly the largest obstacle in accurately predicting propfan flowfields of any 

kind is in the specification of the deflected airfoil shape. Calculations for the SR7 

clearly displayed how even small errors in the blade setting angle could drastically 

change the predicted results. A second obstacle encountered in this study was a lack 

of experimental data for ducted fans which prevented a lengthy verification study for 

a range of geometries and flow conditions.
PRECEDING PAGE BLANK NOT FILMED 

REFERENCES 

[1] Barton, J. M., Yamamoto, 0., and Bober, L. J., "Euler Analysis of Transonic 

Propeller Flows," AIAA Journal, Vol. 3, June 1987, pp. 277-282. 

[2] Usab, W. J., Lee, A. K., and Sullivan, J. P., "A Comparison of Numerical Simulation and Experimental Measurements of Flow through Propellers," AIAA 

Paper 88-0367, 1988. 

[3] Celestina, M. L., Mulac, R. A., and Adamczyk, J. J., "A Numerical Simulation 

of the Inviscid Flow Through a Counterrotating Propeller," ASME Paper 86-

GT-136, 1986. 

[4] Whitfield, D. L., Swafford, T. W., Janus, J. M., Mulac, R. A., and Belk, D. 

M., "Three-Dimensional Unsteady Euler Solutions for Propfans and CounterRotating Propfans in Transonic Flow," AIAA Paper 87-1197, 1987. 

[5] Saito, S., Kobayashi, H., Nakamura, Y., and Matsuo, Y., "Predicted Flow Field 

Around the Advanced Propeller at Take-Off," AIAA Paper 88-3151, 1988. 

[6] Matsuo, Y., Arakawa, C., Saito, S., and Kobayashi, H., "Navier-Stokes Computations for Flowfield of an Advanced Turboprop," AIAA Paper 88-3094, 1988. 

[7] Kobayawa, M., and Hatano, I., "Flow Field around a Propeller by Navier-Stokes 

Equation Analysis," AIAA Paper 88-3150, 1988. 

[8] Sacks, A. H., and Burnell, J. A., "Ducted Propellers—A Critical Review," 

Progress in Aeronautical Sciences Vol III, Pergamon Press Inc., New York, 

1962, pp. 87-135. 

[9] Groeneweg, J. F., and Bober, L. J., "NASA Advanced Propeller Research," 

NASA TM 101361, 1988.
100 

[10]Walatka, P. P., and Buning, P. G., "PLOT3D User's Manual,", rough draft for 

NASA TM, 1988. 

[11] Mulac, R. A., "A Multistage Mesh Generator for Solving the Average-Passage 

Equation System," NASA CR 179539, 1988. 

[12]Adamczyk, J. J., Celestina, M. L., Beach, T. A., and Barnett, M., "Simulation of Three-Dimensional Viscous Flow Within a Multistage Turbine," ASME 

Paper 89-GT-152, 1989. 

[13]Brackbill, J. B., and Saltzman, J. S., "Adaptive Zoning for Singular Problems 

in Two Dimensions," Journal of Computational Physics, 46, 342-368, 1982. 

[14]Thompson, J. F., "Elliptic Grid Generation," Numerical Grid Generation 

(Edited by Joe Thompson), Elsevier Science Publishing Co., 1982. 

[15]Adamczyk, J. J., "Model Equation for Simulating Flows in Multistage Turbomachinery," ASME Paper 85-GT-226, 1985. 

[16]Adamczyk, J. J., Mulac, R. A., and Celestina, M. L., "A Model for Closing the 

Inviscid Form of the Average-Passage Equation System," ASME Paper 86-GT227, 1986. 

[17]Hung, C. M., and Kordulla, W., "A Time-Split Finite Volume Algorithm for 

Three-Dimensional Flow-Field Simulation," AIAA Paper 83-1957, 1983. 

[18]Jameson, A., Schmidt, W., and Turkel, E., "Numerical Solutions of the Euler Equations by Finite Volume Methods Using Runge-Kutta Time-Stepping 

Schemes," AIAA Paper 81-1259, 1981. 

[19]Hollanders, H., Lerat, A., and Peyret, R., "Three-Dimensional Calculation of 

Transonic Viscous Flows by an Implicit Method," AIAA Journal, 23, pp. 1670-

1678, 1985. 

[20]Jameson, A., and Baker, T. J., "Solution of the Euler Equations for Complex 

Configurations," AIAA Paper 83-1929, 1983. 

[21]Jorgensen, P. C. E., and Chima, R. V., "An Unconditionally Stable RungeKutta Method for Unsteady Flows," NASA TM 101347, 1989.
101 

[22] Stefko, C. L., Rose, G. E., and Podboy, G. C., "Wind Tunnel Performance 

Results of an Aeroelstically Scaled 2/9 Model of the PTA Flight Test Propfan," 

NASA TM 89917, 1987. 

[23] Bushnell, P., "Measurement of the Steady Surface Pressure on a Single Rotation 

Large Scale Advanced Prop-Fan Blade at Mach Numbers from 0.03 to 0.78," 

NASA CR 182124, 1988. 

[24] Hall, E. J., "Simulation of Time-Dependent Compressible Viscous Flows Using 

Central and Upwind-Biased Finite-Difference Techniques," Ph.D. dissertation, 

Iowa State University, 1989. 

[25] Williams, M. H., Cho, J., and Dalton, W. N., "Unsteady Aerodynamic Analysis 

of Ducted Fans," (to appear in the AIAA Journal of Propulsion and Power). 

[26] Kruger, W., "On Wind Tunnel Tests and Computations Concerning the Problem of Shrouded Propellers," NACA TM 1202, 1949. 

[27] Steffen, F. W., "Cruise Performance of an Isolated 1.15 Pressure Ratio Turbofan 

Propulsion Simulator at Mach Numbers from 0.6 to 0.85," NASA TM X-3064, 

1974. 

[28] Wesoky, H. L., Abbott, J. M., Albers, J. A., and Dietrich, D. A., "Low-Speed 

Wind Tunnel Tests of a 50.8 Centimeter (20-in.) 1.15 Pressure Ratio Fan Engine 

Model," NASA TM X-3062, 1974. 

[29] Wesoky, F. L., and Steffen, F. W., "Wind Tunnel Tests of a 20 in. Diameter 

1.15 Pressure Ratio Fan Engine Model," NASA TM X-71445, 1973. 

[30] Osborn, W. M., and Steinke, R. J., "Performance of a 1.15 Pressure-Ratio 

Axial-Flow Fan Stage with a Blade Tip Solidity of 0.5," NASA TM X-3052, 

1974.
PRECEDING PAGE BLANK NOT FILMED 

103 

APPENDIX A. GRID GENERATOR OPERATING INSTRUCTIONS

A.1 Introduction 

Grid generation for both the H-grid and C-grid Euler analyses is performed by a 

single program with identical input formats. The grid generation source program is 

written in FORTRAN 77, and has been used successfully on Cray UNICOS and IBM 

VM/CMS mainframe computer systems, as well as Silicon Graphics 4D workstations 

using a UNIX operating system. Array dimensions in the programs are specified by 

a PARAMETER statement in the MAIN routine. The PARAMETER statement 

appears as:

PARAMETER( IJX=100, KX=15, NX=2) 

and appears in every subroutine. Another important parameter, PRECIS, is based 

on the precision of the host computer and is determined internally. 

The individual variables in the PARAMETER statement are defined as: 

IJX Used to determine the maximum number of points in the axial direction 

KX The maximum number of points in the circumferential direction 

NX The maximum number of blade rows 

Array storage is determined through a variable IXT which is defined as IXT = 

(2NX + 1)IJX. The parameter IXT then determines the maximum number of
104 

points in the axial direction. Several rules for storage must be observed for proper 

operation of the code. Storage is assigned under the assumption that IJX > KX. 

The parameter NX must be greater than or equal to the number of blade rows. For 

multiple-block C-grid generation, KX must be greater than or equal to JMAX + 1 

(JMAX is an input parameter described below). The first step in any grid generation 

sequence should be to estimate the size of the grid to be generated and to make sure 

that adequate array storage is provided for the anticipated grid. Unusual results can 

occur when an array size is exceeded, and many unexplained errors can be traced to 

inadequate array dimensions. 

The grid generation code can be used for either ducted or unducted, single- or 

counter-rotation propfans systems in the H-grid mode, and ducted single-rotation 

propellers in the C-grid mode. An option is provided to arbitrarily adjust the 3/4 

radius blade setting angle without modifying the blade geometry parameters. This 

assumes that the pitch change axis is perpendicular to the axis of rotation, and only a 

solid body rotation can occur (no deflections). This procedure is somewhat limited by 

the airfoil cross section coordinate distribution due to the requirement of maintaining 

the index of the true leading and trailing edge points. 

The program produces plotting output suitable for viewing using the PLOT3D 

graphics software package developed at the NASA Ames Research Center. This software permits workstation grid generation and processing before invoking the flow 

solver. PLOT3D output is written to FORTRAN unit 4 for H-grid generation, while 

the C-grid version writes individual block grid information to units 41-45 for the 

five-block architecture. Additional PLOT3D outputs can be generated for debugging
105 

purposes by small changes in the input. Some PostScript output is also available on 

FORTRAN unit 15, although this is primarily for analyzing the quality of the cowl 

C-grid in the multiblock grid generation. 

The input/output files used by the grid generation codes are described in the 

first section below. This section is then followed by a more general section which 

describes the function of the individual subroutines in the computer program. 

Geometric data may be input in either dimensional, or nondimensional form. 

However, each component must be defined in a consistent (either dimensional or 

nondimensional) manner. The final grid coordinates will be nondimensionalized by 

the maximum propeller diameter. 

A.2 Input/Output File Description 

An example input data file is shown in Fig. A.1. A brief description of the 

variables used in the data file is given below. 

VARIABLE DESCRIPTION 

TITLE An 80 character title for the mesh generation. 

NBLROW Number of blade rows: 

for C-grids, this must be 1, 

for H-grids, this must be less than or equal to PARAMETER NX. 

IGEOM Internal/external geometry parameter: 

if = 1, internal flow (e.g. compressor, turbine), 

if = 2, external flow (e.g. propfan, ducted propfan).
106 

+TITLE 

GENERIC DUCTED PROPFAN 

+NBLROW__++IGEOM___++ITHETA__+±IWRITE--++NCOWL---++JMAX++IORTH"-++IDEBUG 

1 2 2 2 1 5 1 

+NINPLS--++NEXPLS--++NBLPTZ--++NBLPTR--++NBLPTT++CITER++ISOLVE++------- -

1 1 15 15 7 50 0 

+ZINLET--++ZEXIT---++CFACLE--++CFACTE--++DFACLE--++DFACTE++--------++--------

-0.50000 1.09000 0.02 0.02 1.0 1.0 

+RTIN---++RATEX---++PATBB---++RATBLZ--++ATBLR--++RATBLT"++RATTOB++--------

1.32000 1.32000 1.20000 1.200000 1.1000 1.05000 1.36000 

****************************************************************************** 

+NHUB----++NTIP----++--------++--------++--------++--------++--------++--------

53 2 

+ZHUB, RHUB+--------++--------++--------++--------++--------++--------++_______ 

-80.00000 -20.00000 -0.45000 -0.40000 -0.37000 -0.35000 -0.33600 -0.330( 

-0.32685 -0.32336 -0.31773 -0.31097 -0.30302 -0.29383 -0.28285 -0.257 

-0.22232 -0.18130 -0.13606 -0.09774 -0.07757 -0.07042 -0.06853 -0.061 

-0.05049 -0.03837 -0.02810 -0.01930 -0.01181 -0.00522 0.00080 0.006 

0.01215 0.01808 0.02469 0.03237 0.04148 0.05238 0.06434 0.073: 

0.07788 0.07965 0.09454 0.11764 0.15665 0.21074 0.27769 0.354 

0.43958 0.58799 0.75179 0.92019 1.09000 

0.00200 0.00200 0.00200 0.00200 0.00260 0.00400 0.00600 0.008( 

0.01000 0.01300 0.02177 0.03076 0.03877 0.04635 0.05436 0.068 

0.08364 0.09408 0.09971 0.10164 0.10204 0.10239 0.10251 0.103( 

0.10421 0.10611 0.10819 0.11031 0.11257 0.11463 011653 0.118: 

0.12065 0.12307 0.12572 0.12874 0.13223 0.13630 0.14061 0.143 

0.14530. 0.14590 0.15081 0.15819 0.16767 0.17458 0.17754 0.174 

0.16806 0.15661 0.15465 0.15465 0.15465 

+ZTIP,RTIP+--------++--------++--------++--------++--------++--------++-------

-10.0000 10.000 

1.7000 1.700 ** *** ******** ******   **  **   **  **   **  **** ** ****** ************ ************* 

+NSHR----++NSHRIN--++NSFIROU--++NPBCAB--++--------++--------++--------++_______ 

75 7 7 3 

+ZSHR, RSHR+--------++--------++--------++--------++--------++--------++_______ 

0.250000 0.249464 0.247860 0.245194 0.241480 0.236734 0.230970 0.2242 

0.216504 0.207864 0.198340 0.187960 0.176774 -0.164834 -0.176774 -0.1879' 

-0.198340 -0.207864 -0.216504 -0.224220 -0.230970 -0.236734 -0.241480 -0.2451 

-0.247860 -0.249464 -0.250000 -0.249464 -0.247860 -0.245194 -0.241480 -0.2367 

-0.230970 -0.224220 -0.216504 -0.207864 -0.198340 -0.187960 -0.176774 -0.1648 

-0.152190 -0.138894 -0.125000 -0.110570 -0.09570 -0.080360 -0.064706 -0.0487 

-0.032630 -0.016350 0.000000 0.016350 0.032630 0.048774 0.064706 0.0803 

0.095670 0.110570 0.125000 0.138894 0.152190 0.164834 0.176774 0.1879 

0.198340 0.207864 0.216504 0.224220 0.230970 0.236734 0.241480 0.2451 

0.247860 0.249464 0.250000 

0.410000 0.410010 0.410040 0.410067 0.410083 0.410075 0.410047 0.4100' 

0.409958 0.409910 0.409857 0.409797 0.409725 0.404410 0.404362 0.4043 

0.404480 0.404638 0.404887 0.405247 0.405617 0.406077 0.406807 0.4074 

0.408125 0.408895 0.410000 0.410750 0.411665 0.413025 0.414677 0.4161 

0.417597 0.419175 0.420660 0.422092 0.423482 0.424795 0.426020 0.4271 

Figure A.1: Sample input data file for grid generation
107

0.428175 0.429077 0.429845 0.430467 0.430922 0.431187 0.431250 0.431113 

0.430775 0.430250 0.429552 0.428700 0.427725 0.426660 0.425538 0.424382 

0.423212 0.422045 0.420887 0.419752 0.418655 0.417597 0.416587 0.415635 

0.414738 0.413900 0.413117 0.412395 0.411742 0.411178 0.410720 0.410382 0.410162 0.410040 0.410000 

*********************** ********************************************************* 

-------------------------------------------------------------------------------+ 

6 8 21 1 0.0 0.0 60.2 

+ZBLA, RBLA, THBLA, TTBLA-------++--------++--------++--------++--------++--------+ 

-0.07042 -0.06259 -0.05477 -0.04699 -0.03924 -0.03153 -0.02386 -0.01628 -0.00878 -0.00130 0.00619 0.01346 0.02074 0.02804 0.03535 0.04269 0.05004 0.05742 0.06481 0.07222 0.07965 

-0.07840 -0.07040 -0.06241 -0.05443 -0.04648 -0.03857 -0.03070 -0.02287 -0.01513 -0.00745 0.00021 0.00786 0.01530 0.02277 0.03026 0.03777 0.04530 0.05285 0.06041 0.06800 0.07560 

-0.08927 -0.08096 -0.07266 -0.06436 -0.05608 -0.04782 -0.03959 -0.03140 -0.02326 -0.01519 -0.00719 0.00081 0.00875 0.01654 0.02436 0.03219 0.04005 0.04792 0.05581 0.06372 0.07164 

-0.10222 -0.09346 -0.08470 -0.07594 -0.06718 -0.05844 -0.04971 -0.04101 

-0.03235 -0.02372 -0.01517 -0.00667 0.00183 0.01025 0.01855 0.02690 0.03526 0.04363 0.05203 0.06043 0.06886 

-0.11092 -0.10175 -0.09257 -0.08339 -0.07421 -0.06504 -0.05588 -0.04674 

-0.03762 -0.02853 -0.01948 -0.01050 -0.00152 0.00741 0.01626 0.02510 0.03396 0.04283 0.05172 0.06062 0.06953 

-0.10567 -0.09655 -0.08743 -0.07831 -0.06919 -0.06008 -0.05097 -0.04189 

-0.03281 -0.02376 -0.01475 -0.00577 0.00321 0.01214 0.02102 0.02992 0.03883 0.04775 0.05668 0.06562 0.07457 

-0.08982 -0.08126 -0.07270 -0.06414 -0.05558 -0.04704 -0.03850 -0.02997 -0.02145 -0.01295 -0.00447 0.00401 0.01246 0.02089 0.02933 0.03777 0.04622 0.05468 0.06314 0.07160 0.08007 

-0.06358 -0.05613 -0.04869 -0.04125 -0.03381 -0.02638 -0.01895 -0.01153 -0.00411 0.00330 0.01070 0.01810 0.02549 0.03289 0.04029 0.04769 0.05509 0.06249 0.06990 0.07731 0.08472 0.10239 0.10294 0.10369 0.10469 0.10595 0.10746 0.10915 0.11119 0.11351 0.11588 0.11825 0.12120 0.12414 0.12704 0.12989 0.13269 0.13544 0.13813 0.14078 0.14337 0.14590 0.11840 0.11876 0.11930 0.12005 0.12104 0.12229 0.12378 0.12550 0.12754 0.12983 0.13216 0.13456 0.13753 0.14041 0.14325 0.14603 0.14877 0.15145 0.15409 0.15666 0.15919 0.14224 0.14233 0.14263 0.14312 0.14379 0.14471 0.14589 0.14732 0.14897 0.15095 0.15319 0.15546 0.15789 0.16079 0.16361 0.16639 0.16910 0.17176 0.17438 0.17693 0.17943 0.17874 0.17901 0.17911 0.17927 0.17967 0.18025 0.18105 0.18213 0.18345 0.18499 0.18689 0.18902 0.19118 0.19357 0.19638 0.19906 0.20169 0.20427 0.20679 0.20926 0.21167 0.22445 0.22475 0.22492 0.22504 0.22525 0.22562 0.22620 0.22701 0.22805 0.22934 0.23084 0.23273 0.23467 0.23678 0.23921 0.24169 0.24409 0.24644 0.24874 0.25098 0.25317 0.27612 0.27636 0.27645 0.27653 0.27679 0.27718 0.27773 0.27849 0.27945 0.28059 0.28199 0.28357 0.28518 0.28704 0.28914 0.29112

Figure A.1: Sample input data file for grid generation (continued)
0.29306 0.29497 

0.32968 0.32972 

0.33312 0.33418 

0.34365 0.34500 

0.38143 0.38163 

0.38513 0.38585 

0.39196 0.39278 

-0.22623 -0.11448 

0.15770 0.17116 

0.18416 0.18144 

-0.22180 -0.14123 

0.07632 0.09372 

0.14188 0.14653 

-0.21743 -0.16069 

0.02138 0.03888 

0.11278 0.12067 

-0.22009 -0.17785 

-0.01496 0.00299 

0.09651 0.10748 

-0.22349 -0.18910 

-0.03311 -0.01442 

0.09562 0.10941 

-0.19509 -0.16591 

-0.02054 -0.00199 

0.11322 0.12781 

-0.15033 -0.12553 

0.00442 0.02140 

0.13073 0.14494 

-0.09282 -0.07252 

0.03766 0.05229 

0.14816 0.16076 

-0.22623 -0.25766 

-0.19041 -0.16474 

0.05918 0.09269 

-0.22180 -0.23798 

-0.13762 -0.11651 

0.06803 0.09214 

-0.21743 -0.22212 

-0.10906 -0.08898 

0.06427 0.08517 

-0.22009 -0.21557 

-0.09505 -0.07598 

0.06456 0.08392 

-0.22349 -0.21227 

-0.08176 -0.06255 

0.07264 0.09166 

-0.19509 -0.18096 

-0.05226 -0.03375 

0.09554 0.11384 

-0.15032 -0.13572

0.29683

0.32988 

0.33533 

0.34632 

0.38192 

0.38668 

0.39357 

-0.06580 

0.17992 

0.17914 

-0.10573 

0.10771 

0.15096 

-0.13115 

0.05401 

0.12867 

-0.15193 

0.01943 

0.11820 

-0.16466 

0.00338

0.12280 

-0.14364 

0.01598

0.14191 

-0.10575 

0.03806 

0.15867 

-0.05588 

0.06671 

0.17300 

-0.25881 

-0.13510 

0.12376 

-0.22938 

-0.09377 

0.11529 

-0.20847 

-0.06796 

0.10544 

-0.19920 

-0.05637 

0.10277 

-0.19376 

-0.04321 

0.11059

-0.16248 

-0.01518 

0.13209

108 

0.29863 

0.33011 

0.33651 

0.34761 

0.38227 

0.38761 

0.39436 

-0.00593 

0.18148 

0.17788 

-0.06270 

0.11843 

0.15519 

-0.09805 

0.06720 

0.13656 

-0.12440 

0.03450 

0.12866 

-0.13966 

0.02025 

0.13576 

-0.12121 

0.03346 

0.15547 

-0.08594 

0.05435 

0.17184 

-0.03926 

0.08089 

0.18478 

-0.25715 

-0.10277 

0.15196 

-0.22032 

-0.06880 

0.13745 

-0.19620 

-0.04625 

0.12510 

-0.18427 

-0.03634 

0.12103 

-0.17620 

-0.02382 

0.12936 

-0.14470 

0.00331 

0.15024

0.30040 

0.33045 

0.33789 

0.34886 

0.38269 

0.38851 

0.39512 

0.03783 

0.18524 

0.17801 

-0.02902 

0.12448 

0.15910 

-0.07011 

0.07870 

0.14430 

-0.09990 

0.04855 

0.13889 

-0.11662 

0.03646 

0.14805 

-0.09998 

0.05037 

0.16812 

-0.06705 

0.07034

0.18402 

-0.02329 

0.09482 

0.19570

-0.25046 

-0.07418 

0.17801

-0.20653 

-0.04048 

0.15910 

-0.18074 

-0.02392 

0.14430 

-0.16738 

-0.01600 

0.13890 

-0.15761 

-0.00435 

0.14805 

-0.12627 

0.02187 

0.16813

0.33092 0.33153 0.3322 

0.33939 0.34084 0.3422 

0.38317 0.38374 0.3844: 

0.38940 0.39028 0.391L 

0.07747 0.11208 0.1384 

0.18762 0.18808 0.1867 

0.00160 0.03009 0.0550 

0.12892 0.13273 0.1371

-0.04451 -0.02039 0.0015: 

0.08821 0.09698 0.1050 

-0.07691 -0.05491 

0.06161 0.07367 

-0.09455 -0.07318 

0.05204 0.06701 

-0.07939 -0.05926 

0.06683 0.08275 

-0.04864 -0.03057 

0.08590 0.10119 

-0.00766 0.00771 

0.10855 0.12203 

-0.24091 -0.22859 

-0.04316 -0.01028 

-0.19103 -0.17463 

-0.01104 0.01729 

-0.16395 -0.14652 

-0.00143 0.02083 

-0.14978 -0.13197 

0.00444 0.02471 

-0.13876 -0.11988 

0.01505 0.03434 

-0.10775 -0.08929 -0.0707 

0.04033 0.05875 0.0771. 

-0.11844 -0.10164 -0.08457 -0.06750 -0.05050 -0.0335 

Figure A.1: Sample input data file for grid generation (continued) 
-0.01658 0.00034 

0.11762 0.13434 

-0.09282 -0.07935 

0.02374 0.03826 

0.13894 0.15326

0.01715 

0.15109 

-0.06428 

0.05271 

0.16758

109 

0.013399 0.05071 

0.16771 0.18402 

-0.04962 -0.03484 

0.06713 0.08149 

0.18182 0.19570

0.06742 0.08415 0.10087 

-0.02008 -0.00543 0.00917 

0.09586 0.11021 0.12457 

Figure A.1: Sample input data file for grid generation (continued)
110 

ITHETA Theta input format for blade definitions: 

if = 1, a mean camber line and tangential thickness are to be input, 

if = 2, surface theta coordinates are to be input. 

IWRITE Mesh write to disk format (not used in this version). 

NCOWL Cowl geometry parameter: 

if = 0 no cowl, no cowl geometry is read, 

if = 1 generate an H-grid for the ducted geometry (cowl geometry is 

read), 

if = 2 generate a multiple-block C-grid for the ducted geometry (cowl 

geometry is read). 

JMAX Defines the number of grid lines from the blade to the outer boundary 

of the axisymmetric projection of the C-grid block. 

(must be <PARAMETER KX). 

IORTH C-grid orthogonality control parameter: 

if = 0, cowl surface points for a C-grid are based on arc length interpolation, 

if = 1, cowl surface points for a C-grid are adjusted for orthogonality. 

IDEBUG Debugging parameter: 

if = 0 no extra debugging plots, 

if = 1 several extra PLOT3D data sets are produced for debugging. 

NINPLS Number of additional constant spacing inlet planes added to grid. 

NEXPLS Number of additional constant spacing exit planes added to grid. 

NBLPTZ Number of points on the airfoil in the axial direction
111 

(this must be an odd number). 

NBLPTR Number of points on the airfoil in the radial direction 

(this must be an odd number). 

NBLPTT Number of points between airfoils in the circumferential 

direction (this must be an odd number < PARAMETER KX). 

NCITER Number of iterations for the C-grid solver. 

if = 0 C-grid is generated using the algebraic solver. 

ISOLVE C-grid generation algorithm control parameter: 

if = 0 C-grid is generated using the Poisson solver. 

if = 1 C-grid is generated using the variational solver. 

ZINLET Axial location of the initial inlet boundary plane 

(before constant spacing inlet planes are added). 

ZEXIT Axial location of the initial exit boundary plane 

(before constant spacing exit planes are added). 

CFACLE Factor which determines the axial extent of the C-grid 

upstream of the cowl leading edge (multiplied by 

cowl axial chord). 

CFACTE Factor which determines the axial extent of the C-grid 

downstream of the cowl trailing edge (multiplied by 

cowl axial chord). 

DFACLE Factor which determines the radial width of the C-grid 

upstream of the cowl leading edge (if > 1, then 

larger than tip gap, if < 1, then smaller than tip
112 

gap). 

DFACTE Factor which determines the radial width of the C-grid 

downstream of the cowl trailing edge (if > 1, then 

larger than tip gap, if < 1, then smaller than tip 

gap). 

RATIN Maximum axial ratio of cells in the inlet region 

(must be> 1). 

RATEX Maximum axial ratio of cells in the exit region 

(must be> 1). 

RATBB Maximum axial ratio of cells between blade rows 

(must be > 1). 

RATBLZ Maximum axial ratio of cells on blades 

(must be > 1). 

RATBLR Maximum radial ratio of cells on blades 

(must be > 1). 

RATBLT Maximum circumferential ratio of cells blade to blade 

(must be > 1). 

RATTOB Maximum radial ratio of cells from tip (unducted H-grids) or 

the cowl upper surface (ducted H-grids) or the C-grid upper 

boundary (ducted C-grids) to the outer boundary 

(must be > 1). 

NHUB Number of X,R hub coordinate input pairs to be read. 

NTIP Number of X,R outer boundary coordinate pairs to be read.
113 

XHUB,RHUB Hub coordinates. 

All of the hub axial coordinates are given first (8 per 

line), followed separately by the radial coordinates. 

The number of points input is determined by NHUB. 

XTIP,RTIP Outer boundary coordinates. 

All of the outer axial coordinates are given first (8 per 

line, followed separately by the radial coordinates. 

The number of points input is determined by NTIP. 

The following 5 variables are used only for a ducted geometry. They should be 

deleted for an unducted case. 

NSHR Number of X,R cowl coordinate input pairs to be read. 

NSHRIN Number of axial grid lines between first blade row 

leading edge and cowl leading edge. 

(this must be an odd number). 

NSHROU Number of axial grid lines between last blade row 

trailing edge and cowl trailing edge. 

(this must be an odd number). 

NPBCAB Number of radial grid lines between lower cowl surface 

and blade tip. 

XSHR,RSHR Cowl coordinates. The coordinates are input in a 

wraparound fashion beginning at the trailing edge 

and proceeding in the clockwise direction. The 

trailing edge point must be duplicated as the last coordinate.
114 

All of the cowl axial coordinates are given first (8 per 

line), followed separately by the radial coordinates. 

The number of points input is determined by NSHR. 

The remaining variables define the propfan blades themselves. The blades are 

defined by radial slices with identical axial and radial coordinates on each side of the 

airfoil. The slices are given from the hub to the tip. The tip cross section must be 

specified in a manner which provides some clearance between the tip of the blade 

and the cowl lower surface, even if the two parts are actually mated. The 8 variables 

which follow are duplicated for each blade row when more than one row is considered. 

NBLA Number of blades in this blade row. 

NBLCRC Number of radial cross section used to define the blade. 

NPPRC Number of axial points per radial cross section used to 

define the blade. 

IPCH Pitch change trigger: 

if = 0 no pitch change is performed, 

if = 1 the pitch change variables ZPCA, THPCA, and BETA 

are read in and the blade setting angle is adjusted 

to the value specified by BETA. 

ZPCA Axial location of the blade pitch change axis. 

The pitch change axis is assumed perpendicular to 

the axis of rotation. 

THPCA Circumferential position of the blade pitch change axis. 

The pitch change axis is assumed perpendicular to
115 

the axis of rotation. 

BETA Absolute 3/4 radius blade setting angle. (Measured positive 

from the plane of rotation. The blade setting angle reference 

plane is the same for either positive or negative rotation 

[see Fig. A.21). 

ZBLA,RBLA, 

THBLA,TTBLA Blade coordinates. 

Each of the blade axial coordinates is given first 

(8 per line) for each cross section, cross sections are 

given from hub to tip) followed by all the radial 

coordinates. The blade circumferential coordinates are given 

next, and may be in two different formats depending 

on the variable ITHETA. The preferred format is 

ITHETA=2, where the actual tangential surface 

coordinates are given for the surfaces ordered 

in the counter-clockwise direction (for a 

positive rotation, this would be the suction 

surface first, then the pressure surface). 

The grid generation system produces two primary types of output files: the 

standard output and the grid or grids (depending on whether a C-grid or H-grid is 

being produced) themselves. A description of these files is given in Table A.1. 

A typical H-grid generated from the sample data set given in Fig. A.1 is shown 

in Fig. A.3. The blade tip gap has been expanded to illustrate the grid in this region.
116

,I YVUl _l I III II 

Figure A.2: Blade setting angle reference plane description
117

Table A.1: Description of input/output flies and UNIX-based filenames for H-grid 

generation 

Name Lu No. Description 

fort.4 4 Unformatted grid data set for first blade row. 

This may be viewed using PLOT3D, or used 

directly as input for the aerodynamic 

analysis. 

fort.14 14 Unformatted grid data set for second blade row. 

fort.24 24 Unformatted grid data set for third blade row (etc.). 

fort.98 98 Geometry and input for grid generator. 

fort.99 99 Standard output from grid generator.

The C-grid generation scheme utilizes the same input and standard output files, 

but produces a number of additional flies which are described in Table A.2. 

A typical C-grid generated from the sample input data file given in Fig. A.1 

modified by the condition NCOWL = 2 is shown in Fig. A.4. 

A.3 Subroutine Description 

A list of the grid generation program subroutines and their functions is given 

below for reference. A skeleton program flowchart is illustrated in Fig. A.5. 

SUBROUTINE DESCRIPTION 

CHGRID Main calling routine. 

ALPHM Plotting routine - specifically for alphanumeric labels, etc. 

000WL Routine for setting up and constructing C-grid boundaries. 

CGRID Routine for determining C-grid interior point distribution. 
118 

Figure A.3: Sample grid based on H-grid generation
119 

Figure A.4: Sample grid based on multiple-block C-grid generation
120 

Table A.2: Description of input/output files and UNIX-based filenames for C-grid 

generation 

Name Lu No. Description 

fort.4 4 Unformatted grid data set of axisymmetric projection 

of coupled grid blocks without C-grid. 

fort.15 41 PostScript plot data of C-grid quality. 

fort .41 41 Block #1 unformatted grid data set. 

fort.42 42 Block #2 unformatted grid data set. 

fort.43 43 Block #3 unformatted grid data set. 

fort.44 44 Block #4 unformatted grid data set. 

fort .45 45 Block #5 unformatted grid data set. 

fort .98 98 Geometry and input for grid generator. 

fort.99 99 Standard output from grid generator. 

fort.17 17 Unformatted grid data set for debugging. 

fort.27 27 Unformatted grid data set for debugging. 

fort.37 37 Unformatted grid data set for debugging. 

fort.47 47 Unformatted grid dataset for debugging. 

fort.57 57 Unformatted grid data set for debugging. 

fort .67 67 Unformatted grid data set for debugging.
121 

COWEX Routine which determines the outer boundary of the exit region. 

COWIN Routine which determines the outer boundary of the inlet region. 

CURPLT Plotting routine - specifically for plotting data as curves. 

END Plotting routine - specifically for ending a plot. 

ERROR Routine for internal evaluation of potential input data 

or grid errors. 

EULER Variational-based interior point solver for C-grids. 

GENBB Routine for generating grid point distribution between blades. 

GENBL Routine for determining blade grid point distribution. 

GENIN Routine for determining inlet region grid point distribution. 

GUESS Routine for generating an initial guess for C-grid interior 

point grid generation scheme. 

INRSCT Spline-based routine for determining the intersection of two 

curves. 

INTERP Polynomial-based interpolation routine. 

METRIC Routine for determining metric terms in C-grid generation. 

NUMBER Plotting routine - specifically for plotting numbers. 

NUMPTS Routine for determining the number of grid lines required for 

a specified grid cell ratio in a given region. 

POISON Poisson-based interior point solver for C-grid. 

P OLFIT Polynomial-based interpolation routine. 

PRINTO Output routine for grid generation. 

SCALE Plotting routine - specifically for generating scales.
122 

SPLINT Cubic spline polynomial interpolation routine. 

START Plotting routine - specifically for initiating a plot. 

WEIGHT Grid point weighting routine for variational grid generation. 

ZCOWL Routine for generating H-grids about cowls.
CHGRID

-- PRINT 123 

-- ERROR 

-- INTERP 

-- GENBL 

-- COWIN 

-- GENIN 

-- GENBB 

-- COWEX 

-- GENEX 

-- ZCOWL 

I -- COWIN 

I -- COWEX 

I -- CCOWL 

I -- CGRID 

I -- GUESS 

I -- WEIGHT 

I -- METRIC 

I -- POISON 

-- EULER 

I -- ANALYZ 

I 

I

-- START 

I

I 

-- SCALE 

-- CURPLT 

-- ALPHM 

--END

NUMPTS 

I --GENBL 

I -- GENIN 

I I 

I -- GENBB 

I I 

I -- GENEX 

I I 

I STOP 

Figure A.5: Program flowchart for ducted propfan analysis grid generation 
PRECEDING PAGE BLANK NOT FILMED

125

APPENDIX B. 3D EULER SOLVER OPERATING INSTRUCTIONS

B.1 Introduction 

This appendix contains the computer program User's Manual for the 31) Euler 

ducted propfan aerodynamic analysis. Again, as in the grid generation, both an Hgrid and a multiple-block C-grid solution scheme have been developed for the analysis 

of ducted propfans. In this case, separate programs are used for each of these analyses, 

although a number of similarities exist for the operation of each code. The H-grid 

flow solver can also be used as is for the analysis of unducted propfans. The H-grid 

version also has the capability of analyzing counter-rotating ducted propfans using 

the average-passage system of equations. 

The flow solver source programs are written in FORTRAN 77, and have been 

used successfully on Cray UNICOS and IBM VM/CMS mainframe computer systems 

as well as Silicon Graphics 4D workstations using a UNIX operating system. Array 

dimensions in the programs are specified by a PARAMETER statement in each code. 

The H-grid analysis parameters are specified in the statement: 

PARAMETER( IMX=125, JMX=25, KMX=16) 

which appears in every subroutine. 

The H-grid Euler analysis PARAMETER variables are defined as:
126 

IMX Number of axial grid lines + 1 

JMX Number of radial grid lines + 1 

KMX Number of circumferential grid lines + 1 

The C-grid Euler analysis code uses several PARAMETER statements all in the 

main routine given as: 

PARAMETER( NBLKS =5) 

PARAMETER( IMX1=125, JMX1=25, KMX1=16) 

PARAMETER( !MX2=125, JMX2=25, KMX2=16) 

PARAMETER( IMX3=125, JMX3=25, KMX3=16) 

PARAMETER( IMX4=125, JMX4=25, KMX4=16) 

PARAMETER( IMX5=125, JMX5=25, KMX5=16) 

The C-grid Euler analysis PARAMETER variables are defined as: 

NBLKS Number of grid blocks (this need not be changed) 

IMX1 Number of axial grid lines in block # 1 + 1 

IMX2 Number of axial grid lines in block # 2 + 1 

IMX3 Number of axial grid lines in block # 3 + 1 

IMX4 Number of axial grid lines in block # 4 + 1 

IMX5 Number of axial grid lines in block # 5 + 1 

JMX1 Number of radial grid lines in block # 1 + 1 

JMX2 Number of radial grid lines in block # 2 + 1 

JMX3 Number of radial grid lines in block # 3 + 1 

JMX4 Number of radial grid lines in block # 4 + 1 

JMX5 Number of radial grid lines in block # 5 + 1
127 

KMX1 Number of circumferential grid lines in block // 1 + 1 

KMX2 Number of circumferential grid lines in block # 2 + 1 

KMX3 Number of circumferential grid lines in block # 3 + 1 

KMX4 Number of circumferential grid lines in block # 4 + 1 

KMX5 Number of circumferential grid lines in block # 5 + 1 

The requirement that the PARAMETER variables be 1 element larger than the 

grid dimensions results from the use of phantom points outside the computational 

domain to impose the numerical boundary conditions. 

An important requirement here is that the first array dimension must be larger 

than the second and third array dimensions. This is based on the dimensioning scheme 

used for certain temporary variables in the code. This condition therefore imposes 

some limitations on the grid sizes that may be run without modifying the code. Under 

certain conditions, it may also be required that the second grid dimension be larger 

than the third. 

Approximate computational storage and CPU requirements can be estimated for 

either the H-grid or C-grid Euler codes from the following formulas: 

CPU sec 4.0X10T5(# grid points)(# iterations) 

Memory Mw 2.7z10 5 (# grid points) 

These formulas are valid for a Cray-II computer operating under the UNICOS 

environment and the cf77 compiler. 

The programs produce plotting output suitable for viewing using the PLOT3D 

graphics software package developed at the NASA Ames Research Center. PLOT3D 

output is written to FORTRAN unit 9 for the H-grid solver, while the C-grid version
128

writes individual block grid information to units 91-95 for the five-block architecture.

B.2 Input/Output Files 

An example input data file for the H-grid and multiple-block C-grid 3D Euler 

solvers is shown in Fig. B.1. A brief description of the form of the data file follows. 

It should be noted that the H-grid Euler solver was developed from an existing 3D 

Euler code with a number of acceleration features including axisymmetric averaging, 

reduced solution domain, grid coarsening and grid refinement. Many of the input 

variables listed below deal with controlling these functions. Due to the nature of the 

analysis for ducted propfans, however, several of these algorithms are not enabled. 

Therefore, a description of the controlling variables is provided for reference only, and 

should not be considered functional in the present analysis. 

VARIABLE DESCRIPTION 

MACH Freestream Mach number. 

GAMMA Specific heat ratio (cp/cv). 

PEXIT Freestream static pressure. This is applied at the

outer exit boundary and integrated radially inward 

for propellers (DUCT=O.0), or applied at the hub and 

integrated radially outward for compressors (DIJCT=1.0) 

to satisfy radial equilibrium. 

OMEGA Nondimensional rotational speed (this value will be 

adjusted to match the advance ratio when ADVR. 54 0.0). 
129 

+__MACH__++_GAMA__++_PEXIT__++_OMEGA__++ADR++DUCT+W 

0.8000 1.4000 0.65602 0.0000 -3.0600 0.0 1.0

+__CFL___++__VIS2__++__VIS4__++QFIL++BCBLDR0T 

-5.5 2.00 1.000 0.0000 0.0000 1.0 0.80 

+_FNCMAX_++__REST__++__SAVE__++_FISTEP_++FNPRNTR0WM 

720.0 1.0 1.0 1.0 1.0 1.0 0.0 

+_REFINI_++_REFINJ_++_REFINK_++EPSX++EPSYEPSZ 

1.0 1.0 1.0 2.00 2.00 2.00 0.0

+-RFINI I_++_RFINJJ_++_RFINKK_++__FAXI__++-FAXIMX++FR3DMX+FTURB 

1.0 1.0 1.0 0.0 99999.0 0.0 0020.0 

+_FAXIS__++_FAXIE__++_FSTART_++FEX1T+JBASE+0GS 1-++-OMEGS2-+ 

-10.0 -9.0 5.00 4.00 1.0 0.0000 0.0000 

Figure B.1: Sample input data ifie for ducted propfan 3D Euler solver
130 

ADVR Advance ratio (J = U/nD): 

if = 0.0, rotational speed is determined by OMEGA. 

DUCT Internal flow duct option: 

if = 0.0, external flow options are utilized (propeller), 

if = 1.0, an internal flow is assumed (compressor). No slip 

boundary conditions are applied at the outer 

boundary rather than the characteristic far field 

condition. 

COWL Cowl geometry trigger: 

if = 0.0, an unducted geometry is assumed, 

if = 1.0, duct boundary conditions are applied 

(only meaningful for H-grid flow solver). 

CFL Time step parameter: 

if < 0 then this is the CFL number used 

to determine the calculation time step for 

local time stepping, 

if> 0.0 then this is the CFL number used 

to determine the calculation time step without local 

time stepping. 

VIS2 Second order damping coefficient ( 2.0) 

(divided by 4 in the code) 

VIS4 Fourth order damping coefficient ( 1.0) 

(divided by 64 in the code)
131 

QFIL Not active in this version (0.0). 

BC Exit boundary condition trigger: 

if = 0.0, characteristic boundary condition based on Ptot is used, 

if = 1.0, characteristic boundary condition based on mass flow is used. 

BLDROW Blade row parameter. This value determines which blade row 

is being calculated during a multiple blade row solution. 

For single blade row calculations, this should be = 1.0. 

DTAXI Not active in this version. 

FNCMAX Maximum number of time steps to be performed. 

REST Restart option parameter: 

if = —1.0 a restart file and body force file are read, 

if = 0.0 no restart file, initialize variables in code, 

if = 1.0 a restart file is read. 

SAVE Save restart file option parameter: 

if = —1.0 restart files and body source terms are output 

at the end of a run, 

if = 0.0 no restart file is output at the end of the run, 

if = 1.0 a restart file is output at the end of the run. 

FISTEP Not active in this version (= 1.0). 

FNPRNT Not active in this version (= 0.0). 

ROWMAX Maximum number of blade rows in the current solution. 

For single rotation, = 1.0. For multiple blade rows, 

this parameter indicates that additional grids and
132 

body source terms are to be read. 

FPRAXI Print axisyinmetric solution, 

if = 0.0 axisymmetric solution is not printed, 

if = 1.0 axisymmetric solution is printed. 

REFINI Grid-coarsening trigger, not active in this version (1.0). 

REFINJ Grid-coarsening trigger, not active in this version (1.0). 

REFINK Grid-coarsening trigger, not active in this version (1.0). 

EPSX Implicit residual smoothing coefficient in the axial 

direction ( 2.0 is a typical value). 

EPSY Implicit residual smoothing coefficient in the radial 

direction ( 2.0 is a typical value). 

EPSZ Implicit residual smoothing coefficient in the circumferential 

direction ( 2.0 is a typical value). 

FPRBF Print trigger for body forces: 

if = 0.0, body forces are not printed, 

if = 1.0, body forces are printed. 

RFINII Fine mesh conversion, not active in this version (1.0). 

RFINJJ Fine mesh conversion, not active in this version (1.0). 

RFINKK Fine mesh conversion, not active in this version (1.0). 

FAXI Number of axisymmetric solutions, not active in this version (0.0).

FAXIMX Number of 3D passes before axisymmetric solutions, not active in this 

version (99999.0). 

FR3DMX Number of reduced domain 3D passes, not active in this version (0.0).
133 

FTURB Not active in this version (0.0).

FAXIS Not active in this version (0.0).

FAXIE Not active in this version (0.0).

FSTART Not active in this version (0.0).

FEXIT Not active in this version (0.0).

JBASE Not active in this version (0.0).

OMEGS1 Not active in this version (0.0).

OMEGS2 Not active in this version (0.0). 

The Euler solver system produces three primary types of output files: the standard output, plot output, and restart files. A list and description of each of these 

files is given in Table B.1 for the H-grid Euler analysis and Table B.2 for the C-grid 

Euler analysis.

B.3 H-grid Euler Solver Subroutine Description 

A list of the 3D H-grid Euler solver program subroutines and their functions is 

given below for reference. A skeleton program flowchart is illustrated in Fig. B.2. 

SUBROUTINE DESCRIPTION 

HPR03D Main calling routine. 

BAXI Calculates axisymmetric flow quantities. 

BAX12D Calculates axisymmetric flow quantities for 2D solution. 

BAXIA Axisymmetric flow averaging routine. 

BCCL Boundary condition routine for cowl lower surface.
134 

Table B.1: Description of input/output files and UNIX-based filenames for H-grid 

Euler solver 

Name Lu No. Description 

fort.1 1 New restart file output for first blade row. 

fort .2 2 Old restart file input for first blade row. 

fort.3 3 Old body force output for first blade row. 

fort .4 4 Unformatted grid data set for first blade row. 

This may be viewed using PLOT3D, or used directly 

as input for the aerodynamic analysis. 

fort.5 5 Main program input data file (may be redirected). 

fort.6 6 Main program output data file (may be redirected). 

fort.9 9 Unformatted flow data set. This may be viewed 

using PLOT3D. 

fort. 11 11 New body force output for current blade row. 

fort.12 12 Old restart file for second blade row. 

fort.22 22 Old restart file for third blade row (etc.). 

fort.13 13 Old body force output for second blade row. 

fort.23 23 Old body force output for third blade row (etc.). 

fort.14 14 Unformatted grid data set for second blade row. 

fort.24 24 Unformatted grid data set for third blade row (etc.).
135 

Table B.2: Description of input/output files and UNIX-based filenames for multiple. 

block C-grid Euler solver 

Name Lu No. Description 

fort.11 11 New restart file output for block #1. 

fort.12 12 New restart file output for block #2. 

fort.13 13 New restart file output for block #3. 

fort.14 14 New restart file output for block #4. 

fort.15 15 New restart file output for block #5. 

fort.21 21 Old restart file input for block #1. 

fort.22 22 Old restart file input for block #2. 

fort.23 23 Old restart ifie input for block #3. 

fort.24 24 Old restart file input for block #4. 

fort.25 25 Old restart file input for block #5. 

fort.41 41 Block #1 unformatted grid data set. 

fort.42 42 Block #2 unformatted grid data set. 

fort.43 43 Block #3 unformatted grid data set. 

fort.44 44. Block #4 unformatted grid data set. 

fort.45 45 Block #5 unformatted grid data set. 

fort.51 51 Main program input data file for block #1. 

fort.52 52 Main program input data file for block #2. 

fort.53 53 Main program input data file for block #3. 

fort.54 54 Main program input data file for block #4. 

fort.55 55 Main program input data file for block #5. 

fort .6 6 Main program output data file (may be redirected). 

fort.61 61 Block #1 printed output file. 

fort.62 62 Block #2 printed output file. 

fort .63 63 Block #3 printed output file. 

fort.64 64 Block #4 printed output file. 

fort.65 65 Block #5 printed output file. 

fort.91 91 Unformatted flow plot data set for block #1. 

fort.92 92 Unformatted flow plot data set for block #2. 

fort.93 93 Unformatted flow plot data set for block #3. 

fort .94 94 Unformatted flow plot data set for block #4. 

fort.95 95 Unformatted flow plot data set for block #5.

All unformatted grid and flow plot data sets may be 

viewed using PLOT3. 
136 

BCCU Boundary condition routine for cowl upper surface. 

BCEXT Boundary condition routine for exit cells. 

BCINL Boundary condition routine for inlet cells. 

BCJ Boundary condition routine for j=constant surfaces (hub and 

outer boundaries). 

BCJ2 Boundary condition routine for j=constant surfaces (hub and 

outer boundaries for 2D solution). 

BCK Boundary condition routine for the blade surfaces. 

BCWAL2 Axisymmetric boundary conditions at the exit for the 

2D solution. 

BCWALL Boundary condition routine for periodic cells. 

BFORCE Routine for calculating body force terms using the 

axisymmetric equations. 

BFREAD Routine to read in body force terms for all adjacent blade 

rows. 

ERROR Convergence-checking routine for Runge-Kutta solver (not used 

- for propfans). 

ERROR2 Convergence-checking routine for 2D Runge-Kutta solver. 

ERROR3 Convergence-checking routine for reduced domain 

3D solution. 

ERRORP Convergence-checking routine for Runge-Kutta solver for 

propfan calculations. 

EULER Runge-Kutta solver.
137 

EULER2 2D Runge-Kutta solver. 

EULER3 Reduced domain Runge-Kutta solver. 

FILTE2 2D artificial dissipation routine. 

FILTE3 Reduced domain artificial dissipation routine. 

FILTER Artificial dissipation routine. 

FIXUP Routine to correct boundary cells before output averaging 

is performed. 

GRIDG Routine to read and set up the grid. 

INITA Routine to initialize flowfield with axisymmetric solution. 

INITB Routine to replace flowfield with latest axisymmetric solution. 

INPUT Routine to read in input data and set up reference values. 

METRI2 Routine to calculate cell volumes and surface normals for 

2D solution. 

METRIC Routine to calculate cell volumes and surface normals. 

OP2D Routine to calculate axisymmetric residual of axisymmetric 

flow variables. 

OUTPUT Routine to print output and save restart files. 

PTCALC Routine to determine inlet total pressure profile. 

REFIN Routine to refine a 3D solution. 

REFIN2 Routine to refine a 2D solution. 

REFIN3 Routine to refine a reduced domain 3D solution. 

RESID Implicit residual smoothing routine. 

RESID2 2D implicit residual smoothing routine.
RiJNGE Runge-Kutta flux calculation routine. 

RUNGE2 2D Runge-Kutta flux calculation routine. 

STEP Routine to determine time step. 

STEP2' Routine to determine time step for 2D solution. 

TEST2D Routine to calculate the axisyminetric average of the 

residuals of the 3D flowfleld. 

VOLTJX Routine to calculate an individual cell volume. 

VOLIJX2 Routine to calculate an individual 2D cell volume. 

B.4 Multiple-Block C-Grid Euler Solver Subroutine Description 

A list of the multiple-block C-grid 3D Euler solver program subroutines and their 

functions is given below for reference. A skeleton program flowchart is illustrated in 

Fig. B.3 for the multiple-block C-grid 3D Euler solvers. 

SUBROUTINE DESCRIPTION 

CPR03D Main calling routine. 

BCI12 Boundary condition routine between block 1 and 2. 

BCI13 Boundary condition routine between blocks 1 and 3. 

BCI14 Boundary condition routine between blocks 1 and 4. 

BC123 Boundary condition routine between blocks 2 and 3. 

BC125 Boundary condition routine between blocks 2 and 5. 

BC133 Boundary condition routine along C-grid branch cut. 

BC134 Boundary condition routine between blocks 3 and 4.
HPRO3D 

-- INPUT 

I -- GRIDG 

I -- PTCALC 

I -- INIT 

I I 

-- INITA 

I I 

I -- BFREAD 

I I 

I -- INITB 

I -- METRIC 

I I I 

I I -- VOLUX 

I -- METRI2 

I -- VOLU2 

-- STEP 

-- EULER 

I -- FILTER 

I -- RtJNGE 

I I I 

I I -- RESID 

I -- BCINL 

I -- BCWALL 

-- BCEXT 

--BCJ 

-- BCK 

-- BCCL 

-- BCCU 

-- FIXUP

139 

-- EULER2 

I I 

I -- BAXI2D 

I I 

I -- TEST2D 

I I 

I --OP2D 

I I I 

I I --BCJ2 

I I I 

I I -- BCWAL2 

I I 

I -- STEP2 

I I 

I -- FILTE2 

I I 

I -- RUNGE2 

I I I 

I I --RESID 

I I 

I -- BCWAL2 

I I 

--BCJ2 

-- ERROR2 

-- BCINL 

-- BCEXT 

--BCJ 

--BCK

-- EULER3 

STEP 

I I 

I -- FILTER 

I I 

I -- RUNGE 

I I I 

I -- RESID 

-- BCINL 

-- BCWALL 

-- BCEXT 

BCJ 

-- BCK 

-- ERROR3 

-- REFIN 

I I 

I -- GRIDG 

I -- REFIN3 

I -- REFIN2 

I -- METRIC 

I I 

I -- VOLUX 

-- OUTPUT 

-- BAXI 

TEST2D 

-- BFORCE 

I -- BCJ2 

-- BCWAL2 

STOP -- ERRORP I 

Figure B.2: Program flowchart for ducted propfan analysis H-grid Euler solver
140 

BC135 Boundary condition routine between blocks 3 and 5. 

BC145 Boundary condition routine between blocks 4 and 5. 

BCINL Boundary condition routine for inlet cells. 

BCJ Boundary condition routine for j=constant surfaces (hub and 

outer boundaries). 

BCK Boundary condition routine for the blade surfaces. 

BCPER Boundary condition routine for periodic cells. 

ERRORP Convergence-checking routine for Runge-Kutta solver for 

propfan calculations. 

EULER Runge-Kutta solver. 

FILTER Artificial dissipation routine. 

FIXUP Routine to correct boundary cells before output averaging 

is performed. 

GRIDG Routine to read and set up the grid. 

INIT Routine to initialize solution. 

INPUT Routine to read in input data and set up reference values. 

METRIC Routine to calculate cell volumes and surface normals. 

OUTPUT Routine to print output and save restart files. 

PTCALC Routine to determine inlet total pressure profile. 

REFIN Routine to refine a 3D solution. 

RESID Implicit residual smoothing routine. 

RUNGE Runge-Kutta flux calculation routine. 

STEP Routine to determine time step.
141 

VOLIJX Routine to calculate an individual cell volume.
CPI tO3D 142 

INPUT 

-- GRIDG 

-- PTCALC 

-- INIT 

-- METRIC 

VOLUX

I -- BC133 

I -- ERRORP 

-- FIXUP 

I I 

I -- OUTPUT 

I I 

I -- P3DOUT 

I STOP 

Figure B.3:

-- STEP I 

-- FILTER I 

-- RtJNGE 

I I I 

I -- RESID I 

-- BCPER I 

-- BCBLA I 

-- BCCOWL 

-- BCHUB I 

-- BCINL I 

- - BCEXT I 

-- BCOUT 

-- BCI12 

-- BCI13 I 

-- BCI14 I 

-- BC123 I 

--BC134 I 

-- BC145 I 

-- BC125 I 

Program flowchart for ducted propfan analysis multiple-block C-grid 

Euler solver
USA Report Documentation Page Aflcs arc 

So" Ao,v.arsol 

1. Report No. 2. Government Accession No. 3. Recipients Catalog No. 

NASA CR-185217 

4. Title and Subtitle S. Report Date 

Investigation of Advanced Counterrotation Blade ConfigApril 1990 

uration Concepts for High Speed Turboprop Systems, Task 

1--Ducted Propfan Analysis

6. Performing Organization Code 

7. Author(s) S. Performing Organization Report No. 

Edward J. Hall EL/E' H&Z2 

Robert A. Delaney 10. Work Unit No. 

James L. Bettner 

9. Performing Organization Name and Address

Allison Gas Turbine Division of 11. Contractor Grant No. 

General Motors Corporation NAS3-25270 

P.O. Box 420

13. Type of Report and Period Covered 

Contractor Report 

Indianapolis, IN 46206-0420 

12. Sponsoring Agency Name and Address 

National Aeronautics and Space Administration Final 

14. Sponsoring Agency Code Lewis Research Center 

Cleveland, OH 44135-3191 

15. Supplementary Notes 

Prepared in cooperation with NASA Project Manager C. J. Miller, NASA Lewis Research 

Center, Cleveland, OH 

16. Abstract This report summarizes the progress of contract NAS3-25270, Task 1--Ducted 

Propfan Analysis. The time-dependent three-dimensional Euler equations of gas dynamics have been solved numerically to study the steady compressible transonic flow 

about ducted propfan propulsion systems. Aerodynamic calculations were based on a 

four-stage Runge-Kutta time-marching finite volumd solution technique with added numerical dissipation. An implicit residual smoothing operator was used to aid convergence. Two calculation grids were employed in this study. The first grid utilized 

an H-type mesh network with a branch cut opening to represent the axisyinmetric cowl. 

The second grid utilized a multiple-block mesh system with a C-type grid about the 

cowl. The individual blocks were numerically coupled in the Euler solver. Grid systems were generated by a combined algebraic/elliptic algorithm developed specifically 

for ducted propfans. The Euler analyses were applied for the calculation of both unducted and ducted propfan flows, and predicted results were compared with experimental data for several cases. The three-dimensional Euler analyses displayed exceptional accuracy, although predicted performance was very sensitive to geometric deflections. Both solution schemes were very robust and demonstrated nearly equal 

efficiency and accuracy, although the multiple-block C-grid formulation provided 

somewhat better resolution of the cowl leading edge region. 

17. Key Words (Suggested by Author(s)) 18. Distribution Statement 

Turbofans Unclassified-Unlimited 

Propfans 

Ducted propfans 

Euler - 

19. Security Classif. (of this report) 20. Security Classif. (of this page) 21. No. of pages 22. Price 

Unclassified Unclassified 159 A08

NASA FORM 1626 OCT 86