cedille

[desynchronisation-controle.git] / exp_controle_asynchrone / simulMWSN.py

diff --git a/exp_controle_asynchrone/simulMWSN.py b/exp_controle_asynchrone/simulMWSN.py
index fbd48aa9323587da7fd791e0049e88def6547823..bf174629478e737918827f3b6bc0c04d2bdc5805 100644 (file)
--- a/exp_controle_asynchrone/simulMWSN.py
+++ b/exp_controle_asynchrone/simulMWSN.py
@@ -6,6 +6,9 @@ import pylab as pb
 from itertools import *
 from scipy import optimize as opt
 from copy import deepcopy 
 from itertools import *
 from scipy import optimize as opt
 from copy import deepcopy 

+import sys as sy
+
+

 error  = 0.1
 epsilon = 1E-10
 vrate = 0.8
 error  = 0.1
 epsilon = 1E-10
 vrate = 0.8


@@ -17,6 +20,7 @@ POS = 1
 POS_NUL = 2
 POSINF1 = 3
 init = []
 POS_NUL = 2
 POSINF1 = 3
 init = []

+fichier_init="config_initiale_default.txt"

 
 
 
 
 
 


@@ -103,6 +107,7 @@ alpha = 0.5
 beta = 1.3E-8
 gamma = 55.54
 delta = 0.2
 beta = 1.3E-8
 gamma = 55.54
 delta = 0.2

+zeta = 0.1

 amplifieur = 1
 sigma2 = 3500
 Bi = 5
 amplifieur = 1
 sigma2 = 3500
 Bi = 5


@@ -240,7 +245,7 @@ def maj(k,maj_theta,mxg,idxexp):
     vp = {}
     for h in V:
         if not ASYNC or  random() < taux_succes:
     vp = {}
     for h in V:
         if not ASYNC or  random() < taux_succes:

-            s = Rh[h]- mt.log(float(sigma2)/D)/(gamma*mt.pow(Ps[h],float(1)/3))
+            s = Rh[h]- mt.log(float(sigma2)/D)/(gamma*mt.pow(Ps[h],float(2)/3))

             if abs(s) > mxg :
                 print "ds calcul v",abs(s),idxexp
                 mxg = abs(s) 
             if abs(s) > mxg :
                 print "ds calcul v",abs(s),idxexp
                 mxg = abs(s) 


@@ -304,14 +309,20 @@ def maj(k,maj_theta,mxg,idxexp):
     Psp={}
     #print "maj des des Psh" 
     def f_Ps(psh,h):
     Psp={}
     #print "maj des des Psh" 
     def f_Ps(psh,h):

-        #print "ds f_ps",psh, v[h]* mt.log(float(sigma2)/D)/(gamma*((psh**2)**(float(1)/3))) +la[h]*psh
-        return v[h]* mt.log(float(sigma2)/D)/(gamma*mt.pow(float(2)/3)) +la[h]*psh
-    for h in V:
-        if not ASYNC or  random() < taux_succes:
-            lah = 0.05 if la[h] == 0 else  la[h]
-            rep = (float(2*v[h]*mt.log(float(sigma2)/D))/mt.pow(3*gamma*lah,float(3)/5))
-            Psp[h] = epsilon if rep <= 0 else rep
-        else :
+        #print "ds f_ps",psh, v[h]* mt.log(float(sigma2)/D)/(gamma*((psh**2)**(float(2)/3))) +la[h]*psh
+         return v[h]* mt.log(float(sigma2)/D)/(gamma*mt.pow(float(2)/3)) +la[h]*psh
+     for h in V:
+         if not ASYNC or  random() < taux_succes:
+             """
+             lah = 0.05 if la[h] == 0 else  la[h]
+             rep = mt.pow(float(2*v[h]*mt.log(float(sigma2)/D))/(3*gamma*lah),float(3)/5)
+             Psp[h] = epsilon if rep <= 0 else rep
+             """
+             t= float(-3*la[h]+mt.sqrt(9*(la[h]**2)+64*zeta*v[h]*mt.log(float(sigma2)/D)/gamma))/(16*zeta)
+             #print t
+             rep = mt.pow(t,float(3)/5)
+             Psp[h]=rep
+         else :

             Psp[h] = Ps[h]
 
 
             Psp[h] = Ps[h]
 
 


@@ -442,19 +453,37 @@ def initialisation():
 
 
 
 
 
 

-def __evalue_maj_theta__():
+def initialisation_():
+    global u, v, la, w, theta , q,  Ps, Rh,  eta, x,init 
+    fd = open(fichier_init,"r")
+    l= fd.readline()
+    init_p = eval(l)
+    print init_p
+    theta = omega
+    (q,Ps,Rh,eta,x,u,v,la,w) = tuple([deepcopy(x) for x in init_p])
+    init = [deepcopy(q),deepcopy(Ps),deepcopy(Rh),deepcopy(eta),
+            deepcopy(x),deepcopy(u),deepcopy(v),deepcopy(la),deepcopy(w)]
+
+
+
+def __evalue_maj_theta__(nbexp,out=False):

     global u, v, la, w, theta , q,  Ps, Rh,  eta, x, valeurFonctionDuale 
     global u, v, la, w, theta , q,  Ps, Rh,  eta, x, valeurFonctionDuale 

