summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
authorJaron Kent-Dobias <jaron@kent-dobias.com>2017-06-22 23:36:35 -0400
committerJaron Kent-Dobias <jaron@kent-dobias.com>2017-06-22 23:36:35 -0400
commitc902f022795157b28ef01cffcc499753e2212c28 (patch)
treed7bc1cf982d36a67894efecec45ee34e2a1100d2
parentefcad567ea9414b807017cbaadd7013e9310ca3e (diff)
downloadc++-c902f022795157b28ef01cffcc499753e2212c28.tar.gz
c++-c902f022795157b28ef01cffcc499753e2212c28.tar.bz2
c++-c902f022795157b28ef01cffcc499753e2212c28.zip
added correct support for negative external fields
-rw-r--r--lib/wolff.h2
-rw-r--r--lib/wolff_tools.c8
-rw-r--r--src/wolff.c4
3 files changed, 10 insertions, 4 deletions
diff --git a/lib/wolff.h b/lib/wolff.h
index 4270f33..428ce85 100644
--- a/lib/wolff.h
+++ b/lib/wolff.h
@@ -31,6 +31,8 @@ typedef struct {
int32_t dHb;
} cluster_t;
+int32_t sign(double x);
+
cluster_t *flip_cluster(const graph_t *g, const double *ps, bool *x,
gsl_rng *r);
diff --git a/lib/wolff_tools.c b/lib/wolff_tools.c
index b3758cc..7886bab 100644
--- a/lib/wolff_tools.c
+++ b/lib/wolff_tools.c
@@ -2,6 +2,10 @@
#include "queue.h"
#include "wolff.h"
+int32_t sign(double x) {
+ return x > 0 ? 1 : -1;
+}
+
graph_t *graph_add_ext(const graph_t *g) {
graph_t *h = (graph_t *)calloc(1, sizeof(graph_t));
@@ -131,8 +135,8 @@ uint32_t wolff_step(double T, double H, ising_state_t *s, gsl_rng *r,
cluster_t *c = flip_cluster(s->g, ps, s->spins, r);
- s->M += -2 * c->dHb;
- s->H += 2 * (c->dJb + H * c->dHb);
+ s->M += - sign(H) * 2 * c->dHb;
+ s->H += 2 * (c->dJb + sign (H) * H * c->dHb);
uint32_t n_flips = c->nv;
diff --git a/src/wolff.c b/src/wolff.c
index a8266ff..1176fde 100644
--- a/src/wolff.c
+++ b/src/wolff.c
@@ -84,8 +84,8 @@ int main(int argc, char *argv[]) {
s->g = graph_add_ext(h);
s->spins = (bool *)calloc(h->nv + 1, sizeof(bool));
- s->M = h->nv;
- s->H = -(1.0 * h->ne) - H * h->nv;
+ s->M = sign(H) * h->nv;
+ s->H = -(1.0 * h->ne) - sign (H) * H * h->nv;
double *bond_probs = (double *)malloc(2 * sizeof(double));
bond_probs[0] = 1 - exp(-2 / T);