Equations and events reupload

[TorkelE]

I have tried to rebase #736 on master. Unfortunately, there have been too many changes, and it is not feasible anymore. This PR builds that one.

Creating events in the DSL

MTK does not currently have a wrapper around the @continuous_events or @discrete_events options. However, we here use the same syntax as for programmatic description, but in a begin ... end block, with one event on each line. If there is only a single event, it can just be given after @X_events.

# When t == 2.5, increase p by 0.2.
rs = @reaction_network begin
    @continuous_events begin
        [t ~ 2.5] => [p ~ p + 0.2]
    end
    (p,d), 0 <--> X
end

# Every 2.0 time unit, increase d by 0.2.
rs = @reaction_network begin
    @discrete_events begin
        2.0 => [d ~ d + 0.2]
    end
    (p,d), 0 <--> X
end

# At time unit 2.0, increase d by 0.2.
rs = @reaction_network begin
    @discrete_events begin
        [2.0] => [d ~ d + 0.2]
    end
    (p,d), 0 <--> X
end

# After every timestep, if X > 0.2, increase d by 0.2.
rs = @reaction_network begin
    @discrete_events begin
        (X > 0.1) => [d ~ d + 0.2]
    end
    (p,d), 0 <--> X
end

Creating equations in the DSL

Equations use the @equations option. This utilises some MTK stuff, but the notation is similar as when creating stuff programmatically. Both differential and algebraic equations are supported.

rs_dsl = @reaction_network hybrid_rs begin
    @variables X2
    @parameters v
    @equations begin
        D(V) ~ X - v*V
        X2  ~ 2*X
    end
(p,d), 0 <--> X
end

Parameters and variables used within equations must be declared as such The exception is if we have a differential equation on the form

D(sym) ~ ...

where sym is a single symbol. Here, sym is inferred to be a variable (it can still be declared in @variables, if you e.g. want to attach metadata). In these cases, D is inferred to be a differential with respect to the independent variable.

It is also possible to create new differentials using the @differentials option:

rs_dsl = @reaction_network hybrid_rs begin
    @differentials begin
        DD = Differential(t)
    end
    @variables V(t)
    @equations begin
        DD(V) ~ -V
    end
    (p,d), 0 <--> X
end

these can then be used in equations. Note that if a non-default differential is used, the V variable have to be decalred.

Working with hybrid equations

These work like before, however, with the following additions:

• Algebraic equations can now be used. These require using structural_simplify on the final XSystem. All relevant problem calls now have a structural_simpligy = false kwarg. Setting this to true ensures that structural_simplify (not compelte) is called on the internally created ODESystem.
• Hybrid systems can now be used to create Nonlinear Systems. When this happens, andy term D(...) (where D is some differential) are automatically set to 0).
• Hybrid systems can be converted into SDEs. Currently, we do not support adding noise equations to variables (which would have to be done at a later point).

There are a couple of things which are broken, awaiting fixes in MTK, but these are rather tangential, and the main functionality works well.

[TorkelE]

Figured I should have another go through and see if I missed anything, and was immediately punished. On the plus side I discovered another MTK bug: https://github.com/SciML/SciMLBase.jl/issues...

[TorkelE]

Ok, I've updated to the new MTK version where this was fixed, all good now.

[isaacsas]

We can leave the reformats I didn't address as is. But please stop rewriting / refactoring perfectly fine existing code as part of PRs that are adding new functionality. It just makes the PR much longer, increasing the review time (it is much easier to find time here and there to review 10 200-line PRs than one 2000 line PR). Modifying existing code in ways that does not change what it is doing or fix a bug should be left for separate PRs.

[isaacsas]

I don't think we should take bad input and fix it, instead we should just error if a user gives something that is defined as an observable as a species or variable. (And maybe this should be a separate function to encapsulate the check.)

[isaacsas]

(Part of the reason for this is that a user might be incorrectly understanding how observables and variables/species work, so it is better to error and let them correct their model I think.)

[TorkelE]

This is needed because if you want to e.g. change whether an observable is a species/variable you do

@species O(t)

in the macro. However, then that observable is considered a species of the system (and have to be removed at the next stage).

If you want we can have an error check in the programmatic ReactionSystem constructor that checks for this and throws an error.

[isaacsas]

Then shouldn't this be part of the DSL? I still don't think that ReactionSystems should be "fixing" bad input but should instead just error.

[isaacsas]

I'd be willing to have this in the two-argument constructor, but I don't think it should be in the inner constructor. The inner constructor implicitly assumes the collection of unknowns and parameters are correct and doesn't do discovery.

[TorkelE]

I will look at revamping where this is done. Generally (in agreement with our previous discussion) I think that stuff like this should be handled outside of the DSL.