Hybrid composed reaction networks

[rveltz]

Hi,

I know you can do this:

prob = DiscreteProblem(10.0,(0.0,3.0))
jump_prob = JumpProblem(prob,Direct(),jump,jump2)

but I would like to do

sir_model1 = @reaction_network begin
           c1, s + i --> 2i
           c2, i --> r
       end c1 c2

sir_model2 = @reaction_network begin
           c1, s + i --> 2i
           c2, i --> r
       end c1 c2
jump_prob = JumpProblem(prob,Direct(),sir_model1,sir_model2)

which returns

ERROR: MethodError: no method matching JumpProblem(::DiscreteProblem{Array{Int64,1},Tuple{Float64,Float64},true,Tuple{Float64,Float64,Float64},DiscreteFunction{true,getfield(DiffEqBase, Symbol("##161#162")),Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}}, ::Direct, ::SIR, ::SIR)

Moreover, it would be good to have sir_model implemented as ConstantRateJump or provide a way to enforce this. Do you have a way to do this please?

Lastly, it would be awesome to be able to use sir_model1 for the ODE part (ODEProblem) of the PDMP and sir_model2 for the jump part. But I am probably dreaming...

[ChrisRackauckas]

I'd wait for a change to ModelingToolkit with this. With a proper IR this is a simple merge operation. In the current form it might be possible but it would be a little hacky.

[isaacsas]

I transferred your issue to DiffEqBio as it is more relevant here I think.

I'm not sure I completely understand what you want to do in general, so I'll give you a few answers to some of your questions:

jumps(sir_model1) will return a tuple of ConstantRateJumps for each reaction. When you call JumpProblem passing in the sir_model, however, we convert any reactions that are in mass action form to a MassActionJump as this provides better performance. So the jump_prob would in this case only have MassActionJumps. Is there an issue for you with having MassActionJumps? They are essentially a more restricted form of ConstantRateJump.

Building a hybrid model is something we'd like to support, but would take some work at this point. It is on the TODO list, but probably won't materialize until we switch to a ModellingToolkit backend, or at least not for a little while.

As for combining networks, I think this is something we could add if you were using min_reaction_networks and had not yet constructed the underlying jumps (i.e. called addjumps!). You could do it manually by constructing a new min_reaction_network, sir_model_merged, looping through sir_model1 and sir_model2 and calling addspecies!, addparam! and addreaction! on the sir_model_merged network.

[rveltz]

My issue is to build a jump problem (here a PDMP) from several reaction_networks choosing which one is used as an ODE or as a Jump part.

[isaacsas]

So we’re part way now on this issue as we now have merge! on generated ReactionSystems. Still working towards hybrid models.

[TorkelE]

With ModelingToolkit you have some decent flexibility with regards to this now. What do you think @isaacsas ?

[isaacsas]

We still need ModelingToolkit support for hybrid systems, and then an interface here for specifying which sub-components are represented in which model type.