Model building reference

Mark a system constructor function as building a component. For example,

@component function AddOne(; name)
    @variables in(t) out(t)
    eqs = [out ~ in + 1]
    return System(eqs, t, [in, out], []; name)
end

ModelingToolkit systems are either components or connectors. Components define dynamics of the model. Connectors are used to connect components together. See the Model building reference section of the documentation for more information.

See also: @connector.

struct LocalScope <: SymScope

The default scope of a variable. It belongs to the system whose equations it is involved in and is namespaced by every level of the hierarchy.

struct ParentScope <: SymScope

Denotes that the variable does not belong to the system whose equations it is involved in. It is not namespaced by this system. In the immediate parent of this system, the scope of this variable is given by parent.

Fields

• parent::SymScope

struct GlobalScope <: SymScope

Denotes that a variable belongs to the root system in the hierarchy, regardless of which equations of subsystems in the hierarchy it is involved in. Variables with this scope are never namespaced and only added to the unknowns/parameters of a system when calling complete or mtkcompile.

apply_to_variables(f, ex) -> Any

Apply function f to each variable in expression ex. f should be a function that takes a variable and returns the replacement to use. A "variable" in this context refers to a symbolic quantity created directly from a variable creation macro such as Symbolics.@variables, @independent_variables, @parameters, @constants or @brownians.

@named y = foo(x)
@named y[1:10] = foo(x)
@named begin
    y[1:10] = foo(x)
    z = foo(x)
end # returns `[y; z]`
@named y 1:10 i -> foo(x*i)  # This is not recommended

Pass the LHS name to the model. When it's calling anything that's not an AbstractSystem, it wraps all keyword arguments in default_to_parentscope so that namespacing works intuitively when passing a symbolic default into a component.

Examples:

julia> using ModelingToolkit

julia> foo(i; name) = (; i, name)
foo (generic function with 1 method)

julia> x = 41
41

julia> @named y = foo(x)
(i = 41, name = :y)

julia> @named y[1:3] = foo(x)
3-element Vector{NamedTuple{(:i, :name), Tuple{Int64, Symbol}}}:
 (i = 41, name = :y_1)
 (i = 41, name = :y_2)
 (i = 41, name = :y_3)

hierarchy(sys::AbstractSystem; describe = false, bold = describe, kwargs...)

Print a tree of a system's hierarchy of subsystems.

Keyword arguments

• describe: Whether to also print the description of each subsystem, if present.
• bold: Whether to print the name of the system in bold font.

connect(
    sys1::ModelingToolkit.AbstractSystem,
    sys2::ModelingToolkit.AbstractSystem,
    syss::ModelingToolkit.AbstractSystem...
) -> Equation

Connect multiple connectors created via @connector. All connected connectors must be unique.

connect(
    var1::Union{Num, SymbolicUtils.BasicSymbolic, Symbolics.Arr},
    var2::Union{Num, SymbolicUtils.BasicSymbolic, Symbolics.Arr},
    vars::Union{Num, SymbolicUtils.BasicSymbolic, Symbolics.Arr}...
) -> Equation

Connect multiple causal variables. The first variable must be an output, and all subsequent variables must be inputs. The statement connect(var1, var2, var3, ...) expands to:

var1 ~ var2
var1 ~ var3
# ...

connect(output_connector, ap_name::Symbol, input_connector; verbose = true)
connect(output_connector, ap::AnalysisPoint, input_connector; verbose = true)

Connect output_connector and input_connector with an AnalysisPoint inbetween. The incoming connection output_connector is expected to be an output connector (for example, ModelingToolkitStandardLibrary.Blocks.RealOutput), and vice versa.

PLEASE NOTE: The connection is assumed to be causal, meaning that

@named P = FirstOrder(k = 1, T = 1)
@named C = Gain(; k = -1)
connect(C.output, :plant_input, P.input)

is correct, whereas

connect(P.input, :plant_input, C.output)

typically is not (unless the model is an inverse model).

Arguments

• output_connector: An output connector
• input_connector: An input connector
• ap: An explicitly created AnalysisPoint
• ap_name: If a name is given, an AnalysisPoint with the given name will be created automatically.

Keyword arguments

