PR to sync with latest breaking SBML change to Reaction

[anandijain]

There is something going on with negative stoichs.

The tests will fail since a new release of SBML hasn't happened.

I also probably messed up the boundary condition handling and something in reversible reactions.

Everything passes locally but 7 failures that all either look like one of the following:

Test Failed at /Users/anand/.julia/dev/SBMLToolkit/test/reactionsystem.jl:24
  Expression: isequal(Catalyst.get_eqs(rs), ModelingToolkit.Reaction[ModelingToolkit.Reaction(k1, nothing, [s1], nothing, [1.0])])
   Evaluated: isequal(Reaction[Reaction{Any, Float64}(k1, Term{Real, Nothing}[s1(t)], Any[], [1.0], Float64[], Pair{Any, Float64}[s1(t) => -1.0], true, nothing)], Reaction[Reaction{Any, Float64}(k1, Any[], Term{Real, Nothing}[s1(t)], Float64[], [1.0], Pair{Any, Float64}[s1(t) => 1.0], false, nothing)])

Test Failed at /Users/anand/.julia/dev/SBMLToolkit/test/reactionsystem.jl:129
  Expression: isequal((nothing, [s1], nothing, [1.0]), SBMLToolkit.getreagents(REACTION1.reactants, REACTION1.products, MODEL1))
   Evaluated: isequal((nothing, Num[s1(t)], nothing, [1.0]), (Num[s1(t)], nothing, [1.0], nothing))

@paulflang understandable if you don't have time to review but if you get the chance

[anandijain]

Now there are two test failures, but I'm not sure what wen't wrong here.

Model to MTK conversions: Test Failed at /Users/anand/.julia/dev/SBMLToolkit/test/reactionsystem.jl:140
  Expression: isequal(([ncs], [ncs], [1.0], [1.0]), SBMLToolkit.getreagents(constreac.reactants, constreac.products, constmod))
   Evaluated: isequal((Term{Real, Nothing}[constspec(t)], Term{Real, Nothing}[constspec(t)], [1.0], [1.0]), (nothing, Num[constspec(t)], nothing, [1.0]))

Model to MTK conversions: Test Failed at /Users/anand/.julia/dev/SBMLToolkit/test/reactionsystem.jl:141
  Expression: isequal((nothing, nothing, nothing, nothing), SBMLToolkit.getreagents(constreac.reactants, constreac.products, constmod))
   Evaluated: isequal((nothing, nothing, nothing, nothing), (nothing, Num[constspec(t)], nothing, [1.0]))

See the tests where I updated reactions, I'd appreciate a double check that I converted the netstoich to reactants and products correctly.

I'm also wondering what I need to do for the boundary condition stuff, since I effectively removed that logic

[paulflang]

I'm also wondering what I need to do for the boundary condition stuff, since I effectively removed that logic

This is the reason why the tests are failing iiuc. If boundary_condition=true the species does not change in concentration when the reaction happens. I discussed with Sam how to best implement this without messing up auto-detection from mass action kinetics. The solution is very simple: wherever boundary_condition=true for a reactant, we add the reactant and its stoichiometry to the products and vice versa. You'll have to restructure your code a little to implement this correctly.

[anandijain]

@paulflang okay yeah I think i see. what is the reason we don't push to reactants if bc is true here

[anandijain]

i guess because that would double things up right

[paulflang]

if bc is true we need to make sure that the stoichiometry of the species in question is the same on the right and left side of the reaction equation. however, we never want to modify the left hand side. this is because the left hand side can be used to determine if a reaction is mass action kinetics (i.e. A + B -> C gives mass action kinetic kAB, A -> B gives mass action kinetic k*A, etc.). This leaves us with the option to modify the product side. But this is already taken care of in line 149. In line 154 we only need to take care of the default case of bc = false, which is that reaction products are actually produced.

[anandijain]

this latest commit only fails on

Model to MTK conversions: Test Failed at /Users/anand/.julia/dev/SBMLToolkit/test/reactionsystem.jl:141
  Expression: isequal((nothing, nothing, nothing, nothing), SBMLToolkit.getreagents(constreac.reactants, constreac.products, constmod))
   Evaluated: isequal((nothing, nothing, nothing, nothing), (Num[constspec(t)], Num[constspec(t)], [1.0], [1.0]))

I don't see why thouhg

[paulflang]

This is because you call SBMLToolkit.getreagents(constreac.reactants, constreac.products, constmod when you should call SBMLToolkit.getreagents(constreac.reactants, constreac.products, constmod; rev=true. The rev=true just swaps reagents and products.

The reaction is X -> nothing, k*X. Then the reversible reaction is nothing -> X, k. But since species is constant we modify to X -> X, k*X and nothing -> nothing, k for fw and rv reaction, respectively. Note that we are never allowed to modify the substrates of a reaction, as we need to be able to determine whether the kineticLaw of a reaction follows mass action kinetics.

[anandijain]

thanks for the help paul #28 is a much cleaner copy of this