RV-8F

N772RW

Building, testing and operating an experimental aircraft

RV-8 Flight Performance Analysis

11 Jan 2017

I had to write a short paper for my propulsion class last semester, so I put together a quick and dirty (incomplete) performance analysis of a garden-variety RV-8. I used the open source Jupyter Notebook and the Python programming language to run the numbers. The code is pretty ugly in some spots, but the primary goal was to quickly generate some pretty graphs for the report. Comments and corrections are welcome.

This notebook is also available on GitHub.

Required Software

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import seaborn as sns
from skaero.atmosphere import coesa as atm
# Plot style with Seaborn
%matplotlib inline
sns.set_context("paper", font_scale=1.8, rc={"lines.linewidth": 3, "figure.figsize": (10,6)})
sns.set_style("ticks")

Van’s Aircraft RV-8 Specifications

Ref: https://www.vansaircraft.com/public/rv8specs.htm

# Wing span (ft)
b = 24

# Wing area (sq ft) 
S = 116

# Wing chord (constant, no taper)
c = S/b

# Aspect ratio
AR = b**2/S

# Gross weight (lbs)
W = 1800

# Wing loading (lb/sq ft)
WL = W/S

# Fuel capacity (US gal)
fuel_capacity = 42

Zero-Lift Drag and Oswald Efficiency

Initially, we cheat a bit here, using data derived by Kevin Horton in his AeroCalc Python package.

# Oswald efficiency factor estimation (ref: D. Raymer [2], Eqn. (12.49), p. 347.)
e = 1.78 * (1 - 0.045 * AR**0.68) - 0.64

# Oswald efficiency factor estimation (ref.: K. Horton)
e = 0.86

# Parasite drag coefficient, zero-lift (ref.: K. Horton)
C_d_0 = 0.0209

Free Stream Velocity

A range of likely velocity values for a light aircraft, from about 60 to 220 KTAS (knots true airspeed).

# Free stream velocity (ft/s)
V = np.arange(100,370)

U.S. Standard Atmosphere

The scikit-aero COESA algorithm is based on the 1976 U.S. Standard Atmosphere. The previous 1962 U.S. Standard Atmosphere is the standard used by the FAA and the flight test community for aircraft performance assessments, however the updated 1976 version differs from this only above 32 km [1].

Compare computed results here with an online Standard Atmosphere Calculator.

# COESA 1976 U.S. Standard Atmosphere (SI units!)

# Altitude to test (ft)
altitude = 23000

# Air density (kg/m^3)
rho = atm.density(altitude * 0.3048)

# Air density (slugs/ft^3)
rho = rho * 0.00194032

# Compare computed results with: https://www.digitaldutch.com/atmoscalc/
print("Air density at {:,.0f} ft = {:,.6f} slugs/ft^3".format(altitude, rho))
Air density at 23,000 ft = 0.001142 slugs/ft^3

The Drag Polar

Equations from Anderson [3].

\begin{equation} \tag{6.17} C_L = \frac{W}{\frac{1}{2} \rho_{\infty} V_{\infty}^2 S} \end{equation}

\begin{equation} \tag{6.1c} C_D = C_{D,0} + \frac{C_L^2}{\pi e AR} \end{equation}

def plot_drag_polar():
    """
    Plot drag polar for altitudes from sea level to Max_Alt in 4,000-ft increments.
    """
    # This range allows the +/- curves to meet in the middle
    V = np.arange(0,2000)
        
    # COESA 1976 U.S. Standard Atmosphere (SI units!)
    # Air density (kg/m^3)
    rho = atm.density(altitude * 0.3048)

    # Air density (slugs/ft^3)
    rho = rho * 0.00194032

    # Lift coefficient
    CL = W / (0.5 * rho * V**2 * S)

    # Drag coefficient
    CD = C_d_0 + CL**2 / (np.pi * e * AR)

    # Add a curve
    plt.plot(CD,CL)
    plt.plot(CD,-CL,color="#4c72b0")
 
    # Locate the axis labels and limits
    ax = plt.gca()
    ax.set_xlim([0,0.14])
    ax.set_ylim([-0.5,1.5])
    ax.xaxis.set_label_coords(0.5, -0.08)
    ax.yaxis.set_label_coords(-0.11, 0.5)
    
    # Add title and axis labels
    #plt.title("\n Drag Polar, RV-8, 1800 lbs \n")
    plt.ylabel(r'$C_L$', rotation=0)
    plt.xlabel(r'$C_D$')
    sns.despine()
    plt.show()