• verbose: Warn if an input is connected to an output (reverse causality). Silence this warning if you are analyzing an inverse model.

domain_connect(
    sys1::ModelingToolkit.AbstractSystem,
    sys2::ModelingToolkit.AbstractSystem,
    syss::ModelingToolkit.AbstractSystem...
) -> Equation

Adds a domain only connection equation, through and across state equations are not generated.

Mark a system constructor function as building a connector. For example,

@connector function ElectricalPin(; name, v = nothing, i = nothing)
    @variables begin
        v(t) = v, [description = "Potential at the pin [V]"]
        i(t) = i, [connect = Flow, description = "Current flowing into the pin [A]"]
    end
    return System(Equation[], t, [v, i], []; name)
end

Since connectors only declare variables, the equivalent shorthand syntax can also be used:

@connector Pin begin
    v(t), [description = "Potential at the pin [V]"]
    i(t), [connect = Flow, description = "Current flowing into the pin [A]"]
end

ModelingToolkit systems are either components or connectors. Components define dynamics of the model. Connectors are used to connect components together. See the Model building reference section of the documentation for more information.

See also: @component.

expand_connections(
    sys::ModelingToolkit.AbstractSystem;
    tol
) -> Any

Given a hierarchical system with connect equations, expand the connection equations and return the new system. tol is the tolerance for handling the singularities in stream connection equations that happen when a flow variable approaches zero.

struct Equality <: ModelingToolkit.AbstractConnectType

Flag which is meant to be passed to the connect metadata of a variable to affect how it behaves when the connector it is in is part of a connect equation. Equality is the default value and such variables when connected are made equal. For example, electric potential is equated at a junction.

For more information, refer to the Connection semantics section of the docs.

See also: connect, @connector, Flow, Stream.

struct Flow <: ModelingToolkit.AbstractConnectType

Flag which is meant to be passed to the connect metadata of a variable to affect how it behaves when the connector it is in is part of a connect equation. Flow denotes that the sum of marked variable in all connectors in the connection set must sum to zero. For example, electric current sums to zero at a junction (assuming appropriate signs are used for current flowing in and out of the function).

For more information, refer to the Connection semantics section of the docs.

See also: connect, @connector, Equality, Stream.

struct Stream <: ModelingToolkit.AbstractConnectType

Flag which is meant to be passed to the connect metadata of a variable to affect how it behaves when the connector it is in is part of a connect equation. Stream denotes that the variable is part of a special stream connector.

For more information, refer to the Connection semantics section of the docs.

See also: connect, @connector, Equality, Flow.

instream(a) -> SymbolicUtils.BasicSymbolic

instream is used when modeling stream connections. It is only allowed to be used on Stream variables.

Refer to the Connection semantics section of the docs for more information.

extend(
    sys::ModelingToolkit.AbstractSystem,
    basesys::ModelingToolkit.AbstractSystem;
    name,
    description,
    gui_metadata
) -> Any

Extend basesys with sys. This can be thought of as the merge operation on systems. Values in sys take priority over duplicates in basesys (for example, defaults).

By default, the resulting system inherits sys's name and description.

The & operator can also be used for this purpose. sys & basesys is equivalent to extend(sys, basesys).

See also compose.

extend(
    sys,
    basesys::Array{T<:ModelingToolkit.AbstractSystem, 1}
) -> Any

Extend sys with all systems in basesys in order.

compose(sys, systems; name)

Compose multiple systems together. This adds all of systems as subsystems of sys. The resulting system inherits the name of sys by default.

The ∘ operator can also be used for this purpose. sys ∘ basesys is equivalent to compose(sys, basesys).

See also extend.

compose(syss...; name) -> Any

Syntactic sugar for adding all systems in syss as the subsystems of first(syss).

substitute_component(
    sys::ModelingToolkit.AbstractSystem,
    rule::Pair{T<:ModelingToolkit.AbstractSystem, T<:ModelingToolkit.AbstractSystem}
) -> Any

Given a hierarchical system sys and a rule lhs => rhs, replace the subsystem lhs in sys by rhs. The lhs must be the namespaced version of a subsystem of sys (e.g. obtained via sys.inner.component). The rhs must be valid as per the following conditions:

