The AbstractSystem Interface
Overview
The AbstractSystem
interface is the core of the system level of ModelingToolkit.jl. It establishes a common set of functionality that is used between systems from ODEs and chemical reactions, allowing users to have a common framework for model manipulation and compilation.
Composition and Accessor Functions
Each AbstractSystem
has lists of variables in context, such as distinguishing parameters vs states. In addition, an AbstractSystem
also can hold other AbstractSystem
types. Direct accessing of the values, such as sys.states
, gives the immediate list, while the accessor functions states(sys)
gives the total set, which includes that of all systems held inside.
The values which are common to all AbstractSystem
s are:
equations(sys)
: All equations that define the system and its subsystems.states(sys)
: All the states in the system and its subsystems.parameters(sys)
: All parameters of the system and its subsystems.nameof(sys)
: The name of the current-level system.get_eqs(sys)
: Equations that define the current-level system.get_states(sys)
: States that are in the current-level system.get_ps(sys)
: Parameters that are in the current-level system.get_systems(sys)
: Subsystems of the current-level system.
Optionally, a system could have:
observed(sys)
: All observed equations of the system and its subsystems.get_observed(sys)
: Observed equations of the current-level system.get_default_u0(sys)
: ADict
that maps states into their default initial condition.get_default_p(sys)
: ADict
that maps parameters into their default value.independent_variable(sys)
: The independent variable of a system.get_noiseeqs(sys)
: Noise equations of the current-level system.
Note that there's get_iv(sys)
, but it is not advised to use, since it errors when the system has no field iv
. independent_variable(sys)
returns nothing
for NonlinearSystem
s.
A system could also have caches:
get_jac(sys)
: The Jacobian of a system.get_tgrad(sys)
: The gradient with respect to time of a system.
Transformations
Transformations are functions which send a valid AbstractSystem
definition to another AbstractSystem
. These are passes, like optimizations (e.g., Block-Lower Triangle transformations), or changes to the representation, which allow for alternative numerical methods to be utilized on the model (e.g., DAE index reduction).
Function Calculation and Generation
The calculation and generation functions allow for calculating additional quantities to enhance the numerical methods applied to the resulting system. The calculations, like calculate_jacobian
, generate ModelingToolkit IR for the Jacobian of the system, while the generations, like generate_jacobian
, generate compiled output for the numerical solvers by applying build_function
to the generated code. Additionally, many systems have function-type outputs, which cobble together the generation functionality for a system, for example, ODEFunction
can be used to generate a DifferentialEquations-based ODEFunction
with compiled version of the ODE itself, the Jacobian, the mass matrix, etc.
Below are the possible calculation and generation functions:
ModelingToolkit.calculate_tgrad
— Functioncalculate_tgrad(sys::AbstractSystem)
Calculate the time gradient of a system.
Returns a vector of Num
instances. The result from the first call will be cached in the system object.
ModelingToolkit.calculate_gradient
— Functioncalculate_gradient(sys::AbstractSystem)
Calculate the gradient of a scalar system.
Returns a vector of Num
instances. The result from the first call will be cached in the system object.
ModelingToolkit.calculate_jacobian
— Functioncalculate_jacobian(sys::AbstractSystem)
Calculate the jacobian matrix of a system.
Returns a matrix of Num
instances. The result from the first call will be cached in the system object.
ModelingToolkit.calculate_factorized_W
— Functioncalculate_factorized_W(sys::AbstractSystem)
Calculate the factorized W-matrix of a system.
Returns a matrix of Num
instances. The result from the first call will be cached in the system object.
ModelingToolkit.calculate_hessian
— Functioncalculate_hessian(sys::AbstractSystem)
Calculate the hessian matrix of a scalar system.
Returns a matrix of Num
instances. The result from the first call will be cached in the system object.
ModelingToolkit.generate_tgrad
— Functiongenerate_tgrad(sys::AbstractSystem, dvs = states(sys), ps = parameters(sys), expression = Val{true}; kwargs...)
Generates a function for the time gradient of a system. Extra arguments control the arguments to the internal build_function
call.
ModelingToolkit.generate_gradient
— Functiongenerate_gradient(sys::AbstractSystem, dvs = states(sys), ps = parameters(sys), expression = Val{true}; kwargs...)
Generates a function for the gradient of a system. Extra arguments control the arguments to the internal build_function
call.
ModelingToolkit.generate_jacobian
— Functiongenerate_jacobian(sys::AbstractSystem, dvs = states(sys), ps = parameters(sys), expression = Val{true}; sparse = false, kwargs...)
Generates a function for the jacobian matrix matrix of a system. Extra arguments control the arguments to the internal build_function
call.
ModelingToolkit.generate_factorized_W
— Functiongenerate_factorized_W(sys::AbstractSystem, dvs = states(sys), ps = parameters(sys), expression = Val{true}; sparse = false, kwargs...)
Generates a function for the factorized W-matrix matrix of a system. Extra arguments control the arguments to the internal build_function
call.
ModelingToolkit.generate_hessian
— Functiongenerate_hessian(sys::AbstractSystem, dvs = states(sys), ps = parameters(sys), expression = Val{true}; sparse = false, kwargs...)
Generates a function for the hessian matrix matrix of a system. Extra arguments control the arguments to the internal build_function
call.
Additionally, jacobian_sparsity(sys)
and hessian_sparsity(sys)
exist on the appropriate systems for fast generation of the sparsity patterns via an abstract interpretation without requiring differentiation.
Problem Constructors
At the end, the system types have DEProblem
constructors, like ODEProblem
, which allow for directly generating the problem types required for numerical methods. The first argument is always the AbstractSystem
, and the proceeding arguments match the argument order of their original constructors. Whenever an array would normally be provided, such as u0
the initial condition of an ODEProblem
, it is instead replaced with a variable map, i.e., an array of pairs var=>value
, which allows the user to designate the values without having to know the order that ModelingToolkit is internally using.