Index: Ising2D-MP.py
===================================================================
--- Ising2D-MP.py	(revision 87)
+++ Ising2D-MP.py	(revision 88)
@@ -33,7 +33,7 @@
         # Random access coordonate
         X,Y=numpy.random.randint(SizeX),numpy.random.randint(SizeY)
         
-        DeltaE=J*sigma[X,Y]*(2*(sigma[X,(Y+1)%SizeY]+
+        DeltaE=sigma[X,Y]*(2*J*(sigma[X,(Y+1)%SizeY]+
                                 sigma[X,(Y-1)%SizeY]+
                                 sigma[(X-1)%SizeX,Y]+
                                 sigma[(X+1)%SizeX,Y])+B)
@@ -46,13 +46,13 @@
 def Magnetization(sigma,M):
     return(numpy.sum(sigma)/(sigma.shape[0]*sigma.shape[1]*1.0))
 
-def Energy(sigma,J):
+def Energy(sigma,J,B):
     # Copier et caster 
     E=numpy.copy(sigma).astype(numpy.float32)
     
     # Appel par slice
-    E[1:-1,1:-1]=-J*E[1:-1,1:-1]*(E[:-2,1:-1]+E[2:,1:-1]+
-                                  E[1:-1,:-2]+E[1:-1,2:])
+    E[1:-1,1:-1]=E[1:-1,1:-1]*(2.*J*(E[:-2,1:-1]+E[2:,1:-1]+
+                                     E[1:-1,:-2]+E[1:-1,2:])+B)
     
     # Bien nettoyer la peripherie
     E[:,0]=0
@@ -64,6 +64,13 @@
 
     return(Energy/(E.shape[0]*E.shape[1]*1.0))
 
+def CriticalT(T,E):
+
+    Epoly=numpy.poly1d(numpy.polyfit(T,E,T.size/3))
+    dEpoly=numpy.diff(Epoly(T))
+    dEpoly=numpy.insert(dEpoly,0,0)
+    return(T[numpy.argmin(dEpoly)])
+
 def DisplayCurves(T,E,M,J,B):
 
     plt.xlabel("Temperature")
@@ -99,12 +106,12 @@
     try:
         opts, args = getopt.getopt(sys.argv[1:],"hcj:b:z:i:s:e:p:u:",["coupling=","magneticfield=","size=","iterations=","tempstart=","tempend=","tempstep=","units"])
     except getopt.GetoptError:
-        print '%s -j <Coupling Factor> -b <Magnetic Field> -z <Size of Lattice> -i <Iterations> -s <Minimum Temperature> -e <Maximum Temperature> -p <steP Temperature> -u <Unit of process> -c (Print Curves)' % sys.argv[0]
+        print '%s -j <Coupling Factor> -b <Magnetic Field> -z <Size of Square Lattice> -i <Iterations> -s <Minimum Temperature> -e <Maximum Temperature> -p <steP Temperature> -u <Unit of process> -c (Print Curves)' % sys.argv[0]
         sys.exit(2)
     
     for opt, arg in opts:
         if opt == '-h':
-            print '%s -j <Coupling Factor> -b <Magnetic Field> -z <Size of Lattice> -i <Iterations> -s <Minimum Temperature> -e <Maximum Temperature> -p <steP Temperature> -c (Print Curves)' % sys.argv[0]
+            print '%s -j <Coupling Factor> -b <Magnetic Field> -z <Size of Square Lattice> -i <Iterations> -s <Minimum Temperature> -e <Maximum Temperature> -p <steP Temperature> -c (Print Curves)' % sys.argv[0]
             sys.exit()
         elif opt == '-c':
             Curves=True
@@ -142,19 +149,17 @@
 
     Trange=numpy.arange(Tmin,Tmax+Tstep,Tstep)
 
-    E=[]
-    M=[]
-
     def MetropolisStrip(T):
         # Indispensable d'utiliser copy : [:] ne fonctionne pas avec numpy !
         sigma=numpy.copy(sigmaIn)
         duration=Metropolis(sigma,J,B,T,Iterations)
-        ImageOutput(sigma,"Ising2D_Serial_%i_%1.1f_Final" % (Size,T))
+        ImageOutput(sigma,"Ising2D_MP_%i_%1.1f_Final" % (Size,T))
         print "CPU Time : %f" % (duration)
-        E,M=Energy(sigma,J),Magnetization(sigma,B)
+        indice=numpy.where(Trange==T)[0][0]
+        E,M=Energy(sigma,J,B),Magnetization(sigma,B)
         print "Total Energy at Temperature %f : %f" % (T,E)
         print "Total Magnetization at Temperature %f : %f" % (T,M)
-        return(T,E,M)
+        return([T,E,M])
 
     pool=Pool(processes=Procs)
     print "Start on %i processors..." % Procs
@@ -162,9 +167,15 @@
     # MAP: distribution of usecases T to be applied to MetropolisStrip 
     Results=pool.map(MetropolisStrip,Trange)
 
+    E=numpy.array(Results)[:,1]
+    M=numpy.array(Results)[:,2]
+    
     if Curves:
         DisplayCurves(Trange,E,M,J,B)
 
     # Save output
     numpy.savez("Ising2D_MP_%i_%.8i" % (Size,Iterations),(Trange,E,M))
-    
+      
+    # Estimate Critical temperature
+    print "The critical temperature on %ix%i lattice with J=%f, B=%f is %f " % (Size,Size,J,B,CriticalT(Trange,E))
+  
