SIS

#coding=utf-8
def writeLine(file, i1, i2):
    file.write("%s,%s\n" % (i1, i2))

if __name__ == "__main__":
    file = open("test","w")
    i1 = 25.0
    i2 = 20.0
    i12 = 5.0
    N = 1000.0
    b1 = 0.0025
    b2 = 0.0020
    sigma1 = 3.0
    sigma2 = 2.0
    a1 = b1/sigma1
    a2 = b2/sigma2
    writeLine(file, i1, i2)
    e = 0
    for i in range(1, 100):
        s = N - i1 -i2 - i12
        
        q1 = (N if b1*s*(i1+i12)>N else b1*s*(i1+i12))
        q2 = (i12 if a2*i12>i12 else a2*i12)
        q3 = (i1 if a1*i1>i1 else a1*i1)
        q4 = (i1 if e*b2*i1*(i2+i12)>i1 else e*b2*i1*(i2+i12))
        print "%s,%s,%s,%s" % (q1,q2,q3,q4)
        delta1 = q1+q2-q3-q4
        i1_new = i1+delta1
        delta2 = (N if b2*s*(i2+i12)>N else b2*s*(i2+i12))\
            +(i12 if a1*i12>i12 else a1*i12)\
            -(i2 if a2*i2>i2 else a2*i2)\
            -(i2 if e*b1*i2*(i1+i12)>i2 else e*b1*i2*(i1+i12))
        i2_new = i2+delta2
        delta12 = (i2 if e*b1*i2*(i1+i12)>i2 else e*b1*i2*(i1+i12))\
            +(i1 if e*b2*i1*(i2+i12)>i1 else e*b2*i1*(i2+i12))\
            -(i12 if (a1+a2)*i12>i12 else (a1+a2)*i12)
        i12_new = i12+delta12
        if i1<0 or i2<0:
          break
        i1 = i1_new
        i2 = i2_new
        i12 = i12_new
        writeLine(file,(i1+i12)/N,(i2+i12)/N)
    file.close()