CkTable.hpp

#pragma once
#include "Angular/Wigner369j.hpp"
#include <algorithm>
#include <iostream>
#include <vector>
 
namespace Angular {
 
//==============================================================================
// "Helper" functions
constexpr int twoj(int jindex) { return 2 * jindex + 1; }
constexpr int jindex(int twoj) { return (twoj - 1) / 2; }
constexpr int jindex_kappa(int ka) { return (ka > 0) ? ka - 1 : -ka - 1; }
 
//==============================================================================
class CkTable {
 
public:
  CkTable(const int in_max_twoj = 0) { fill(in_max_twoj); }
 
public:
 
  void fill(const int in_max_twoj);
 
  double get_Ckab_mutable(int k, int ka, int kb);
  double get_tildeCkab_mutable(int k, int ka, int kb);
  double get_3jkab_mutable(int k, int ka, int kb);
 
  double get_Ckab(int k, int ka, int kb) const;
  double get_tildeCkab(int k, int ka, int kb) const;
  double get_3jkab(int k, int ka, int kb) const;
 
  double operator()(int k, int ka, int kb) const { return get_Ckab(k, ka, kb); }
 
  double get_Lambdakab(int k, int ka, int kb) const;
 
  int max_tj() const { return twoj(m_max_jindex_sofar); }
  int max_k() const { return m_max_k_sofar; }
 
private:
  std::vector<std::vector<std::vector<double>>> m_3j_k_a_b = {};
  std::vector<std::vector<double>> m_Rjab_a_b = {}; // Sqrt([ja][jb])
  int m_max_jindex_sofar = -1;
  int m_max_k_sofar = -1;
};
 
} // namespace Angular