1. rhs must not be namespaced.
2. The name of rhs must be the same as the unnamespaced name of lhs.
3. Neither one of lhs or rhs can be marked as complete.
4. Both lhs and rhs must share the same independent variable.
5. rhs must contain at least all of the unknowns and parameters present in lhs.
6. Corresponding unknowns in rhs must share the same connection and causality (input/output) metadata as their counterparts in lhs.
7. For each subsystem of lhs, there must be an identically named subsystem of rhs. These two corresponding subsystems must satisfy conditions 3, 4, 5, 6, 7. If the subsystem of lhs is a connector, the corresponding subsystem of rhs must also be a connector of the same type.

sys also cannot be marked as complete.

flatten(sys::System) -> System
flatten(sys::System, noeqs) -> System

Flatten the hierarchical structure of a system, collecting all equations, unknowns, etc. into one top-level system after namespacing appropriately.

mtkcompile(
    sys;
    additional_passes,
    simplify,
    split,
    allow_symbolic,
    allow_parameter,
    conservative,
    fully_determined,
    inputs,
    outputs,
    disturbance_inputs,
    kwargs...
)

Compile the given system into a form that ModelingToolkit can generate code for. Also performs a variety of symbolic-numeric enhancements. For ODEs, this includes processes such as order reduction, index reduction, alias elimination and tearing. A subset of the unknowns of the system may be eliminated as observables, eliminating the need for the numerical solver to solve for these variables.

Does not rely on metadata to identify variables/parameters/brownians. Instead, queries the system for which symbolic quantites belong to which category. Any variables not present in the equations of the system will be removed in this process.

Keyword Arguments

• When simplify=true, the simplify function will be applied during the tearing process.
• allow_symbolic=false, allow_parameter=true, and conservative=false limit the coefficient types during tearing. In particular, conservative=true limits tearing to only solve for trivial linear systems where the coefficient has the absolute value of $1$.
• fully_determined=true controls whether or not an error will be thrown if the number of equations don't match the number of inputs, outputs, and equations.
• inputs, outputs and disturbance_inputs are passed as keyword arguments.All inputs get converted to parameters and are allowed to be unconnected, allowing models where n_unknowns = n_equations - n_inputs.
• sort_eqs=true controls whether equations are sorted lexicographically before simplification or not.

Macro shorthand for building and compiling a system in one step.

@mtkcompile sys = Constructor(args...; kwargs....)

Is shorthand for

@named sys = Constructor(args...; kwargs...)
sys = mtkcompile(sys)

complete(
    sys::ModelingToolkit.AbstractSystem;
    split,
    flatten,
    add_initial_parameters
) -> Any

Mark a system as completed. A completed system is a system which is done being defined/modified and is ready for structural analysis or other transformations. This allows for analyses and optimizations to be performed which require knowing the global structure of the system.

One property to note is that if a system is complete, the system will no longer namespace its subsystems or variables, i.e. isequal(complete(sys).v.i, v.i).

This namespacing functionality can also be toggled independently of complete using toggle_namespacing.

substitute_observed(
    sys::ModelingToolkit.AbstractSystem,
    expr;
    simplify
) -> Any

Recursively substitute the observed equations of sys into expr. If simplify, call Symbolics.simplify on the resultant expression.

empty_substitutions(sys) -> Any

Check if any variables were eliminated from the system as part of mtkcompile.

get_substitutions(sys) -> Dict

Get a dictionary mapping variables eliminated from the system during mtkcompile to the expressions used to calculate them.

If lifting converts (nested) if statements into a series of continuous events + a logically equivalent if statement + parameters.

Lifting proceeds through the following process:

• rewrite comparisons to be of the form eqn [op] 0; subtract the RHS from the LHS
• replace comparisons with generated parameters; for each comparison eqn [op] 0, generate an event (dependent on op) that sets the parameter

SymbolicContinuousCallback(eqs::Vector{Equation}, affect = nothing, iv = nothing; 
                           affect_neg = affect, initialize = nothing, finalize = nothing, rootfind = SciMLBase.LeftRootFind, alg_eqs = Equation[])

A ContinuousCallback specified symbolically. Takes a vector of equations eq as well as the positive-edge affect and negative-edge affect_neg that apply when any of eq are satisfied. By default affect_neg = affect; to only get rising edges specify affect_neg = nothing.

