SciML / ModelingToolkit.jl

An acausal modeling framework for automatically parallelized scientific machine learning (SciML) in Julia. A computer algebra system for integrated symbolics for physics-informed machine learning and automated transformations of differential equations
https://mtk.sciml.ai/dev/
Other
1.41k stars 203 forks source link

MTK allows you to create faulty events without throwing an error. #2612

Open TorkelE opened 5 months ago

TorkelE commented 5 months ago

I have gone through the events in MTK. There are several cases where MTK allow your to build a system with a misformated event, and an error is not actually thrown until you create a problem or attempts to simulate it. Especially since it is a fit trick how to format the events, the errors should be thrown at the first possible time (when a system is created).

using ModelingToolkit, OrdinaryDiffEq

@parameters p d
@variables t X(t)
eqs = [Differential(t)(X) ~ p - d*X]

# Example 1: Discrete events with dosage given as a Tuple
discrete_events = (2.0, 5.0) => [X ~ X + 10.0]
@mtkbuild osys = ODESystem(eqs, t; discrete_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.

# Example 2: Discrete events where single expression condition is given in a vector.
discrete_events = [X < 5.0] => [X ~ X + 10.0]
@mtkbuild osys = ODESystem(eqs, t; discrete_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5])
sol = solve(oprob, Tsit5()) # Error is thrown here.

# Example 3: Discrete events where multiple expression conditions are given in a vector.
discrete_events = [X < 5.0, X > 10.0] => [X ~ X + 10.0]
@mtkbuild osys = ODESystem(eqs, t; discrete_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.

# Example 4: Discrete events that affects parameters and variables.
discrete_events = 5.0 => [X ~ X + 10.0, p ~ 2.0]
@mtkbuild osys = ODESystem(eqs, t; discrete_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.

# Example 5: Discrete events with non-trivial affect right-hand side.
discrete_events = 5.0 => [X + 1.0 ~ X + 10.0]
@mtkbuild osys = ODESystem(eqs, t; discrete_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.

# Example 6: Continuous events that affects parameters and variables.
continuous_events = [X ~ 5.0] => [X ~ X + 10.0, p ~ 2.0]
@mtkbuild osys = ODESystem(eqs, t; continuous_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.

# Example 7: Continuous events with non-trivial affect right-hand side.
continuous_events = [X ~ 5.0] => [X + 1.0 ~ X + 10.0]
@mtkbuild osys = ODESystem(eqs, t; continuous_events)
oprob = ODEProblem(osys, [X => 10], (0.0, 10.), [p => 2.0, d => 0.5]) # Error is thrown here.
baggepinnen commented 5 months ago

If we used the shift-index notation to describe the affect, we could actually handle cases like

discrete_events = 5.0 => [X + 1.0 ~ X + 10.0]

since it would be well-defined which X to solve for, and easy to solve

discrete_events = 5.0 => [X(k) + 1.0 ~ X(k-1) + 10.0]

for X(k).

ChrisRackauckas commented 5 months ago

That's a pretty interesting proposal. We'd have to make validate that only a single shift is given but that would make that implicitness become explicit.