nzm

c-dichteG/a-dichteG.py

@@ -4,6 +4,16 @@ OCTO = 8
 N=5
 E=2**N
 
+T=1000 # number of iterations to reach the stationary distrib. on pi
+epsilon=0.05 #step for moving with the gradient
+
+# ---- iteration on phi -------
+NSTEPS_ITER=20
+LAMBDA=0.4 # coefficient for weighting the derivative
+
+# starting point
+phiFuks = [0.,1/3,1/3,2/3,1/3,2/3,2/3,1.]
+
 #Randomisation Function which return a vector
 def randomization(n):
         """
@@ -58,23 +68,23 @@ def tindex(config,pos):
 def Ptrans(phiT,config1,config2):
     assert(len(config1)==len(config2))
     L=len(config1)
-    w12=1.
+    p12=1.
     for pos in range(L):
+        assert( config2[pos]==0 or config2[pos]==1 )
         itr=tindex(config1,pos)
         if config2[pos]==0 :
-            w12 *= 1. - phiT[itr]
-        elif config2[pos]==1 :
-            w12 *= phiT[itr]
-        else:
-            raise ValueError("binary inputs only")
-        m=phiT[itr] if config2[pos]==1 else 1-phiT[itr]
+            p12 *= 1. - phiT[itr]  # product
+        else :
+            p12 *= phiT[itr]       # product
+        #m=phiT[itr] if config2[pos]==1 else 1-phiT[itr]
         #print("pos:%d indx:%d pres:%f m1:%f"%(pos,itr,w12,m))
-    return w12
+    return p12
 
 # uses Ptrans
+# can be optimised with shifts
 def calculateTransitionMatrix(phiT):
         """
-        from phi to transition matrix
+        from phi to transition matrix, 2^n transitions to perform
         """
         tMatrix= np.zeros((E, E))
         for j in range(E):
@@ -83,7 +93,7 @@ def calculateTransitionMatrix(phiT):
         return tMatrix
 
 
-# prints for all the configurations pi and calculates Q
+# prints pi for all the configurations and calculates Q
 def printDistribStats(PIst):
         """
         prints the distribution of pi
@@ -136,20 +146,23 @@ def iterateDistrib(PI,k,transMatrix):
         p=0.
         for i in range(1,E): # remove 0 because it is zero ! 
                 p+= transMatrix[k,i]*PI[i]
-        delta=p-PI[k]
+        #delta=p-PI[k]
         #print (" %s p=%f p'=%f delta=%f"%(configT[k],PI[k],p,delta))
         PI[k]=p
 
+# procedure for calculating the stationary distrtibution by iterating
+# T times ; fully asynchronous updating
 def iterateOnPi(PI, transMatrix):
         """
         T iterations, changes PI
         """
         for t in range(T):
-                # draw between 1 and E-2 (included)
-               k=np.random.randint(1,E-1)
-               #print("updated:::",k)
+               # draws RANDOMLY A CONFIGURATION between 1 and E-2 (included)
+               k= np.random.randint(1,E-1)
+               # print("updated:::",k)
                iterateDistrib(PI,k,transMatrix)
 
+# stationary ditrib of the proba to go to a 1
 def findStationaryDistrib(phi):
     """
     finds the stationary distribution from a random init distribution
@@ -190,14 +203,10 @@ def evaluateGradient(phiZ):
 
 ############
 print("----- starting test ------")
-T=5000 # number of iterations to reach the stationary distrib. on pi
-epsilon=0.05
 
-phiFuks = [0,1/3,1/3,2/3,1/3,2/3,2/3,1.]
 #phiG=randomization(OCTO)
 #print("PHI:%s"%phiG)
 
-
 #fills the config array
 configT=[]
 for i in range(E):
@@ -205,22 +214,21 @@ for i in range(E):
 
 
 # ---- for next phi -------
-def fadd(x,y):
-    coef= 1;#np.random.random() # [0,1] random weight for each coordinate
+def fWeightedAddition(x,y):
+    coef= np.random.random() # [0,1] random weight for each coordinate
     val= x + coef*LAMBDA*y
     if val < 0.: return 0.
     if val > 1.: return 1.
     return val
 
-# ---- iteration on phi -------
-NSTEPS_ITER=20
-LAMBDA=0.4 # coefficient for weighting the derivative
+
 
 phiIter = phiFuks.copy()  # deep copy
+#phiIter = randomization(8)
 for stepG in range(NSTEPS_ITER):
     print("----- step %d -----"%stepG)
-    grad= evaluateGradient(phiIter)
-    phiNew= list( map(fadd,phiIter,grad) )
+    grad= evaluateGradient(phiIter)  #vector with eight dimensions
+    phiNew= list( map(fWeightedAddition,phiIter,grad) ) # "moving" according to grad
     phiIter= phiNew.copy()
     #T=+2000
     LAMBDA*= 0.98