diff options
| -rw-r--r-- | space_wolff.hpp | 2 | ||||
| -rw-r--r-- | transformation.hpp | 56 |
2 files changed, 46 insertions, 12 deletions
diff --git a/space_wolff.hpp b/space_wolff.hpp index 9d529a3..141396a 100644 --- a/space_wolff.hpp +++ b/space_wolff.hpp @@ -62,7 +62,7 @@ public: for (Spin<U, D, S>* s : t->toConsider()) { double ΔE = t->ΔE(s); if (ΔE > 0 && 1.0 - exp(-ΔE / T) > rng.uniform<double>(0, 1)) { - queue.push(t->createNew(s)); + queue.push(t->createNew(T, visited, s, rng)); } } diff --git a/transformation.hpp b/transformation.hpp index 0731ec9..fc1833e 100644 --- a/transformation.hpp +++ b/transformation.hpp @@ -12,7 +12,7 @@ public: virtual std::set<Spin<U, D, S>*> current() const { return {NULL}; } virtual std::set<Spin<U, D, S>*> toConsider() const { return {}; } virtual double ΔE(const Spin<U, D, S>*) const { return 0; } - virtual Transformation* createNew(Spin<U, D, S>*) const { return new Transformation(); } + virtual Transformation* createNew(double T, const std::set<Spin<U, D, S>*>&, Spin<U, D, S>*, Rng& rng) const { return new Transformation(); } virtual void apply(){}; }; @@ -56,7 +56,7 @@ public: } } - Transformation<U, D, R, S>* createNew(Spin<U, D, S>* s) const override; + Transformation<U, D, R, S>* createNew(double, const std::set<Spin<U, D, S>*>&, Spin<U, D, S>* s, Rng&) const override; void apply() override { M.dict.remove(sOld); @@ -112,8 +112,8 @@ public: } } - Transformation<U, D, R, S>* createNew(Spin<double, D, S>* s) const; - Transformation<U, D, R, S>* createNew(Spin<signed, D, S>* s) const; + Transformation<U, D, R, S>* createNew(double T, const std::set<Spin<U, D, S>*>& v, Spin<double, D, S>* s, Rng&) const; + Transformation<U, D, R, S>* createNew(double T, const std::set<Spin<U, D, S>*>& v, Spin<signed, D, S>* s, Rng&) const; void apply() override { M.dict.remove(s1Old); @@ -150,11 +150,11 @@ public: return M.B(s0s_new) - M.B(s0s_old); } - Transformation<double, D, R, S>* createNew(Spin<double, D, S>* s) const { + Transformation<double, D, R, S>* createNew(double T, const std::set<Spin<U, D, S>*>& v, Spin<double, D, S>* s, Rng&) const { return new SpinFlip<U, D, R, S>(M, r, s); } - Transformation<signed, D, R, S>* createNew(Spin<signed, D, S>* s) const { + Transformation<signed, D, R, S>* createNew(double, const std::set<Spin<U, D, S>*>&, Spin<signed, D, S>* s, Rng&) const { Vector<signed, 2> v = r.act(*s).x; std::set<Spin<U, D, S>*> on_site = M.dict.at(v); if (on_site.empty()) { @@ -168,25 +168,59 @@ public: }; template <class U, int D, class R, class S> -Transformation<U, D, R, S>* SpinFlip<U, D, R, S>::createNew(Spin<U, D, S>* s) const { +Transformation<U, D, R, S>* SpinFlip<U, D, R, S>::createNew(double T, const std::set<Spin<U, D, S>*>& v, Spin<U, D, S>* s, Rng& rng) const { if (s == NULL) { return new FieldFlip<U, D, R, S>(M, r); } else { - return new SpinFlip(M, r, s); + SpinFlip<U, D, R, S>* t = new SpinFlip<U, D, R, S>(M, r, s); + Spin<U, D, S>* sMax = NULL; + double ΔEMax = 0.0; + for (Spin<U, D, S>* ss : t->toConsider()) { + if (ss != NULL && ss != s && !v.contains(ss)) { + double ΔE = t->ΔE(ss); + if (ΔE > ΔEMax) { + sMax = ss; + ΔEMax = ΔE; + } + } + } + if (sMax == NULL || 1.0 - exp(-ΔEMax / T) < rng.uniform<double>(0, 1)) { + return t; + } else { + delete t; + return new PairFlip<U, D, R, S>(M, r, s, sMax); + } } } template <class U, int D, class R, class S> -Transformation<U, D, R, S>* PairFlip<U, D, R, S>::createNew(Spin<double, D, S>* s) const { +Transformation<U, D, R, S>* PairFlip<U, D, R, S>::createNew(double T, const std::set<Spin<U, D, S>*>& v, Spin<double, D, S>* s, Rng& rng) const { if (s == NULL) { return new FieldFlip<U, D, R, S>(M, r); } else { - return new SpinFlip<U, D, R, S>(M, r, s); + SpinFlip<U, D, R, S>* t = new SpinFlip<U, D, R, S>(M, r, s); + Spin<U, D, S>* sMax = NULL; + double ΔEMax = 0.0; + for (Spin<U, D, S>* ss : t->toConsider()) { + if (ss != NULL && ss != s && !v.contains(ss)) { + double ΔE = t->ΔE(ss); + if (ΔE > ΔEMax) { + sMax = ss; + ΔEMax = ΔE; + } + } + } + if (sMax == NULL || 1.0 - exp(-ΔEMax / T) < rng.uniform<double>(0, 1)) { + return t; + } else { + delete t; + return new PairFlip<U, D, R, S>(M, r, s, sMax); + } } } template <class U, int D, class R, class S> -Transformation<U, D, R, S>* PairFlip<U, D, R, S>::createNew(Spin<signed, D, S>* s) const { +Transformation<U, D, R, S>* PairFlip<U, D, R, S>::createNew(double, const std::set<Spin<U, D, S>*>&, Spin<signed, D, S>* s, Rng&) const { if (s == NULL) { return new FieldFlip<U, D, R, S>(M, r); } else { |