Assume without loss of generality that the equation is of the form c(u,p,t) ~ 0; we denote the integrator state as i.u. For compactness, we define prev_sign = sign(c(u[t-1], p[t-1], t-1)) and cur_sign = sign(c(u[t], p[t], t)). A condition edge will be detected and the callback will be invoked iff prev_sign * cur_sign <= 0. The positive edge affect will be triggered iff an edge is detected and if prev_sign < 0; similarly, affect_neg will be triggered iff an edge is detected and prev_sign > 0.

Inter-sample condition activation is not guaranteed; for example if we use the dirac delta function as c to insert a sharp discontinuity between integrator steps (which in this example would not normally be identified by adaptivity) then the condition is not guaranteed to be triggered.

Once detected the integrator will "wind back" through a root-finding process to identify the point when the condition became active; the method used is specified by rootfind from SciMLBase.RootfindOpt. If we denote the time when the condition becomes active as tc, the value in the integrator after windback will be:

• u[tc-epsilon], p[tc-epsilon], tc if LeftRootFind is used,
• u[tc+epsilon], p[tc+epsilon], tc if RightRootFind is used,
• or u[t], p[t], t if NoRootFind is used.

For example, if we want to detect when an unknown variable x satisfies x > 0 using the condition x ~ 0 on a positive edge (that is, D(x) > 0), then left root finding will get us x=-epsilon, right root finding x=epsilon and no root finding will produce whatever the next step of the integrator was after it passed through 0.

Multiple callbacks in the same system with different rootfind operations will be grouped by their rootfind value into separate VectorContinuousCallbacks in the enumeration order of SciMLBase.RootfindOpt. This may cause some callbacks to not fire if several become active at the same instant. See the SciMLBase documentation for more information on the semantic rules.

Affects (i.e. affect and affect_neg) can be specified as either:

• A list of equations that should be applied when the callback is triggered (e.g. x ~ 3, y ~ 7) which must be of the form unknown ~ observed value where each unknown appears only once. Equations will be applied in the order that they appear in the vector; parameters and state updates will become immediately visible to following equations.
• A tuple (f!, unknowns, read_parameters, modified_parameters, ctx), where:
  • f! is a function with signature (integ, u, p, ctx) that is called with the integrator, a state index vector u derived from unknowns, a parameter index vector p derived from read_parameters, and the ctx that was given at construction time. Note that ctx is aliased between instances.
  • unknowns is a vector of symbolic unknown variables and optionally their aliases (e.g. if the model was defined with @variables x(t) then a valid value for unknowns would be [x]). A variable can be aliased with a pair x => :y. The indices of these unknowns will be passed to f! in u in a named tuple; in the earlier example, if we pass [x] as unknowns then f! can access x as integ.u[u.x]. If no alias is specified the name of the index will be the symbol version of the variable name.
  • read_parameters is a vector of the parameters that are used by f!. Their indices are passed to f in p similarly to the indices of unknowns passed in u.
  • modified_parameters is a vector of the parameters that are modified by f!. Note that a parameter will not appear in p if it only appears in modified_parameters; it must appear in both parameters and modified_parameters if it is used in the affect definition.
  • ctx is a user-defined context object passed to f! when invoked. This value is aliased for each problem.
• A ImperativeAffect; refer to its documentation for details.

reinitializealg is used to set how the system will be reinitialized after the callback.

• Symbolic affects have reinitialization built in. In this case the algorithm will default to SciMLBase.NoInit(), and should not be provided.
• Functional and imperative affects will default to SciMLBase.CheckInit(), which will error if the system is not properly reinitialized after the callback. If your system is a DAE, pass in an algorithm like SciMLBase.BrownBasicFullInit() to properly re-initialize.

Initial and final affects can also be specified identically to positive and negative edge affects. Initialization affects will run as soon as the solver starts, while finalization affects will be executed after termination.

SymbolicDiscreteCallback(conditions::Vector{Equation}, affect = nothing, iv = nothing;
                         initialize = nothing, finalize = nothing, alg_eqs = Equation[])

A callback that triggers at the first timestep that the conditions are satisfied.

The condition can be one of:

