[Tutor] Optimizing a simple radioactive decay simulation

Alfred Milgrom fredm@smartypantsco.com
Thu Dec 12 19:45:02 2002


Hi Tim:

The most time-consuming part of your program appears to be the calls to the 
random routine. 2 million calls to random.randint(0,1) (which is what you 
need with an initial population of 1 million) takes about 5 minutes on my 
computer.

The rest of the time is due to your use of classes.
Here is a much simpler program which accomplishes similar result to yours 
without the use of classes. It still uses the random.randint routine, but 
it's about 1/3 faster than yours:

def main2():
     """Run the decay simulator."""
    numcoins = int(raw_input("How many coins? "))
    halflives = 0
     while numcoins > 0:
         print "Undecayed (heads): %s" % numcoins
         ds = [1 for coin in range(numcoins) if random.randint(0,1)]
         halflives += 1
         numcoins = len(ds)
     print "\n%s half-lives required." % halflives

if __name__ == '__main__':
     main2()

I'm sure that there is an even faster solution out there.

HTH, Best regards,
Alfred Milgrom

At 05:22 PM 12/12/02 -0600, Tim Wilson wrote:
>Hi everyone,
>
>My students and I did a fun little simulation of C-14 decay today using
>pennies to represent the atoms. We flipped the coins repeatedly until all
>the pennies had "decayed" (by turning up heads) and used the results to
>make a graph of the decay curve.
>
>I thought it would be fun to replicate the demo on the computer so I could
>expand the simulation to many more pennies. My first approach, which I will
>include below, was to create Coin and DecaySim classes that contain the
>needed methods. This worked perfectly, but was quite slow.
>-Tim
>
>
>Here's the output on my 1-GHz Athlon with 512 MB RAM.
>
>wilsont@galileo:~/Documents> time python decay.py
>How many coins? 1000000
>...
>21 half-lives required.
>
>real    12m6.267s
>user    8m26.820s
>sys     0m1.190s
>
>--
>Tim Wilson