06 Mar

Introducing GEVfoil

I have decided to rename my little aerofoil project. I found no fewer than four other projects using the previous name, which I’ll avoid stating to reduce confusing search engines any further.

I’ve settled on GEVfoil (Ground Effect Vehicle Aerofoil) which may seem a little obvious, but clarity should at least be in my favour. The objectives remain the same;

  • To create a set of tools to aid in the analysis of aerofoils in ground effect.
  • To primarily focus on 2D effects, but not exclude 3D effects where appropriate and relevant.
  • To utilise open-source tools as much as possible, for example links into the OpenFOAM CFD (Computational Fluid Dynamics) project.
  • To maintain a blog of the development, and open-source the whole project when it has content enough to be worthy.

I have renamed previous blog posts with this in mind:

18 Feb

Ramer–Douglas–Peucker algorithm

Generally speaking, for low-speed flow around an aerofoil, the largest gradients in the flow-properties coincide with the greatest curvature of the aerofoil, and thus are the region that requires the greatest concentration of cells (or nodes) in the computational domain.

Many of the simplified treatments for spacing points along an aerofoil use a logarithmic spacing biased towards the leading edge. This of course can create issues if spaced in the chord-wise direction (x/c). The surface-normal at the leading edge of a symmetric aerofoil is perpendicular to the stream-wise direction, and thus can become under-resolved, and even faceted, with a simple linear or logarithmic chord-wise spacing.

One of the more significant issues with linear or logarithmic spacing is that it simply does not discriminate on the basis of curvature, thus the growth factor becomes a qualitative parameter controlled by the aerodynamicist. Any other curvature, for example the upper-surface curvature caused by the s-shaped chord of a DHTMU aerofoil, is neglected.

S-shaped chord of the DHMTU aerofoils. Source: JavaFoil Manual

Some computational treatments, though able to nicely capture curvature, are a little too complex for this simplified tool intended for 2D aerofoil analyses. I believe that the Ramer-Douglas-Peucker (RDP) algorithm is a simple yet powerful tool ideally suited to be applied to spacing points/nodes along a computational grid for for 2D aerofoil analyses.

RDP is a a decimation (down-sampling) algorithm, which although cannot increase resolution (create points) in areas of high-curvature, it can remove remove points in areas where curvature is low. This approach is quite simple and elegant. The algorithm works by creating a line between the start and end points, then recursively divides this line whilst each time finding the furthest point, marking these points to be kept if they are greater than a user-specified distance (epsilon).

Ramer-Douglas-Peucker Animation. Source: Mysid

As this requires a set of points, rather than a spline or other analytical curve, there are upsides and downsides. In the case of of analytic aerofoils, such as the DHTMU, or the NACA 4 that I originally scripted into GEVgoil, these need to be discretised before being down-sampled. However, for aerofoils that are already being supplied as a set of points, for example using the CSV method that is now also in GEVfoil, RDP can be applied directly.

To simplify matters further, RDP is available as a Python module via PIP: https://pypi.org/project/rdp/ Using this excellent ready-made library, I have constructed perhaps the simplest Python class, implementing this RDP module into GEVfoil:

from rdp import rdp
import numpy as np

class RDP(object):
    """A class to decimate aerofoil points"""
    def __init__(self, x, y, e=0.0001):
        xy = np.vstack((x, y)).T
        self.reducedXY = rdp(xy, epsilon=e)
        self.reducedX = self.reducedXY[:,0]
        self.reducedY = self.reducedXY[:,1]

And of course, the proof is in the pudding. Some examples of the RDP module in action in GEVfoil with varying epsilon:

Original NACA 2412 as generated by GEVfoil with 800 points each on the upper and lower surfaces
Decimated with RDP using an epsilon of 1e-05
Decimated with RDP using an epsilon of 1e-04
Decimated with RDP using an epsilon of 1e-03
21 Oct

GEVfoil – NACA 4 Series Part 2 and CSV Reading

Previously I had made a simple function that will generate a number of points/coordinates for any specified NACA 4 series aerofoil. The limitation with this approach was that it did a fairly poor job at capturing the leading-edge. The leading edge is of critical importance to the aerodynamic performance and characteristics of the aerofoil. After some consideration, I have decided not to modify this technique. I am designing the GEVfoil tools to be as modular, generic and extendable as possible. This isn’t to say I am neglecting the importance of the leading edge, or that the method is flawed.

