Index

Public API

Creating an aerodynamic derivatives solver

AeroBeams.FDType
FD (Finite Differences) DerivationMethod composite type

Fields

  • method::FiniteDifferenceMethod
source

Creating an aerodynamic solver

Creating an airfoil

AeroBeams.create_AirfoilFunction
create_Airfoil(; kwargs...)

Airfoil constructor

Arguments

  • name::String: name of the airfoil
  • Re::Real: Reynolds number
  • Ma::Real: Mach number
  • U::Real: relative airspeed
  • b::Real: semichord
source
AeroBeams.create_flapped_AirfoilFunction
create_flapped_Airfoil(; kwargs...)

Airfoil constructor (with a trailing-edge flap)

Arguments

  • name::String: name of the airfoil
  • flapSiteID::Int64: flap site ID
  • Re::Real: Reynolds number
  • Ma::Real: Mach number
  • U::Real: relative airspeed
  • b::Real: semichord
source
AeroBeams.update_Airfoil_params!Function
update_Airfoil_params!(airfoil::Airfoil; kwargs...)

Updates the airfoil parameters

Arguments

  • airfoil::Airfoil

Keyword arguments

  • Re::Real: Reynolds number
  • Ma::Real: Mach number
  • U::Real: reference relative airspeed
  • b::Real: reference semichord
source

Creating an aerodynamic surface

AeroBeams.create_AeroSurfaceFunction
create_AeroSurface(; kwargs...)

Aerodynamic surface constructor

Keyword arguments

  • solver::AeroSolver: aerodynamic solver for pitch-plunge-induced loads
  • flapLoadsSolver::FlapAeroSolver: aerodynamic solver for flap-induced loads
  • gustLoadsSolver::GustAeroSolver: aerodynamic solver for gust-induced loads
  • derivationMethod::DerivationMethod: method for calculation of aerodynamic derivatives
  • airfoil::Airfoil: airfoil section
  • c::Union{<:Function,Real}: chord
  • Λ::Union{<:Function,Real}: sweep angle
  • φ::Union{<:Function,Real}: twist angle
  • normSparPos::Union{<:Function,Float64}: normalized position of the spar (beam reference line) on the chord
  • normFlapSpan::Union{Nothing,Vector{<:Real}}: normalized position of the trailing-edge flap along the span (beam arclength)
  • normFlapPos::Union{Nothing,Float64}: normalized position of the trailing-edge flap hinge on the chord
  • δIsTrimVariable::Bool: flag for trailing-edge deflection being a trim variable
  • δ::Union{Nothing,<:Function,Real}: trailing-edge deflection [rad]
  • flapSiteID::Union{Nothing,Int64}: trailing-edge flap site ID
  • updateAirfoilParameters::Bool: flag to update airfoil parameters with local airspeed
  • hasTipCorrection::Bool: flag to employ a tip correction on aerodynamic coefficients
  • tipLossFunction::String: respective tip loss function
  • tipLossDecayFactor::Float64: respective tip loss factor
  • smallAngles::Bool: flag to employ small angles approximation on the calculation of the angle of attack
  • hasInducedDrag::Bool: flag to include induced drag generated by the lifting surface
  • hasSymmetricCounterpart::Bool: flag to indicate the surface has a symmetric counterpart (for the computation of aspect ratio)
source

Creating an atmosphere

AeroBeams.standard_atmosphereFunction
standard_atmosphere(altitude::Real)

