evalPRNG/compareFonctionMixingTime.py

   1 import random as rn
   2 import numpy as np
   3 from math import *
   4 from optparse import OptionParser
   5
   6
   7 def bin(elem,n):
   8     """Convertit un nombre en binaire"""
   9     q = -1
  10     res =[0 for i in range(n)]
  11     i = 1
  12     while q != 0:
  13         q = elem // 2
  14         r = elem % 2
  15         res[n-i] =  r
  16         elem = q
  17         i+=1
  18     return res
  19
  20 def dec(ch,n):
  21     l = len(ch)
  22     acc = 0
  23     for i in range(l):
  24         if ch[i]==1 :
  25             acc = acc + 2**(n-i-1)
  26     return acc
  27
  28
  29 def MarkovMatrixUnPas(fbin,n,p2n,lp2nm1):
  30     #print p2n
  31     MM = np.zeros((p2n,p2n))
  32     # diagonal
  33     for i in lp2nm1 :
  34         ib = bin(i,n)
  35         MM[i,i] += 0.5
  36         image = fbin[tuple(ib)]
  37         for j in range(n):
  38             ipb =[_ for _ in ib]
  39             ipb[j] = image[j]
  40             MM[i, dec(ipb,n)] +=float(1)/(2*n)
  41     return MM
  42 """
  43
  44 def MarkovMatrixUnPas(fbin,n,p2n,lp2nm1):
  45     MM = np.zeros((p2n,p2n))
  46     # diagonal
  47     for i in lp2nm1 :
  48         ib = bin(i,n)
  49         #MM[i,i] += 0.5
  50         image = fbin[tuple(ib)]
  51         for j in range(n):
  52             ipb =[_ for _ in ib]
  53             ipb[j] = image[j]
  54             MM[i, dec(ipb,n)] +=float(1)/(n)
  55     return MM
  56 """
  57
  58
  59
  60
  61 def MarkovMatrixSaut(fbin,n,p2n,lp2nm1):
  62     #print p2n
  63     MM = np.zeros((p2n,p2n))
  64     # diagonal
  65     for i in lp2nm1 :
  66         ib = bin(i,n)
  67         MM[i,i] += 0.5
  68         image = fbin[tuple(ib)]
  69         for j in lp2nm1:
  70             st = bin(j,n)
  71             ipb =[_ for _ in ib]
  72             for k in range(n):
  73                 if st[k] == 1 :
  74                     ipb[k] = image[k]
  75             MM[i, dec(ipb,n)] +=float(1)/(2*2**n)
  76     return MM
  77
  78 """
  79
  80 def MarkovMatrixSaut(fbin,n,p2n,lp2nm1):
  81     MM = np.zeros((p2n,p2n))
  82     # diagonal
  83     for i in lp2nm1 :
  84         ib = bin(i,n)
  85         #MM[i,i] += 0.5
  86         image = fbin[tuple(ib)]
  87         for j in lp2nm1:
  88             st = bin(j,n)
  89             ipb =[_ for _ in ib]
  90             for k in range(n):
  91                 if st[k] == 1 :
  92                     ipb[k] = image[k]
  93             MM[i, dec(ipb,n)] +=float(1)/(2**n)
  94     return MM
  95
  96
  97 """
  98
  99
 100
 101 ## n'est pas un circuit hmiltonien
 102 lf= []
 103 #lf +=[ [13, 10, 9, 14, 3, 11, 1, 12, 15, 4, 7, 5, 2, 6, 0, 8]] #4
 104 #lf += [[29, 22, 25, 30, 19, 27, 24, 16, 21, 6, 5, 28, 23, 26, 1, 17, 31, 12, 15, 8, 10, 14, 13, 9, 3, 2, 7, 20, 11, 18, 0, 4]] #5
 105 lf +=  [[55, 60, 45, 44, 58, 62, 61, 48, 53, 50, 52, 36, 59, 34, 33, 49, 15, 42, 47, 46, 35, 10, 57, 56, 7, 54, 39, 37, 51, 2, 1, 40, 63, 26, 25, 30, 19, 27, 17, 28, 31, 20, 23, 21, 18, 22, 16, 24, 13, 12, 29, 8, 43, 14, 41, 0, 5, 38, 4, 6, 11, 3, 9, 32]] #6
 106 #lf+= [[47, 58, 57, 44, 51, 42, 41, 60, 63, 22, 37, 53, 19, 54, 48, 32, 45, 14, 13, 46, 11, 43, 33, 8, 7, 36, 39, 4, 34, 50, 1, 40, 29, 62, 61, 30, 59, 27, 17, 56, 55, 20, 23, 52, 18, 2, 49, 24, 31, 10, 9, 28, 3, 26, 25, 12, 15, 38, 21, 5, 35, 6, 0, 16]]
 107
 108 #lf += [[111, 94, 93, 116, 122, 90, 125, 88, 115, 126, 119, 84, 123, 98, 81, 120, 109, 106, 105, 110, 99, 107, 104, 72, 71, 118, 117, 96, 103, 102, 113, 64, 79, 86, 95, 124, 83, 91, 121, 24, 85, 22, 69, 20, 19, 114, 17, 112, 77, 76, 13, 108, 74, 10, 9, 73, 67, 66, 101, 100, 75, 82, 97, 0, 127, 54, 57, 62, 51, 59, 56, 48, 53, 38, 37, 60, 55, 58, 33, 49, 63, 44, 47, 40, 42, 46, 45, 41, 35, 34, 39, 52, 43, 50, 32, 36, 29, 28, 61, 92, 26, 18, 89, 25, 87, 30, 23, 4, 27, 2, 16, 80, 31, 78, 15, 14, 3, 11, 8, 12, 5, 70, 21, 68, 7, 6, 65, 1]] #7
 109
 110 #lf += [[223, 190, 249, 254, 187, 251, 233, 232, 183, 230, 247, 180, 227, 178, 240, 248, 237, 236, 253, 172, 203, 170, 201, 168, 229, 166, 165, 244, 163, 242, 241, 192, 215, 220, 205, 216, 218, 222, 221, 208, 213, 210, 212, 214, 219, 211, 217, 209, 239, 202, 207, 140, 139, 234, 193, 204, 135, 196, 199, 132, 194, 130, 225, 200, 159, 62, 185, 252, 59, 250, 169, 56, 191, 246, 245, 52, 243, 50, 176, 48, 173, 238, 189, 44, 235, 42, 137, 184, 231, 38, 37, 228, 35, 226, 177, 224, 151, 156, 141, 152, 154, 158, 157, 144, 149, 146, 148, 150, 155, 147, 153, 145, 175, 206, 143, 136, 11, 142, 129, 8, 7, 198, 197, 4, 195, 2, 161, 160, 255, 124, 109, 108, 122, 126, 125, 112, 117, 114, 116, 100, 123, 98, 97, 113, 79, 106, 111, 110, 99, 74, 121, 120, 71, 118, 103, 101, 115, 66, 65, 104, 127, 90, 89, 94, 83, 91, 81, 92, 95, 84, 87, 85, 82, 86, 80, 88, 77, 76, 93, 72, 107, 78, 105, 64, 69, 102, 68, 70, 75, 67, 73, 96, 55, 58, 45, 188, 51, 186, 61, 40, 119, 182, 181, 53, 179, 54, 33, 49, 15, 174, 47, 60, 171, 46, 57, 32, 167, 6, 36, 164, 43, 162, 1, 0, 63, 26, 25, 30, 19, 27, 17, 28, 31, 20, 23, 21, 18, 22, 16, 24, 13, 10, 29, 14, 3, 138, 41, 12, 39, 134, 133, 5, 131, 34, 9, 128]]#8 pas totally
 111
 112 #lf +=[[223, 250, 249, 254, 187, 251, 233, 232, 183, 230, 247, 180, 227, 178, 240, 248, 237, 236, 173, 172, 171, 238, 201, 168, 229, 166, 228, 244, 235, 242, 241, 192, 215, 220, 205, 216, 218, 222, 153, 152, 151, 210, 212, 214, 219, 146, 217, 209, 239, 142, 141, 206, 195, 234, 193, 136, 231, 196, 199, 197, 194, 226, 225, 200, 63, 188, 253, 252, 59, 190, 189, 176, 191, 246, 245, 164, 243, 162, 161, 177, 143, 170, 45, 44, 43, 138, 185, 184, 135, 38, 167, 165, 179, 34, 129, 224, 31, 154, 221, 158, 147, 26, 25, 156, 159, 22, 213, 149, 211, 150, 144, 208, 207, 14, 13, 204, 203, 202, 169, 8, 133, 198, 132, 4, 139, 131, 1, 160, 255, 124, 109, 108, 122, 126, 125, 112, 117, 114, 116, 100, 123, 98, 97, 113, 79, 106, 111, 110, 99, 74, 121, 120, 71, 118, 103, 101, 115, 66, 65, 104, 127, 90, 89, 94, 83, 91, 81, 92, 95, 84, 87, 85, 82, 86, 80, 88, 77, 76, 93, 72, 107, 78, 105, 64, 69, 102, 68, 70, 75, 67, 73, 96, 55, 58, 57, 62, 51, 186, 41, 40, 119, 182, 181, 53, 35, 54, 48, 56, 175, 174, 61, 60, 11, 46, 9, 32, 37, 6, 36, 52, 163, 50, 49, 0, 23, 28, 157, 24, 155, 30, 29, 16, 21, 18, 20, 148, 27, 19, 145, 17, 47, 10, 15, 140, 3, 42, 137, 12, 39, 134, 7, 5, 2, 130, 33, 128]] #8 totally
 113
 114
 115 #lf +=[
 116 #    [29, 22, 21, 30, 19, 27, 24, 28, 7, 20, 5, 4, 23, 26, 25, 17, 31, 12, 15, 8#, 10, 14, 13, 9, 3, 2, 1, 6, 11, 18, 0, 16],#sylvain 1
 117 #    [29, 22, 25, 30, 19, 27, 24, 16, 21, 6, 5, 28, 23, 26, 1, 17, 31, 12, 15, 8, 10, 14, 13, 9, 3, 2, 7, 20, 11, 18, 0, 4],#sylvain 2
 118 #    [29, 22, 21, 24, 31, 30, 27, 26, 7, 23, 5, 28, 3, 19, 25, 17, 13, 12, 9, 8, 15, 2, 1, 10, 6, 14, 4, 20, 11, 18, 0, 16], #mini avec 7
 119 #    [29, 22, 25, 30, 19, 27, 24, 16, 21, 6, 5, 28, 23, 26, 1, 17, 31, 12, 15, 8, 10, 14, 13, 9, 3, 2, 7, 20, 11, 18, 0, 4]] #5
 120
 121
 122 # 14 = (1,1,1,0) = f(0,0,0,0)
 123 #  6 = (0,1,1,0) = f(0,0,0,1)
 124 #....
 125 # 8  = (1,0,0,0) = f(1,1,1,1)
 126
 127
 128 def traite_f(f,dev):
 129     n = int(log(len(f))/log(2))
 130     #nbre d'elements
 131     p2n = int(pow(2,n))
 132     lp2nm1 = range(p2n)
 133     # pour eviter de le calculer a chaque fois
 134     # fbin est la fonction (representee comme un dico)
 135     # qui au tuple binaire (x3,x2,x1,x0) associe
 136     # le nombre f(x3,x2,x1,x0)
 137     fbin={}
 138     for j in range(len(f)):
 139         fbin[tuple(bin(j,n))] = bin(f[j],n)
 140
 141     MM = MarkovMatrixUnPas(fbin,n,p2n,lp2nm1)
 142     M = np.zeros((p2n,p2n))
 143     np.copyto(M,MM)
 144
 145     error = 1
 146     cpt = 2
 147     while error > dev and cpt < 1000:
 148         M = np.dot(M,MM)
 149         error =max([sqrt(sum([(M[i,j] - float(1)/p2n)**2  for i in range(p2n)])) for j in range(p2n)])
 150         #print cpt, error, M
 151         cpt +=1
 152
 153
 154     """
 155     MMs = MarkovMatrixSaut(fbin,n,p2n,lp2nm1)
 156     Ms = np.zeros((p2n,p2n))
 157     np.copyto(Ms,MMs)
 158
 159
 160     error = 1
 161     cpts = 2
 162     while error > dev :
 163         if n==4 and cpts==2:
 164             print Ms
 165         Ms = np.dot(Ms,MMs)
 166         error =max([sqrt(sum([(Ms[i,j] - float(1)/p2n)**2  for i in range(p2n)])) for j in range(p2n)])
 167         #print cpt, error, M
 168         cpts +=1
 169
 170     """
 171     cpts=1
 172     return n, cpt-1, cpts-1
 173
 174
 175 def get_treatedlist(rf):
 176     treatedlist = []
 177     fd = open(rf,"r")
 178     ll = fd.readlines()
 179     compteur = 0
 180     for l in ll:
 181         compteur +=1
 182         # suppr [ et ]
 183         #print compteur, l
 184
 185         l = l[1:len(l)-2]
 186         #print l
 187         lkm = [x for x in l.split(",")]
 188         #print lkm
 189         treatedlist.append([int(x) for x in lkm])
 190     return treatedlist
 191
 192 def lcompare(l1,l2):
 193     taille= len(l1)
 194     cpt = 0
 195     ok = True
 196     while ok and cpt < taille:
 197         if l1[cpt] == l2[cpt]:
 198             cpt +=1
 199         else:
 200             ok=False
 201     if ok :
 202         return 0
 203     else :
 204         return cpt
 205
 206
 207
 208
 209
 210 def main():
 211     global lf
 212     if rf != False:
 213         lf =get_treatedlist(rf)
 214     def compareMarche(x,y):
 215         return cmp(x[1],y[1])
 216     def compareSaute(x,y):
 217         return cmp(x[2],y[2])
 218
 219     brut=[]
 220     cpt = 0
 221     for f in lf:
 222         #print "cpt",cpt
 223         dev = 1E-14
 224         n,cptM,cptS = traite_f(f,dev)
 225         brut +=[(f,cptM,cptS)]
 226         cpt+=1
 227         """
 228         print f, n,cptM,cptS
 229         #print 8*n*n + 4*n*log(n+1)#,8*n*n + (n+2)*(log(n)+2)
 230         #print "Pour "+str(n)+" bits et pour eps="+str(dev)+", appels a rand moy pr 1 bit genere.",
 231         print cpt
 232
 233         brut.sort(cmp=compareMarche)
 234         f,cptM,cptS = brut[0]
 235
 236         print f
 237         print "En marchant:", cptM, cptM*log(n)/(log(2)*n),
 238         print "En sautant:", cptS
 239
 240         brut.sort(cmp=compareSaute)
 241         fp,cptM,cptS = brut[0]
 242         print fp
 243         print "En sautant:", cptS
 244
 245         cp = lcompare(f,fp)
 246         if cp > 0:
 247             print "diff",f[cp],fp[cp]
 248         else :
 249             print "egales"
 250         """
 251         #print brut
 252     brut.sort(cmp=compareMarche)
 253     print brut
 254     (meilleure,k,kk) = brut[0]
 255     n,cpM,cPS = traite_f(f,1E-6)
 256     print "Meilleure fonction :\n", f
 257     print "CPM :", cpM
 258
 259
 260
 261 rf = False
 262 def options():
 263     global rf
 264     parser = OptionParser()
 265     parser.add_option("-i", "--input", dest="i",
 266                       help="file of sequences")
 267
 268     (options, args) = parser.parse_args()
 269     if (options.i != None):
 270         rf = options.i
 271
 272
 273 if __name__ == "__main__":
 274     options()
 275     print "------"
 276     main()
 277