• Δt::Real - periodic events with period Δt
• ts::Vector{Real} - events trigger at these preset times given by ts
• eqs::Vector{Symbolic} - events trigger when the condition evaluates to true

Arguments:

• iv: The independent variable of the system. This must be specified if the independent variable appears in one of the equations explicitly, as in x ~ t + 1.
• alg_eqs: Algebraic equations of the system that must be satisfied after the callback occurs.

ImperativeAffect(f::Function; modified::NamedTuple, observed::NamedTuple, ctx)

ImperativeAffect is a helper for writing affect functions that will compute observed values and ensure that modified values are correctly written back into the system. The affect function f needs to have the signature

    f(modified::NamedTuple, observed::NamedTuple, ctx, integrator)::NamedTuple

The function f will be called with observed and modifiedNamedTuples that are derived from their respective NamedTuple definitions. Each declarationNamedTuple should map an expression to a symbol; for example if we pass observed=(; x = a + b) this will alias the result of executing a+b in the system as x so the value of a + b will be accessible as observed.x in f. modified currently restricts symbolic expressions to only bare variables, so only tuples of the form (; x = y) or (; x) (which aliases x as itself) are allowed.

The argument NamedTuples (for instance (;x=y)) will be populated with the declared values on function entry; if we require (;x=y) in observed and y=2, for example, then the NamedTuple (;x=2) will be passed as observed to the affect function f.

The NamedTuple returned from f includes the values to be written back to the system after f returns. For example, if we want to update the value of x to be the result of x + y we could write

ImperativeAffect(observed=(; x_plus_y = x + y), modified=(; x)) do m, o
    @set! m.x = o.x_plus_y
end

Where we use Setfield to copy the tuple m with a new value for x, then return the modified value of m. All values updated by the tuple must have names originally declared in modified; a runtime error will be produced if a value is written that does not appear in modified. The user can dynamically decide not to write a value back by not including it in the returned tuple, in which case the associated field will not be updated.

modelingtoolkitize(
    prob::ODEProblem;
    u_names,
    p_names,
    return_symbolic_u0_p,
    kwargs...
) -> Any

Convert an ODEProblem to a ModelingToolkit.System.

Keyword arguments

• u_names: An array of names of the same size as prob.u0 to use as the names of the unknowns of the system. The names should be given as Symbols.
• p_names: A collection of names to use for parameters of the system. The collection should have keys corresponding to indexes of prob.p. For example, if prob.p is an associative container like NamedTuple, then p_names should map keys of prob.p to the name that the corresponding parameter should have in the returned system. The names should be given as Symbols.
• INTERNAL return_symbolic_u0_p: Also return the symbolic state and parameter objects.

All other keyword arguments are forwarded to the created System.

modelingtoolkitize(
    prob::SDEProblem;
    u_names,
    p_names,
    kwargs...
) -> Any

Convert an SDEProblem to a ModelingToolkit.System.

Keyword arguments

• u_names: an array of names of the same size as prob.u0 to use as the names of the unknowns of the system. The names should be given as Symbols.
• p_names: a collection of names to use for parameters of the system. The collection should have keys corresponding to indexes of prob.p. For example, if prob.p is an associative container like NamedTuple, then p_names should map keys of prob.p to the name that the corresponding parameter should have in the returned system. The names should be given as Symbols.

All other keyword arguments are forwarded to the created System.

modelingtoolkitize(
    prob::OptimizationProblem;
    u_names,
    p_names,
    kwargs...
) -> System

Convert an OptimizationProblem to a ModelingToolkit.System.

Keyword arguments

• u_names: An array of names of the same size as prob.u0 to use as the names of the unknowns of the system. The names should be given as Symbols.
• p_names: A collection of names to use for parameters of the system. The collection should have keys corresponding to indexes of prob.p. For example, if prob.p is an associative container like NamedTuple, then p_names should map keys of prob.p to the name that the corresponding parameter should have in the returned system. The names should be given as Symbols.

All other keyword arguments are forwarded to the created System.

modelingtoolkitize(
    prob::Union{SciMLBase.NonlinearLeastSquaresProblem, NonlinearProblem};
    u_names,
    p_names,
    kwargs...
) -> System

Convert a NonlinearProblem or NonlinearLeastSquaresProblem to a ModelingToolkit.System.

