Allow input types commonly used in Catlayst/MTK

[TorkelE]

I was looking into writing a short tutorial on how to use FiniteStateProjection.jl in Catalyst, and then link to the package. While I was looking through the documentation I had some thoughts:

Current workflow for simulating a model using FSP:

using FiniteStateProjection
using OrdinaryDiffEq

rn = @reaction_network begin
    σ, 0 --> A
    d, A --> 0
end

sys = FSPSystem(rn)

# Parameters for our system
ps = [ 10.0, 1.0 ]

# Initial distribution (over 1 species)
# Here we start with 0 copies of A
u0 = zeros(50)
u0[1] = 1.0

prob = convert(ODEProblem, sys, u0, (0, 10.0), ps)
sol = solve(prob, Vern7())

Comment 1: Right now, using Catalyst and MTK, parameter values are no longer given as vectors, but rather maps from parameters to values, e.g:

ps = [:σ => 10.0, :d => 1.0]

This has the advantage that one does not need to consider the order with which the parameters are stored within the model (which, after a recent MTK update, is no longer fully reliable).

Comment 2: Technically the same holds for initial conditions. Now, since the initial condition for a FSP simulation is a multidimensional array, this does not directly work. However, it might be possible to enable FSPSystem to take an optional argument which specifies the order with which the species occur? (it is possible to manually check using species(rn), but something like this might make sense)

Comment 3: Mostly we have a dispatch on ODEProblem for special system types, rather than using convert. E.g. something like

ODEProblem(sys, u0, (0, 10.0), ps)

instead of the current

prob = convert(ODEProblem, sys, u0, (0, 10.0), ps)

[kaandocal]

Good points! This should be an easy change using the MTK API, I hope to do this soon.

[TorkelE]

Sounds good! I was thinking of adding an example of how to use FiniteStateProjection.jl in the Catalyst docs, and then linking this repo for people who comes from Catalyst and want more information.

[kaandocal]

I've updated the package to include comments 1 and 3, although some of the tests fail on GitHub because of issues in SimpleNonLinearsolve.jl. I am not sure what the best way to proceed with comment 2 is, but will think about it.

[TorkelE]

Sounds great!

I haven't thought super much about it, I think the simplest might be to simply provide the species order as an additional argument somewhere? E.g. for

rn = @reaction_network begin
    (k1,k2), X + Y <--> XT
end

one would probably provide a vector [:X, :Y, :XY] (or use symbolic form [X, Y, XY] or [rn.X, rn.Y, rn.XY]), but one could select any order. It would also help the user be sure which part of the solution corresponds to which species.