plot_drag_polar()

png

Thrust Required for Level Flight

Equations from Anderson [3].

\begin{equation} \tag{6.13} T = D = q_{\infty} S C_D \end{equation}

\begin{equation} \tag{6.14} L = W = q_{\infty} S C_D \end{equation}

\begin{equation} \tag{6.15} \frac{T}{W} = \frac{C_D}{C_L} \end{equation}

\begin{equation} \tag{6.16} T_R = \frac{W}{C_L/C_D} \end{equation}

def plot_thrust_required(V, Max_Alt=0, min_color="#4c72b0", drag_curves='N'):
    """
    Plot thrust required curves for altitudes from sea level to Max_Alt in 4,000-ft increments.
    
    Keyword arguments:
    V -- a 1D array of velocity values in ft/s
    Max_Alt -- an altitude, in feet, at which to stop plotting curves (default 0, sea level)
    min_color -- a hex color code for the min value marker
    drag_curves -- include drag curves? 'Y' or 'N'
    """
    legend = ['Sea Level']
    for altitude in range(0, Max_Alt + 1, 4000):
        
        # COESA 1976 U.S. Standard Atmosphere (SI units!)
        # Air density (kg/m^3)
        rho = atm.density(altitude * 0.3048)
        
        # Air density (slugs/ft^3)
        rho = rho * 0.00194032
        
        # Lift coefficient
        CL = W / (0.5 * rho * V**2 * S)
        
        # Drag coefficient
        CD = C_d_0 + CL**2 / (np.pi * e * AR)
        
        # Lift force
        L = CL * 0.5 * rho * V**2 * S
        
        # Parasite drag
        Dp = 0.5 * C_d_0 * rho * V**2 * S
        
        # Induced drag
        Di = (2 * L**2) / (rho * V**2 * S * np.pi * e * AR)
        
        # Thrust required (lbs)
        Tr = W / (CL / CD)
        
        # Add a curve
        plt.plot(V,Tr)
        
        # Plot drag curves
        if drag_curves == "Y":
            plt.plot(V,Dp)
            plt.plot(V,Di)
            legend.append("Parasite Drag")
            legend.append("Induced Drag")
        
        # Update the legend
        if altitude > 0:
            legend.append("{:,.0f} ft".format(altitude))
 
    # Locate the axis labels and limits
    ax = plt.gca()
    ax.set_xlim([75,400])
    ax.set_ylim([0,500])
    ax.xaxis.set_label_coords(0.5, -0.15)
    ax.yaxis.set_label_coords(-0.15, 0.5)
    ax.text(138, 107, 'Best Range', size=13)
    
    # Add title and axis labels
    plt.legend(legend,loc=(0.22, 0.45))
    #plt.title("\n Thrust Required for Level Flight \n")
    plt.ylabel('Thrust\nRequired\n(lb)', rotation=0)
    plt.xlabel('Free Stream Velocity (ft/s)')
    
    # Plot (L/D)max point
    plt.plot(V[np.argmin(Tr)], np.amin(Tr), marker='o', ms=9, color=min_color)
    sns.despine(trim=True)
    plt.show()

#plot_thrust_required(V,0)
plot_thrust_required(V,0,"#4c72b0","Y")
#plot_thrust_required(V,16000,"#64b5cd")

png

Power Required for Level Flight

def plot_power_required(V, Max_Alt=0, min_color="#4c72b0"):
    """
    Plot power required curves for altitudes from sea level to Max_Alt in 4,000-ft increments.
    
    Keyword arguments:
    V -- a 1D array of velocity values in ft/s
    Max_Alt -- an altitude, in feet, at which to stop plotting curves (default 0, sea level)
    min_color -- a hex color code for the min value marker
    """
    legend = ['Sea Level']
    for altitude in range(0, Max_Alt + 1, 4000):
        
        # COESA 1976 U.S. Standard Atmosphere (SI units!)
        # Air density (kg/m^3)
        rho = atm.density(altitude * 0.3048)
        
        # Air density (slugs/ft^3)
        rho = rho * 0.00194032
        
        # Lift coefficient
        CL = W / (0.5 * rho * V**2 * S)
        
        # Drag coefficient
        CD = C_d_0 + CL**2 / (np.pi * e * AR)
        
        # Thrust required (lbs)
        Tr = W / (CL / CD)
        
        # Power required (hp)
        Pr = (Tr * V) / 550
        
        # Add a curve (1 ft/s = 0.592484 knots)
        plt.plot(V*0.592484,Pr)
        
        # Update the legend
        if altitude > 0:
            legend.append("{:,.0f} ft".format(altitude))
 
    # Locate the axis labels and limits
    ax = plt.gca()
    ax.set_xlim([25,250])
    ax.set_ylim([0,300])
    ax.xaxis.set_label_coords(0.5, -0.15)
    ax.yaxis.set_label_coords(-0.15, 0.5)
    #ax.text(75, 15, 'Best Endurance', size=13)
    
    # Add title and axis labels
    plt.legend(legend,loc=(0.12, 0.4))
    #plt.title("\n Power Required for Level Flight \n")
    plt.ylabel('Power\nRequired\n(hp)', rotation=0)
    plt.xlabel('Free Stream Velocity (knots)')
    
    # Plot (L/D)max point
    plt.plot(V[np.argmin(Pr)]*0.592484, np.amin(Pr), marker='o', ms=9, color=min_color)
    sns.despine()
    plt.show()

