Skip to content

Commit

Permalink
Merge branch 'develop' into feature-climb_segment_bug
Browse files Browse the repository at this point in the history
  • Loading branch information
mclarke2 authored Dec 20, 2022
2 parents 76b23c7 + 24f153e commit ea2fbe1
Show file tree
Hide file tree
Showing 3 changed files with 186 additions and 8 deletions.
23 changes: 17 additions & 6 deletions regression/scripts/airfoil_import/airfoil_import_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
# Created:
# Modified: Sep 2020, M. Clarke
# May 2021, R. Erhard
# Dec 2022, J. Smart

#----------------------------------------------------------------------
# Imports
Expand All @@ -12,10 +13,8 @@
from SUAVE.Core import Units, Data
from SUAVE.Plots.Geometry import plot_airfoil
import matplotlib.pyplot as plt
from SUAVE.Methods.Geometry.Two_Dimensional.Cross_Section.Airfoil.import_airfoil_geometry\
import import_airfoil_geometry
from SUAVE.Methods.Geometry.Two_Dimensional.Cross_Section.Airfoil.compute_airfoil_properties \
import compute_airfoil_properties
from SUAVE.Methods.Geometry.Two_Dimensional.Cross_Section.Airfoil\
import import_airfoil_geometry, compute_airfoil_properties, convert_airfoil_to_meshgrid
from SUAVE.Plots.Performance.Airfoil_Plots import *
import os
import numpy as np
Expand Down Expand Up @@ -68,8 +67,20 @@ def main():
for j in range(0, len(airfoil_geometry_3.camber_coordinates)):
assert( np.abs(airfoil_geometry_3.camber_coordinates[j] - airfoil_geometry_4.camber_coordinates[j]) < 1E-8 )

plot_airfoil(airfoil_geometry_with_selig[0])

# Multiple meshes use too much memory on AppVeyor

A_MASK_1 = convert_airfoil_to_meshgrid(airfoil_geometry_1)
# A_MASK_2 = convert_airfoil_to_meshgrid(airfoil_geometry_2)
# A_MASK_3 = convert_airfoil_to_meshgrid(airfoil_geometry_3)
# A_MASK_4 = convert_airfoil_to_meshgrid(airfoil_geometry_4)

assert (len(np.where(A_MASK_1)[0]) == 32313)
# assert (len(np.where(A_MASK_2)[0]) == 32313)
# assert (len(np.where(A_MASK_3)[0]) == 122051849)
# assert (len(np.where(A_MASK_4)[0]) == 122051849)

plot_airfoil(airfoil_geometry_with_selig[1])

return

if __name__ == '__main__':
Expand Down
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
## @defgroup Methods-Geometry-Two_Dimensional-Cross_Section-Airfoil Airfoil
# Geometry functions for two dimensional airfoils.
# Geometry functions for two-dimensional airfoils.
# @ingroup Methods-Geometry-Two_Dimensional-Cross_Section

from .compute_naca_4series import compute_naca_4series
from .compute_airfoil_properties import compute_airfoil_properties
from .import_airfoil_dat import import_airfoil_dat
from .import_airfoil_geometry import import_airfoil_geometry
from .import_airfoil_polars import import_airfoil_polars
from .import_airfoil_polars import import_airfoil_polars
from .convert_airfoil_to_meshgrid import convert_airfoil_to_meshgrid
Original file line number Diff line number Diff line change
@@ -0,0 +1,166 @@
## @ingroup Methods-Geometry-Two_Dimensional-Cross_Section-Airfoil
# convert_airfoil_to_meshgrid.py
#
# Created: Dec 2022, J. Smart
# Modified:

# ----------------------------------------------------------------------
# Imports
# ----------------------------------------------------------------------

from SUAVE.Core import Data
import numpy as np

## @ingroup Methods-Geometry-Two_Dimensional-Cross_Section-Airfoil
def convert_airfoil_to_meshgrid(airfoil_geometry, *args, **kwargs):
"""Converts a SUAVE airfoil geometry representation to a Numpy meshgrid
array mask of boolean values.
Assumptions:
None
Source:
None
Inputs:
airfoil_geometry [SUAVE Data Structure]
.x_lower_surface [Numpy Array, float32]
.y_lower_surface [Numpy Array, float32]
.y_upper_surface [Numpy Array, float32]
Outputs:
airfoil_meshgrid [Numpy Array, bool]
Properties Used:
N/A
"""

# Unpack Values

x_lower_surface = airfoil_geometry.x_lower_surface
y_lower_surface = airfoil_geometry.y_lower_surface
y_upper_surface = airfoil_geometry.y_upper_surface