Atmosphere constructor based on the International Standard Atmosphere (https://en.wikipedia.org/wiki/InternationalStandardAtmosphere)

Arguments

  • altitude::Real: altitude
source

Creating boundary conditions

AeroBeams.create_BCFunction
create_BC(; kwargs...)

BC constructor

Keyword arguments

  • name::String: name of the BC
  • beam::Beam: beam at which the BC is applied
  • node::Int64: node of the beam at which the BC is applied
  • types::Vector{String}= types of BCs applied to the node (generalized forces and displacements)
  • values: corresponding values of the applied BCs (constants or functions of time)
  • toBeTrimmed::Union{BitVector,Vector{Bool}}: TF on whether the BC is to be trimmed
source

Creating a beam

AeroBeams.create_BeamFunction
create_Beam(; kwargs...)

Beam constructor

Keyword Arguments

  • name::String: name of the beam
  • length::Real: (arc)length of the beam
  • rotationParametrization::String: type of rotation parametrization to define basis b
  • p0::Vector{<:Real}: rotation parameters from basis A to basis b
  • k::Union{Vector{<:Real},<:Function}: undeformed beam's curvatures per unit length
  • initialPosition::Vector{<:Real}: initial position of the beam's first node relative to the beam's origin (which may be another beam's node)
  • nElements::Int64: number of elements for discretization
  • normalizedNodalPositions::Vector{Float64}: normalized nodal positions of beam elements
  • S::Vector{<:Matrix{<:Real}}: array of sectional stiffness matrices
  • I::Vector{<:Matrix{<:Real}}: array of sectional inertia matrices
  • connectedBeams::Union{Nothing,Vector{Beam}}: array of beams to which this beam is connected (a non-recursive property)
  • connectedNodesThis::Vector{Int64}: nodes of this beam which are connected to other beams' nodes
  • connectedNodesOther::Vector{Int64}: respective nodes of the other beams
  • pointInertias::Vector{PointInertia}: attached point inertias
  • hingedNodes::Vector{Int64}: nodes with a hinge
  • hingedNodesDoF::Union{Vector{Vector{Bool}},Vector{BitVector}}: respective hinged degrees-of-freedom
  • u0_of_x1::Union{Vector{<:Real},<:Function,Nothing}: initial displacement (resolved in the undeformed beam basis, b) of the beam as a function of its arclength coordinate (x1)
  • p0_of_x1::Union{Vector{<:Real},<:Function,Nothing}: initial rotation parameters (resolved in the undeformed beam basis, b) of the beam as a function of its arclength coordinate (x1)
  • udot0_of_x1::Union{Vector{<:Real},<:Function,Nothing}: initial displacement's rates (resolved in the undeformed beam basis, b) of the beam as a function of its arclength coordinate (x1)
  • pdot0_of_x1::Union{Vector{<:Real},<:Function,Nothing}: initial rotation parameters' rates (resolved in the undeformed beam basis, b) of the beam as a function of its arclength coordinate (x1)
  • f_A_of_x1t::Union{Nothing,<:Function}: distributed dead forces initially resolved in basis A, as a function of the beam arclength coordinate (x1) and time (t)
  • m_A_of_x1t::Union{Nothing,<:Function}: distributed dead moments initially resolved in basis A, as a function of the beam arclength coordinate (x1) and time (t)
  • f_b_of_x1t::Union{Nothing,<:Function}: distributed dead forces initially resolved in basis b, as a function of the beam arclength coordinate (x1) and time (t)
  • m_b_of_x1t::Union{Nothing,<:Function}: distributed dead moments initially resolved in basis b, as a function of the beam arclength coordinate (x1) and time (t)
  • ff_A_of_x1t::Union{Nothing,<:Function}: distributed follower forces initially resolved in basis A, as a function of the beam arclength coordinate (x1) and time (t)
  • mf_A_of_x1t::Union{Nothing,<:Function}: distributed moments forces initially resolved in basis A, as a function of the beam arclength coordinate (x1) and time (t)
  • ff_b_of_x1t::Union{Nothing,<:Function}: distributed follower forces initially resolved in basis b, as a function of the beam arclength coordinate (x1) and time (t)
  • mf_b_of_x1t::Union{Nothing,<:Function}: distributed follower moments initially resolved in basis b, as a function of the beam arclength coordinate (x1) and time (t)
  • aeroSurface::Union{Nothing,AeroSurface}: attached aerodynamic surface
  • springs::Vector{Spring}: array of attached springs
source
AeroBeams.add_point_inertias_to_beam!Function
add_point_inertias_to_beam!(beam::Beam; inertias::Vector{PointInertia})

Adds point inertias to the beam

Arguments

  • beam::Beam

Keyword arguments

  • inertias::Vector{PointInertia}
source
AeroBeams.add_loads_to_beam!Function
add_loads_to_beam!(beam::Beam; loadTypes::Vector{String},loadFuns::Vector{<:Function})

Adds loads to the beam

Arguments

  • beam::Beam

Keyword arguments

  • loadTypes::Vector{String}
  • loadFuns::Vector{<:Function}
source
AeroBeams.add_initial_displacements_and_velocities_to_beam!Function
add_initial_displacements_and_velocities_to_beam!(beam::Beam; conditionTypes::Vector{String},conditionFuns::Vector{<:Function})

Adds initial generalized displacements and velocities to the beam

Arguments

  • beam::Beam

Keyword arguments

  • conditionTypes::Vector{String}
  • conditionFuns::Vector{<:Function}
source
AeroBeams.add_springs_to_beam!Function
add_springs_to_beam!(; beam::Beam,springs::Vector{Spring})

Adds simply-attached springs to a beam

Keyword arguments

  • beam::Beam
  • springs::Vector{Spring}
source
AeroBeams.add_spring_to_beams!Function
add_spring_to_beams!(; beams::Vector{Beam},spring::Spring)

Adds a doubly-attached spring to the beams

Keyword arguments

  • beams::Vector{Beam}
  • spring::Spring
source

Creating a gust

AeroBeams.create_SharpEdgedGustFunction
create_SharpEdgedGust(; kwargs...)

Creates a sharp-edged gust

Keyword arguments

  • initialTime::Real: time when the gust begins
  • duration::Real: duration of the gust
  • convectiveVelocity::Real: convective velocity of the gust (in longitudinal direction)
  • verticalVelocity::Real: peak vertical velocity of the gust (in lift direction)
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
source
AeroBeams.create_OneMinusCosineGustFunction
create_OneMinusCosineGust(; kwargs...)

Creates a one-minus-cosine gust

Keyword arguments

  • initialTime::Real: time when the gust begins
  • duration::Real: duration of the gust
  • convectiveVelocity::Real: convective velocity of the gust (in longitudinal direction)
  • verticalVelocity::Real: peak vertical velocity of the gust (in lift direction)
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
source
AeroBeams.create_Continuous1DGustFunction
create_Continuous1DGust(; kwargs...)

Creates a continuous 1D gust

Keyword arguments

  • spectrum::String: spectrum of the PSD: von Kármán ("vK") or Dryden ("Dryden")
  • generationMethod::String: method for generation of gust velocity ("sinusoids" or "whiteNoise")
  • initialTime::Real: time when the gust begins
  • duration::Real: duration of the gust
  • generationDuration::Real: duration of the gust considered for its generation
  • components::Vector{Int}: components of gust velocity to consider
  • L::Real: turbulence length scale [m]
  • ωmin::Real: minimum frequency of the PSD [rad/s]
  • ωmax::Real: maximum frequency of the PSD [rad/s]
  • Uref::Real: reference airspeed
  • convectiveVelocity::Real: convective velocity of the gust (in longitudinal direction)
  • σ::Real: turbulence intensity (RMS) of gust velocity component
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
  • seed::Int64: seed for random numbers generation (reproducibility)
  • plotPSD::Bool: flag to plot the PSD
source
AeroBeams.create_DiscreteSpaceGustFunction
create_DiscreteSpaceGust(; kwargs...)

Creates a discrete space gust

Keyword arguments

  • type::String: type of gust
  • gustLength::Real: length of the gust (in longitudinal direction)
  • gustWidth::Real: width of the gust (in spanwise direction)
  • convectiveVelocity::Real: convective velocity of the gust (in longitudinal direction)
  • verticalVelocity::Real: peak vertical velocity of the gust (in lift direction)
  • c0::Vector{<:Real}: position vector of the front of the gust, resolved in the inertial basis
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
source
AeroBeams.create_Continuous1DSpaceGustFunction
create_Continuous1DSpaceGust(; kwargs...)

Creates a continuous 1D space gust

Keyword arguments

  • spectrum::String: spectrum of the PSD: von Kármán ("vK") or Dryden ("Dryden")
  • gustLength::Real: length of the gust (in longitudinal direction)
  • N::Int64: number of nodes for length discretization
  • L::Real: turbulence length scale
  • σ::Real: turbulence intensity (RMS) of gust velocity components
  • c0::Vector{<:Real}: position vector of the front of the gust, resolved in the inertial basis
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
  • seed::Int64: seed for random numbers generation (reproducibility)
  • plotPSD::Bool: flag to plot the PSD
source
AeroBeams.create_Continuous2DSpaceGustFunction
create_Continuous2DSpaceGust(; kwargs...)

Creates a continuous 2D space gust

Keyword arguments

  • spectrum::String: spectrum of the PSD: von Kármán ("vK") or Dryden ("Dryden")
  • gustLength::Real: length of the gust (in longitudinal direction)
  • gustWidth::Real: width of the gust (in lateral direction)
  • Nx::Int64: number of nodes for length discretization
  • Ny::Int64: number of nodes for width discretization
  • L::Real: turbulence length scale
  • σ::Real: turbulence intensity (RMS) of gust velocity component
  • c0::Vector{<:Real}: position vector of the front of the gust, resolved in the inertial basis
  • p::Vector{<:Real}: Euler rotation parameters (3-2-1 sequence) from the inertial basis to the gust basis
  • seed::Int64: seed for random numbers generation (reproducibility)
source
AeroBeams.create_TrimLoadsLinkFunction
create_TrimLoadsLink(; masterBC::BC,slaveBCs::Vector{BC})

Trim loads link constructor (a link between trim loads so that they are equal)

Keyword arguments

  • masterBC::BC: master BC
  • slaveBCs::Vector{BC}: slave BCs
source
AeroBeams.create_FlapLinkFunction
create_FlapLink(; kwargs...)

Flap link constructor (a link between flapped surfaces)

Keyword arguments

  • masterBeam::Beam: beam of the master surface
  • slaveBeams::Vector{Beam}: beams of the slave surfaces
  • δMultipliers::Vector{<:Real}: multiplication factors of flap deflection in slave surfaces relative to the master surface
source

Creating a model

AeroBeams.create_ModelFunction
create_Model(; kwargs...)

Creates a model

Keyword arguments

  • name::String: name of the model
  • units::create_UnitsSystem: units system
  • beams::Vector{Beam}: beams in the assembly
  • initialPosition::Vector{<:Real}: initial position vector of the first node of the first beam, resolved in the inertial frame I
  • gravityVector::Vector{<:Real}: gravity vector
  • BCs::Vector{BC}: boundary condtions
  • p_A0::Vector{Float64}: initial rotation parameters that bring basis I to basis A
  • u_A::Union{Vector{<:Real},<:Function,Nothing}: displacement of basis A relative to basis I
  • v_A::Union{Vector{<:Real},<:Function,Nothing}: velocity of basis A relative to basis I
  • ω_A::Union{Vector{<:Real},<:Function,Nothing}: angular velocity of basis A relative to basis I
  • vdot_A::Union{Vector{<:Real},<:Function,Nothing}: acceleration of basis A relative to basis I
  • ωdot_A::Union{Vector{<:Real},<:Function,Nothing}: angular acceleration of basis A relative to basis I
  • altitude::Union{Nothing,Real}: altidude
  • atmosphere::Union{Nothing,Atmosphere}: atmosphere
  • gust::Union{Nothing,Gust}: gust
  • trimLoadsLinks::Vector{TrimLoadsLink}: links between trim loads
  • flapLinks::Vector{FlapLink}: links between flapped surfaces
  • hingeAxisConstraints::Vector{HingeAxisConstraint}: hinge axis constraints
source
AeroBeams.set_motion_basis_A!Function
set_motion_basis_A!(; kwargs...)

Sets the motion of basis A into the model

Keyword arguments

  • model::Model: model
  • u_A::Union{Vector{<:Real},<:Function,Nothing}: displacement of basis A relative to basis I
  • v_A::Union{Vector{<:Real},<:Function,Nothing}: velocity of basis A relative to basis I
  • ω_A::Union{Vector{<:Real},<:Function,Nothing}: angular velocity of basis A relative to basis I
  • vdot_A::Union{Vector{<:Real},<:Function,Nothing}: acceleration of basis A relative to basis I
  • ωdot_A::Union{Vector{<:Real},<:Function,Nothing}: angular acceleration of basis A relative to basis I
source

Creating a point inertia

AeroBeams.PointInertiaType
PointInertia composite type

Fields

  • elementID::Int64: local element ID to which the point inertia is attached
  • mass::Real: mass
  • η::Vector{Real}: position relative to element's midpoint's reference line
  • Ixx::Real: mass moment of inertia about the x1-axis
  • Iyy::Real: mass moment of inertia about the x2-axis
  • Izz::Real: mass moment of inertia about the x3-axis
  • Ixy::Real: mass product of inertia with respect to the x1-axis and x2-axis
  • Ixz::Real: mass product of inertia with respect to the x1-axis and x3-axis
  • Iyz::Real: mass product of inertia with respect to the x2-axis and x3-axis
  • inertiaMatrix::Union{Matrix{Real},Matrix{Nothing}}: mass moment of inertia matrix (tensor)
source

Creating a hinge axis constraint

AeroBeams.create_HingeAxisConstraintFunction
HingeAxisConstraint(; kwargs...)

Hinge axis constraint constructor

Keyword arguments

  • beam::Beam: beam
  • localHingeAxis::Vector{<:Real}: three-element vector defining the hinge axis, resolved in the undeformed (b) basis of the beam
  • pHValue::Union{Real,Nothing}: constrained value for the rotation at the hinge
  • solutionMethod::String: solution method for the constraint ("addedResidual" or "appliedMoment")
source

Creating and solving a problem

AeroBeams.InitialVelocitiesUpdateOptionsType
InitialVelocitiesUpdateOptions composite type

Defines variables for the update of the initial velocities states

Fields

  • Δt::Real: time step
  • maxIter::Int64: maximum number of iterations
  • relaxFactor::Float64: relaxation factor
  • tol::Float64: convergence tolerance
  • displayProgress::Bool: flag to display progress
source
AeroBeams.create_SteadyProblemFunction
create_SteadyProblem(; kwargs...)

Steady problem constructor

Keyword arguments

  • model::Model: model
  • systemSolver::systemSolver: nonlinear system solver
  • getLinearSolution::Bool: flag to solve for linear structural solution
  • x0::Vector{Float64}: initial states
source
AeroBeams.create_TrimProblemFunction
create_TrimProblem(; kwargs...)

Trim problem constructor

Keyword arguments

  • model::Model: model
  • systemSolver::systemSolver: nonlinear system solver
  • getLinearSolution::Bool: flag to solve for linear structural solution
  • getInertiaMatrix::Bool: flag to compute inertia matrix
  • x0::Vector{Float64}: initial states
source
AeroBeams.create_EigenProblemFunction
create_EigenProblem(; kwargs...)

Eigen problem constructor

Keyword arguments

  • model::Model: model
  • systemSolver::systemSolver: nonlinear system solver
  • getLinearSolution::Bool: flag to solve for linear structural solution
  • nModes::Int64=Inf64: number of modes to be computed
  • frequencyFilterLimits::Vector{Float64}: limits of the frequency filter
  • normalizeModeShapes::Bool: flag to normalize mode shapes
  • x0::Vector{Float64}: initial states
  • refTrimProblem::Union{Nothing,TrimProblem}: reference trim problem containing steady solution
source
AeroBeams.create_DynamicProblemFunction
create_DynamicProblem(; kwargs...)

Dynamic problem constructor

Keyword arguments

  • model::Model: model
  • systemSolver::systemSolver: nonlinear system solver
  • getLinearSolution::Bool: flag to solve for linear structural solution
  • initialTime::Real: initial time
  • Δt::Union{Nothing,Real}: time step
  • finalTime::Union{Nothing,Real}: final time
  • timeVector::Union{Nothing,Vector{Float64}}: time vector
  • adaptableΔt::Bool: flag for adaptable time step
  • minΔt::Union{Nothing,Real}: minimum time step (when adaptable)
  • maxΔt::Union{Nothing,Real}: maximum time step (when adaptable)
  • δb::Float64: discontinuities boundary convergence norm
  • skipInitialStatesUpdate::Bool: flag to skip update of initial states
  • initialVelocitiesUpdateOptions::InitialVelocitiesUpdateOptions: options for the initial velocities update
  • trackingTimeSteps::Bool: flag to track time step solutions
  • trackingFrequency::Int64: frequency of time steps in which to track solution
  • displayProgress::Bool: flag to display progress
  • displayFrequency::Int64: frequency of time steps in which to display progress
  • x0::Vector{Float64}: initial states
source

Creating a sample model

AeroBeams.tip_loss_function_PazyFunction
tip_loss_function_Pazy(type::String,θ::Real,U::Real,Λ::Real)

Computes the tip loss function for the Pazy wing's tip correction function

Arguments

  • type::String: tip loss function type
  • θ::Real: root pitch angle [rad]
  • U::Real: airspeed [m/s]
  • Λ::Real: sweep angle [rad]
source
AeroBeams.typical_section_dataFunction
typical_section_data(name::String)

Returns the data of the typical section of given name

Arguments

  • name::String: name of the typical section

Outputs

  • a: position of elastic axis (semichords aft of midchord)
  • e: position of mass axis (semichords aft of midchord)
  • μ: mass ratio
  • rα²: squared radius of gyration in pitch
  • σ: ratio of plunge to pitch frequency in vacuo
  • ωα: pitch frequency in vacuo [rad/s]
  • c: chord [m]
source
AeroBeams.create_PazyFunction
create_Pazy(; kwargs...)

Creates the Pazy wing

Returns the wing model and geometrical properties

Arguments

  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • withSkin::Bool: flag for skin on
  • sweepStructuralCorrections::Bool: flag to apply ad hoc structural corrections on sectional stiffness matrix with sweep angle
  • updateAirfoilParameters::Bool: flag to update airfoil parameters with airspeed
  • upright::Bool: flag to set the wing in the upright position
  • airfoil::Airfoil: airfoil section
  • θ::Real: pitch angle at the root [rad]
  • Λ::Real: sweep angle [rad]
  • hasTipCorrection::Bool: flag for aerodynamic tip correction
  • GAy::Real: shear stiffness in the x2 direction
  • GAz::Real: shear stiffness in the x3 direction
  • altitude::Real: altitude
  • g::Real: acceleration of gravity
  • airspeed::Real: airspeed
  • smallAngles::Bool: flag for small angles approximation
  • tipLossType::String: tip loss function type
  • tipMass::Real: mass of a point inertia added to the tip of the wing
  • ηtipMass::Vector{<:Real}: position vector of the tip mass relative to the tip of the spar, resolved in the local basis
  • tipMassInertia::Matrix{<:Real}: mass moment of inertia matrix of the tip mass
  • additionalBCs::Vector{BC}: additional BCs (beyond the clamp)
  • gust::Union{Nothing,Gust}: gust
source
AeroBeams.create_PazyFFWTFunction
create_PazyFFWT(; kwargs...)

Creates a version of the Pazy wing with flared folding wingtip (FFWT)

Arguments

  • solutionMethod::String: solution method for the constraint ("addedResidual" or "appliedMoment")
  • airfoil::Airfoil: airfoil section
  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • withSkin::Bool = flag for skin on
  • sweepStructuralCorrections::Bool: flag to apply ad hoc structural corrections on sectional stiffness matrix with sweep angle
  • hasTipCorrection::Bool: flag for aerodynamic tip correction
  • GAy::Real: shear stiffness in the x2 direction
  • GAz::Real: shear stiffness in the x3 direction
  • hingeNode::Int64: hinge node
  • pitchAngle::Real: pitch angle [rad]
  • Λ::Real: sweep angle [rad]
  • foldAngle::Union{Real,Nothing}: fold angle [rad]
  • flareAngle::Real: flare angle [rad]
  • kSpring::Real: stiffness of the spring around the hinge for folding (combination of OOP bending and twist)
  • kIPBendingHinge::Real: stiffness of the rotational spring around the hinge for in-plane bending
  • g::Real: local acceleration of gravity
  • altitude::Real: altitude
  • airspeed::Real: local airspeed
  • tipLossType::String: tip loss function type
  • flightDirection::Vector{<:Real}: flight direction vector, resolved in basis A
  • tipMass::Real: tip mass for passive flutter control
  • ηtipMass::Vector{<:Real}: position vector of the tip mass, relative to the spar
  • gust::Union{Nothing,Gust}=nothing: gust
source
AeroBeams.create_SMWFunction
create_SMW(; kwargs...)

Creates the 16 meters wing (SMW) of the conventional HALE aircraft described by Patil, Hodges and Cesnik in: Nonlinear Aeroelasticity and Flight Dynamics of HALE (2001)

Returns the wing model and its span

Arguments

  • aeroSolver::AeroSolver: aerodynamic solver
  • flapLoadsSolver::FlapAeroSolver: aerodynamic solver for flap loads
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • airfoil::Airfoil: airfoil section
  • θ::Real: pitch angle [rad]
  • k1::Real: twisting curvature
  • k2::Real: flapwise bending curvature
  • nElem::Int64: number of elements for discretization
  • altitude::Real: altitude
  • airspeed::Real: airspeed
  • g::Real: acceleration of gravity
  • stiffnessFactor::Real: stiffness factor for beam structural properties
  • ∞::Real: value of rigid structural properties
  • Ψ::Real: torsion-IP-bending stiffness coupling factor
  • tipF3::Union{Real,<:Function}: tip dead transverse force applied at the tip
  • cd0::Real: parasite drag coefficient for the wing
  • cnδ::Real: cn vs δ slope for the wing
  • cmδ::Real: cm vs δ slope for the wing
  • cdδ::Real: cd vs δ slope for the wing
  • hasInducedDrag::Bool: flag to include induced drag
  • δAil::Union{Nothing,Real,<:Function}: aileron deflection
  • additionalBCs::Vector{BC}: additional BCs (beyond the clamp and tip load)
source
AeroBeams.create_TDWingFunction
create_TDWing(; kwargs...)

Creates the wing described by Tang and Dowell in: Experimental and Theoretical Study on Aeroelastic Response of High-Aspect-Ratio Wings (2001)

Returns the wing model

Arguments

  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • updateAirfoilParameters::Bool: flag to update airfoil parameters with airspeed
  • airfoil::Airfoil: airfoil section
  • θ::Real: pitch angle
  • nElem::Int64: number of elements for discretization
  • altitude::Real: altitude
  • airspeed::Real: airspeed
  • g::Real: acceleration of gravity
  • ∞::Real: value of rigid structural properties
source
AeroBeams.create_HeliosFunction
create_Helios(; kwargs...)

Creates a model based on the flying-wing aircraft described by Patil and Hodges in: Flight Dynamics of Highly Flexible Flying Wings (2006)

Keyword arguments

  • altitude::Real: altitude
  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • g::Real: local acceleration of gravity
  • wingAirfoil::Airfoil: airfoil section of the wing
  • podAirfoil::Airfoil: airfoil section of the pods
  • beamPods::Bool: flag to include pods
  • stiffnessFactor::Real: stiffness factor for the wing structure
  • ∞::Real: value for rigid structural properties
  • nElemStraightSemispan::Int64: number of elements in the straight section of the semispan
  • nElemDihedralSemispan::Int64: number of elements in the dihedral section of the semispan
  • nElemPod::Int64: number of elements in the pods
  • payloadPounds::Real: payload, in pounds
  • airspeed::Real: local initial/trim airspeed
  • δIsTrimVariable::Bool: flag for flap deflection being a trim variable
  • thrustIsTrimVariable::Bool: flag for motors' thrust being a trim variable
  • δ::Union{Nothing,Real,<:Function}: flap deflection
  • thrust::Union{Real,<:Function}: motors' thrust
  • reducedChord::Bool: flag to employ a reduced (7 ft) chord
  • payloadOnWing::Bool: flag to set the payload on the wing's reference line
  • wingRootAoA::Real: root angle of attack for clamped wing model [rad]
  • gust::Union{Nothing,Gust}: gust
source
AeroBeams.create_conventional_HALEFunction
create_conventional_HALE(; kwargs...)

Creates a model based on the conventional HALE aircraft described by Patil, Hodges and Cesnik in: Nonlinear Aeroelasticity and Flight Dynamics of HALE (2001)

Keyword arguments

  • altitude::Real: altitude
  • aeroSolver::AeroSolver: aerodynamic solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • flapLoadsSolver::FlapAeroSolver: aerodynamic solver for flap loads
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • stabilizersAero::Bool: flag for stabilizers with aerodynamic surfaces
  • includeVS::Bool: flag to include a vertical stabilizer in the model
  • wAirfoil::Airfoil: wing airfoil section
  • sAirfoil::Airfoil: stabilizers airfoil section
  • nElemWing::Int64: number of elements of the full wing
  • nElemHorzStabilizer::Int64: number of elements of the horizontal stabilizer
  • nElemTailBoom::Int64: number of elements of the tail boom
  • nElemVertStabilizer::Int64: number of elements of the vertical stabilizer
  • ∞::Real=1e12: value of rigid structural properties
  • stiffnessFactor::Real: stiffness factor for the wing structure
  • k1::Real: undeformed wing torsional curvature
  • k2::Real: undeformed wing flapwise bending curvature
  • airspeed::Real: local initial/trim airspeed
  • δElevIsTrimVariable::Bool: flag for elevator deflection being a trim variable
  • δAilIsTrimVariable::Bool: flag for aileron deflection being a trim variable
  • δRuddIsTrimVariable::Bool: flag for rudder deflection being a trim variable
  • thrustIsTrimVariable::Bool: flag for motors' thrust being a trim variable
  • δElev::Union{Nothing,Real,<:Function}: elevator deflection
  • δAil::Union{Nothing,Real,<:Function}: aileron deflection
  • δRudd::Union{Nothing,Real,<:Function}: rudder deflection
  • thrust::Union{Real,<:Function}: motors' thrust
  • g::Real: local acceleration of gravity
  • wingCd0::Real: parasite drag coefficient for the wing
  • wingcnδ::Real: cn vs δ slope for the wing
  • wingcmδ::Real: cm vs δ slope for the wing
  • wingcdδ::Real: cd vs δ slope for the wing
  • stabsCd0::Real: parasite drag coefficient for the stabilizers
  • stabscnδ::Real: cn vs δ slope for the stabilizers
  • stabscmδ::Real: cm vs δ slope for the stabilizers
  • stabscdδ::Real: cd vs δ slope for the stabilizers
  • hasInducedDrag::Bool: flag to include induced drag
source
AeroBeams.create_BWBFunction
create_BWB(; kwargs...)

Creates a model based on the blended-wing-body described by Weihua Su's PhD thesis

Keyword arguments

  • altitude::Real: altitude
  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • ∞::Real=1e12: value of rigid structural properties
  • stiffnessFactor::Real: stiffness factor for the wing structure
  • airspeed::Real = local initial/trim airspeed
  • δElevIsTrimVariable::Bool: flag for elevator deflection being a trim variable
  • thrustIsTrimVariable::Bool: flag for motors' thrust being a trim variable
  • δElev::Union{Nothing,Real,<:Function}: elevator deflection
  • thrust::Union{Real,<:Function}: motors' thrust
  • g::Real: local acceleration of gravity
  • updateAirfoilParameters::Bool: flag to update airfoil parameters with airspeed
  • hasTipCorrection::Bool: flag to employ aerodynamic tip correction
  • tipLossDecayFactor::Real: tip loss decay factor
source
AeroBeams.create_HealyFFWTFunction
create_HealyFFWT(; kwargs...)

Creates a version of Healy's wing with flared folding wingtip (FFWT). See Healy's PhD thesis, chapter 7.

Arguments

  • solutionMethod::String: solution method for the constraint ("addedResidual" or "appliedMoment")
  • airfoil::Airfoil: airfoil section
  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • hasTipCorrection::Bool: flag for aerodynamic tip correction
  • tipLossDecayFactor::Real: tip loss decay factor, in case hasTipCorrection is true
  • pitchAngle::Real: root pitch angle [rad]
  • wingtipTwist::Real: wingtip twist angle [rad]
  • foldAngle::Union{Real,Nothing}: fold angle [rad]
  • flareAngle::Real: flare angle [rad]
  • kSpring::Real: stiffness of the spring around the hinge for folding (combination of OOP bending and twist)
  • kIPBendingHinge::Real: stiffness of the rotational spring around the hinge for in-plane bending
  • g::Real: local acceleration of gravity
  • altitude::Real: altitude
  • airspeed::Real: local airspeed
  • flightDirection::Vector{<:Real}: flight direction vector, resolved on basis A
  • nElementsInner::Int64: number of elements for discretization of inner part of the wing
  • nElementsFFWT::Int64: number of elements for discretization of the wingtip
  • gust::Union{Nothing,Gust}=nothing: gust
source
AeroBeams.create_HealyBaselineFFWTFunction
create_HealyBaselineFFWT(; kwargs...)

Creates a version of Healy's baseline wing with flared folding wingtip (FFWT). See Healy's PhD thesis, chapter 3.

Arguments

  • solutionMethod::String: solution method for the constraint ("addedResidual" or "appliedMoment")
  • airfoil::Airfoil: airfoil section
  • aeroSolver::AeroSolver: aerodynamic solver
  • gustLoadsSolver::GustAeroSolver: indicial gust loads solver
  • derivationMethod::DerivationMethod: method for aerodynamic derivatives
  • hasTipCorrection::Bool: flag for aerodynamic tip correction
  • tipLossDecayFactor::Real: tip loss decay factor, in case hasTipCorrection is true
  • hingeConfiguration::String: hinge hingeConfiguration ("free" or "locked")
  • flareAngle::Real: flare angle [rad]
  • pitchAngle::Real: root pitch angle [rad]
  • foldAngle::Union{Real,Nothing}: fold angle [rad]
  • kSpringHinge::Real: stiffness of the rotational spring around the hinge for folding (combination of OOP bending and twist)
  • kIPBendingHinge::Real: stiffness of the rotational spring around the hinge for in-plane bending
  • altitude::Real: altitude
  • airspeed::Real: local airspeed
  • g::Real: local acceleration of gravity
  • nElementsInner::Int64: number of elements for discretization of inner part of the wing
  • nElementsFFWT::Int64: number of elements for discretization of the wingtip
  • gust::Union{Nothing,Gust}=nothing: gust
source

Creating a spring

AeroBeams.create_SpringFunction
create_Spring(; kwargs...)

Creates a spring

Keyword arguments

  • basis::String: basis on which stiffnesses are defined ("b" or "A")
  • elementsIDs::Vector{Int64}: local IDs of the element(s)' node(s) to which the spring is attached
  • nodesSides::Vector{Int64}: sides (1 or 2) of the node(s) to which the spring is attached
  • ku::Vector{<:Real}: translational stiffness vector, resolved in basis A
  • kp::Vector{<:Real}: rotational stiffness vector, resolved in basis A
source

Creating a system solver

AeroBeams.create_NewtonRaphsonFunction
create_NewtonRaphson(; kwargs...)

Newton-Raphson nonlinear system solver constructor

Keyword arguments

  • absoluteTolerance::Float64: absolute convergence tolerance
  • relativeTolerance::Float64: relative convergence tolerance
  • maximumIterations::Int64: maximum number of iterations
  • desiredIterations::Int64: desired number of iterations
  • maximumAbsoluteError::Real: maximum absolute error for divergence detection
  • maximumRelativeError::Real: maximum relative error for divergence detection
  • initialLoadFactor::Real: initial load factor
  • minimumLoadFactor::Float64: minimum load factor
  • maximumLoadFactorStep::Float64: maximum load factor step
  • minimumLoadFactorStep::Float64: minimum load factor step
  • ρ::Real: relaxation factor for trim variables
  • trackingLoadSteps::Bool: flag to track partial load steps solutions
  • displayStatus::Bool: flag to display status
  • minConvRateAeroJacUpdate::Real: minimum convergence rate to skip computation of aerodynamic Jacobians
  • minConvRateJacUpdate::Real: minimum convergence rate to skip computation of structural Jacobians
  • alwaysUpdateJacobian::Bool: flag to update Jacobians on every iteration
  • allowAdvanceThroughUnconvergedAeroStates::Bool: flag to allow the advancement of the algorithm through unconverged aerodynamic states
  • aeroStatesResidualRatioThreshold::Float64: threshold ratio of aerodynamic/total residuals in order to allow the above flag
  • aeroStatesRelativeErrorThreshold::Float64: threshold relative error in order to allow the above flag
source

Creating a units system

AeroBeams.create_UnitsSystemFunction
create_UnitsSystem(; kwargs...)

Creates a system composed of length, force, angle, frequency and mass units (this is only for plotting purposes and does not influence calculations)

Keyword arguments

  • length::String: length unit
  • force::String: force unit
  • angle::String: angle unit
  • frequency::String: frequency unit
  • mass::String: mass unit
source

Utilities

AeroBeams.round_off!Function
round_off!(x)

Rounds the array or number to input tolerance (defaults to machine epsilon)

Arguments

  • x: array or number
  • tol: tolerance
source
AeroBeams.tildeFunction
tilde(v)

Computes the skew-symmetric matrix associated with a three-element vector

Arguments

  • v: vector
source
AeroBeams.gauss_legendre7Function
gauss_legendre7(f, a::Real, b::Real)

Computes the Gauss-Legendre quadrature (integral) for the vector-valued function f in the interval from a to b, using 7 points

Arguments

  • f: function
  • a::Real: lower limit
  • b::Real: upper limit
source
AeroBeams.isotropic_stiffness_matrixFunction
isotropic_stiffness_matrix(; kwargs...)

Creates a 6x6 sectional stiffness matrix for a cross-section made of isotropic material

Arguments

  • : value for rigid properties
  • EA: axial stiffness
  • GAy: shear stiffness in the x2 direction
  • GAz: shear stiffness in the x3 direction
  • GJ: torsional stiffness
  • EIy: bending stiffness in the x2 direction
  • EIz: bending stiffness in the x3 direction
  • t2: offset from reference line to tension center in the x2 direction
  • t3: offset from reference line to tension center in the x3 direction
  • s2: offset from reference line to shear center in the x2 direction
  • s3: offset from reference line to shear center in the x3 direction
source
AeroBeams.inertia_matrixFunction
inertia_matrix(; kwargs...)

Creates a 6x6 sectional inertia matrix

Arguments

  • ρA: mass per unit length
  • ρIy: mass moment of inertia per unit length about the x2-axis
  • ρIz: mass moment of inertia per unit length about the x3-axis
  • ρIs: mass moment of inertia per unit length about the x1-axis
  • ρIyz: mass product of inertia per unit length
  • e2: offset from reference line to center of gravity in the x2 direction
  • e3: offset from reference line to center of gravity in the x3 direction
source
AeroBeams.rotation_angleFunction
rotation_angle(p)

Computes the rotation angle given the Wiener-Milenkovic rotation parameters

Arguments

  • p: rotation parameters
source
AeroBeams.rotation_angle_limitedFunction
rotation_angle_limited(p)

Computes the rotation angle (in the range -360 to 360 degrees) given the Wiener-Milenkovic rotation parameters

Arguments

  • p: rotation parameters
source
AeroBeams.yrp_from_rotation_tensorFunction
yrp_from_rotation_tensor(R; ϵround=1e-10,ϵsingularity=1e-3,assumeNullYawInSingularity=true)

Computes the Euler angles from the sequence 3-2-1 (yaw, roll, pitch) given the rotation tensor

Arguments

  • R: rotation tensor

Keyword arguments

  • ϵround: tolerance for rounding off elements of R to zero
  • ϵsingularity: tolerance to consider singular case
  • assumeNullYawInSingularity: flag to assume zero yaw angle in singularity
source
AeroBeams.WM_to_yrpFunction
WM_to_yrp(p)

Transforms Wiener-Milenkovic parameters to Euler parameters of sequence 3-2-1 (yaw, roll, pitch)

Arguments

  • p: rotation parameters
source
AeroBeams.yrp_to_WMFunction
yrp_to_WM(p)

Transforms Euler parameters of sequence 3-2-1 (yaw, roll, pitch) to Wiener-Milenkovic parameters

Arguments

  • p: rotation parameters
source
AeroBeams.rotation_between_WMFunction
rotation_between_WM(p1,p2)

Computes the Wiener-Milenkovic parameters describing the rotation from p1 to p2, i.e., p12 such that R2(p2) = R12(p12)*R1(p1)

Arguments

  • p1: initial rotation parameters
  • p2: final rotation parameters
source
AeroBeams.mode_trackingFunction
mode_tracking(controlParam::AbstractVector{<:Real},untrackedFreqs::Array{Vector{Float64}},untrackedDamps::Array{Vector{Float64}},untrackedEigenvectors::Array{Matrix{ComplexF64}})

Applies mode tracking based on eigenvectors match and a greedy search

Arguments

  • controlParam::AbstractVector{<:Real}: vector of control parameter
  • untrackedFreqs::Array{Vector{Float64}}: frequencies vector
  • untrackedDamps::Array{Vector{Float64}}: dampings vector
  • untrackedEigenvectors::Array{Matrix{ComplexF64}}: complex-valued eigenvectors
source
AeroBeams.mode_tracking_hungarianFunction
mode_tracking_hungarian(controlParam::AbstractVector{<:Real},untrackedFreqs::Array{Vector{Float64}},untrackedDamps::Array{Vector{Float64}},untrackedEigenvectors::Array{Matrix{ComplexF64}}; showProgress::Bool=false,α::Real=1.0,β::Real=0.01,γ::Real=0.01)

Applies mode tracking based on eigenvectors match and frequency and damping continuity, with a Hungarian algorithm

Arguments

  • controlParam::AbstractVector{<:Real}: vector of control parameter
  • untrackedFreqs::Array{Vector{Float64}}: frequencies vector
  • untrackedDamps::Array{Vector{Float64}}: dampings vector
  • untrackedEigenvectors::Array{Matrix{ComplexF64}}: complex-valued eigenvectors

Keyword arguments:

  • showProgress::Bool: flag to show progress over control parameter
  • α::Real: weight for eigenvector similarity
  • β::Real: weight for frequency continuity
  • γ::Real: weight for damping continuity
source
AeroBeams.moving_averageFunction
moving_average(v::AbstractVector{<:Real}, n::Int)

Computes the moving average of a series

Arguments

  • v::AbstractVector{<:Real}: series
  • n::Int: moving average window size
source
AeroBeams.get_FFT_and_PSDFunction
get_FFT_and_PSD(t::AbstractVector{<:Real},y::Vector{<:Real}; tol::Float=1e3*eps())

Computes the FFT and PSD of signal y(t)

Arguments

  • t::AbstractVector{<:Real}: time signal
  • y::Vector{<:Real}: quantity signal

Keyword arguments

  • tol::AbstractFloat: tolerance for time signal being equally spaced
source
AeroBeams.Newton_solverFunction
Newton_solver(f::Function, x0::AbstractArray; absTol::Real=1e-9, relTol::Real=1e-9, maxIter::Int=50)

Solves a nonlinear algebraic system of equations $( f(x) = 0 $) using the Newton-Raphson method.

Arguments

  • f::Function: the nonlinear system of equations to solve. Must return a vector of residuals for a given x.
  • x0::AbstractArray: initial guess for the solution.

Keyword arguments

  • absTol::Real: absolute tolerance for convergence (default: 1e-9).
  • relTol::Real: relative tolerance for convergence (default: 1e-9).
  • maxIter::Int: maximum number of iterations allowed (default: 50).

Returns

  • x::AbstractArray: the solution vector, if convergence is achieved.
  • converged::Bool: true if the solution converged within the given tolerances, false otherwise.
source
AeroBeams.backward_extrapolationFunction
backward_extrapolation(v::AbstractVector{<:Real}; order::Int=5)

Estimates the values of an array from backward finite difference extrapolation

Arguments

  • v::AbstractVector{<:Real}: vector

Keyword arguments

  • order::Int: order of the approximation
source
AeroBeams.count_of_exceedanceFunction
count_of_exceedance(time::AbstractVector{<:Real}, series::Vector{<:Real}, datum::Real, marker::Vector{<:Real})

Computes the count of exceedance of the quantities in the marker vector above the datum for the given time series

Keyword arguments

  • time::AbstractVector{<:Real}: time vector
  • series::Vector{<:Real}: time series
  • datum::Real: datum value for the time series
  • marker::Vector{<:Real}: markers for the count of exceedance computation
source
AeroBeams.frequency_of_exceedanceFunction
frequency_of_exceedance(timeLength::Real, count::Vector{<:Real})

Computes the frequency of exceedance

Keyword arguments

  • count::Vector{<:Real}: count of exceedance
  • timeLength::Real: length of time
source
AeroBeams.NASA_frames_loaderFunction
NASA_frames_loader(frame::Union{Int,String})

Loads dynamic stall data from a frame of the NASA report (https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19830003778.pdf)

Arguments

  • frame::Union{Int,String}: frame number
source
AeroBeams.GU_frames_loaderFunction
GU_frames_loader(frame::Union{Int,String})

Loads dynamic stall data from a frame of the University of Glasgow (DOI:10.5525/gla.researchdata.464)

Arguments

  • frame::Union{Int,String}: frame number
source

Visualizing the results

AeroBeams.plot_undeformed_assemblyFunction
plot_undeformed_assembly(model::Model; kwars...)

Plots the undeformed assembly of beams in the model

Arguments

  • model::Model

Keyword arguments

  • view::Tuple{Real,Real}: view angles
  • equalAspectRatio::Bool: flag to set equal aspect ratio plot
  • plotNodes::Bool: flag to plot nodes
  • nodesColor=:black: color of nodes
  • linesColor=:black: color of lines (beams)
source
AeroBeams.plot_steady_deformationFunction
plot_steady_deformation(problem::Problem; kwargs...)

Plots the initial and final deformed states for the model in the given problem

Arguments

  • problem::Problem: problem

Keyword arguments

  • interactive::Bool: flag for interactive plot
  • plotBCs::Bool: flag to plot BCs
  • plotDistLoads::Bool: flag to plot distributed loads (includes gravitational and aerodynamic loads)
  • view::Union{Nothing,Tuple{Real,Real}}: view angles
  • scale::Real: displacements and rotations scale
  • lw::Real: linewidth
  • colorUndef: color for undeformed assembly
  • colorDef: color for deformed assembly
  • plotGrid::Bool: flag to plot grid
  • legendPos: legend position
  • tolPlane::Real: displacement tolerance to plot as plane
  • plotAeroSurf::Bool: flag to plot aerodynamic surfaces
  • surfα::Float64: transparency factor of aerodynamic surfaces
  • ΔuDef::Vector{<:Real}: displacement vector for first node of deformed assembly relative to the undeformed one
  • σ2plot::Vector{<:Real}: (tentative) load factors to plot
  • partialLoadStepsLineStyle::Symbol: linestyle for partial load steps
  • plotLimits::Union{Nothing,Tuple{Vector{<:Real},Vector{<:Real},Vector{<:Real}}}: plot axis limits
  • showScale::Bool: flag to show scale on plot
  • scalePos::Vector{<:Real}: position of scale on plot
  • save::Bool: flag to save the figure
  • savePath::String: relative path on which to save the figure
source
AeroBeams.plot_steady_outputsFunction
plot_steady_outputs(problem::Problem; kwargs...)

Plots outputs of a steady problem

Arguments

  • problem::Problem: problem

Keyword arguments

  • outputs::Vector{String}: list of outputs
  • beamGroups: list of beams in each group for which arclengths are concatenated
  • lw::Real: linewidth
  • colorScheme: color scheme
  • legendPos: legend position
  • save::Bool: flag to save figures
  • saveFolder::String: relative path of folder where to save figures
  • figureExtension::String: figure extension
source
AeroBeams.plot_mode_shapesFunction
plot_mode_shapes(problem::EigenProblem; kwargs...)

Plots the mode shapes of the model in the given eigenproblem

Arguments

  • problem::EigenProblem: problem

Keyword arguments

  • interactive::Bool: flag for interactive plot
  • plotBCs::Bool: flag to plot BCs
  • view::Union{Nothing,Tuple{Real,Real}}: view angles
  • nModes::Union{Nothing,Int64}: number of modes to plot
  • scale::Real: displacements and rotations scale
  • ΔuDef::Vector{<:Real}: displacement vector for first node of deformed assembly relative to the undeformed one
  • frequencyLabel::String: option for frequency label
  • lw::Real: linewidth
  • colorSteady: color for steadily deformed assembly
  • modalColorScheme: color scheme for mode shapes
  • plotAxes::Bool: flag to plot axes
  • plotGrid::Bool: flag to plot grid
  • plotLimits::Union{Nothing,Tuple{Vector{Int64},Vector{Int64}}} = plot limits
  • legendPos: legend position
  • modeLabels::Union{Nothing,Vector{String}} = labels for each mode
  • tolPlane::Real: displacement tolerance to plot as plane
  • plotAeroSurf: flag to plot aerodynamic surfaces
  • surfα::Float64: transparency factor of aerodynamic surfaces
  • save::Bool: flag to save the figure
  • savePath::String: relative path on which to save the figure
source
AeroBeams.plot_dynamic_deformationFunction
plot_dynamic_deformation(problem::DynamicProblem; kwargs...)

Plots the animated deformation of the model in the given problem

Arguments

  • problem::DynamicProblem: problem

Keyword arguments

  • refBasis::String: reference observer basis for plot
  • plotFrequency::Int64: frequency of time steps to plot
  • plotUndeformed::Bool: flag to plot undeformed assembly
  • plotBCs::Bool: flag to plot BCs
  • plotDistLoads::Bool: flag to plot distributed loads (includes gravitational and aerodynamic loads)
  • view::Union{Nothing,Tuple{Real,Real}}: view angles
  • fps::Real: frame rate for gif
  • scale::Real: displacements and rotations scale
  • frequencyLabel::String: option for frequency label (only frequency or frequency and damping)
  • lw::Real: linewidth
  • colorUndef: color for undeformed assembly
  • colorDef: color for deformed assembly
  • plotGrid::Bool: flag to plot grid
  • legendPos: legend position
  • tolPlane::Real: displacement tolerance to plot as plane
  • plotAeroSurf: flag to plot aerodynamic surfaces
  • surfα::Float64: transparency factor of aerodynamic surfaces
  • plotLimits::Union{Nothing,Tuple{Vector{<:Real},Vector{<:Real},Vector{<:Real}}}: plot axis limits
  • followAssembly::Bool: flag to adjust the plot limits in order to "follow" the assembly of beams with the movement of basis A
  • save::Bool: flag to save the figure
  • savePath::String: relative path on which to save the figure
  • showScale::Bool: flag to show scale on plot
  • showTimeStamp::Bool: flag to show time stamp on plot
  • scalePos::Vector{<:Real}: position of scale on plot
  • timeStampPos::Vector{<:Real}: position of time stamp on plot
  • displayProgress::Bool: flag to display progress of gif creation
source
AeroBeams.plot_time_outputsFunction
plot_time_outputs(problem::Problem; kwargs...)

Plots outputs of a dynamic problem

Arguments

  • problem::Problem: problem

Keyword arguments

  • nodes::Vector{Tuple{Int64,Int64}} = global IDs of nodes for which to plot
  • elements::Vector{Int64} = global IDs of elements for which to plot
  • nodalOutputs::Vector{String} = nodal outputs to plot
  • elementalOutputs::Vector{String}: elemental outputs to plot
  • lw::Real: linewidth
  • colorScheme: color scheme
  • showLegend::Bool: flag to show legend
  • legendPos: legend position
  • save::Bool: flag to save figures
  • saveFolder::String: relative path of folder where to save figures
  • figureExtension::String: figure extension
source
AeroBeams.plot_snapshotsFunction
plot_snapshots(problem::DynamicProblem; kwargs...)

Plots the deformed states for the model in the given dynamic problem at specified time instances (snapshots)

Arguments

  • problem::DynamicProblem: problem

Keyword arguments

  • refBasis::String: reference observer basis for plot
  • plotBCs::Bool: flag to plot BCs
  • plotDistLoads::Bool: flag to plot distributed loads (includes gravitational and aerodynamic loads)
  • view::Union{Nothing,Tuple{Real,Real}}: view angles
  • scale::Real: displacements and rotations scale
  • lw::Real: linewidth
  • ls::Symbol= line style
  • color: color
  • plotAxes::Bool: flag to plot axes
  • plotGrid::Bool: flag to plot grid
  • tolPlane::Real: displacement tolerance to plot as plane
  • plotAeroSurf::Bool: flag to plot aerodynamic surfaces
  • surfα::Float64: transparency factor of aerodynamic surfaces
  • snapshots::Vector{<:Real}: (tentative) time instances to plot
  • plotLimits::Union{Nothing,Tuple{Vector{<:Real},Vector{<:Real},Vector{<:Real}}}: plot axis limits
  • showScale::Bool: flag to show scale on plot
  • scalePos::Vector{<:Real}: position of scale on plot
  • save::Bool: flag to save the figure
  • savePath::String: relative path on which to save the figure
source