Keyword arguments

• u_names: An array of names of the same size as prob.u0 to use as the names of the unknowns of the system. The names should be given as Symbols.
• p_names: A collection of names to use for parameters of the system. The collection should have keys corresponding to indexes of prob.p. For example, if prob.p is an associative container like NamedTuple, then p_names should map keys of prob.p to the name that the corresponding parameter should have in the returned system. The names should be given as Symbols.

All other keyword arguments are forwarded to the created System.

FMIComponent(
    ::Val{Ver};
    fmu,
    tolerance,
    communication_step_size,
    reinitializealg,
    type,
    name
)

A component that wraps an FMU loaded via FMI.jl. The FMI version (2 or 3) should be provided as a Val to the function. Supports Model Exchange and CoSimulation FMUs. All inputs, continuous variables and outputs must be FMI.fmi2Real or FMI.fmi3Float64. Does not support events or discrete variables in the FMU. Does not support automatic differentiation. Parameters of the FMU will have defaults corresponding to their initial values in the FMU specification. All other variables will not have a default. Hierarchical names in the FMU of the form namespace.variable are transformed into symbolic variables with the name namespace__variable.

Keyword Arguments

• fmu: The FMU loaded via FMI.loadFMU.
• tolerance: The tolerance to provide to the FMU. Not used for v3 FMUs since it is not supported by FMI.jl.
• communication_step_size: The periodic interval at which communication with CoSimulation FMUs will occur. Must be provided for CoSimulation FMU components.
• reinitializealg: The DAE initialization algorithm to use for the callback managing the FMU. For CoSimulation FMUs whose states/outputs are used in algebraic equations of the system, this needs to be an algorithm that will solve for the new algebraic variables. For example, OrdinaryDiffEqCore.BrownFullBasicInit().
• type: Either :ME or :CS depending on whether fmu is a Model Exchange or CoSimulation FMU respectively.
• name: The name of the system.

liouville_transform(sys::System; kwargs...) -> System

Generates the Liouville transformed set of ODEs, which is the original ODE system with a new variable trJ appended, corresponding to the -tr(Jacobian). This variable is used for properties like uncertainty propagation from a given initial distribution density.

For example, if $u'=p*u$ and p follows a probability distribution $f(p)$, then the probability density of a future value with a given choice of $p$ is computed by setting the initial trJ = f(p), and the final value of trJ is the probability of $u(t)$.

Example:

using ModelingToolkit, OrdinaryDiffEq

@independent_variables t
@parameters α β γ δ
@variables x(t) y(t)
D = Differential(t)
eqs = [D(x) ~ α*x - β*x*y, D(y) ~ -δ*y + γ*x*y]
@named sys = System(eqs, t)

sys2 = liouville_transform(sys)
sys2 = complete(sys2)
u0 = [x => 1.0, y => 1.0, sys2.trJ => 1.0]
prob = ODEProblem(sys2, u0, tspan, p)
sol = solve(prob, Tsit5())

Where sol[3,:] is the evolution of trJ over time.

Sources:

Probabilistic Robustness Analysis of F-16 Controller Performance: An Optimal Transport Approach

Abhishek Halder, Kooktae Lee, and Raktim Bhattacharya https://abhishekhalder.bitbucket.io/F16ACC2013Final.pdf

change_of_variables(
    sys::System,
    iv,
    forward_subs,
    backward_subs;
    simplify,
    t0,
    isSDE
) -> Any

Generates the set of ODEs after change of variables.

Example:

using ModelingToolkit, OrdinaryDiffEq, Test

# Change of variables: z = log(x)
# (this implies that x = exp(z) is automatically non-negative)

@independent_variables t
@parameters α
@variables x(t)
D = Differential(t)
eqs = [D(x) ~ α*x]

tspan = (0., 1.)
def = [x => 1.0, α => -0.5]

@mtkcompile sys = System(eqs, t;defaults=def)
prob = ODEProblem(sys, [], tspan)
sol = solve(prob, Tsit5())

@variables z(t)
forward_subs  = [log(x) => z]
backward_subs = [x => exp(z)]
new_sys = change_of_variables(sys, t, forward_subs, backward_subs)
@test equations(new_sys)[1] == (D(z) ~ α)