plot_power_required(V,0)
plot_power_required(V,16000,"#64b5cd")

png

png

Propeller Efficiency

def prop_eta(D=82, eta_guess=0.9, P=300, V=50, altitude=0, es=0.0001):
    """
    Estimate propeller efficiency using a numerical method (ref: Solies [9]).
    
    Keyword arguments:
    D -- propeller diameter, inches
    eta_guess -- an initial guess of propeller efficiency, percent
    P -- engine brake horsepower available at selected altitude
    V -- a velocity or array of velocities in knots
    altitude -- operating altitude, feet
    es -- percent error at which solution converges (smaller = more iterations)
    """
    
    # Propulsive disk area, ft^2
    A = (np.pi * (D/12)**2) / 4
    
    # COESA 1976 U.S. Standard Atmosphere (SI units!)
    # Air density (kg/m^3)
    rho = atm.density(altitude * 0.3048)

    # Air density (slugs/ft^3)
    rho = rho * 0.00194032
    
    # Convert V from knots to ft/s
    V = V * 1.68781
    
    # Convert P from hp to ft*lb/s
    P = P * 550
    
    # Initial error, %
    ea = 100
    
    # Initial guess
    eta_0 = eta_guess
    
    i = 0
    while np.any(ea) > es:
    
        # (2) Propeller thrust, lbs
        T = (eta_0 * P) / V

        # (3) Velocity at the propeller disk is V + dV/2
        dV_by_2 = -(V/2) + np.sqrt((V**2)/4 + T/(2 * rho * A))

        # (4) New efficiency
        eta_i = V / (V + dV_by_2)

        # (1) New eta
        eta = eta_guess * eta_i

        # Error
        ea = np.abs((eta - eta_0)/eta) * 100
        
        #print("\n", i, eta)
        
        if (np.any(ea) > es) or (i>20):
            eta_0 = eta
    
    return eta
# From Lycoming O-360 manual, Fig. 3-26, "Sea Level and Altitude Performance, IO-360-M1B"
io360_bhp ={0:180, 4000:160, 8000:140, 12000:120, 16000:100, 20000:90}

# Hartzell propeller diameter, inches
D = 74

Estimation of Piston Engine BHP with Altitude

Equation from Raymer [2].

\begin{equation} \tag{13.10} P = P_{SL} \bigg(\frac{\rho}{\rho_0} - \frac{1 - \rho/\rho_0}{7.55} \bigg) \end{equation}

rho_0 = atm.density(0 * 0.3048)
rho_0 = rho_0 * 0.00194032
rho_20 = atm.density(20000 * 0.3048)
rho_20 = rho_20 * 0.00194032

P_20 = 180 * ((rho_20/rho_0) - ((1 - (rho_20/rho_0))/7.55))

Power Available and Maximum Velocity

# Velocity range to cross power req'd curve
V_range_0 = np.arange(173,193,5)
V_range_8 = np.arange(170,190,5)

# Propeller efficiency at Sea Level
eta_0 = prop_eta(D, 0.9, io360_bhp[0], V_range_0, 0)

