utils.py
- openaerostruct.geometry.utils.build_multi_spline(out_name, num_sections, control_points)[source]
This function returns an OpenMDAO Independent Variable Component with an output vector appropriately named and sized to function as an unified set of B-spline control poitns that join multiple sections by construction.
- Parameters:
- out_name: string
Name of the output to assign to the B-spline
- num_sectionsint
Number of sections
- control_points: list
List of B-spline control point arrays corresponding to each section
- Returns:
- spline_controlOpenMDAO component object
The unified B-spline control point indpendent variable component
- openaerostruct.geometry.utils.build_section_dicts(surface)[source]
This utility function takes a multi-section surface dictionary and outputs a list of individual section surface dictionaries so the geometry group for each individual section can be initialized.
- Parameters:
- surface: dict
OpenAeroStruct multi-section surface dictionary
- Returns:
- section_surfaceslist
List of OpenAeroStruct surface dictionaries for each individual surface
- openaerostruct.geometry.utils.connect_multi_spline(prob, section_surfaces, sec_cp, out_name, comp_name, geom_name, return_bind_inds=False)[source]
This function connects the the unified B-spline component with the individual B-splines of each section. There is a point of overlap at each section so that each edge control point control the edge controls points of each section’s B-spline. This is how section joining by consturction is acheived. An issue occurs however when a B-spline in a particular section only has one control point. In this case the one section control point is bound to the left edge B-spline component control point. As result, there is nothing to maintain C0 continuity with the next section. As result a constraint will need to be manually set. To facilitate this, the array bind_inds will contain a list of the B-spline control point indicies that will need to be manually constrained to their previous sections.
- Parameters:
- probOpenMDAO problem object
The OpenAeroStruct problem object with the unified B-spline component added.
- section_surfaceslist
List of the surface dictionaries for each section.
- sec_cplist
List of B-spline control point arrays for each section.
- out_name: string
Name of the unified B-spline component output to connect from
- comp_name: string
Name of the unified B-spline component added to the problem object
- geom_namestring
Name of the multi-section geometry group
- return_bind_inds: bool
Return list of unjoined unified B-spline inidices. Default is False.
- Returns:
- bind_indslist
List of unified B-spline control point indicies not connected due to the presence of a single control point section.(Only if return bind_inds specified)
- openaerostruct.geometry.utils.dihedral(mesh, dihedral_angle, symmetry)[source]
Apply dihedral angle. Positive angles up.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- dihedral_anglefloat
Dihedral angle in degrees.
- symmetryboolean
Flag set to true if surface is reflected about y=0 plane.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the aerodynamic surface with dihedral angle.
- openaerostruct.geometry.utils.gen_crm_mesh(num_x, num_y, span_cos_spacing=0.0, chord_cos_spacing=0.0, wing_type='CRM:jig')[source]
- openaerostruct.geometry.utils.gen_rect_mesh(num_x, num_y, span, chord, span_cos_spacing=0.0, chord_cos_spacing=0.0)[source]
- openaerostruct.geometry.utils.generate_vsp_surfaces(vsp_file, symmetry=False, include=None, scale=1.0)[source]
Generate a series of VLM surfaces based on geometries in an OpenVSP model.
- Parameters:
- vsp_filestr
OpenVSP file to generate meshes from.
- symmetrybool
Flag specifying if the full model should be read in (False) or only half (True). Half model only reads in right side surfaces. Defaults to full model.
- includelist[str]
List of body names defined in OpenVSP model that should be included in VLM mesh output. Defaults to all bodies found in model.
- scale: float
A global scale factor from the OpenVSP geometry to incoming VLM mesh geometry. For example, if the OpenVSP model is in inches, and the VLM in meters, scale=0.0254. Defaults to 1.0.
- Returns:
- surfaceslist[dict]
List of surfaces dictionaries, one (two if symmetry==False) for each body requested in include. This is a relatively empty surface dictionary that contains only basic information about the VLM mesh (i.e. name, symmetry, mesh).
- openaerostruct.geometry.utils.plot3D_meshes(file_name, zero_tol=0)[source]
Reads in multi-surface meshes from a Plot3D mesh file for VLM analysis.
- Parameters:
- fileNamestr
Plot3D file name to be read in.
- zero_tolfloat
If a node location read in the file is below this magnitude we will just make it zero. This is useful for getting rid of noise in the surface that may be due to the meshing tools geometry tolerance.
- Returns:
- mesh_dictdict
Dictionary holding the mesh of every surface included in the plot3D sorted by surface name.
- openaerostruct.geometry.utils.rotate(mesh, theta_y, symmetry, rotate_x=True, ref_axis_pos=0.25)[source]
Compute rotation matrices given mesh and rotation angles in degrees.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- theta_y[ny]numpy array
1-D array of rotation angles about y-axis for each wing slice in degrees.
- symmetryboolean
Flag set to True if surface is reflected about y=0 plane.
- rotate_xboolean
Flag set to True if the user desires the twist variable to always be applied perpendicular to the wing. To clarify, this is to ensure that non-planar surfaces such as winglets are twisted correctly about their axis. The winglets themselves have to be created with zshear or a user created mesh.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the twisted aerodynamic surface.
- openaerostruct.geometry.utils.scale_x(mesh, chord_dist, ref_axis_pos=0.25)[source]
Modify the chords along the span of the wing by scaling only the x-coord.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- chord_dist[ny]numpy array
Spanwise distribution of the chord scaler.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh with the new chord lengths.
- openaerostruct.geometry.utils.shear_x(mesh, xshear)[source]
Shear the wing in the x direction (distributed sweep).
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- xshear[ny]numpy array
Distance to translate wing in x direction.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh with the new chord lengths.
- openaerostruct.geometry.utils.shear_y(mesh, yshear)[source]
Shear the wing in the y direction (distributed span).
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- yshear[ny]numpy array
Distance to translate wing in y direction.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh with the new span widths.
- openaerostruct.geometry.utils.shear_z(mesh, zshear)[source]
Shear the wing in the z direction (distributed dihedral).
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- zshear[ny]numpy array
Distance to translate wing in z direction.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh with the new chord lengths.
- openaerostruct.geometry.utils.stretch(mesh, span, symmetry, ref_axis_pos=0.25)[source]
Stretch mesh in spanwise direction to reach specified span.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- spanfloat
Relative stetch ratio in the spanwise direction.
- symmetryboolean
Flag set to true if surface is reflected about y=0 plane.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the stretched aerodynamic surface.
- openaerostruct.geometry.utils.sweep(mesh, sweep_angle, symmetry)[source]
Apply shearing sweep. Positive sweeps back.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- sweep_anglefloat
Shearing sweep angle in degrees.
- symmetryboolean
Flag set to true if surface is reflected about y=0 plane.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the swept aerodynamic surface.
- openaerostruct.geometry.utils.taper(mesh, taper_ratio, symmetry, ref_axis_pos=0.25)[source]
Alter the spanwise chord linearly to produce a tapered wing. Note that we apply taper around the quarter-chord line.
- Parameters:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the initial aerodynamic surface.
- taper_ratiofloat
Taper ratio for the wing; 1 is untapered, 0 goes to a point.
- symmetryboolean
Flag set to true if surface is reflected about y=0 plane.
- Returns:
- mesh[nx, ny, 3]numpy array
Nodal mesh defining the tapered aerodynamic surface.
- openaerostruct.geometry.utils.unify_mesh(sections, shift_uni_mesh=True)[source]
Function that produces a unified mesh from all the individual wing section meshes.
- Parameters:
- sectionslist
List of section OpenAeroStruct surface dictionaries
- shift_uni_meshbool
Flag that shifts sections so that their leading edges are coincident. Intended to keep sections from seperating or intersecting during scalar span or sweep operations without the use of the constraint component.
- Returns:
- uni_meshnumpy array
Unfied surface mesh in OAS format