From 0c6efee6b0f099d71ab4daa9f79715e9dc2dc902 Mon Sep 17 00:00:00 2001
From: Jaron Kent-Dobias <jaron@kent-dobias.com>
Date: Tue, 2 Apr 2024 20:20:48 +0200
Subject: Simplied implementation, specific for the quadractic models.

---
 tensor.hpp | 108 -------------------------------------------------------------
 1 file changed, 108 deletions(-)
 delete mode 100644 tensor.hpp

(limited to 'tensor.hpp')

diff --git a/tensor.hpp b/tensor.hpp
deleted file mode 100644
index 7f98d02..0000000
--- a/tensor.hpp
+++ /dev/null
@@ -1,108 +0,0 @@
-#pragma once
-
-#include <array>
-#include <functional>
-
-#include <eigen3/unsupported/Eigen/CXX11/Tensor>
-
-#include "types.hpp"
-#include "factorial.hpp"
-
-using Vector = Eigen::Matrix<Real, Eigen::Dynamic, 1>;
-
-using Matrix = Eigen::Matrix<Real, Eigen::Dynamic, Eigen::Dynamic>;
-
-template <int p>
-using Tensor = Eigen::Tensor<Real, p>;
-
-template <int p>
-using constTensor = Eigen::Tensor<const Real, p>;
-
-template <int p, std::size_t... Indices>
-Tensor<p> initializeJHelper(unsigned N, unsigned M, std::index_sequence<Indices...>) {
-  std::array<unsigned, p> Ns;
-  std::fill_n(Ns.begin(), p, N);
-  Ns[p-1] = M;
-  return Tensor<p>(std::get<Indices>(Ns)...);
-}
-
-template <int p>
-Tensor<p> initializeJ(unsigned N, unsigned M) {
-  return initializeJHelper<p>(N, M, std::make_index_sequence<p>());
-}
-
-template <int p, std::size_t... Indices>
-void setJHelper(Tensor<p>& J, const std::array<unsigned, p>& ind, Real val, std::index_sequence<Indices...>) {
-  J(std::get<Indices>(ind)...) = val;
-}
-
-template <int p>
-void setJ(Tensor<p>& J, std::array<unsigned, p> ind, Real val) {
-  setJHelper<p>(J, ind, val, std::make_index_sequence<p>());
-}
-
-template <int p, std::size_t... Indices>
-Real getJHelper(const Tensor<p>& J, const std::array<unsigned, p>& ind, std::index_sequence<Indices...>) {
-  return J(std::get<Indices>(ind)...);
-}
-
-template <int p>
-Real getJ(const Tensor<p>& J, const std::array<unsigned, p>& ind) {
-  return getJHelper<p>(J, ind, std::make_index_sequence<p>());
-}
-
-template <int p>
-void iterateOverHelper(Tensor<p>& J,
-    std::function<void(Tensor<p>&, std::array<unsigned, p>)>& f,
-    unsigned l, std::array<unsigned, p> is) {
-  if (l == 0) {
-    f(J, is);
-  } else {
-    for (unsigned i = 0; i < J.dimension(p - l); i++) {
-      std::array<unsigned, p> js = is;
-      js[p - l] = i;
-      iterateOverHelper<p>(J, f, l - 1, js);
-    }
-  }
-}
-
-template <int p>
-void iterateOver(Tensor<p>& J, std::function<void(Tensor<p>&, std::array<unsigned, p>)>& f) {
-  std::array<unsigned, p> is;
-  iterateOverHelper<p>(J, f, p, is);
-}
-
-template <int p, class Distribution, class Generator>
-Tensor<p> generateCouplings(unsigned N, unsigned M, Distribution d, Generator& r) {
-  Tensor<p> J = initializeJ<p>(N, M);
-
-  std::function<void(Tensor<p>&, std::array<unsigned, p>)> setRandom =
-    [&d, &r] (Tensor<p>& JJ, std::array<unsigned, p> ind) {
-      setJ<p>(JJ, ind, d(r));
-    };
-
-  iterateOver<p>(J, setRandom);
-
-  return J;
-}
-
-template <int p, class Generator>
-Tensor<p> generateRealPSpinCouplings(unsigned N, unsigned M, Generator& r) {
-  Real σp = sqrt(factorial(p-1) / ((Real)pow(N, p - 1)));
-
-  return generateCouplings<p>(N, M, std::normal_distribution<Real>(0, σp), r);
-}
-
-Tensor<3> contractDown(const Tensor<3>& J, const Vector& z) {
-  return J;
-}
-
-const std::array<Eigen::IndexPair<int>, 1> ip00 = {Eigen::IndexPair<int>(0, 0)};
-const std::array<Eigen::IndexPair<int>, 1> ip20 = {Eigen::IndexPair<int>(2, 0)};
-
-template <int r>
-Tensor<3> contractDown(const Tensor<r>& J, const Vector& z) {
-  Tensor<1> zT = Eigen::TensorMap<constTensor<1>>(z.data(), {z.size()});
-  Tensor<r - 1> Jz = J.contract(zT, ip00);
-  return contractDown(Jz, z);
-}
-- 
cgit v1.2.3-70-g09d2