Implemented nor algorithm.

[phys/heatbath.git] / sim-nor.hpp

diff --git a/sim-nor.hpp b/sim-nor.hpp
index aa305ff71abcaaec1e0e1a7f976e671ffe13d74d..7bb4a7c8fa8adb868c5c9755096ca5f48c07f158 100644 (file)
--- a/sim-nor.hpp
+++ b/sim-nor.hpp
@@ -7,6 +7,8 @@
 
 #include "latlib/neigh.h"
 
+#define EPSILONPHI 0.5
+
 class sim : public o815::sim {
 public:
   struct siteconf {
@@ -24,7 +26,12 @@ private:
 
   neigh *nb;
 
-  void updatePhi (const int& x);
+  int updatePhi (const int& x);
+
+  double M;
+
+  double SofPhi(const int& x, const complex<double>& phix);
+  double rhoPhi(const int& x, const complex<double>& phixCandidate);
 };
 
 
@@ -42,17 +49,41 @@ sim::sim(o815 *_O815) : o815::sim( _O815,
   nb = new neigh(2, _O815->comargs.lsize[0], _O815->comargs.lsize[1]);
 }
 
-void sim::updatePhi (const int& x) 
-{
-  const double r = gsl_rng_uniform(rangsl);
-  const double theta = gsl_rng_uniform(rangsl) * 2*M_PI;
-}
-
 void sim::_makeSweep() {  
+  for (int ix = 0; ix < LSIZE2; ix++)
+    updatePhi(ix);
 }
 
 void sim::_newParas() {
-  /* reset variables */
+  M = pow( (*O815->paraQ)["mass"], 2 ) + 4;
+  *log << "SIM: Resetting fields." << endl << flush;
+  memset(conf, 0, sizeof(siteconf)*LSIZE2);
+}
+
+int sim::updatePhi(const int& x) {
+  complex<double> phixCandidate = conf[x].phi +
+    complex<double> ( 2*EPSILONPHI*( 0.5 - gsl_rng_uniform(rangsl) ), 
+                     2*EPSILONPHI*( 0.5 - gsl_rng_uniform(rangsl) ) );
+  
+  if ( gsl_rng_uniform(rangsl) < rhoPhi(x, phixCandidate) ) {
+    conf[x].phi = phixCandidate;
+    return 1;
+  }
+  
+  return 0;
+}
+
+double sim::rhoPhi(const int& x, const complex<double>& phixCandidate) {
+  return exp( SofPhi(x, conf[x].phi) - SofPhi(x, phixCandidate) );
+}
+
+double sim::SofPhi(const int& x, const complex<double>& phix) {
+  double sofphi = M * norm(phix);
+  
+  for (int nu=0; nu<4; nu++)
+    sofphi -= 2 * real( conj(phix) * conf[ (*nb)[x*4+nu] ].phi );
+  
+  return sofphi;
 }
 
 #endif