# Propeller efficiency, 16000
eta_8 = prop_eta(D, 0.9, io360_bhp[8000], V_range_8, 8000)
def plot_max_speed(V, altitude, V_range, BHPa, prop_eta):
    """
    Plot 
    
    Keyword arguments:
    V -- 
    """
    legend = ['Sea Level', '8000 ft', 'THP Available']
    
    # COESA 1976 U.S. Standard Atmosphere (SI units!)
    # Air density (kg/m^3)
    rho = atm.density(0 * 0.3048)

    # Air density (slugs/ft^3)
    rho = rho * 0.00194032

    # Lift coefficient
    CL = W / (0.5 * rho * V**2 * S)

    # Drag coefficient
    CD = C_d_0 + CL**2 / (np.pi * e * AR)

    # Thrust required (lbs)
    Tr = W / (CL / CD)

    # Power required (hp)
    Pr = (Tr * V) / 550
    
    # Thrust horsepower available
    THPa = BHPa * prop_eta

    # Add a curve (1 ft/s = 0.592484 knots)
    plt.plot(V*0.592484, Pr)
    
    
    # COESA 1976 U.S. Standard Atmosphere (SI units!)
    # Air density (kg/m^3)
    rho = atm.density(8000 * 0.3048)

    # Air density (slugs/ft^3)
    rho = rho * 0.00194032

    # Lift coefficient
    CL = W / (0.5 * rho * V**2 * S)

    # Drag coefficient
    CD = C_d_0 + CL**2 / (np.pi * e * AR)

    # Thrust required (lbs)
    Tr = W / (CL / CD)

    # Power required (hp)
    Pr = (Tr * V) / 550
    
    # Thrust horsepower available
    THPa_0 = io360_bhp[0] * eta_0
    THPa_8 = io360_bhp[8000] * eta_8

    # Add a curve (1 ft/s = 0.592484 knots)
    plt.plot(V*0.592484, Pr)
    plt.plot(V_range_0, THPa_0, color="#c44e52")
    plt.plot(V_range_8, THPa_8, color="#c44e52")

    # Locate the axis labels and limits
    ax = plt.gca()
    ax.set_xlim([176,182])
    ax.set_ylim([110,170])
    ax.xaxis.set_label_coords(0.5, -0.15)
    ax.yaxis.set_label_coords(-0.15, 0.5)
    
    ax.xaxis.set_major_locator(ticker.MultipleLocator(0.5))
    
    # Add title and axis labels
    plt.legend(legend,loc=(0.15, 0.57))
    #plt.title("\n Power Required for Level Flight \n")
    plt.ylabel('Power\nRequired\n(hp)', rotation=0)
    plt.xlabel('Free Stream Velocity (knots)')
    sns.despine(trim=True)
    plt.show()

sns.set_style("whitegrid")
plot_max_speed(V, 8000, V_range_8, io360_bhp[8000], eta_8)

png

The intersection of the Thrust Horsepower Available (THP) curves and the Power Required curves should give a rough indication of maximum velocity at a given altitude. Here, we show a max speed of 180.5 knots at sea level, and 177 knots at 8,000 feet. These numbers compare favorably to Van’s published data of 213 mph (185 knots) at sea level and 203 mph (176 knots) at 8,000 feet.

References

[ 1 ] R. D. Kimberlin, Flight Testing of Fixed-Wing Aircraft. Reston, Va.: AIAA, 2003.

[ 2 ] D. P. Raymer, Aircraft Design: A Conceptual Approach. Reston, Va: AIAA, 4th ed., 2006.

[ 3 ] J. D. Anderson, Introduction to Flight. New York: McGraw-Hill, 5th ed., 2005.

[ 4 ] J. D. Mattingly and K. M. Boyer, Elements of Propulsion: Gas Turbines and Rockets. Reston, Va.: AIAA, 2nd ed., 2016.

[ 5 ] M. Nita and D. Scholz, “Estimating the Oswald Factor from Basic Aircraft Geometrical Parameters,” in German Aerospace Congress, pp. 1–19, Dec. 2012.

[ 6 ] C. E. Jobe, “Prediction of Aerodynamic Drag,” Tech. Rep. AFWAL-TM-84-203, Air Force Wright Aeronautical Laboratories, July 1984.

[ 7 ] C. M. Jackson Jr., “Estimation of Flight Performance with Closed-Form Approximations to the Equations of Motion,” pp. 1–33, Mar. 2001.

[ 8 ] I. H. Abbot, A. E. Von Doenhoff, and L. S. Stivers Jr., Summary of Airfoil Data. NACA: National Advisory Committee for Aeronautics, 1945.

[ 9 ] U. P. Solies, “Numerical method for estimation of propeller efficiencies,” Journal of Aircraft, vol. 31, pp. 996–998, July 1994.

comments powered by Disqus