From d3b3e39a525d0c3aa9663f824ba96e0c429a8313 Mon Sep 17 00:00:00 2001
From: Jaron Kent-Dobias <jaron@kent-dobias.com>
Date: Thu, 26 Jul 2018 00:32:38 -0400
Subject: partially class-ified, ising and On work but potts and height do not

---
 lib/orthogonal.h | 27 ++++++++++-----------------
 1 file changed, 10 insertions(+), 17 deletions(-)

(limited to 'lib/orthogonal.h')

diff --git a/lib/orthogonal.h b/lib/orthogonal.h
index 2a391ae..c0958f2 100644
--- a/lib/orthogonal.h
+++ b/lib/orthogonal.h
@@ -53,20 +53,19 @@ void free_spin (orthogonal_t <q, T> m) {
 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));
 
   if (m.is_reflection) {
     double prod = 0;
     for (q_t i = 0; i < q; i++) {
-      prod += v.x[i] * m.x[i];
+      prod += v[i] * m.x[i];
     }
     for (q_t i = 0; i < q; i++) {
-      v_rot.x[i] = v.x[i] - 2 * prod * m.x[i];
+      v_rot[i] = v[i] - 2 * prod * m.x[i];
     }
   } else {
     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];
+        v_rot[i] += m.x[q * i + j] * v[j];
       }
     }
   }
@@ -112,11 +111,10 @@ vector_t <q, T> act_inverse (orthogonal_t <q, T> m, vector_t <q, T> v) {
     return act(m, v); // reflections are their own inverse
   } else {
     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];
+        v_rot[i] += m.x[q * j + i] * v[j];
       }
     }
 
@@ -176,33 +174,32 @@ orthogonal_t <q, double> generate_rotation_perturbation (gsl_rng *r, vector_t <q
 
   if (n > 1) {
     unsigned int rotation = gsl_rng_uniform_int(r, n);
-    v.x = (double *)malloc(q * sizeof(double));
 
     double cosr = cos(2 * M_PI * rotation / (double)n / 2.0);
     double sinr = sin(2 * M_PI * rotation / (double)n / 2.0);
 
-    v.x[0] = v0.x[0] * cosr - v0.x[1] * sinr;
-    v.x[1] = v0.x[1] * cosr + v0.x[0] * sinr;
+    v[0] = v0[0] * cosr - v0[1] * sinr;
+    v[1] = v0[1] * cosr + v0[0] * sinr;
 
     for (q_t i = 2; i < q; i++) {
-      v.x[i] = v0.x[i];
+      v[i] = v0[i];
     }
   } else {
-    v.x = v0.x;
+    v = v0;
   }
 
   double m_dot_v = 0;
 
   for (q_t i = 0; i < q; i++) {
     m.x[i] = gsl_ran_ugaussian(r);
-    m_dot_v += m.x[i] * v.x[i];
+    m_dot_v += m.x[i] * v[i];
   }
 
   double v2 = 0;
   double factor = epsilon * gsl_ran_ugaussian(r);
 
   for (q_t i = 0; i < q; i++) {
-    m.x[i] = m.x[i] - m_dot_v * v.x[i] + factor * v.x[i];
+    m.x[i] = m.x[i] - m_dot_v * v[i] + factor * v[i];
     v2 += pow(m.x[i], 2);
   }
 
@@ -212,10 +209,6 @@ orthogonal_t <q, double> generate_rotation_perturbation (gsl_rng *r, vector_t <q
     m.x[i] /= mag_v;
   }
 
-  if (n > 1) {
-    free(v.x);
-  }
-
   return m;
 }
 
-- 
cgit v1.2.3-70-g09d2