# Determine necessary resolution of the meshgrid. We do this by dividing the
# x-length of the lower surface by the minimum separation between
# any two x-coordinates of the geometry (ceil-rounded to an int). Later
# we'll instantiate the meshgrid with this number of x-indices, so that the
# separation between any two points in the meshgrid is equal to the minimum
# separation between any two x-coordinates.

x_length = (
np.max(x_lower_surface)
)

Nx = np.ceil(
x_length / np.abs(np.min(np.diff(x_lower_surface)))
).astype(int)

# We determine the necessary number of y-coordinate points by taking the
# maximum separation between the highest point of the upper surface and the
# lowest point of the lower surface and multiplying that by the number of
# x-points in order to re-normalize to our future meshgrid coordinates,
# then ciel-rounding to an int.

Ny = np.ceil(
Nx * ( np.max(y_upper_surface) - np.min(y_lower_surface) )
).astype(int)

# Instantiate the meshgrid, using ij-indexing so that X[i,j] returns i
# for all points, and Y[i,j] returns j for all coordinates.

X, Y = np.meshgrid(np.arange(Nx), np.arange(Ny), indexing="ij")

# Create the indexing arrays for the meshgrid. These convert the airfoil
# geometry coordinates into meshgrid array indices. The X_INDICES are found
# just by multplying/stretching the x_lower_surface coordinates across the
# number of x-coodinates in the meshgrid.

X_INDICES = np.ceil(
Nx / x_length * x_lower_surface
).astype(int)

# The Y_INDICES are similarly stretched, but first are offset by the
# minimum of the lower surface to bring them to a relative zero

Y_LOWER_INDICES = np.floor(
Nx / x_length * (
y_lower_surface - np.min(y_lower_surface)
)
).astype(int)

Y_UPPER_INDICES = np.ceil(
Nx /x_length * (
y_upper_surface - np.min(y_lower_surface)
)
).astype(int)

# We then repeat the elements of the Y_INDICES by the number of gridpoints
# between each x-coordinate, essentially treating the y-surface as flat
# between those points. We trim the final point by telling it to repeat 0
# times

REPEATS = np.append(
np.diff(X_INDICES),
0
)

# Need to hand the case where the X_INDICES aren't sorted, and swap
# some elements around to allow the masks to be created

if np.any(REPEATS<0):

REPEAT_FLAG = True

NEG_REPEATS = np.where(REPEATS<0)[0]

if np.any(np.diff(NEG_REPEATS) == 1):
print("Airfoil geometry contains sequential negative x-steps. Meshing Failed.")
return None

(X_INDICES[NEG_REPEATS],
X_INDICES[NEG_REPEATS + 1]) = (X_INDICES[NEG_REPEATS + 1],
X_INDICES[NEG_REPEATS])

(Y_LOWER_INDICES[NEG_REPEATS],
Y_LOWER_INDICES[NEG_REPEATS + 1]) = (Y_LOWER_INDICES[NEG_REPEATS + 1],
Y_LOWER_INDICES[NEG_REPEATS])

(Y_UPPER_INDICES[NEG_REPEATS],
Y_UPPER_INDICES[NEG_REPEATS + 1]) = (Y_UPPER_INDICES[NEG_REPEATS + 1],
Y_UPPER_INDICES[NEG_REPEATS])

REPEATS = np.append(
np.diff(X_INDICES),
0
)

Nx = np.sum(REPEATS)

Ny = np.ceil(
Nx * (np.max(y_upper_surface) - np.min(y_lower_surface))
).astype(int)

X, Y = np.meshgrid(np.arange(Nx), np.arange(Ny), indexing="ij")

Y_LOWER_INDICES = np.repeat(Y_LOWER_INDICES, REPEATS)
Y_UPPER_INDICES = np.repeat(Y_UPPER_INDICES, REPEATS)

# We then create masks for the upper and lower surfaces by tiling the
# indices over the meshgrid (taking a transpose to comport with our earlier
# indexing style).

Y_LOWER_GRID = np.tile(Y_LOWER_INDICES, (Ny,1)).T
Y_UPPER_GRID = np.tile(Y_UPPER_INDICES, (Ny,1)).T

# We then create our airfoil meshgrid mask by comparing our Y coordinates
# from the meshgrid to our upper and lower grids, intermediately treating
# them as ints to simplify the multi-condition comparison

Y_LOWER = (Y > Y_LOWER_GRID).astype(int)
Y_UPPER = (Y < Y_UPPER_GRID).astype(int)

AIRFOIL_MASK = (Y_LOWER + Y_UPPER) > 1

return AIRFOIL_MASK

0 comments on commit ea2fbe1

Please sign in to comment.