High Level API
The high-level API allows modelers to interactively build models in a symbolic manner. It is designed as a semi-DSL for easily building large complex models and manipulating the models to generate optimal forms to be used in numerical methods.
High-Level API Documentation
ModelingToolkit.@parameters
— MacroDefine one or more known variables.
ModelingToolkit.@variables
— MacroDefine one or more unknown variables.
@parameters t α σ(..) β[1:2]
@variables w(..) x(t) y() z(t, α, x)
expr = β₁* x + y^α + σ(3) * (z - t) - β₂ * w(t - 1)
Note that @parameters
and @variables
implicitly add ()
to values that are not given a call. The former specifies the values as known, while the latter specifies it as unknown. (..)
signifies that the value should be left uncalled.
Sometimes it is convenient to define arrays of variables to model things like x₁,…,x₃
. The @variables
and @parameters
macros support this with the following syntax:
@variables x[1:3];
x
3-element Array{Operation,1}:
x₁()
x₂()
x₃()
# support for arbitrary ranges and tensors
@variables y[2:3,1:5:6];
y
2×2 Array{Operation,2}:
y₂̒₁() y₂̒₆()
y₃̒₁() y₃̒₆()
# also works for dependent variables
@parameters t; @variables z[1:3](t);
z
3-element Array{Operation,1}:
z₁(t())
z₂(t())
z₃(t())
Missing docstring for Differential
. Check Documenter's build log for details.
Base.:~
— Method~(lhs::Num, rhs::Num) -> Equation
Create an Equation
out of two Num
instances, or an Num
and a Number
.
Examples
julia> using ModelingToolkit
julia> @variables x y;
julia> x ~ y
Equation(x(), y())
julia> x - y ~ 0
Equation(x() - y(), 0)
ModelingToolkit.modelingtoolkitize
— Functionmodelingtoolkitize(prob::ODEProblem) -> Union{Tuple{Any,Any,Any}, ODESystem}
Generate ODESystem
, dependent variables, and parameters from an ODEProblem
.
modelingtoolkitize(prob::SDEProblem) -> Union{Tuple{Any,Any,Any}, SDESystem}
Generate SDESystem
, dependent variables, and parameters from an SDEProblem
.
modelingtoolkitize(prob::OptimizationProblem) -> OptimizationSystem
Generate OptimizationSystem
, dependent variables, and parameters from an OptimizationProblem
.
Differentiation Functions
Missing docstring for ModelingToolkit.derivative
. Check Documenter's build log for details.
Missing docstring for ModelingToolkit.gradient
. Check Documenter's build log for details.
Missing docstring for ModelingToolkit.jacobian
. Check Documenter's build log for details.
ModelingToolkit.sparsejacobian
— Functionsparsejacobian(ops::AbstractVector, vars::AbstractVector; simplify=false)
A helper function for computing the sparse Jacobian of an array of expressions with respect to an array of variable expressions.
Missing docstring for ModelingToolkit.hessian
. Check Documenter's build log for details.
ModelingToolkit.sparsehessian
— Functionsparsehessian(O, vars::AbstractVector; simplify=false)
A helper function for computing the sparse Hessian of an expression with respect to an array of variable expressions.
Sparsity Detection
ModelingToolkit.jacobian_sparsity
— Functionjacobian_sparsity(ops::AbstractVector, vars::AbstractVector)
Return the sparsity pattern of the Jacobian of an array of expressions with respect to an array of variable expressions.
ModelingToolkit.hessian_sparsity
— Functionhessian_sparsity(ops::AbstractVector, vars::AbstractVector)
Return the sparsity pattern of the Hessian of an array of expressions with respect to an array of variable expressions.
Latexification
ModelingToolkit.jl's expressions support Latexify.jl, and thus
using Latexify
latexify(ex)
will produce LaTeX output from ModelingToolkit models and expressions. This works on basics like Term
all the way to higher primitives like ODESystem
and ReactionSystem
.
The Auto-Detecting System Constructors
For the high-level interface, the system constructors, such as ODESystem
, have high-level constructors, which just take in the required equations and automatically parse the expressions to figure out the states and parameters of the system. The following high-level constructors exist:
ODESystem(eqs)
NonlinearSystem(eqs)
Direct Tracing
Because ModelingToolkit expressions respect Julia semantics, one way to generate symbolic expressions is to simply place ModelingToolkit variables as inputs into existing Julia code. For example, the following uses the standard Julia function for the Lorenz equations to generate the symbolic expression for the Lorenz equations:
function lorenz(du,u,p,t)
du[1] = 10.0(u[2]-u[1])
du[2] = u[1]*(28.0-u[3]) - u[2]
du[3] = u[1]*u[2] - (8/3)*u[3]
end
@variables t u[1:3](t) du[1:3](t)
@parameters p[1:3]
lorenz(du,u,p,t)
du
3-element Array{Num,1}:
10.0 * (u₂(t) - u₁(t))
u₁(t) * (28.0 - u₃(t)) - u₂(t)
u₁(t) * u₂(t) - 2.6666666666666665 * u₃(t)
Or similarly:
@variables t x(t) y(t) z(t) dx(t) dy(t) dz(t)
@parameters σ ρ β
du = [dx,dy,dz]
u = [x,y,z]
p = [σ,ρ,β]
lorenz(du,u,p,t)
du
3-element Array{Num,1}:
10.0 * (y(t) - x(t))
x(t) * (28.0 - z(t)) - y(t)
x(t) * y(t) - 2.6666666666666665 * z(t)
Intermediate Calculations
The system building functions can handle intermediate calculations by simply defining and using an Term
of Sym
s. For example:
@variables x y z
@parameters σ ρ β
a = y - x
eqs = [0 ~ σ*a,
0 ~ x*(ρ-z)-y,
0 ~ x*y - β*z]
ns = NonlinearSystem(eqs, [x,y,z], [σ,ρ,β])
nlsys_func = generate_function(ns)[2] # second is the inplace version
expands to:
:((var"##MTIIPVar#368", var"##MTKArg#365", var"##MTKArg#366")->begin
@inbounds begin
let (x, y, z, σ, ρ, β) = (var"##MTKArg#365"[1], var"##MTKArg#365"[2], var"##MTKArg#365"[3], var"##MTKArg#366"[1], var"##MTKArg#366"[2], var"##MTKArg#366"[3])
var"##MTIIPVar#368"[1] = (*)(σ, (-)(y, x))
var"##MTIIPVar#368"[2] = (-)((*)(x, (-)(ρ, z)), y)
var"##MTIIPVar#368"[3] = (-)((*)(x, y), (*)(β, z))
end
end
nothing
end)
In addition, the Jacobian calculations take into account intermediate variables to appropriately handle them.
I/O and Saving
Note that Julia's standard I/O functionality can be used to save ModelingToolkit expressions out to files. For example, here we will generate an in-place version of f
and save the anonymous function to a .jl
file:
using ModelingToolkit
@variables u[1:3]
function f(u)
[u[1]-u[3],u[1]^2-u[2],u[3]+u[2]]
end
ex1, ex2 = build_function(f(u),u)
write("function.jl", string(ex2))
Now we can do something like:
f = include("function.jl")
and that will load the function back in. Note that this can be done to save the transformation results of ModelingToolkit.jl so that they can be stored and used in a precompiled Julia package.