Создание модели Gillesp ie

Question

Я работаю в Python и пытаюсь воспроизвести Gillesp ie Stochasti c Моделирование с помощью циклов. Я пока закодировал ниже. Мне нужно добавить al oop, который проходит через Schemati c моделирования стохастики c. Кроме того, мне нужно добавить реакции на точку. Я не могу найти код для графа. Заранее спасибо.

import numpy as np
from scipy.special import factorial
import matplotlib.pyplot as plt

#reaction constrants
c1= 0.5
print ("c1: ", c1)

#initial molecular population numbers of S1 and S2
X0=1000
print ("X0: ", X0)
X1=999
print ("X1: ", X1)

#initialize time variable t and reaction n
t=0
n=0

print ("t: ", t)
print ("n: ",n)



  #calculate a1=h1*c1 where h1 is the total number of all possible # combos of pairs of s1

    a1= X0*X1*c1
    a = np.array([a1])
    a0=sum(a)

    #print debug
    print ("a1: ", a1)
    print ("a: ", a)
    print ("a0: ", a0)

     #generating random #s r1 for normal distribution using np.random.random
    r = np.random.random(2)
    r1= r[0]
    r2= r[1]
    #find the incrememnt T in time as (1/a0)*ln(1/r1)
    T= (1/a0)*np.log(1/r1)
    print ("T: ", T)

    #choose next reaction
    mu= 0
    N= r2*a0 - a[mu

]
while N>0:
    mu= mu+1
    N = N- a[mu]
next_r = mu
print ("next_r: ", next_r)

 #define the time of the next reaction
t1= t + T
print ("t1: ", t1)

#add one to the n
#number of reactions count n
n1 = n+1
print ("n1: ", n1)
#update the system according to the reaction that was chosen
X2=X1+1
X3= X1-1
print ("X2: ", X2)
print ("X3: ", X3)

while N>0:
    mu= mu+1
    N = N- a[mu]
next_r = mu
print ("next_r: ", next_r)

 #define the time of the next reaction
t2=t1+ T
print ("t2: ", t2)

#add one to the n
#number of reactions count n
n2 = n1+1
print ("n2: ", n2)
#update the system according to the reaction that was chosen


X1_all= np.array([X0,X1])
t_all= np.array([0,10])
ax=plt.subplot()
ax.plot(t_all,X1_all)
plt.show()