Initially I had considered returning the aerofoil as a mathematical function, rather than a set of points. This concept may well be revisited in the future. However, most of the other input aerofoils will arrive as a set of coordinates of varying fidelity, from sources such as a the UIUC Airfoil Coordinates Database or from Comma Separated Variable (CSV) files. Therefore, the classes which will receive the aerofoil coordinate objects will deal with interpolating and appropriately converting coordinates into functions anyway. Consequently, as long as the linear-spacing is arbitrarily high enough a simple piece-wise cubic interpolation will suffice.

Note that the CSV reader itself is an abstraction of the AerofoilReader class, which I plan to later extend to read the DAT file from UIUC,

To that end, the CSV reader currently looks as follows:

import os

class AerofoilReader(object):
    """ An abstract class to be inherited by classes for specific
        aerofoil file types """
    def __init__(self, filename):
        super(AerofoilReader, self).__init__()
        self.coord = []
        self.title = []  # Column titles
        self.label = []  # Aerofoil Label

        foil = open(filename, 'r')
        self.lines = foil.readlines()


class CSVReader(AerofoilReader):
    """ Read in CSV. Assumes a header is present."""
    def __init__(self, filepath, xcol, ycol, zcol=False):
        super(CSVReader, self).__init__(filepath)

        self.label = os.path.basename(filepath).split('.')[0]

        inpcols = 3        # Two or three columns
        if not zcol:
            inpcols = 2

        for i, line in enumerate(self.lines):
            parts = line.rstrip().replace('"', '').split(',')
            additional = []

            if i > 0:
                # Not the header
                if len(parts) > inpcols:
                    # More than x, y & z - store other columns
                    for j in range(len(parts)):
                        if j != xcol and j != ycol and j != zcol:
                            additional.append(parts[j])

                xyz = [parts[xcol], parts[ycol], parts[zcol]]
                self.coord.append(xyz + additional)
            else:
                # It's the header
                self.title = parts
13 Feb

GEVfoil – NACA 4 Series

The intention behind GEVfoil is to establish an open-source set of Python classes to handle aerofoils (American English: airfoils). The intention is that these can be used for anything from generating simple point clouds, analytical representations all the way to STL files for CFD meshing and OpenFOAM BlockMesh output.

The first step is starting with something basic, in this case a simple class to generate NACA 4 series aerofoils:

import numpy as np
import matplotlib.pyplot as plt


class NACA4(object):
    """A class for generating  NACA 4 series aerofoils
        TO-DO - linspace not good at LE capture """
    def __init__(self, foil, points=200, chord=1):
        super(NACA4, self).__init__()
        m = float(foil[0])/100      # max camber
        p = float(foil[1])/10       # chordwise position of max camber
        t = float(foil[2:])/100     # thickness
        x = np.linspace(0, chord, points)
        self.x = x
        yt = self.thickness(x, t)
        yc = []
        for coord in x:
            if m:
                yc.append(self.camber(m, p, coord))
            else:
                # No camber
                yc.append(0)
        y1 = yc + yt
        y2 = yc - yt

        self.y = y1

        plt.plot(self.x, yc)
        plt.plot(self.x, y1, '.')
        plt.plot(self.x, y2, '.')
        plt.axis('equal')
        plt.show()

    def thickness(self, x, t):
        # Y thickness at given x point
        return t / 0.2 * \
            (0.2969*np.sqrt(x)-0.126*x-0.3516*x**2+0.2843*x**3-0.1015*x**4)

    def camber(self, m, p, x):
        # Return the camber of the aerofoil
        if x <= p:
            return m/p**2 * (2*p*x-x**2)
        return m/(1-p)**2*((1-2*p)+2*p*x-x**2)

if __name__ == "__main__":
    test = NACA4('2412')

For now, Matplotlib has been built into the class for visual reference. For example, 2412 generates:

And 0010 creates:

At this stage there is a clear issue around the leading edge, the function was interrogated with a linearly spaced x/chord-wise interval, this must be improved for the future to cluster points around the greatest curvature.