import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

def f(x,t):
    k1, k2, k3 =0.04, 3*10**7, 10**4
    return np.array([-k1*x[0]+k3*x[1]*x[2], k1*x[0]-k3*x[1]*x[2]-k2*x[1]**2,
        k2*x[1]**2])

y0=np.array([1. ,0., 0.])
tpt=np.linspace(0,1,1000)
Y=odeint(f,y0,tpt)

plt.plot(tpt,Y[:,1])
plt.show()