import math,random
import matplotlib.pyplot as plt

def interact(p1,p2):
    interaction=0
    A=random.getrandbits(3)
    daughters=[]
    B=0
    if A<=3:
        B=random.getrandbits(1)
        interaction=1


    if interaction==1:
        if (p1=='E' or p2=='E') and (p1=='N' or p2=='N') and B==1:
            #print'ABSORBTION'
            daughters.append('H')
        elif (p1=='E' or p2=='E') and (p1=='N' or p2=='N') and B==0:
            #print'NEUTRON TRANSMUTATION'
            daughters.append('H')
            daughters.append('E')
        elif (p1=='NOBLE' and p2=='E') or (p1=='E' and p2=='NOBLE'):
            daughters.append('NOBLE')
        else:
            #print'NO INTERACTION'
            daughters.append(p1)
            daughters.append(p2)
    else:
        #print'NO INTERACTION'
        daughters.append(p1)
        daughters.append(p2)
    return daughters

def fusion(fuel,t,p):
    products=[]
    if fuel=='H':
        if t>=2000 and p>=50:
            products.append('N')
            products.append('P')
            products.append('E')
            products.append('PLASMA')
            products.append('NOBLE')
        else:
            products.append(fuel)

    elif fuel=='NOBLE':
        if t>=5000 and p>=100:
            products.append('N')
            products.append('P')
            products.append('E')
            products.append('PLASMA')
            products.append('CO2')
        else:
            products.append(fuel)

    elif fuel=='CO2':
        if t>=9500 and p>=200:
            products.append('O2')
            products.append('P')
            products.append('N')
            products.append('PLASMA')
            products.append('E')
        else:
            products.append(fuel)
    elif fuel=='O2':
        if t>=9700 and p>=250:
            products.append('BRMT')
            products.append('P')
            products.append('N')
            products.append('PLASMA')
            products.append('E')
        else:
            products.append(fuel)
    else:
        products.append(fuel)
    return products

popt=[]
ppt=[]
intt=[]
fust=[]
end=0
zoo=['E','N']
t=0
temp=float(3500)
pres=float(100)
time=0

while t<=25:
    if temp>9725:
        temp=float(9725)
    if pres>256:
        pres=float(256)
    c=0
    pp=0
    intc=0
    fusc=0
    print''
    print''
    print'### - ITERATION',t
    print'TEMP',temp
    print'PRES',pres
    print'NEUTRONS',zoo.count('N')
    print'ELECTRONS',zoo.count('E')
    print'PHOTONS',zoo.count('P')
    print'HYDROGEN',zoo.count('H')
    print'NOBLE',zoo.count('NOBLE')
    print'CO2',zoo.count('CO2')
    print'OXYGEN',zoo.count('O2')
    print'ZOO SIZE',len(zoo)
    for i in zoo:
        prt=zoo[c]
        if prt=='E' or prt=='N':
            prt=zoo.pop(c)
            find=0
            for j in range(len(zoo)):
                p2=zoo[j]
                if (prt=='E' and p2=='N') or (prt=='N' and p2=='E'):
                    find=1
                    this=j
                    break
            if find==1:
                prt2=zoo.pop(this)
                products=interact(prt,prt2)
                if products!=[prt,prt2] and products!=[prt2,prt]: #we have a winner!
                    intc+=1
                for j in products:
                    zoo.append(j)
            else:
                zoo.append(prt)
        elif prt=='PLASMA':
            if random.getrandbits(1)==1:
                zoo.pop(c)
        c+=1
    print''
    print'PRE-FUSION'
    print'NEUTRONS',zoo.count('N')
    print'ELECTRONS',zoo.count('E')
    print'PHOTONS',zoo.count('P')
    print'HYDROGEN',zoo.count('H')
    print'NOBLE',zoo.count('NOBLE')
    print'CO2',zoo.count('CO2')
    print'OXYGEN',zoo.count('O2')
    print'INTERACTIONS',intc
    print'ZOO SIZE',len(zoo)
    presS=pres
    for i in zoo:
        A=zoo.pop(0)
        B=fusion(A,temp,presS)
        if len(B)>1:
            fusc+=1
            pres+=float(0.2)
        for j in B:
            zoo.append(j)
            if j=='P' and len(B)>1:
                pp+=1
    ppt.append(pp)
    popt.append(len(zoo))
    intt.append(intc)
    fust.append(fusc)
    print''
    print'FINAL'
    print'NEUTRONS',zoo.count('N')
    print'ELECTRONS',zoo.count('E')
    print'PHOTONS',zoo.count('P')
    print'HYDROGEN',zoo.count('H')
    print'NOBLE',zoo.count('NOBLE')
    print'CO2',zoo.count('CO2')
    print'OXYGEN',zoo.count('O2')
    print'FUSIONS',fusc
    print'PHOTONS PRODUCED',pp
    print'ZOO SIZE',len(zoo)
    t+=1
    temp+=100
print'END SIZE',len(zoo)
plt.plot(range(t),popt)
plt.xlabel('TIME')
plt.ylabel('POPULATION')
plt.show()

plt.plot(range(t),intt)
plt.xlabel('TIME')
plt.ylabel('INTERACTIONS')
plt.show()

plt.plot(range(t),fust)
plt.xlabel('TIME')
plt.ylabel('FUSION EVENTS')
plt.show()

plt.plot(range(t),ppt)
plt.xlabel('TIME')
plt.ylabel('PHOTONS PRODUCED')
plt.show()