diff options
Diffstat (limited to 'lib/cluster.h')
-rw-r--r-- | lib/cluster.h | 177 |
1 files changed, 1 insertions, 176 deletions
diff --git a/lib/cluster.h b/lib/cluster.h index 29dd0cb..a20912e 100644 --- a/lib/cluster.h +++ b/lib/cluster.h @@ -20,6 +20,7 @@ #include "graph.h" #include "tree.h" #include "measurement.h" +#include "vector.h" #include "orthogonal.h" #include "dihedral.h" #include "dihinf.h" @@ -95,182 +96,6 @@ class state_t { } }; -template <q_t q, class T> -struct vector_t { T *x; }; - -template <q_t q, class T> -void init(vector_t <q, T> *ptr) { - ptr->x = (T *)calloc(q, sizeof(T)); - - ptr->x[0] = (T)1; -} - -template <q_t q, class T> -vector_t <q, T> copy (vector_t <q, T> v) { - vector_t <q, T> v_copy; - - v_copy.x = (T *)calloc(q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - v_copy.x[i] = v.x[i]; - } - - return v_copy; -} - -template <q_t q, class T> -void add (vector_t <q, T> v1, vector_t <q, T> v2) { - for (q_t i = 0; i < q; i++) { - v1.x[i] += v2.x[i]; - } -} - -template <q_t q, class T> -void subtract (vector_t <q, T> v1, vector_t <q, T> v2) { - for (q_t i = 0; i < q; i++) { - v1.x[i] -= v2.x[i]; - } -} - -template <q_t q, class T> -vector_t <q, T> scalar_multiple(v_t a, vector_t <q, T> v) { - vector_t <q, T> multiple; - multiple.x = (T *)malloc(q * sizeof(T)); - for (q_t i = 0; i < q; i++) { - multiple.x[i] = a * v.x[i]; - } - - return multiple; -} - -template <q_t q, class T> -T dot(vector_t <q, T> v1, vector_t <q, T> v2) { - T prod = 0; - - for (q_t i = 0; i < q; i++) { - prod += v1.x[i] * v2.x[i]; - } - - return prod; -} - -template <q_t q, class T> -void free_spin (vector_t <q, T> v) { - free(v.x); -} - -template <q_t q, class T> -struct orthogonal_t { T *x; }; - -template <q_t q, class T> -void init(orthogonal_t <q, T> *ptr) { - ptr->x = (T *)calloc(q * q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - ptr->x[q * i + i] = (T)1; - } -} - -template <q_t q, class T> -orthogonal_t <q, T> copy (orthogonal_t <q, T> m) { - orthogonal_t <q, T> m_copy; - m_copy.x = (T *)calloc(q * q, sizeof(T)); - - for (q_t i = 0; i < q * q; i++) { - m_copy.x[i] = m.x[i]; - } - - return m_copy; -} - -template <q_t q, class T> -void free_spin (orthogonal_t <q, T> m) { - free(m.x); -} - -template <q_t q, class T> -vector_t <q, T> act (orthogonal_t <q, T> m, vector_t <q, T> v) { - vector_t <q, T> v_rot; - v_rot.x = (T *)calloc(q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - for (q_t j = 0; j < q; j++) { - v_rot.x[i] += m.x[q * i + j] * v.x[j]; - } - } - - return v_rot; -} - - -template <q_t q, class T> -orthogonal_t <q, T> act (orthogonal_t <q, T> m1, orthogonal_t <q, T> m2) { - orthogonal_t <q, T> m2_rot; - m2_rot.x = (T *)calloc(q * q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - for (q_t j = 0; j < q; j++) { - for (q_t k = 0; k < q; k++) { - m2_rot.x[i * q + j] += m1.x[i * q + j] * m2.x[j * q + k]; - } - } - } - - return m2_rot; -} - -template <q_t q, class T> -vector_t <q, T> act_inverse (orthogonal_t <q, T> m, vector_t <q, T> v) { - vector_t <q, T> v_rot; - v_rot.x = (T *)calloc(q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - for (q_t j = 0; j < q; j++) { - v_rot.x[i] += m.x[q * j + i] * v.x[j]; - } - } - - return v_rot; -} - -template <q_t q, class T> -orthogonal_t <q, T> act_inverse (orthogonal_t <q, T> m1, orthogonal_t <q, T> m2) { - orthogonal_t <q, T> m2_rot; - m2_rot.x = (T *)calloc(q * q, sizeof(T)); - - for (q_t i = 0; i < q; i++) { - for (q_t j = 0; j < q; j++) { - for (q_t k = 0; k < q; k++) { - m2_rot.x[i * q + j] += m1.x[j * q + i] * m2.x[j * q + k]; - } - } - } - - return m2_rot; -} - -template <q_t q> -void generate_rotation (gsl_rng *r, orthogonal_t <q, double> *ptr) { - double *v = (double *)malloc(q * sizeof(double)); - double v2 = 0; - - for (q_t i = 0; i < q; i++) { - v[i] = gsl_ran_ugaussian(r); - v2 += v[i] * v[i]; - } - - ptr->x = (double *)calloc(q * q, sizeof(double)); - - for (q_t i = 0; i < q; i++) { - ptr->x[q * i + i] = 1.0; - for (q_t j = 0; j < q; j++) { - ptr->x[q * i + j] -= 2 * v[i] * v[j] / v2; - } - } - - free(v); -} - template <class R_t, class X_t> v_t flip_cluster(state_t <R_t, X_t> *state, v_t v0, R_t r, gsl_rng *rand) { v_t nv = 0; |