Different trajectories for equivalent models using gillespie integrator

[leon9812]

The following two models should produce the same trajectories, but only r1 yields the expected results.

import tellurium as te
r1 = te.loada("""
   S1 => 2 S2; c1*S1
   S1 = 100; c1 = 1.0
""")
r2 = te.loada("""
   S1 => S2 + S2; c1*S1
   S1 = 100; c1 = 1.0
""")

def simulate(r):
    r.integrator = "gillespie"
    r.integrator.seed = 0
    return r.simulate(0, 2)

s1 = simulate(r1)
s2 = simulate(r2)

import matplotlib.pyplot as plt
plt.plot(s1["time"], s1["[S1]"], label="r1-S1")
plt.plot(s2["time"], s2["[S1]"], label="r2-S1", linestyle="dashed")
plt.plot(s1["time"], s1["[S2]"], label="r1-S2")
plt.plot(s2["time"], s2["[S2]"], label="r2-S2", linestyle="dashed")
plt.legend()
plt.show()

trajectories:

[luciansmith]

Thanks for the report! I'll investigate.

I will note that the problem is not, as I first thought, that the seed does different things in the two models. The problem is that the gillespie integrator doesn't handle 'S1 + S1' correctly.

[luciansmith]

This is now fixed! It should be included in the next roadrunner release. Thank you again for the report!