new_prob = ODEProblem(new_sys, [], tspan)
new_sol = solve(new_prob, Tsit5())

@test isapprox(new_sol[x][end], sol[x][end], atol=1e-4)

stochastic_integral_transform(
    sys::System,
    correction_factor
) -> Any

Choose correction_factor=-1//2 (1//2) to convert Ito -> Stratonovich (Stratonovich->Ito).

Girsanov_transform(sys::System, u; θ0) -> Any

Measure transformation method that allows for a reduction in the variance of an estimator Exp(g(X_t)). Input: Original SDE system and symbolic function u(t,x) with scalar output that defines the adjustable parameters d in the Girsanov transformation. Optional: initial condition for θ0. Output: Modified SDE System with additional component θ_t and initial value θ0, as well as the weight θ_t/θ0 as observed equation, such that the estimator Exp(g(X_t)θ_t/θ0) has a smaller variance.

Reference: Kloeden, P. E., Platen, E., & Schurz, H. (2012). Numerical solution of SDE through computer experiments. Springer Science & Business Media.

Example

using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D

@parameters α β
@variables x(t) y(t) z(t)

eqs = [D(x) ~ α*x]
noiseeqs = [β*x]

@named de = System(eqs,t,[x],[α,β]; noise_eqs = noiseeqs)

# define u (user choice)
u = x
θ0 = 0.1
g(x) = x[1]^2
demod = ModelingToolkit.Girsanov_transform(de, u; θ0=0.1)

u0modmap = [
    x => x0
]

parammap = [
    α => 1.5,
    β => 1.0
]

probmod = SDEProblem(complete(demod),u0modmap,(0.0,1.0),parammap)
ensemble_probmod = EnsembleProblem(probmod;
          output_func = (sol,i) -> (g(sol[x,end])*sol[demod.weight,end],false),
          )

simmod = solve(ensemble_probmod,EM(),dt=dt,trajectories=numtraj)

change_independent_variable(
    sys::System, iv, eqs = [];
    add_old_diff = false, simplify = true, fold = false
)

Transform the independent variable (e.g. $t$) of the ODE system sys to a dependent variable iv (e.g. $u(t)$). The transformation is well-defined when the mapping between the new and old independent variables are one-to-one. This is satisfied if one is a strictly increasing function of the other (e.g. $du(t)/dt > 0$ or $du(t)/dt < 0$).

Any extra equations eqs involving the new and old independent variables will be taken into account in the transformation.

Keyword arguments

• add_old_diff: Whether to add a differential equation for the old independent variable in terms of the new one using the inverse function rule $dt/du = 1/(du/dt)$.
• simplify: Whether expanded derivative expressions are simplified. This can give a tidier transformation.
• fold: Whether internal substitutions will evaluate numerical expressions.

Usage before structural simplification

The variable change must take place before structural simplification. In following calls to mtkcompile, consider passing allow_symbolic = true to avoid undesired constraint equations between between dummy variables.

Usage with non-autonomous systems

If sys is non-autonomous (i.e. $t$ appears explicitly in its equations), consider passing an algebraic equation relating the new and old independent variables (e.g. $t = f(u(t))$). Otherwise the transformed system can be underdetermined. If an algebraic relation is not known, consider using add_old_diff instead.

Usage with hierarchical systems

It is recommended that iv is a non-namespaced variable in sys. This means it can belong to the top-level system or be a variable in a subsystem declared with GlobalScope.

Example

Consider a free fall with constant horizontal velocity. Physics naturally describes position as a function of time. By changing the independent variable, it can be reformulated for vertical position as a function of horizontal position:

julia> @variables x(t) y(t);

julia> @named M = System([D(D(y)) ~ -9.81, D(D(x)) ~ 0.0], t);

julia> M = change_independent_variable(M, x);

julia> M = mtkcompile(M; allow_symbolic = true);

julia> unknowns(M)
3-element Vector{SymbolicUtils.BasicSymbolic{Real}}:
 xˍt(x)
 y(x)
 yˍx(x)

add_accumulations(sys::System) -> Any
add_accumulations(sys::System, vars) -> Any

Add accumulation variables for vars. For every unknown x in vars, add D(accumulation_x) ~ x as an equation.

