hilbert_space.hpp

Go to the documentation of this file.

// Copyright (c) 2015-2018 Commissariat à l'énergie atomique et aux énergies alternatives (CEA)
// Copyright (c) 2015-2018 Centre national de la recherche scientifique (CNRS)
// Copyright (c) 2018-2022 Simons Foundation
// Copyright (c) 2016 Igor Krivenko
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You may obtain a copy of the License at
//     https://www.gnu.org/licenses/gpl-3.0.txt
//
// Authors: Michel Ferrero, Igor Krivenko, Olivier Parcollet, Nils Wentzell

 
#pragma once
 
#include "./fundamental_operator_set.hpp"
#include "../utility/exceptions.hpp"
 
#include <boost/container/flat_map.hpp>
#include <h5/h5.hpp>
 
#include <cstdint>
#include <set>
#include <string>
#include <vector>
 
namespace triqs::hilbert_space {


  using fock_state_t = uint64_t;

  class hilbert_space {
    public:
    hilbert_space() : dim_(0) {}

    hilbert_space(fundamental_operator_set const &fops) : dim_(static_cast<int>(1ull << fops.size())) {}

    [[nodiscard]] auto size() const { return dim_; }

    bool operator==(hilbert_space const &hs) const { return dim_ == hs.dim_; }

    [[nodiscard]] bool has_state(fock_state_t f) const { return f < dim_; }

    [[nodiscard]] auto get_state_index(fock_state_t f) const {
      if (f >= dim_) TRIQS_RUNTIME_ERROR << "Error in hilbert_space::get_state_index: Fock state is not in H: f = " << f;
      return static_cast<int>(f);
    }

    [[nodiscard]] fock_state_t get_fock_state(int f) const {
      if (f >= dim_) TRIQS_RUNTIME_ERROR << "Error in hilbert_space::get_fock_state: Index is >= dim(H): " << f;
      return f;
    }

    [[nodiscard]] fock_state_t get_fock_state(fundamental_operator_set const &fops,
                                              std::set<fundamental_operator_set::indices_t> const &alphas) const {
      fock_state_t f = 0;
      for (auto const &alpha : alphas) f += 1 << fops[alpha];
      return f;
    }

    [[nodiscard]] static std::string hdf5_format() { return "hilbert_space"; }

    friend void h5_write(h5::group g, std::string const &name, hilbert_space const &hs) {
      auto gr = g.create_group(name);
      h5::write(gr, "dim", hs.dim_);
    }

    friend void h5_read(h5::group g, std::string const &name, hilbert_space &hs) {
      auto gr = g.open_group(name);
      h5::read(gr, "dim", hs.dim_);
    }
 
    private:
    int dim_;
  };

  class sub_hilbert_space {
    public:
    sub_hilbert_space(int m = -1) : m_(m) {}
 
#ifdef TRIQS_WORKAROUND_INTEL_COMPILER_BUGS
    // Workaround needed for icc, checked with 17.0.1 20161005)
    sub_hilbert_space(sub_hilbert_space const &) = default;
    sub_hilbert_space(sub_hilbert_space &&)      = default;
    sub_hilbert_space &operator=(sub_hilbert_space const &x) {
      index         = x.index;
      fock_states   = x.fock_states;
      fock_to_index = x.fock_to_index;
      return *this;
    }
    sub_hilbert_space &operator=(sub_hilbert_space &&) = default;
#endif

    void add_fock_state(fock_state_t f) {
      int const M = static_cast<int>(fock_states_.size());
      fock_states_.push_back(f);
      fock_to_index_.insert(std::make_pair(f, M));
    }

    [[nodiscard]] auto size() const { return static_cast<int>(fock_states_.size()); }

    bool operator==(sub_hilbert_space const &hs) const { return m_ == hs.m_ && fock_states_ == hs.fock_states_; }

    [[nodiscard]] auto get_state_index(fock_state_t f) const { return fock_to_index_.find(f)->second; }

    [[nodiscard]] bool has_state(fock_state_t f) const { return fock_to_index_.count(f) == 1; }

    [[nodiscard]] auto get_fock_state(int i) const { return fock_states_[i]; }

    [[nodiscard]] auto const &get_all_fock_states() const { return fock_states_; }

    [[nodiscard]] auto get_index() const { return m_; };

    void set_index(int m) { m_ = m; }

    [[nodiscard]] static std::string hdf5_format() { return "sub_hilbert_space"; }

    friend void h5_write(h5::group g, std::string const &name, sub_hilbert_space const &hs) {
      auto gr = g.create_group(name);
      h5::write(gr, "index", hs.m_);
      h5::write(gr, "fock_states", hs.fock_states_);
    }

    friend void h5_read(h5::group g, std::string const &name, sub_hilbert_space &hs) {
      auto gr = g.open_group(name);
      h5::read(gr, "index", hs.m_);
      h5::read(gr, "fock_states", hs.fock_states_);
      hs.fock_to_index_.clear();
      for (auto f : hs.fock_states_) hs.fock_to_index_.insert(std::make_pair(f, static_cast<int>(hs.fock_to_index_.size())));
    }
 
    private:
    int m_;
    std::vector<fock_state_t> fock_states_;
    boost::container::flat_map<fock_state_t, int> fock_to_index_;
  };

 
} // namespace triqs::hilbert_space