@meta
CurrentModule = SymbolicAWEModels
API Reference
This page provides a detailed reference for all public functions exported by the SymbolicAWEModels.jl
package.
High-Level Simulation Functions
These functions provide convenient wrappers for running common simulation scenarios.
SymbolicAWEModels.sim!
— Functionsim!(sam, set_values; dt, total_time, vsm_interval, prn, lin_model)
Run a generic simulation for a given AWE model and a matrix of control inputs. Optionally, also simulate a provided linear model, returning both logs.
Arguments
sam::SymbolicAWEModel
: Initialized AWE model.set_values::Matrix{Float64}
: A matrix of external control torques [Nm] to be applied at each time step. The number of rows must equal the number of simulation steps, and the number of columns must equal the number of winches.
Keywords
dt::Float64
: Time step [s]. Defaults to1/sam.set.sample_freq
.total_time::Float64
: Total simulation duration [s]. Defaults to 10.0.vsm_interval::Int
: Interval for the value state machine updates. Defaults to 3.prn::Bool
: If true, prints a performance summary upon completion. Defaults to true.lin_model
: (optional) a continuous-timeStateSpace
object fromControlSystemsBase
. If provided, the linear model is simulated in parallel and a second log is returned.
Returns
- If
lin_model
is not provided:(SysLog, Nothing)
(nonlinear log, nothing) - If
lin_model
is provided:(SysLog, SysLog)
(nonlinear, linear logs)
SymbolicAWEModels.sim_oscillate!
— Functionsim_oscillate!(sam; dt, total_time, steering_freq, steering_magnitude, vsm_interval,
bias, prn, lin_model)
Run a simulation with oscillating steering input on the given AWE model. Optionally also simulate a provided linear model.
Keywords (see sim!)
lin_model
: (optional) a continuous-timeStateSpace
object fromControlSystemsBase
.
Returns
- If
lin_model
is not provided:(SysLog, Nothing)
(nonlinear log, nothing) - If
lin_model
is provided:(SysLog, SysLog)
(nonlinear, linear logs)
SymbolicAWEModels.sim_turn!
— Functionsim_turn!(sam; dt, total_time, steering_time, steering_magnitude, vsm_interval, prn,
lin_model)
Run a simulation with a constant steering input for a specified duration. Optionally also simulate a provided linear model.
Keywords (see sim!)
lin_model
: (optional) a continuous-timeStateSpace
object fromControlSystemsBase
.
Returns
- If
lin_model
is not provided:(SysLog, Nothing)
(nonlinear log, nothing) - If
lin_model
is provided:(SysLog, SysLog)
(nonlinear, linear logs)
SymbolicAWEModels.sim_reposition!
— Functionsim_reposition!(sam; dt, total_time, reposition_interval_s, target_elevation_deg,
target_azimuth_deg, prn)
Run a simulation that periodically resets the kite's elevation and azimuth.
This function simulates the AWE model and, at a specified time interval, calls reposition!
to reposition the kite to a target elevation and azimuth. It logs the entire simulation and returns a SysLog
.
Arguments
sam::SymbolicAWEModel
: Initialized AWE model.
Keywords
dt::Float64
: Time step [s]. Defaults to1/sam.set.sample_freq
.total_time::Float64
: Total simulation duration [s]. Defaults to 20.0.reposition_interval_s::Float64
: The interval in seconds at which to reset the pose. Defaults to 5.0.target_elevation::Float64
: The target elevation in rad for repositioning. Defaults to deg2rad(45.0).target_azimuth::Float64
: The target azimuth in rad for repositioning. Defaults to 0.0.target_heading::Float64
: The target heading in rad for repositioning. Defaults to 0.0.prn::Bool
: If true, prints status messages during the simulation. Defaults to true.
Returns
SysLog
: A log of the complete simulation.
Low-Level Simulation and Analysis
These functions provide direct control over the simulation and tools for model analysis.
KiteUtils.init!
— Functioninit!(sam::SymbolicAWEModel; ...)
Initialize the SymbolicAWEModel
.
This is the main entry point for setting up the model. It handles:
- Loading or building the symbolic model (
full_sys
). - Creating the
ODEProblem
,LinearizationProblem
, and control functions as needed. - Serializing the model to disk if it was newly built.
- Initializing the ODE integrator.
Keyword Arguments
solver
: The ODE solver to use. Ifnothing
, a default is chosen based on settings.adaptive::Bool
: Enable adaptive time-stepping for the solver.prn::Bool
: Enable printing of progress messages.remake::Bool
: Force a full rebuild of the symbolic model, ignoring any cached versions.reload::Bool
: Force reloading of the serialized model from disk.outputs
: A vector of variables to be treated as system outputs.create_prob::Bool
: Whether to create theODEProblem
.create_lin_prob::Bool
: Whether to create theLinearizationProblem
.create_control_func::Bool
: Whether to generate the control functions.lin_vsm::Bool
: Whether to linearize the aerodynamics using the Vortex Step Method (VSM) after initialization.
Returns
- The initialized
ODEIntegrator
.
KiteUtils.next_step!
— Functionnext_step!(s::SymbolicAWEModel, integrator::ODEIntegrator; set_values, dt, vsm_interval)
Take a simulation step, using the provided integrator.
This is a convenience method that calls the main next_step!
function.
next_step!(s::SymbolicAWEModel; set_values, dt, vsm_interval)
Take a simulation step forward in time.
This function advances the simulation by one time step, optionally updating control inputs and re-linearizing the VSM model. It then updates the SystemStructure
with the new state from the ODE integrator.
Arguments
s::SymbolicAWEModel
: The kite power system state object.
Keyword Arguments
set_values=nothing
: New values for the control inputs. Ifnothing
, the current values are used.dt=1/s.set.sample_freq
: Time step size [s].vsm_interval=1
: The interval (in steps) to re-linearize the VSM model. If 0, it is not re-linearized.
SymbolicAWEModels.find_steady_state!
— Functionfind_steady_state!(s::SymbolicAWEModel, integ=s.integrator; t=1.0, dt=1/s.set.sample_freq)
Run the simulation for a short period to allow the system to settle.
During this period, the winches are braked and the wing's elevation and azimuth angles are fixed, but it is free to move radially (in distance). This allows the dynamic components of the bridle and tethers to settle into a stable, steady-state equilibrium before starting a maneuver or analysis.
Arguments
s::SymbolicAWEModel
: The model to be stabilized.integ
: The integrator to use. Defaults tos.integrator
.
Keywords
t::Float64=1.0
: The duration [s] for which to run the settling simulation.dt::Float64
: The time step [s] for the settling simulation.
SymbolicAWEModels.linearize!
— Functionlinearize!(s::SymbolicAWEModel; set_values=s.get_set_values(s.integrator)) -> LinType
Compute the full state-space linearization of the model around the current operating point.
This function uses the LinearizationProblem
generated by ModelingToolkit.jl
to calculate the A, B, C, and D matrices for the complete, high-order system.
Arguments
s::SymbolicAWEModel
: The model to linearize.
Keywords
set_values
: The control input vectoru
around which to linearize.
Returns
LinType
: A NamedTuple(A, B, C, D)
containing the state-space matrices.
SymbolicAWEModels.copy_to_simple!
— Functioncopy_to_simple!(sam, tether_sam, simple_sam; prn=true)
Simplify a detailed AWE model into a 1-segment tether model.
This high-level function orchestrates the simplification process:
- Calculates the equivalent axial stiffness and damping of the detailed model (
tether_sam
) by analyzing its step response. - Assigns these calculated properties to the single-segment tethers of the simple model (
simple_sam
). - Copies the dynamic state (wing position, orientation, tether attachment points, etc.) from the detailed model (
sam
) to the simple model (simple_sam
). - Reinitializes the simple model's integrator to apply the new state.
Arguments
sam::SymbolicAWEModel
: The detailed source model, used as a reference for state.tether_sam::SymbolicAWEModel
: A copy of the detailed model, used to perform the step response test.simple_sam::SymbolicAWEModel
: The destination simple model to be updated.
Keywords
prn::Bool=true
: If true, enables printing during the process.
Returns
SymbolicAWEModel
: The updated simple modelsimple_sam
.
copy_to_simple!(sys::SystemStructure, ssys::SystemStructure)
Copy the dynamic state from a detailed SystemStructure
to a simplified one.
This is a low-level utility that maps the state of a complex model (e.g., "ram" with 4 groups and bridle pulleys) to a simpler model (e.g., "simple_ram" with 2 groups and direct connections). It ensures the simplified model matches the key dynamic properties of the detailed one.
Arguments
sys::SystemStructure
: The sourceram
model structure.ssys::SystemStructure
: The destinationsimple_ram
model structure.
SymbolicAWEModels.simple_linearize!
— Functionsimple_linearize!(s::SymbolicAWEModel; tstab=10.0) -> LinType
Compute a simplified, low-order state-space model by numerically linearizing the full simulation.
This function performs system identification by perturbing the states and inputs of the full nonlinear model and observing the effect on the state derivatives and outputs. It runs the simulation for a short duration (tstab
) after each perturbation to find the steady-state response.
Arguments
s::SymbolicAWEModel
: The model to linearize.
Keywords
tstab::Float64=10.0
: The simulation time [s] to run after each perturbation to reach a steady state.
Returns
LinType
: A NamedTuple(A, B, C, D)
containing the identified low-order state-space matrices.
System Structure Constructors
These functions are used to procedurally generate predefined SystemStructure
topologies.
SymbolicAWEModels.create_ram_sys_struct
— Functioncreate_ram_sys_struct(set::Settings)
Create a SystemStructure
for the primary "ram" model with a stability-enhancing bridle.
This is the main, detailed model configuration. It differs from the 4_attach
version by having each of its four Group
sections defined by three deforming bridle points and one statically attached (non-deforming) point. This design improves stability by ensuring the kite's z-axis remains aligned with the bridle system.
The model features:
- A flexible wing with 4 deformable groups (3 deforming points + 1 static point each).
- A complex bridle system with pulleys.
- Four main tethers and three winches.
Arguments
set::Settings
: Configuration parameters defining the kite's geometry and properties.
Returns
SystemStructure
: A newSystemStructure
object representing the "ram" model.
SymbolicAWEModels.create_tether_sys_struct
— Functioncreate_tether_sys_struct(set::Settings; axial_stiffness, axial_damping)
Create a simplified SystemStructure
for testing tether dynamics.
This model consists of only four independent tethers, each represented by a dynamic point mass connected to a fixed ground anchor. It does not include a wing or bridle system, making it ideal for isolating and analyzing the behavior of the tethers themselves.
Arguments
set::Settings
: Configuration parameters.
Keywords
axial_stiffness::Vector{Float64}
: Predefined axial stiffness [N] for each tether.axial_damping::Vector{Float64}
: Predefined axial damping [Ns] for each tether.
Returns
SystemStructure
: A newSystemStructure
representing the 4-tether test system.
SymbolicAWEModels.create_simple_ram_sys_struct
— Functioncreate_simple_ram_sys_struct(set::Settings; axial_stiffness, axial_damping)
Create a simplified SystemStructure
for a ram-air kite with direct tether connections.
This model represents a kite with a flexible wing (2 deformable groups) but simplifies the bridle by connecting the four main tethers directly to the wing attachment points, omitting the complex pulley system. Each tether is modeled as a single segment.
Arguments
set::Settings
: Configuration parameters.
Keywords
axial_stiffness::Vector{Float64}
: Predefined axial stiffness [N] for each tether.axial_damping::Vector{Float64}
: Predefined axial damping [Ns] for each tether.
Returns
SystemStructure
: A newSystemStructure
representing the simplified model.
State Accessor Functions (Getters)
Use these functions to retrieve state information and calculated values from a model instance.
SymbolicAWEModels.winch_force
— Functionwinch_force(s::SymbolicAWEModel)
Returns the winch force [N] for each winch.
SymbolicAWEModels.unstretched_length
— Functionunstretched_length(s::SymbolicAWEModel)
Returns the unstretched tether length [m] for each winch.
SymbolicAWEModels.tether_length
— Functiontether_length(s::SymbolicAWEModel)
Returns the current stretched tether length [m] for each tether.
Visualization Functions
SymbolicAWEModels provides plotting functionality through a package extension that automatically loads when you import GLMakie. The extension provides the following functions:
3D System Visualization
Plot the 3D structure of the kite system with interactive features:
using GLMakie
plot(sys::SystemStructure; kwargs...)
Keyword Arguments:
size::Tuple=(1200, 800)
: Figure size in pixelsmargin::Float64=10.0
: Margin around the system in world unitssegment_color=:black
: Default color for segmentshighlight_color=:red
: Color for highlighted segmentsshow_points::Bool=true
: Show point markersshow_segments::Bool=true
: Show tether segmentsshow_orient::Bool=true
: Show wing orientation axes
Interactive Features:
- Hover over segments to highlight them
- Click on a segment to zoom in
- Click in empty space to zoom out
- Rotate, pan, and zoom with mouse
Time-Series Visualization
Plot simulation results as multi-panel time-series:
plot(sys::SystemStructure, log::SysLog; kwargs...)
Keyword Arguments:
plot_default::Bool=true
: Enable default plot panelsplot_reelout::Bool=plot_default
: Show reel-out velocitiesplot_aero_force::Bool=plot_default
: Show aerodynamic forcesplot_twist::Bool=plot_default
: Show wing twist anglesplot_aoa::Bool=plot_default
: Show angle of attackplot_heading::Bool=plot_default
: Show heading angleplot_winch_force::Bool=plot_default
: Show winch forcesplot_aero_moment::Bool=false
: Show aerodynamic momentsplot_turn_rates::Bool=false
: Show angular velocitiesplot_elevation::Bool=false
: Show elevation angleplot_azimuth::Bool=false
: Show azimuth angleplot_tether_moment::Bool=false
: Show tether-induced momentsplot_set_values::Bool=false
: Show set torque valuessuffix::String=" - " * sys.name
: Suffix for plot labelssize::Tuple=(1200, 800)
: Figure size in pixels
Example:
using SymbolicAWEModels, GLMakie
set = Settings("system.yaml")
sam = SymbolicAWEModel(set, "ram")
init!(sam)
(log, _) = sim_oscillate!(sam)
# Plot only heading and angle of attack
plot(sam.sys_struct, log;
plot_default=false,
plot_heading=true,
plot_aoa=true)
You don't need to explicitly load the plotting extension. Simply using GLMakie
after loading SymbolicAWEModels will automatically make the plot
functions available.
Utility and Helper Functions
General helper functions for package management and setup.
SymbolicAWEModels.init_module
— Functioninit_module(; force=false, add_pkg=true)
Initialize the module in the current working directory.
This function performs the following actions:
- Copies all files from the module's
data
directory to the current working directory'sdata
folder (pwd()/data
). Existing files in the destination are NOT overwritten unlessforce=true
. - Copies all example scripts from the module to the current working directory's
examples
folder (pwd()/examples
). The folder is created if it does not exist. Existing files are NOT overwritten unlessforce=true
. - Installs all required packages if they are not already installed. This occurs only if
add_pkg=true
(default). The packages are automatically determined fromexamples/Project.toml
.
Keyword Arguments
force::Bool=false
: Iftrue
, existing files in the destination directories will be overwritten. Iffalse
(default), existing files will be preserved.add_pkg::Bool=true
: Iftrue
(default), installs required packages if they are not already present. Iffalse
, package installation is skipped.