4249.sagews

# define rate constants
from scipy.integrate import odeint
import numpy as np
from pylab import *
from matplotlib import ticker


def Pd100_LH_model():
    #set some rate constants
    kadsO2 = 1000000
    kdesO2 = 55000

    kadsCO = 2000000
    kdesCO = 70000
    kreact = 200000

    # define the right-hand side of your rate equations
    def dThetadt((O2, CO), t):
        return [kadsO2*(1-O2-CO)*(1-O2-CO)-kdesO2*O2-kreact*O2*CO,
                kadsCO*(1-O2-CO)-kdesCO*CO-kreact*O2*CO]
    return dThetadt, kreact

t = np.linspace(0, 5.e-5, 1000000)
model, kreact = Pd100_LH_model()
O2, CO = odeint(model, (0,0), t).transpose()
figure()
subplot(211)
plot(t, O2, 'g-',label='$\Theta_{O_{2}}$')
plot(t, CO, 'b-', label='$\Theta_{CO}$')
ylabel('$\mathrm{Coverage}$')
ylim(0,1)
legend(loc=2)
subplot(212)
ylabel('$\mathrm{TOF}$')
plot(t, kreact*O2*CO, 'r-')
xlabel('$t\ \mathrm{\ in\ s}$')
savefig('TOF CO oxidation')

show()