AstrodynamicalModels.jl
Common models within astrodynamics!
Overview
This package extends ModelingToolkit to represent common astrodynamical models. All available models are shown on the Docstrings page. Consult the Models pages for more detail about each model in this package!
Usage
If you're familiar with ModelingToolkit.jl, then you'll be able to use this package! Some AstrodynamicalModels-specific usage instructions are provided here. Please don't be shy about making Discourse posts, or filing issues on GitHub!
Installation & Setup
This package can be installed just like any other registered Julia package.
# To install wherever Julia code runs...
import Pkg
Pkg.add("AstrodynamicalModels") # or ]add AstrodynamicalModels in Julia's REPLTo load the package, simply enter using AstrodynamicalModels.
julia> using AstrodynamicalModels
Retrieving a Model
Each model within this package is implemented with a function – each function returns some AbstractSystem from ModelingToolkit.jl. Typically, this will be an ODESystem. If you're worried about overhead from calling each function every time you need a particular model, don't! Each function is implemented with @memoize, so all results are cached the first time you call a model's function with a particular function signature.
julia> R2BPModel = R2BP() # Restricted Two-body Problem dynamicsERROR: UndefVarError: `R2BP` not definedjulia> CR3BPModel = CR3BP() # Circular Restricted Three-body Problem dynamicsERROR: UndefVarError: `CR3BP` not definedjulia> CR3BPModelWithSTM = CR3BP(; stm=true) # Optionally include state transition matrix dynamicsERROR: UndefVarError: `CR3BP` not defined
Using a Model
To actually use each model, you probably also want to load ModelingToolkit (and any other SciML packages of your choice).
julia> using ModelingToolkit
Now you can use any method defined for ModelingToolkit.AbstractSystem instances. Once again, the ModelingToolkit Documentation are the best place to learn how to interact with AbstractSystem instances! Some quick examples are shown below.
Check the Equations of Motion
julia> eqs = equations(R2BSystem())6-element Vector{Equation}: Differential(t)(x(t)) ~ ẋ(t) Differential(t)(y(t)) ~ ẏ(t) Differential(t)(z(t)) ~ ż(t) Differential(t)(ẋ(t)) ~ (-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) Differential(t)(ẏ(t)) ~ (-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) Differential(t)(ż(t)) ~ (-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3)
List the States and Parameters
julia> x = states(R2BSystem())ERROR: UndefVarError: `states` not definedjulia> p = parameters(R2BSystem())WARNING: both ModelingToolkit and AstrodynamicalModels export "parameters"; uses of it in module Main must be qualified ERROR: UndefVarError: `parameters` not defined
Calculate the Jacobian
julia> J = calculate_jacobian(R2BSystem())6×6 Matrix{Num}: 0 … 1 0 0 0 0 1 0 0 0 0 1 (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3x(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0 -3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*x(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0 -3x(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) … 0 0 0
Generate Code to Replicate the Model
julia> print(build_function(R2BSystem()))function () #= /juliateam/.julia/packages/SymbolicUtils/r1pzW/src/code.jl:373 =# #= /juliateam/.julia/packages/SymbolicUtils/r1pzW/src/code.jl:374 =# #= /juliateam/.julia/packages/SymbolicUtils/r1pzW/src/code.jl:375 =# ODESystem(0x0000000000000012, Equation[Differential(t)(x(t)) ~ ẋ(t), Differential(t)(y(t)) ~ ẏ(t), Differential(t)(z(t)) ~ ż(t), Differential(t)(ẋ(t)) ~ (-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ẏ(t)) ~ (-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ż(t)) ~ (-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3)], t, SymbolicUtils.BasicSymbolic{Real}[x(t), y(t), z(t), ẋ(t), ẏ(t), ż(t)], SymbolicUtils.BasicSymbolic{Real}[μ], nothing, Dict{Any, Any}(:ż => ż(t), :μ => μ, :y => y(t), :ẏ => ẏ(t), :ẋ => ẋ(t), :z => z(t), :x => x(t)), Any[], Equation[], Base.RefValue{Vector{Num}}(Num[]), Base.RefValue{Any}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Any}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), :R2B, ODESystem[], Dict{Any, Any}(), Dict{Any, Any}(), nothing, nothing, Equation[], nothing, nothing, nothing, ModelingToolkit.SymbolicContinuousCallback[ModelingToolkit.SymbolicContinuousCallback(Equation[], Equation[])], ModelingToolkit.SymbolicDiscreteCallback[], nothing, nothing, nothing, nothing, nothing, false, nothing, nothing, nothing, nothing, nothing) end
Generate Code which Implements the Dynamics
julia> print(R2BFunction())ODEFunction{true, SciMLBase.FullSpecialize, ModelingToolkit.var"#f#684"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x85184888, 0xb17a535f, 0x82fddb8f, 0x46fde24e, 0xe0ba6752), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x4feca1d2, 0xcfce7c7e, 0x9e2fea20, 0xdcc4aa6a, 0xca1fa926), Nothing}}, LinearAlgebra.UniformScaling{Bool}, Nothing, Nothing, ModelingToolkit.var"#_jac#689"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x67420626, 0x1a4595aa, 0x37e5bbe7, 0xf10f7fc6, 0xf0e5b87a), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x0923dfa2, 0x8c524e4a, 0x4b8fb075, 0x723709b6, 0xf6c592f1), Nothing}}, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, Nothing, ModelingToolkit.var"#51211#generated_observed#693"{Bool, ODESystem, Dict{Any, Any}, Vector{SymbolicUtils.BasicSymbolic{Real}}}, Nothing, ODESystem, Nothing, Nothing}(ModelingToolkit.var"#f#684"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x85184888, 0xb17a535f, 0x82fddb8f, 0x46fde24e, 0xe0ba6752), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x4feca1d2, 0xcfce7c7e, 0x9e2fea20, 0xdcc4aa6a, 0xca1fa926), Nothing}}(RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x85184888, 0xb17a535f, 0x82fddb8f, 0x46fde24e, 0xe0ba6752), Nothing}(nothing), RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x4feca1d2, 0xcfce7c7e, 0x9e2fea20, 0xdcc4aa6a, 0xca1fa926), Nothing}(nothing)), LinearAlgebra.UniformScaling{Bool}(true), nothing, nothing, ModelingToolkit.var"#_jac#689"{RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x67420626, 0x1a4595aa, 0x37e5bbe7, 0xf10f7fc6, 0xf0e5b87a), Nothing}, RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x0923dfa2, 0x8c524e4a, 0x4b8fb075, 0x723709b6, 0xf6c592f1), Nothing}}(RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x67420626, 0x1a4595aa, 0x37e5bbe7, 0xf10f7fc6, 0xf0e5b87a), Nothing}(nothing), RuntimeGeneratedFunctions.RuntimeGeneratedFunction{(:ˍ₋out, :ˍ₋arg1, :ˍ₋arg2, :t), ModelingToolkit.var"#_RGF_ModTag", ModelingToolkit.var"#_RGF_ModTag", (0x0923dfa2, 0x8c524e4a, 0x4b8fb075, 0x723709b6, 0xf6c592f1), Nothing}(nothing)), nothing, nothing, nothing, nothing, nothing, nothing, nothing, nothing, ModelingToolkit.var"#51211#generated_observed#693"{Bool, ODESystem, Dict{Any, Any}, Vector{SymbolicUtils.BasicSymbolic{Real}}}(Core.Box(nothing), Core.Box(nothing), false, Core.Box(Equation[]), ODESystem(0x000000000000000d, Equation[Differential(t)(x(t)) ~ ẋ(t), Differential(t)(y(t)) ~ ẏ(t), Differential(t)(z(t)) ~ ż(t), Differential(t)(ẋ(t)) ~ (-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ẏ(t)) ~ (-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ż(t)) ~ (-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3)], t, SymbolicUtils.BasicSymbolic{Real}[x(t), y(t), z(t), ẋ(t), ẏ(t), ż(t)], SymbolicUtils.BasicSymbolic{Real}[μ], nothing, Dict{Any, Any}(:ż => ż(t), :μ => μ, :y => y(t), :ẏ => ẏ(t), :ẋ => ẋ(t), :z => z(t), :x => x(t)), Any[], Equation[], Base.RefValue{Vector{Num}}(Num[]), Base.RefValue{Any}((Num[0 0 0 1 0 0; 0 0 0 0 1 0; 0 0 0 0 0 1; (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3x(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3y(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3z(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0; -3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*x(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*y(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*z(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0; -3x(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3y(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3z(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0], (false, false))), Base.RefValue{Any}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), :R2B, ODESystem[], Dict{Any, Any}(), Dict{Any, Any}(), nothing, nothing, Equation[], nothing, nothing, nothing, ModelingToolkit.SymbolicContinuousCallback[ModelingToolkit.SymbolicContinuousCallback(Equation[], Equation[])], ModelingToolkit.SymbolicDiscreteCallback[], nothing, nothing, nothing, nothing, nothing, true, nothing, nothing, nothing, nothing, nothing), Dict{Any, Any}(), SymbolicUtils.BasicSymbolic{Real}[μ]), nothing, ODESystem(0x000000000000000d, Equation[Differential(t)(x(t)) ~ ẋ(t), Differential(t)(y(t)) ~ ẏ(t), Differential(t)(z(t)) ~ ż(t), Differential(t)(ẋ(t)) ~ (-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ẏ(t)) ~ (-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3), Differential(t)(ż(t)) ~ (-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3)], t, SymbolicUtils.BasicSymbolic{Real}[x(t), y(t), z(t), ẋ(t), ẏ(t), ż(t)], SymbolicUtils.BasicSymbolic{Real}[μ], nothing, Dict{Any, Any}(:ż => ż(t), :μ => μ, :y => y(t), :ẏ => ẏ(t), :ẋ => ẋ(t), :z => z(t), :x => x(t)), Any[], Equation[], Base.RefValue{Vector{Num}}(Num[]), Base.RefValue{Any}((Num[0 0 0 1 0 0; 0 0 0 0 1 0; 0 0 0 0 0 1; (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3x(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3y(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3z(t)*((-x(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0; -3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*x(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*y(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3((-y(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*z(t)*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0; -3x(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) -3y(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) (-μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^3) - 3z(t)*((-z(t)*μ) / (sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t)))^6))*sqrt(abs2(y(t)) + abs2(x(t)) + abs2(z(t))) 0 0 0], (false, false))), Base.RefValue{Any}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), Base.RefValue{Matrix{Num}}(Matrix{Num}(undef, 0, 0)), :R2B, ODESystem[], Dict{Any, Any}(), Dict{Any, Any}(), nothing, nothing, Equation[], nothing, nothing, nothing, ModelingToolkit.SymbolicContinuousCallback[ModelingToolkit.SymbolicContinuousCallback(Equation[], Equation[])], ModelingToolkit.SymbolicDiscreteCallback[], nothing, nothing, nothing, nothing, nothing, true, nothing, nothing, nothing, nothing, nothing), nothing, nothing)
Generate C/C++ and MATLAB Code
julia> print(build_function([eq.rhs for eq in equations(R2BSystem())], states(R2BSystem()), parameters(R2BSystem()); target=Symbolics.CTarget()))ERROR: UndefVarError: `states` not definedjulia> print(build_function([eq.rhs for eq in equations(R2BSystem())], states(R2BSystem()), parameters(R2BSystem()); target=Symbolics.MATLABTarget()))ERROR: UndefVarError: `states` not defined
If you're interested in learning a bit about each astrodynamical model, or you'd like more specific examples which show how to use each model, consult the Models pages!