1 # This file is part of Marius Peter's airfoil analysis package (this program). 2 # 3 # This program is free software: you can redistribute it and/or modify 4 # it under the terms of the GNU General Public License as published by 5 # the Free Software Foundation, either version 3 of the License, or 6 # (at your option) any later version. 7 # 8 # This program is distributed in the hope that it will be useful, 9 # but WITHOUT ANY WARRANTY; without even the implied warranty of 10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11 # GNU General Public License for more details. 12 # 13 # You should have received a copy of the GNU General Public License 14 # along with this program. If not, see <https://www.gnu.org/licenses/>. 15 """ 16 The creator.py module contains class definitions for coordinates 17 and various components we add to an airfoil (spars, stringers, and ribs). 18 19 Classes: 20 Airfoil: instantiated with class method to provide coordinates to heirs. 21 Spar: inherits from Airfoil. 22 Stringer: also inherits from Airfoil. 23 24 Functions: 25 plot_geom(airfoil): generates a 2D plot of the airfoil & any components. 26 """ 27 28 import sys 29 import os.path 30 import numpy as np 31 from math import sin, cos, atan 32 import bisect as bi 33 import matplotlib.pyplot as plt 34 35 36 class Airfoil: 37 """This class represents a single NACA airfoil. 38 39 The coordinates are saved as two lists 40 for the x- and z-coordinates. The coordinates start at 41 the leading edge, travel over the airfoil's upper edge, 42 then loop back to the leading edge via the lower edge. 43 44 This method was chosen for easier future exports 45 to 3D CAD packages like SolidWorks, which can import such 46 geometry as coordinates written in a CSV file. 47 """ 48 49 # Defaults 50 chord = 100 51 semi_span = 200 52 53 def __init__(self): 54 # mass and area 55 self.mass = float() 56 self.area = float() 57 # Component material 58 self.material = str() 59 # Coordinate 60 self.x = [] 61 self.z = [] 62 63 @classmethod 64 def from_dimensions(cls, chord, semi_span): 65 """Create airfoil from its chord and semi-span.""" 66 if chord > 20: 67 cls.chord = chord 68 else: 69 cls.chord = 20 70 print('Chord too small, using minimum value of 20.') 71 cls.semi_span = semi_span 72 return Airfoil() 73 74 def __str__(self): 75 return type(self).__name__ 76 77 def add_naca(self, naca_num): 78 """Generate surface geometry for a NACA airfoil. 79 80 The nested functions perform the required steps to generate geometry, 81 and can be called to solve the geometry y-coordinate for any 'x' input. 82 Equation coefficients were retrieved from Wikipedia.org. 83 84 Parameters: 85 naca_num: 4-digit NACA wing 86 87 Return: 88 None 89 """ 90 # Variables extracted from 'naca_num' argument passed to the function 91 self.naca_num = naca_num 92 m = int(str(naca_num)[0]) / 100 93 p = int(str(naca_num)[1]) / 10 94 t = int(str(naca_num)[2:]) / 100 95 # x-coordinate of maximum camber 96 p_c = p * self.chord 97 98 def get_camber(x): 99 """ 100 Returns camber z-coordinate from 1 'x' along the airfoil chord. 101 """ 102 z_c = float() 103 if 0 <= x < p_c: 104 z_c = (m / (p**2)) * (2 * p * (x / self.chord) - 105 (x / self.chord)**2) 106 elif p_c <= x <= self.chord: 107 z_c = (m / 108 ((1 - p)**2)) * ((1 - 2 * p) + 2 * p * 109 (x / self.chord) - (x / self.chord)**2) 110 return (z_c * self.chord) 111 112 def get_thickness(x): 113 """Return thickness from 1 'x' along the airfoil chord.""" 114 x = 0 if x < 0 else x 115 z_t = 5 * t * self.chord * (+0.2969 * 116 (x / self.chord)**0.5 - 0.1260 * 117 (x / self.chord)**1 - 0.3516 * 118 (x / self.chord)**2 + 0.2843 * 119 (x / self.chord)**3 - 0.1015 * 120 (x / self.chord)**4) 121 return z_t 122 123 def get_theta(x): 124 dz_c = float() 125 if 0 <= x < p_c: 126 dz_c = ((2 * m) / p**2) * (p - x / self.chord) 127 elif p_c <= x <= self.chord: 128 dz_c = (2 * m) / ((1 - p)**2) * (p - x / self.chord) 129 130 theta = atan(dz_c) 131 return theta 132 133 def get_upper_coord(x): 134 x = x - get_thickness(x) * sin(get_theta(x)) 135 z = get_camber(x) + get_thickness(x) * cos(get_theta(x)) 136 return (x, z) 137 138 def get_lower_coord(x): 139 x = x + get_thickness(x) * sin(get_theta(x)) 140 z = get_camber(x) - get_thickness(x) * cos(get_theta(x)) 141 return (x, z) 142 143 # Densify x-coordinates 10 times for first 1/4 chord length 144 x_chord_25_percent = round(self.chord / 4) 145 146 x_chord = [i / 10 for i in range(x_chord_25_percent * 10)] 147 x_chord.extend(i for i in range(x_chord_25_percent, self.chord + 1)) 148 # Reversed list for our lower airfoil coordinate densification 149 x_chord_rev = [i for i in range(self.chord, x_chord_25_percent, -1)] 150 extend = [i / 10 for i in range(x_chord_25_percent * 10, -1, -1)] 151 x_chord_rev.extend(extend) 152 153 # Generate our airfoil geometry from previous sub-functions. 154 self.x_c = [] 155 self.z_c = [] 156 for x in x_chord: 157 self.x_c.append(x) 158 self.z_c.append(get_camber(x)) 159 self.x.append(get_upper_coord(x)[0]) 160 self.z.append(get_upper_coord(x)[1]) 161 for x in x_chord_rev: 162 self.x.append(get_lower_coord(x)[0]) 163 self.z.append(get_lower_coord(x)[1]) 164 return None 165 166 def add_mass(self, mass): 167 self.mass = mass 168 169 def info_print(self, round): 170 """Print all the component's coordinates to the terminal.""" 171 name = ' CREATOR DATA FOR {} '.format(str(self).upper()) 172 num_of_dashes = len(name) 173 print(num_of_dashes * '-') 174 print(name) 175 for k, v in self.__dict__.items(): 176 if type(v) != list: 177 print('{}:\n'.format(k), v) 178 print(num_of_dashes * '-') 179 for k, v in self.__dict__.items(): 180 if type(v) == list: 181 print('{}:\n'.format(k), np.around(v, round)) 182 return None 183 184 def info_save(self, save_path, number): 185 """Save all the object's coordinates (must be full path).""" 186 file_name = '{}_{}.txt'.format(str(self).lower(), number) 187 full_path = os.path.join(save_path, file_name) 188 try: 189 with open(full_path, 'w') as sys.stdout: 190 self.info_print(6) 191 # This line required to reset behavior of sys.stdout 192 sys.stdout = sys.__stdout__ 193 print('Successfully wrote to file {}'.format(full_path)) 194 except IOError: 195 print( 196 'Unable to write {} to specified directory.\n'.format( 197 file_name), 'Was the full path passed to the function?') 198 return None 199 200 201 class Spar(Airfoil): 202 """Contains a single spar's location.""" 203 def __init__(self): 204 super().__init__() 205 self.x_start = [] 206 self.x_end = [] 207 self.thickness = float() 208 self.z_start = [] 209 self.z_end = [] 210 211 def add_coord(self, airfoil, x_loc_percent): 212 """Add a single spar at the % chord location given to function. 213 214 Parameters: 215 airfoil: gives the spar access to airfoil's coordinates. 216 x_loc_percent: spar's location as a % of total chord length. 217 218 Return: 219 None 220 """ 221 222 # Scaled spar location with regards to chord 223 loc = x_loc_percent * self.chord 224 # bi.bisect_left: returns index of first value in airfoil.x > loc 225 # This ensures that spar geom intersects with airfoil geom. 226 # Spar upper coordinates 227 spar_x = bi.bisect_left(airfoil.x, loc) - 1 228 x = [airfoil.x[spar_x]] 229 z = [airfoil.z[spar_x]] 230 # Spar lower coordinates 231 spar_x = bi.bisect_left(airfoil.x[::-1], loc) 232 x += [airfoil.x[-spar_x]] 233 z += [airfoil.z[-spar_x]] 234 self.x.append(x) 235 self.z.append(z) 236 return None 237 238 def add_spar_caps(self, spar_cap_area): 239 self.cap_area = spar_cap_area 240 return None 241 242 def add_mass(self, mass): 243 self.mass = len(self.x) * mass 244 return None 245 246 def add_webs(self, thickness): 247 """Add webs to spars.""" 248 for _ in range(len(self.x)): 249 self.x_start.append(self.x[_][0]) 250 self.x_end.append(self.x[_][1]) 251 self.z_start.append(self.z[_][0]) 252 self.z_end.append(self.z[_][1]) 253 self.thickness = thickness 254 return None 255 256 257 class Stringer(Airfoil): 258 """Contains the coordinates of all stringers.""" 259 def __init__(self): 260 super().__init__() 261 self.x_start = [] 262 self.x_end = [] 263 self.thickness = float() 264 self.z_start = [] 265 self.z_end = [] 266 self.area = float() 267 268 def add_coord(self, airfoil, stringer_u_1, stringer_u_2, stringer_l_1, 269 stringer_l_2): 270 """Add equally distributed stringers to four airfoil locations 271 (upper nose, lower nose, upper surface, lower surface). 272 273 Parameters: 274 airfoil_coord: packed airfoil coordinates 275 spar_coord: packed spar coordinates 276 stringer_u_1: upper nose number of stringers 277 stringer_u_2: upper surface number of stringers 278 stringer_l_1: lower nose number of stringers 279 stringer_l_2: lower surface number of stringers 280 281 Returns: 282 None 283 """ 284 285 # Find distance between leading edge and first upper stringer 286 interval = airfoil.spar.x[0][0] / (stringer_u_1 + 1) 287 # initialise first self.stringer_x at first interval 288 x = interval 289 # Add upper stringers from leading edge until first spar. 290 for _ in range(0, stringer_u_1): 291 # Index of the first value of airfoil.x > x 292 i = bi.bisect_left(airfoil.x, x) 293 self.x.append(airfoil.x[i]) 294 self.z.append(airfoil.z[i]) 295 x += interval 296 # Add upper stringers from first spar until last spar 297 # TODO: stringer placement if only one spar is created 298 interval = (airfoil.spar.x[-1][0] - 299 airfoil.spar.x[0][0]) / (stringer_u_2 + 1) 300 x = interval + airfoil.spar.x[0][0] 301 for _ in range(0, stringer_u_2): 302 i = bi.bisect_left(airfoil.x, x) 303 self.x.append(airfoil.x[i]) 304 self.z.append(airfoil.z[i]) 305 x += interval 306 307 # Find distance between leading edge and first lower stringer 308 interval = airfoil.spar.x[0][1] / (stringer_l_1 + 1) 309 x = interval 310 # Add lower stringers from leading edge until first spar. 311 for _ in range(0, stringer_l_1): 312 i = bi.bisect_left(airfoil.x[::-1], x) 313 self.x.append(airfoil.x[-i]) 314 self.z.append(airfoil.z[-i]) 315 x += interval 316 # Add lower stringers from first spar until last spar 317 interval = (airfoil.spar.x[-1][1] - 318 airfoil.spar.x[0][1]) / (stringer_l_2 + 1) 319 x = interval + airfoil.spar.x[0][1] 320 for _ in range(0, stringer_l_2): 321 i = bi.bisect_left(airfoil.x[::-1], x) 322 self.x.append(airfoil.x[-i]) 323 self.z.append(airfoil.z[-i]) 324 x += interval 325 return None 326 327 def add_area(self, area): 328 self.area = area 329 return None 330 331 def add_mass(self, mass): 332 self.mass = len(self.x) * mass + len(self.x) * mass 333 return None 334 335 def add_webs(self, thickness): 336 """Add webs to stringers.""" 337 for _ in range(len(self.x) // 2): 338 self.x_start.append(self.x[_]) 339 self.x_end.append(self.x[_ + 1]) 340 self.z_start.append(self.z[_]) 341 self.z_end.append(self.z[_ + 1]) 342 self.thickness = thickness 343 return None 344 345 def info_print(self, round): 346 super().info_print(round) 347 print('Stringer Area:\n', np.around(self.area, round)) 348 return None 349 350 351 def plot_geom(airfoil, view: False): 352 """This function plots the airfoil's + sub-components' geometry.""" 353 fig, ax = plt.subplots() 354 355 # Plot chord 356 x = [0, airfoil.chord] 357 y = [0, 0] 358 ax.plot(x, y, linewidth='1') 359 # Plot quarter chord 360 ax.plot(airfoil.chord / 4, 361 0, 362 '.', 363 color='g', 364 markersize=24, 365 label='Quarter-chord') 366 # Plot mean camber line 367 ax.plot(airfoil.x_c, 368 airfoil.z_c, 369 '-.', 370 color='r', 371 linewidth='2', 372 label='Mean camber line') 373 # Plot airfoil surfaces 374 ax.plot(airfoil.x, airfoil.z, color='b', linewidth='1') 375 376 # Plot spars 377 try: 378 for _ in range(len(airfoil.spar.x)): 379 x = (airfoil.spar.x[_]) 380 y = (airfoil.spar.z[_]) 381 ax.plot(x, y, '-', color='y', linewidth='4') 382 except AttributeError: 383 print('No spars to plot.') 384 # Plot stringers 385 try: 386 for _ in range(0, len(airfoil.stringer.x)): 387 x = airfoil.stringer.x[_] 388 y = airfoil.stringer.z[_] 389 ax.plot(x, y, '.', color='y', markersize=12) 390 except AttributeError: 391 print('No stringers to plot.') 392 393 # Graph formatting 394 plot_bound = max(airfoil.x) 395 ax.set(title='NACA ' + str(airfoil.naca_num) + ' airfoil', 396 xlabel='X axis', 397 xlim=[-0.10 * plot_bound, 1.10 * plot_bound], 398 ylabel='Z axis', 399 ylim=[-(1.10 * plot_bound / 2), (1.10 * plot_bound / 2)]) 400 401 plt.grid(axis='both', linestyle=':', linewidth=1) 402 plt.gca().set_aspect('equal', adjustable='box') 403 plt.gca().legend(bbox_to_anchor=(1, 1), 404 bbox_transform=plt.gcf().transFigure) 405 406 if view is True: 407 plt.show() 408 else: 409 pass 410 return fig, ax 411 412 413 def main(): 414 return None 415 416 417 if __name__ == '__main__': 418 main() 419