summaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--lib/wolff.h3
-rw-r--r--lib/wolff_tools.c2
-rw-r--r--src/wolff_vector.c31
3 files changed, 29 insertions, 7 deletions
diff --git a/lib/wolff.h b/lib/wolff.h
index 29286ba..18ae6dd 100644
--- a/lib/wolff.h
+++ b/lib/wolff.h
@@ -42,7 +42,8 @@ typedef struct {
double *spins;
double T;
double (*J)(double);
- double (*H)(double *);
+ double (*H)(q_t, double *, double *);
+ double *H_info;
double *R;
double E;
double *M;
diff --git a/lib/wolff_tools.c b/lib/wolff_tools.c
index e63623c..7a52546 100644
--- a/lib/wolff_tools.c
+++ b/lib/wolff_tools.c
@@ -131,7 +131,7 @@ v_t flip_cluster_vector(vector_state_t *s, v_t v0, double *rot, gsl_rng *r) {
rs_old = vector_rotate_inverse(s->n, s->R, sn);
rs_new = vector_rotate_inverse(s->n, R_tmp, sn);
}
- double dE = s->H(rs_old) - s->H(rs_new);
+ double dE = s->H(s->n, s->H_info, rs_old) - s->H(s->n, s->H_info, rs_new);
prob = 1.0 - exp(-dE / s->T);
vector_subtract(s->n, s->M, rs_old);
vector_add(s->n, s->M, rs_new);
diff --git a/src/wolff_vector.c b/src/wolff_vector.c
index c6ab695..904fcfd 100644
--- a/src/wolff_vector.c
+++ b/src/wolff_vector.c
@@ -5,12 +5,32 @@ double identity(double x) {
return x;
}
-double zero(double *x) {
+double zero(q_t n, double *H, double *x) {
return 0.0;
}
-double test(double *x) {
- return -x[1];
+double dot(q_t n, double *H, double *x) {
+ double total = 0;
+ for (q_t i = 0; i < n; i++) {
+ total += H[i] * x[i];
+ }
+ return total;
+}
+
+double theta(double x, double y) {
+ double val = atan(y / x);
+
+ if (x < 0.0 && y > 0.0) {
+ return M_PI + val;
+ } else if ( x < 0.0 && y < 0.0 ) {
+ return - M_PI + val;
+ } else {
+ return val;
+ }
+}
+
+double modulated(q_t n, double *H_info, double *x) {
+ return H_info[0] * cos(H_info[1] * theta(x[0], x[1]));
}
int main(int argc, char *argv[]) {
@@ -83,8 +103,9 @@ int main(int argc, char *argv[]) {
s->spins[q * i] = 1.0;
}
+ s->H_info = H;
s->T = T;
- s->H = test;
+ s->H = dot;
s->J = identity;
s->R = (double *)calloc(q * q, sizeof(double));
@@ -95,7 +116,7 @@ int main(int argc, char *argv[]) {
s->M = (double *)calloc(q, sizeof(double));
s->M[0] = 1.0 * (double)h->nv;
- s->E = - ((double)h->ne) * s->J(1.0) - s->H(s->M);
+ s->E = - ((double)h->ne) * s->J(1.0) - s->H(s->n, s->H_info, s->M);
printf("%g %g\n", s->E, s->M[0]);