add_accumulations(sys::System, vars::Vector{<:Pair}) -> Any

Add accumulation variables for vars. vars is a vector of pairs in the form of

[cumulative_var1 => x + y, cumulative_var2 => x^2]

Then, cumulative variables cumulative_var1 and cumulative_var2 that computes the cumulative x + y and x^2 would be added to sys.

All accumulation variables have a default of zero.

noise_to_brownians(sys::System; names) -> Any

Given a system with noise in the form of noise equation (get_noise_eqs(sys) !== nothing) return an equivalent system which represents the noise using brownian variables.

Keyword Arguments

• names: The name(s) to use for the brownian variables. If this is a Symbol, variables with the given name and successive numeric _i suffixes will be used. If a Vector, this must have appropriate length for the noise equations of the system. The corresponding number of brownian variables are created with the given names.

convert_system_indepvar(sys::System, t; name) -> Any

Function which takes a system sys and an independent variable t and changes the independent variable of sys to t. This is different from change_independent_variable since this function only does a symbolic substitution of the independent variable. sys must not be a reduced system (observed(sys) must be empty). If sys is time-independent, this can be used to turn it into a time-dependent system.

Keyword arguments

• name: The name of the returned system.

struct Sample <: Symbolics.Operator

Represents a sample operator. A discrete-time signal is created by sampling a continuous-time signal.

Constructors

Sample(clock::Union{TimeDomain, InferredTimeDomain} = InferredDiscrete())Sample(dt::Real)

Sample(x::Num), with a single argument, is shorthand for Sample()(x).

Fields

• clock

Examples

julia> using Symbolics

julia> t = ModelingToolkit.t_nounits

julia> Δ = Sample(0.01)
(::Sample) (generic function with 2 methods)

struct Hold <: Symbolics.Operator

Represents a hold operator. A continuous-time signal is produced by holding a discrete-time signal x with zero-order hold.

cont_x = Hold()(disc_x)

function SampleTime()

SampleTime() can be used in the equations of a hybrid system to represent time sampled at the inferred clock for that equation.

transition(from, to, cond; immediate::Bool = true, reset::Bool = true, synchronize::Bool = false, priority::Int = 1)

Create a transition from state from to state to that is enabled when transitioncondition cond evaluates to true.

Arguments:

• from: The source state of the transition.
• to: The target state of the transition.
• cond: A transition condition that evaluates to a Bool, such as ticksInState() >= 2.
• immediate: If true, the transition will fire at the same tick as it becomes true, if false, the actions of the state are evaluated first, and the transition fires during the next tick.
• reset: If true, the destination state to is reset to its initial condition when the transition fires.
• synchronize: If true, the transition will only fire if all sub-state machines in the source state are in their final (terminal) state. A final state is one that has no outgoing transitions.
• priority: If a state has more than one outgoing transition, all outgoing transitions must have a unique priority. The transitions are evaluated in priority order, i.e., the transition with priority 1 is evaluated first.

activeState(state)

When used in a finite state machine, this operator returns true if the queried state is active and false otherwise.

entry()
entry(state)

When used in a finite-state machine, this operator returns true at the first tick when the state is active, and false otherwise.

When used to query the entry of the enclosing state, the method without arguments is used, when used to query the entry of another state, the state is passed as an argument.

This can be used to perform a unique action when entering a state.

ticksInState()
ticksInState(state)

Get the number of ticks spent in a state in a finite state machine.

When used to query the number of ticks spent in the enclosing state, the method without arguments is used, i.e.,

@mtkmodel FSM begin
    ...
    @equations begin
        var(k+1) ~ ticksInState() >= 2 ? 0.0 : var(k)
    end
end

If used to query the number of ticks in another state, the state is passed as an argument.

This operator can be used in both equations and transition conditions.

See also timeInState and entry

timeInState()
timeInState(state)

Get the time (in seconds) spent in a state in a finite state machine.

When used to query the time spent in the enclosing state, the method without arguments is used, i.e.,

@mtkmodel FSM begin
    ...
    @equations begin
        var(k+1) ~ timeInState() >= 2 ? 0.0 : var(k)
    end
end

If used to query the residence time of another state, the state is passed as an argument.

This operator can be used in both equations and transition conditions.

See also ticksInState and entry