Radio active decay-Monte Carlo Method

#Radio active decay
from pylab import*
from random import *
n,Lambda,t,tm=1000000,0.2,0.0,10
N=[n]
T=[t]
while n>0 and t<tm:
    for i in range(n):
        if random()<=Lambda:
            n-=1
    t+=1
    N.append(n)
    T.append(t)
plot(T,N,'r')
T=array(T)
Y=N[0]*exp(-Lambda*T)
plot(T,Y,'b')
legend(['simulated','exponential'])
show()

output