Author: Vilgot Lötberg, vilgotl@kth.se, 0725079097 Upon finding issues or if you have questions, please call me.
NOTE! This is the MatRocket source-code, but does not have to be downloaded for getting started with MatRocket. Proper installation of MatRocket is done via pip, follow the guide provided below.
For feedback, updates and possibility to chat with others using MatRocket, feel free to join the Discord: https://discord.gg/6xzcjur4fX
- MatRocket installation
- Project setup
- Functionality & documentation
- Changelog
Below is provided a MATLAB-script that, if MatRocket is not installed, installs MatRocket, and adds the path:
pip_libraries = {"MatRocket"};
for index = 1:numel(pip_libraries)
library = pip_libraries{index};
[failiure, library_status] = system("pip show "+library);
if failiure
system("pip install "+library);
[failiure, library_status] = system("pip show "+library);
end
library_path = strip(extractBetween(string(library_status), "Location: ", "Requires:"))+"\"+library;
addpath(genpath(library_path));
disp(library +" loaded successfully.")
end
A suggested use case would be to put this in a startup-script in MATLAB.
MatRocket is installed via pip. In powershell, cmd or similar:
pip install MatRocketFor installation of pip, see https://pip.pypa.io/en/stable/installation/
NOTE! MATLAB does not automatically recognize libraries installed via pip, and as such, in any scripts that uses MatRocket need to add the MatRocket path.
pip show MatRocketThis gives the resulting output (or similar):
Name: MatRocket
Version: 0.10.15
Summary: MatRocket is a library for simulating rockets in MATLAB.
Home-page: https://github.com/spiggen/MatRocket
Author: Vilgot Lötberg
Author-email: vilgotl@kth.se
License:
Location: C:\Users\jonas\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.13_qbz5n2kfra8p0\LocalCache\local-packages\Python313\site-packages
Requires: MatlabBlenderIO
Required-by:This path can then be added in MATLAB using addpath:
addpath(genpath("C:\Users\jonas\AppData\Local\Packages\PythonSoftwareFoundation.Python.3.13_qbz5n2kfra8p0\LocalCache\local-packages\Python313\site-packages"))For a template on how to set up a project, a project example is provided which can be cloned and modified:
https://github.com/spiggen/MatRocket-Single-Stage-Rocket-Example-
base_atmosphere_model.mbase_atmosphere_from_dataset_model.mbase_aerodynamics_model.mbase_equations_of_motion_model.mbase_gravity_model.mbase_propulsion_model.mbase_UI_loading_bar.m
run_simulation.msimulation_logger.mquery_historian.mrealtime_ode.mrocket2state_vector.mstate_vector2rocket.mderivative2rocket.mflatten.msystem_equations.mmesh2aerodynamics.mmoment_vector.mforce_vector.mbody.mcellsum.mlimit.mnormalize.morthonormalize.mvector2orthonormal_basis.mrandom_wind_dataset.m
All models work around the concept of a state-update. It takes in the rocket, and updates it in some way. It then feeds that same rocket as output, but with an updated state.
Due to the fact that the models feed into each other and are dependent on parameters set by each other, their order is important.
See input & output parameters below for more information on which parameters are required by which model.
[NOTE!] If the parameter is not updated by a previous model, the parameter used will be the one set in the initiation script.
Suggested order:
rocket.models = {@base_atmosphere_model, ...
@base_equations_of_motion_model, ...
@base_propulsion_model, ...
@base_aerodynamics_model, ...
@base_gravity_model, ...
@base_equations_of_motion_model };rocket = base_atmosphere_model(rocket)Recursive: No
File: Simulations/Models/base_atmosphere_model.m
rocket.position| 3x1 Double
rocket.atmosphere.pressure| Double | Atmospheric pressurerocket.atmosphere.temperature| Double | Atmospheric temperaturerocket.atmosphere.density| Double | Atmospheric densityrocket.atmosphere.speed_of_sound| Double | Atmospheric speed of soundrocket.atmosphere.wind_velocity| 3x1 Double | Atmospheric wind velocity
This model queries MATLAB's atmoscoesa model, included in the MATLAB aerospace toolbox based on the rockets altitude, and thus determines the above properties.
rocket = base_atmosphere_from_dataset_model(rocket)Introduced in MatRocket 0.10.18
Recursive: No
File: Simulations/Models/base_atmosphere_from_dataset_model.m
rocket.position| 3x1 Doublerocket.atmosphere.dataset| string | Either the url to a valid dataset-csv file download page, or the filepath to a dataset-csv file.
rocket.atmosphere.pressure| Double | Atmospheric pressurerocket.atmosphere.temperature| Double | Atmospheric temperaturerocket.atmosphere.density| Double | Atmospheric densityrocket.atmosphere.speed_of_sound| Double | Atmospheric speed of soundrocket.atmosphere.wind_velocity| 3x1 Double | Atmospheric wind velocity
This model uses a dataset from a weather balloon or similar to calculate the atmospheric conditions. For an example of where to get data, see: https://weather.uwyo.edu/upperair/sounding.shtml
For formatting of the data pointed to by the url/csv the rocket.atmosphere.dataset, see: https://weather.uwyo.edu/wsgi/sounding?datetime=2018-12-20%2009:00:00&id=02185&src=UNKNOWN&type=TEXT:CSV
To get random data (for example, for monte-carlo style simulations), see the helper-function random_wind_dataset.
rocket = base_aerodynamics_model(rocket)Recursive: Yes
File: Simulations/Models/base_aerodynamics_model.m
componenthasaerodynamicspropertycomponenthasis_bodypropertycomponent.is_bodyis set totrue
rocket.atmosphere.density| 3x1 Double (top level only).component.attitude| 3x3 Double | Rotation matrix from component to parent.component.wind_velocity_absolute| 3x1 Double | Relative wind in parent basis.component.rotation_rate_absolute| 3x1 Double | Absolute rotation rate in parent basis.component.aerodynamics.pressure_coefficient| 3x1 Double | Coefficient of drag for each surface.component.aerodynamics.moment_of_area| 3x3x5 Double | Area-moment tensor.component.aerodynamics.length_scale| Double | Characteristic length used in tensor formulation
.component.moments.LiftMoment| struct | Moment induced by lift/drag in parent basis.component.forces.LiftForce| struct | Lift/drag force in parent basis
.component.aerodynamics.wind_velocity_absolute| 3x1 Double | Wind vector rotated into component’s own frame.component.lift_force_tensor| 3x3 Double | Component-wise lift tensor, for diagnostics/plotting
This model does not use center-of-pressure. It instead computes aerodynamic moments using a series expansion of area-moment tensors, relative wind, and angular velocity — all in the parent frame. The tensor approach blends higher-order spatial information to resolve lift-induced moments.
See:
Detailed technical explanation
The inputs for the aerodynamics-model can be set manually, but recommendation is to let mesh2aerodynamics.m do it for you. See documentation for mesh2aerodynamics.m.
rocket = equations_of_motion(rocket)Recursive: Yes
File: Simulations/Models/base_equations_of_motion_model.m
componenthasis_bodyproperty set totrue
rocket.rigid_body.moment_of_inertia| 3x3 Double (top level only)rocket.atmosphere.wind_velocity| 3x1 Double (top level only)rocket.velocity| 3x1 Double | Velocity of rocket (top level only)rocket.angular_momentum| 3x1 Double | Angular momentum of rocket (top level only).component.mass| Double | Mass of component.component.attitude| 3x3 Double | Rotation matrix from component to parent.component.position| 3x1 Double | Position of component.component.forces| struct | Forces acting on the component.component.moments| struct | Force-moments acting on the component
.component.rotation_rate| 3x1 Double | Angular velocity in parent basis.component.rotation_rate_absolute| 3x1 Double | Accumulated from parent + local.component.wind_velocity_absolute| 3x1 Double | Relative wind in parent basis.component.acceleration| 3x1 Double | Linear acceleration.component.derivative("position"),.derivative("velocity"),.derivative("attitude"),.derivative("angular_momentum")for integration of ODE by ODE-solver
.component.center_of_mass_absolute| 3x1 Double | Recursive center of mass.component.force_absolute,.moment_absolute| struct | Aggregated forces and moments in hierarchy.component.mass_absolute| Double | Total mass of the rocket
This model walks the rocket hierarchy, computing each body’s motion using Newton-Euler equations. Communication is two-way:
- Parent to child: passes down
rotation_rate_absoluteandwind_velocity_absolute - Child to parent: aggregates
mass,force,moment, andcenter_of_mass
Basis handling is explicit — all component-frame vectors are transformed via .attitude or its transpose.
rocket = base_gravity_model(rocket)Recursive: No
File: Simulations/Models/base_gravity_model.m
rocket.attitude| 3x3 Double | Rotation matrix from component to parentrocket.enviroment.g| Double | gravitational constantrocket.mass_absolute| Double | Total mass of the rocketrocket.center_of_mass_absolute| 3x1 Double | Total center of mass of rocket
rocket.forces.Gravity| struct | Gravitational force vector
rocket = base_propulsion_model(rocket)Recursive: No
File: Simulations/Models/base_propulsion_model.m
rocket.t| Double | Simulation-timerocket.engine.burn_time| Double | Length of engine-burnrocket.engine.thrust_force| Double | Magnitude of thrust-force
rocket.engine.forces.thrust_force| struct | Thrust-forcerocket.engine.active_burn| Boolean | Whether or not the burn is ongoing or not
[historian, job] = run_simulation(rocket, job)File: Simulations/run_simulation.m
rocket| structjob| struct | contains fields with options for how the simulation is to be run, for example the simulation time, ode-solver, etc.
historian| struct | a carbon copy of the rocket, but with an extra dimension for each parameter pertaining to the time-steps taken during simulation. Seesimulation_logger.mjob| struct | same as input, but with additional fields pertaining to how the simulation was executed, for example elapsed simulaiton time
A shorthand/wrapper for the MATLAB ODE-solvers. By default, run_simulation uses MATLAB's ode_45. For realtime ode-solving, see MatRocket's realtime_ode.m.
To change out the ode-solver, assign job.ode_solver with the function handle to the ode-solver of your choise:
job.ode_solver = @ode23t;To change the simulation-time, assign job.t_max:
job.t_max = 300 %[s];% Simulation-logging-mode:
simulation_logger(rocket);
% Query/extraction-mode:
historian = simulation_logger();
File: Simulations/lib/simulation_logger.m
rocket| struct
historian| struct | a carbon copy of the rocket, but with an extra dimension for each parameter pertaining to the time-steps taken during simulation
simulation_logger is meant to run in two different modes depending on if it's in the simulation loop or not. It is the main simulation logger, and keeps track of the simulation state by using persistent variables. That way, when the ode-solver is done and the simulation scope is cleaned up by MATLAB, the simulation_logger scope is still untouched, allowing us to extract the simulation data. It is primarily used by run_simulation.m.
The historian-struct is a copy of the rocket-struct, but with one extra dimension for each double/logical type property. This extra dimension contains said property over every single timestep, and is a log of how it changes over the course of the simulation. For example:
>> rocket
rocket =
struct with fields:
is_body: 1
position: [3×1 double]
attitude: [3×3 double]
mesh: ".\Rocket_source\rocket.stl"
engine: [1×1 struct]
models: {1×10 cell}
derivative: [4×1 containers.Map]
enviroment: [1×1 struct]
atmosphere: [1×1 struct]
rigid_body: [1×1 struct]
angular_momentum: [3×1 double]
rotation_rate: [3×1 double]
velocity: [3×1 double]
aerodynamics: [1×1 struct]
chute: [1×1 struct]
tank: [1×1 struct]
>> historian
historian =
struct with fields:
is_body: [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 … ] (1×6001 double)
position: [3×1×6001 double]
attitude: [3×3×6001 double]
mesh: ".\Rocket_source\rocket.stl"
engine: [1×1 struct]
enviroment: [1×1 struct]
atmosphere: [1×1 struct]
rigid_body: [1×1 struct]
angular_momentum: [3×1×6001 double]
rotation_rate: [3×1×6001 double]
velocity: [3×1×6001 double]
aerodynamics: [1×1 struct]
chute: [1×1 struct]
tank: [1×1 struct]
t: [0 0.0167 0.0333 0.0500 0.0667 0.0833 0.1000 0.1167 0.1333 0.1500 0.1667 0.1833 0.2000 0.2167 … ] (1×6001 double)
wind_velocity_absolute: [3×1×6001 double]
mass: [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 … ] (1×6001 double)
rotation_rate_absolute: [3×1×6001 double]
forces: [1×1 struct]
moments: [1×1 struct]
assigned_body_parameters: [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 … ] (1×6001 double)
center_of_mass_absolute: [3×1×6001 double]
mass_absolute: [68 67.9700 67.9400 67.9100 67.8800 67.8500 67.8200 67.7900 67.7600 67.7300 67.7000 67.6700 … ] (1×6001 double)
force_absolute: [1×1 struct]
moment_absolute: [1×1 struct]
acceleration: [3×1×6001 double]
rotation_rate_component_basis: [3×1×6001 double]
lift_force_tensor: [3×3×6001 double]
(In the above case the historian contains extra parameters, since it also contains parameters assigned during the simulation, whereas the rocket only contains the initial conditions).
[NOTE!] If querying simulation_logger directly, remember to clear simulation_logger before starting a new simulation, otherwise it will still contain the old simulation data and try to append the new simulation data to the old. To do this, run:
clear simulation_logger% Query at single timestep
rocket = query_historian(historian, t)
% Query at multiple timesteps
historian = query_historian(historian, t)File: Simulations/lib/query_historian.m
historian| struct | seesimulation_logger.mt| 1x1 double OR 1xN double | time at which historian is to be queried
rocket| struct | iftis a 1x1 doublehistorian| struct | iftis a 1xN double, seesimulation_logger.m
query_historian is a helper function that, given a time or a range of times t, interpolates historian and returns a rocket/historian at that given time or range of times. This is especially useful for animation, as MATLAB's ode-solvers don't take a constant timestep, but instead vary the timestep. By using query_historian, we can interpolate it at a range of times that is linearly spaced to get a fixed framerate.
solution = realtime_ode(dvdt, t_span, v_init)File: Simulations/lib/realtime_ode.m
dvdt| function_handle | ordinary differential equation to be solved on the form:dvdt = @(t,v) ...t_span| 1x2 double | simulation time-spant_v_init| Nx1 double | initial value of state vector, where N is the number of elements in the state vector
solution| struct | struct with two fields:solution.x| 1xM | timesteps taken during simulation, where M is the number of timestepssolution.y| NxM | state-vector values during simulation, where M is the number of timesteps, and N is the number of elements in the state vector
realtime_ode is written on the same form as the standard MATLAB ode-solvers, see MATLAB's documentation. realtime_ode is a second-order Adam-Bashworth ODE-solver, that clocks the time it takes to run the simulation computations, and uses that as it's actual timstep stepsize. This way, simulation time is kept in sync with realtime. This is especially important for hardware-in-the-loop simulations, as realtime systems can only run in real time, and thus the simulation has to keep in sync with the hardware.
url = random_wind_dataset(station_number, start_date, end_date)File: Simulations/lib/random_wind_dataset.m
station_number| string | station number of the station you want to query wind data fromstart_date| datetime | start-date on the interval you want to find wind-dataend_date| datetime | end-date on the interval you want to find wind-data
url| string | url to a valid wind-dataset
random_wind_dataset looks for a random wind dataset on the web in between the start date and end date provided, and, once it finds one that is valid, returns the download-url as a string. The data is provided by the university of wyoming: https://weather.uwyo.edu/upperair/sounding.shtml
For formatting of the data pointed to by the url, see: https://weather.uwyo.edu/wsgi/sounding?datetime=2018-12-20%2009:00:00&id=02185&src=UNKNOWN&type=TEXT:CSV
component = propulsion_model(component)Recursive: Yes
File: Simulations/Models/lib/mesh2aerodynamics.m
componenthasmeshproperty
.component.mesh| struct | Filepath to the stl-file describing the body's geometry
.component.aerodynamics| struct | property containing all aerodynamics information.component.aerodynamics.moment_of_area| 3x3x5 Double | Area-moment tensor.component.aerodynamics.surface_area| 3x1 Double | Surface-area of each face.component.aerodynamics.length_scale| 3x1 Double | Characteristic length of each axis of the component.component.aerodynamics.pressure_coefficient| 3x1 Double | Coefficient of drag for each surface (pressure).component.aerodynamics.friction_coefficient| 3x1 Double | Coefficient of drag for each surface (skin-friction).component.aerodynamics.is_body| Boolean, true | Tells equations of motion to include forces and moments specified under aerodynamics
See:
Detailed technical explanation
How MatRocket is structured.
See https://github.com/spiggen/MatRocket/blob/master/CHANGELOG.md

