...

diff --git a/sim-1mr+.hpp b/sim-1mr+.hpp
index 2b141498f03d8e29e610f0939f6536da3ba1ea57..74fc18e604c40615fd7f7e86073d45a66ddc5aea 100644 (file)
--- a/sim-1mr+.hpp
+++ b/sim-1mr+.hpp
@@ -24,6 +24,9 @@ private:
   void _newParas();
   gsl_rng* rangsl;
   void updatePhi (const int& x);
   void _newParas();
   gsl_rng* rangsl;
   void updatePhi (const int& x);

+  void updatePhi_angle_relax (const int& x, const complex<double>& V);
+  void updatePhi_angle_rand (const int& x, const complex<double>& V);
+  void updatePhi_r (const int& x, const complex<double>& V); 

 };
 
 
 };
 
 


@@ -41,18 +44,8 @@ sim::sim(o815 *_O815) : o815::sim( _O815,
   nb = new neigh(2, _O815->comargs.lsize[0], _O815->comargs.lsize[1]);
 }
 
   nb = new neigh(2, _O815->comargs.lsize[0], _O815->comargs.lsize[1]);
 }
 

-void sim::updatePhi (const int& x) 
+void sim::updatePhi_r (const int& x, const complex<double>& V) 

 {
 {

-  complex<double> V=0;
-
-  for (int nu=0; nu<4; nu++)
-    V += conf[ (*nb)[x*4+nu] ].phi;
-
-  const double V2diff = pow(real(V), 2) - pow(imag(V), 2);
-  const double Vprod = real(V)*imag(V);
-  conf[x].phi = complex<double> ( + real(conf[x].phi) * V2diff + 2 * imag(conf[x].phi) * Vprod,
-                                 - imag(conf[x].phi) * V2diff + 2 * real(conf[x].phi) * Vprod ) / norm(V);
-  

   const double r = exp( - M * norm(conf[x].phi)
                        - 1./M * norm(V)
                        + 2 * real( conf[x].phi * conj(V) ) );
   const double r = exp( - M * norm(conf[x].phi)
                        - 1./M * norm(V)
                        + 2 * real( conf[x].phi * conj(V) ) );


@@ -62,6 +55,34 @@ void sim::updatePhi (const int& x)
     + V / M;
 }
 
     + V / M;
 }
 

+void sim::updatePhi_angle_relax (const int& x, const complex<double>& V)
+{
+  const double V2diff = pow(real(V), 2) - pow(imag(V), 2);
+  const double Vprod = real(V)*imag(V);
+  conf[x].phi = complex<double> ( + real(conf[x].phi) * V2diff + 2 * imag(conf[x].phi) * Vprod,
+                                 - imag(conf[x].phi) * V2diff + 2 * real(conf[x].phi) * Vprod ) / norm(V);
+}
+
+void sim::updatePhi_angle_rand (const int& x, const complex<double>& V)
+{
+  conf[x].phi = polar( abs(conf[x].phi), gsl_rng_uniform(rangsl) * 2*M_PI );
+}
+
+void sim::updatePhi (const int& x) 
+{
+  complex<double> V=0;
+
+  for (int nu=0; nu<4; nu++)
+    V += conf[ (*nb)[x*4+nu] ].phi;
+  
+  updatePhi_r(x,V);
+
+  if ( gsl_rng_uniform(rangsl) < 0.5 )
+    updatePhi_angle_relax(x,V);
+  else
+    updatePhi_angle_rand(x,V);
+}
+

 void sim::_makeSweep() {  
   for (int ichecker=0; ichecker<2; ichecker++)
     for (int it=0; it<O815->comargs.lsize[0]; it++)
 void sim::_makeSweep() {  
   for (int ichecker=0; ichecker<2; ichecker++)
     for (int it=0; it<O815->comargs.lsize[0]; it++)

diff --git a/sim-1mr-.hpp b/sim-1mr-.hpp
index 75c41d97183fd1f9d1f0185b81550cb41d15f93a..01701b959f5651b3cfe7d9d25a7ccbb55b484c03 100644 (file)
--- a/sim-1mr-.hpp
+++ b/sim-1mr-.hpp
@@ -45,7 +45,6 @@ void sim::updatePhi (const int& x)
 {
   const double theta = gsl_rng_uniform(rangsl) * 2*M_PI;
   complex<double> V=0;
 {
   const double theta = gsl_rng_uniform(rangsl) * 2*M_PI;
   complex<double> V=0;

-  const double oldarg = arg(conf[x].phi);

 
   for (int nu=0; nu<4; nu++)
     V += conf[ (*nb)[x*4+nu] ].phi;
 
   for (int nu=0; nu<4; nu++)
     V += conf[ (*nb)[x*4+nu] ].phi;


@@ -55,7 +54,7 @@ void sim::updatePhi (const int& x)
                        + 2 * real( conf[x].phi * conj(V) ) );
 
   conf[x].phi = sqrt(std::log( 1./(1-r) )) / sqrt(M)
                        + 2 * real( conf[x].phi * conj(V) ) );
 
   conf[x].phi = sqrt(std::log( 1./(1-r) )) / sqrt(M)

-    * polar(1.0, oldarg)
+    * polar( 1.0, arg(conf[x].phi) )

     + V / M;
 }
 
     + V / M;
 }
 

diff --git a/sim-r-.hpp b/sim-r-.hpp
index f16d97ec20ddd69f3936de15db7b8f7021e2a10f..f00b873eccdda6ad76461982443a7eecaf43ba77 100644 (file)
--- a/sim-r-.hpp
+++ b/sim-r-.hpp
@@ -44,15 +44,13 @@ sim::sim(o815 *_O815) : o815::sim( _O815,
 void sim::updatePhi (const int& x) 
 {
   const double r = gsl_rng_uniform(rangsl);
 void sim::updatePhi (const int& x) 
 {
   const double r = gsl_rng_uniform(rangsl);

-  const double theta = gsl_rng_uniform(rangsl) * 2*M_PI;
-  const double oldarg = arg(conf[x].phi);

   complex<double> V=0;
 
   for (int nu=0; nu<4; nu++)
     V += conf[ (*nb)[x*4+nu] ].phi;
 
   conf[x].phi = sqrt(std::log( 1./(1-r) )) / sqrt(M)
   complex<double> V=0;
 
   for (int nu=0; nu<4; nu++)
     V += conf[ (*nb)[x*4+nu] ].phi;
 
   conf[x].phi = sqrt(std::log( 1./(1-r) )) / sqrt(M)

-    * polar(1.0, oldarg)
+    * polar( 1.0, arg(conf[x].phi) )

     + V / M;
 }
 
     + V / M;
 }