RCAIDE.Library.Methods.Aerodynamics.Vortex_Lattice_Method.VLM
VLM#
- VLM(conditions, settings, geometry)[source]#
Uses the vortex lattice method to compute the lift, induced drag and moment coefficients.
The user has the option to discretize control surfaces using the boolean settings.discretize_control_surfaces. The user should be forwarned that this will cause very slight differences in results for 0 deflection due to the slightly different discretization.
The user has the option to use the boundary conditions and induced velocities from either RCAIDE or VORLAX. See build_RHS in compute_RHS_matrix.py for more details.
By default in Vortex_Lattice, VLM performs calculations based on panel coordinates with float32 precision. The user may also choose to use float16 or float64, but be warned that the latter can be memory intensive.
The user should note that fully capitalized variables correspond to a VORLAX variable of the same name
Assumptions: The user provides either global discretezation (number_spanwise/chordwise_vortices) or separate discretization (wing/fuselage_spanwise/chordwise_vortices) in settings, not both. The set of settings not being used should be set to None.
The VLM requires that the user provide a non-zero velocity that matches mach number. For surrogate training cases at mach 0, VLM uses a velocity of 1e-6 m/s
Source: 1. Miranda, Luis R., Robert D. Elliot, and William M. Baker. “A generalized vortex lattice method for subsonic and supersonic flow applications.” (1977). (NASA CR)
VORLAX Source Code
Inputs: geometry.
reference_area [m^2] wing.
spans.projected [m] chords.root [m] chords.tip [m] sweeps.quarter_chord [radians] taper [Unitless] twists.root [radians] twists.tip [radians] symmetric [Boolean] aspect_ratio [Unitless] areas.reference [m^2] vertical [Boolean] origin [m]
- fuselage.
origin [m] width [m] heights.maximum [m] lengths.nose [m] lengths.tail [m] lengths.total [m] lengths.cabin [m] fineness.nose [Unitless] fineness.tail [Unitless]
settings.number_of_spanwise_vortices [Unitless] <—| settings.number_of_chordwise_vortices [Unitless] <—|
|–Either/or; see generate_vortex_distribution() for more details
settings.wing_spanwise_vortices [Unitless] <—| settings.wing_chordwise_vortices [Unitless] <—| settings.fuselage_spanwise_vortices [Unitless] <—| settings.fuselage_chordwise_vortices [Unitless] <—|
settings.use_surrogate [Unitless] settings.propeller_wake_model [Unitless] settings.discretize_control_surfaces [Boolean], set to True to generate control surface panels settings.use_VORLAX_matrix_calculation [boolean] settings.floating_point_precision [float16/32/64]
conditions.aerodynamics.angles.alpha [radians] conditions.aerodynamics.angles.beta [radians] conditions.freestream.mach_number [Unitless] conditions.freestream.velocity [m/s] conditions.static_stability.pitch_rate [radians/s] conditions.static_stability.roll_rate [radians/s] conditions.static_stability.yaw_rate [radians/s]
Outputs: results.
CL [Unitless], CLTOT in VORLAX CDi [Unitless], CDTOT in VORLAX CM [Unitless], CMTOT in VORLAX CY [Unitless], Total y force coeff CRTOT [Unitless], Rolling moment coeff (unscaled) CL_mom [Unitless], Rolling moment coeff (scaled by b_ref) CNTOT [Unitless], Yawing moment coeff (unscaled) CN [Unitless], Yawing moment coeff (scaled by b_ref) CL_wing [Unitless], CL of each wing CDi_wing [Unitless], CDi of each wing cl_y [Unitless], CL of each strip cdi_y [Unitless], CDi of each strip alpha_i [radians] , Induced angle of each strip in each wing (array of numpy arrays) CP [Unitless], Pressure coefficient of each panel gamma [Unitless], Vortex strengths of each panel
Properties Used: N/A
- compute_rotation_effects(VD, settings, EW_small, GAMMA, len_mach, X, CHORD, XLE, XBAR, rhs, COSINP, SINALF, COSCOS, PITCH, ROLL, YAW, STB, RNMAX)[source]#
This computes the effects of the freestream and aircraft rotation rate on CLE, the induced flow at the leading edge
Assumptions: Several of the values needed in this calculation have been computed earlier and stored in VD
Normally, VORLAX skips the calculation implemented in this function for linear chordwise spacing (the if statement below). However, since the trends are correct, albeit underestimated, this calculation is being forced here.