From: couchot Date: Thu, 28 Nov 2013 07:37:24 +0000 (+0100) Subject: t X-Git-Url: https://bilbo.iut-bm.univ-fcomte.fr/and/gitweb/desynchronisation-controle.git/commitdiff_plain/90e98a89726a12ad0e1d68d39c993aa73b3c4230 t --- diff --git a/exp_controle_asynchrone/simulMWSN.py b/exp_controle_asynchrone/simulMWSN.py index fbd48aa..4844593 100644 --- 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 +import sys as sy + + error = 0.1 epsilon = 1E-10 vrate = 0.8 @@ -17,6 +20,7 @@ POS = 1 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 +zeta = 0.1 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: - 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) @@ -304,14 +309,20 @@ def maj(k,maj_theta,mxg,idxexp): 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] @@ -442,19 +453,36 @@ 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 - nbexp = 5 res = {} m = [] - itermax = 10000 + itermax = 100000 def __maj_theta(k): + mem = [] 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 @@ -469,16 +497,25 @@ def __evalue_maj_theta__(): 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 + sy.exit(0) if k == itermax: print "nbre d'iteration trop grand" print "init" print init - exit + sy.exit(1) print "###############" print k @@ -492,7 +529,10 @@ def __evalue_maj_theta__(): ASYNC = False -__evalue_maj_theta__() +__evalue_maj_theta__(1,True) + + + #ASYNC = True #taux_succes = 0.9 #__evalue_maj_theta__() diff --git a/exp_controle_asynchrone/simulMWSN.py~ b/exp_controle_asynchrone/simulMWSN.py~ index 171c682..863b2ae 100644 --- a/exp_controle_asynchrone/simulMWSN.py~ +++ b/exp_controle_asynchrone/simulMWSN.py~ @@ -5,8 +5,11 @@ import numpy as np import pylab as pb from itertools import * from scipy import optimize as opt +from copy import deepcopy +import sys as sy -error = 0.01 + +error = 0.1 epsilon = 1E-10 vrate = 0.8 p = 0.7 @@ -16,7 +19,8 @@ nbiter = 1000 POS = 1 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 +zeta = 0.1 amplifieur = 1 sigma2 = 3500 Bi = 5 @@ -215,18 +220,23 @@ def maj_theta(k): -def maj(k,maj_theta): +def maj(k,maj_theta,mxg,idxexp): # initialisation des operateurs lagrangiens - global u, v, la, w, theta , q, Ps, Rh, eta, x, valeurFonctionDuale - - + global u, v, la, w, theta , q, Ps, Rh, eta, x, valeurFonctionDuale + + smax = 0 #print "iteration",k up = {} for h in V: for i in N: - if not ASYNC or random() < taux_perte: - up[(h,i)] = u[(h,i)]-theta*(eta[(h,i)]-sum([a[i][l]*x[(h,l)] for l in L])) + if not ASYNC or random() < taux_succes: + s = eta[(h,i)]-sum([a[i][l]*x[(h,l)] for l in L]) + if abs(s) > mxg : + print "ds calcul u",abs(s),idxexp + mxg = abs(s) + smax = max(smax,abs(s)) + up[(h,i)] = u[(h,i)]-theta*s else : up[(h,i)] = u[(h,i)] @@ -234,10 +244,13 @@ def maj(k,maj_theta): vp = {} for h in V: - #print "vp",gamma*mt.pow(Ps[h],float(1)/3) - #vp[h] = max(0,v[h]-theta*(Rh[h]- mt.log(float(sigma2)/D)/(gamma*mt.pow(Ps[h],float(1)/3)))) - if not ASYNC or random() < taux_perte: - vp[h] = max(0,v[h]-theta*(Rh[h]- mt.log(float(sigma2)/D)/(gamma*mt.pow(Ps[h],float(1)/3)))) + if not ASYNC or random() < taux_succes: + s = Rh[h]- mt.log(float(sigma2)/D)/(gamma*mt.pow(Ps[h],float(1)/3)) + if abs(s) > mxg : + print "ds calcul v",abs(s),idxexp + mxg = abs(s) + smax = max(smax,abs(s)) + vp[h] = max(0,v[h]-theta*s) else : vp[h] = v[h] @@ -245,8 +258,13 @@ def maj(k,maj_theta): lap={} for i in N: - if not ASYNC or random() < taux_perte: - resla = la[i]-theta*((q[i]*Bi -sum([aplus[i][l]*cs[l]*sum([x[(h,l)] for h in V]) for l in L])-cr*sum([amoins[i][l]*sum([x[(h,l)] for h in V]) for l in L])-psi(Ps,i))) + if not ASYNC or random() < taux_succes: + s = q[i]*Bi -sum([aplus[i][l]*cs[l]*sum([x[(h,l)] for h in V]) for l in L])-cr*sum([amoins[i][l]*sum([x[(h,l)] for h in V]) for l in L])-psi(Ps,i) + if abs(s) > mxg : + print "ds calcul la",abs(s),idxexp,i + mxg = abs(s) + smax = max(smax,abs(s)) + resla = la[i]-theta*s lap[i] = max(0,resla) else : lap[i] = la[i] @@ -256,8 +274,13 @@ def maj(k,maj_theta): wp={} for l in L: - if not ASYNC or random() < taux_perte: - wp[l] = w[l] + theta*(sum([a[i][l]*q[i] for i in N])) + if not ASYNC or random() < taux_succes: + s = sum([a[i][l]*q[i] for i in N]) + if abs(s) > mxg : + print "ds calcul w",abs(s),idxexp + mxg = abs(s) + smax = max(smax,abs(s)) + wp[l] = w[l] + theta*s else : wp[l] = w[l] @@ -265,8 +288,6 @@ def maj(k,maj_theta): thetap = maj_theta(k) - #print "k",k,"maj theta", thetap - qp={} def f_q(qi,i): @@ -274,7 +295,7 @@ def maj(k,maj_theta): return qi*qi + qi*fa for i in N: - if not ASYNC or random() < taux_perte: + if not ASYNC or random() < taux_succes: c = -float(sum([a[i][l]*w[l] for l in L]) - la[i]*Bi)/(2*amplifieur) rep = epsilon if c <= 0 else c qp[i] = rep @@ -291,10 +312,19 @@ def maj(k,maj_theta): #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_perte: + 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)) + 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))/float(16*zeta) + #print t + rep = mt.pow(t,float(3)/5) + Psp[h]=rep + """ +############reprendre ici + else : Psp[h] = Ps[h] @@ -313,7 +343,7 @@ def maj(k,maj_theta): return delta*rh*rh-v[h]*rh-sum([u[(h,i)]*eta[(h,i)] for i in N]) for h in V: - if not ASYNC or random() < taux_perte: + if not ASYNC or random() < taux_succes: rep = float(v[h])/(2*delta) Rhp[h] = 0 if rep < 0 else rep else : @@ -327,7 +357,7 @@ def maj(k,maj_theta): for h in V: for l in L: - if not ASYNC or random() < taux_perte: + if not ASYNC or random() < taux_succes: rep = -float(cs[l]*sum([la[i]*aplus[i][l] for i in N]) +cr*sum([la[i]*amoins[i][l] for i in N])+sum([u[(h,i)]*a[i][l] for i in N]))/(2*delta) xp[(h,l)] = 0 if rep < 0 else rep else : @@ -351,7 +381,7 @@ def maj(k,maj_theta): arret = abs(valeurFonctionDuale-valeurFonctionDualep) < error - return (up,vp,lap,wp,thetap,etap,qp,Psp,Rhp,xp,valeurFonctionDualep,arret) + return (up,vp,lap,wp,thetap,etap,qp,Psp,Rhp,xp,valeurFonctionDualep,arret,mxg,smax) @@ -373,48 +403,32 @@ print "x",x def initialisation(): - global u, v, la, w, theta , q, Ps, Rh, eta, x, valeurFonctionDuale + global u, v, la, w, theta , q, Ps, Rh, eta, x,init theta = omega - - # TODO Ahmed : initialiser avec qi >0 et somme (ail.qi) = 0 Ei/i = Ti - #q = [0 for i in N] q = {} for i in N : q[i] = 0.15 + random()*0.05 - - - - #TODO Ahmed doit donner les valeurs initiales pour Ps Ps= {} for h in V : Ps[h] = 0.2+random()*0.3 - - - - #TODO Ahmed doit donner les valeurs initiales pour Rh Rh = {} for vi in V: Rh[vi] = 0.1 + random()*0.1 - - eta = {} for h in V : for i in N: eta[(h,i)]= etahi(h,i,Rh) - - - # Ok pour 0 x={} for h in V : for l in L: - x[(h,l)]=0#random() + x[(h,l)]=0 @@ -437,43 +451,94 @@ def initialisation(): for l in L: w[l] = random() + init = [deepcopy(q),deepcopy(Ps),deepcopy(Rh),deepcopy(eta), + deepcopy(x),deepcopy(u),deepcopy(v),deepcopy(la),deepcopy(w)] -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 - nbexp = 50 res = {} m = [] + itermax = 100000 + def __maj_theta(k): + mem = [] om = omega/(mt.pow(k,0.75)) return om for idxexp in range(nbexp): - - initialisation() + mxg = 0 + if not(out): + initialisation() + else : + initialisation_() + k = 1 arret = False - while k < 1000000 and not arret : - (u,v,la,w,theta,eta,q,Ps,Rh,x,valeurFonctionDuale,ar)=maj(k,__maj_theta) - errorq = abs(min(q.values())-max(q.values())) + 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 - if k%5000 ==0 : + variation = "+" if smax > sm else "-" + sm = smax + 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 + sy.exit(0) + 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 -#ASYNC = False -#__evalue_maj_theta() -ASYNC = True -taux_perte = 0.9 -__evalue_maj_theta__() +ASYNC = False +__evalue_maj_theta__(1,True) + + + +#ASYNC = True +#taux_succes = 0.9 +#__evalue_maj_theta__() """