-    nbexp = 5
+    nbexp = 10

     res = {}
     m = []
     res = {}
     m = []

-    itermax = 10000
+    itermax = 100000

  
     def __maj_theta(k):
  
     def __maj_theta(k):

+        mem = []

         om = omega/(mt.pow(k,0.75))
         return om
     for idxexp in range(nbexp):
         mxg = 0
         om = omega/(mt.pow(k,0.75))
         return om
     for idxexp in range(nbexp):
         mxg = 0

-        initialisation()
+        if not(out):
+            initialisation()
+        else :
+            initialisation_()
+            

         k = 1
         arret = False
         sm = 0
         k = 1
         arret = False
         sm = 0


@@ -464,7 +493,71 @@ def __evalue_maj_theta__():
             arret = errorq <  error
             k+=1
             variation = "+" if smax > sm else "-"
             arret = errorq <  error
             k+=1
             variation = "+" if smax > sm else "-"

+            print variation,
+            if k%100 ==0 :
+                print "k:",k,"erreur sur q", errorq, "et q:",q
+                print "maxg=", mxg
+                mem = [deepcopy(q),deepcopy(Ps),deepcopy(Rh),deepcopy(eta),
+                       deepcopy(x),deepcopy(u),deepcopy(v),deepcopy(la),deepcopy(w)]
+            if k%4500 == 0 :
+                print "#########\n",mem,"\#########\n"
+            if k%4600 == 0 :
+                print "#########\n",mem,"\#########\n"
+
+
+
+            if smax - sm  > 500:
+                print "variation trop grande"
+                print "init"
+                print init
+                exit 

             sm = smax
             sm = smax

+
+        if k == itermax:
+            print "nbre d'iteration trop grand"
+            print "init"
+            print init
+            sy.exit(1)
+
+        print "###############"
+        print k
+        print "###############"
+        m += [k]
+
+    print (min(m),max(m),float(sum(m))/nbexp,m),m
+
+
+
+def __une_seule_exp__(fichier_donees):
+    global u, v, la, w, theta , q,  Ps, Rh,  eta, x, valeurFonctionDuale 
+    initialisation()
+        
+
+    fichier = open(fichier_donees, "r")
+    instructions ={}
+    for line in fichier:    
+         l = line.split("=")
+         instructions[l[0]] = eval(l[1])
+    u, v, la, w,  q,  Ps, Rh,  eta, x, = instructions['u'],    instructions['v'],    instructions['la'],    instructions['w'],    instructions['q'],    instructions['Ps'],    instructions['Rh'],    instructions['eta'],    instructions['x']
+    nbexp = 1
+    res = {}
+    m = []
+    itermax = 100000
+ 
+    def __maj_theta(k):
+        om = omega/(mt.pow(k,0.75))
+        return om
+    for idxexp in range(nbexp):
+        mxg = 0
+        k = 1
+        arret = False
+        sm = 0
+        while k < itermax  and not arret : 
+            (u,v,la,w,theta,eta,q,Ps,Rh,x,valeurFonctionDuale,ar,mxg,smax)=maj(k,__maj_theta,mxg,idxexp)
+            errorq =  (max(q.values()) - min(q.values()))/min(q.values())
+            arret = errorq <  error
+            k+=1
+            variation = "+" if smax > sm else "-"

             print variation,
             if k%100 ==0 :
                 print "k:",k,"erreur sur q", errorq, "et q:",q
             print variation,
             if k%100 ==0 :
                 print "k:",k,"erreur sur q", errorq, "et q:",q


@@ -474,6 +567,8 @@ def __evalue_maj_theta__():
                 print "init"
                 print init
                 exit 
                 print "init"
                 print init
                 exit 

+            sm = smax
+

         if k == itermax:
             print "nbre d'iteration trop grand"
             print "init"
         if k == itermax:
             print "nbre d'iteration trop grand"
             print "init"


@@ -490,45 +585,14 @@ def __evalue_maj_theta__():
 
 
 
 
 
 

-

 ASYNC = False
 ASYNC = False

-__evalue_maj_theta__()
+__une_seule_exp__("config_initiale.py")
+#__evalue_maj_theta__()

 #ASYNC = True
 #taux_succes = 0.9
 #__evalue_maj_theta__()
 
 
 #ASYNC = True
 #taux_succes = 0.9
 #__evalue_maj_theta__()
 
 

-"""
-initialisation()
-k = 1
-arret = False
-while k < 10000  and not arret : 
-    (u,v,la,w,theta,eta,q,Ps,Rh,x,valeurFonctionDuale,ar)=maj(k,maj_theta)
-    arret = ar
-    k+=1
-    errorq =  abs(min(q.values())-max(q.values()))
-    print "errorq",errorq
-    arret = errorq <  error
-"""
-
-"""
-    print "u",u
-    print "w",w
-    print "theta",theta
-    print "eta", eta
-    print "q",q
-    print "v",v
-    print "lambda",la
-    print "Ps",Ps
-    print "Rh",Rh
-    print "x",x
-"""
-
-
-"""
-    k +=1
-
-"""