mc_move_set.hpp

Go to the documentation of this file.

// Copyright (c) 2013-2018 Commissariat à l'énergie atomique et aux énergies alternatives (CEA)
// Copyright (c) 2013-2018 Centre national de la recherche scientifique (CNRS)
// Copyright (c) 2018-2022 Simons Foundation
//
// 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, Olivier Parcollet, Nils Wentzell

 
#pragma once
 
#include "./concepts.hpp"
#include "./mc_move.hpp"
#include "./random_generator.hpp"
 
#include <fmt/format.h>
#include <h5/h5.hpp>
 
#include <cassert>
#include <cstddef>
#include <complex>
#include <functional>
#include <map>
#include <stdexcept>
#include <string>
#include <utility>
#include <vector>
 
namespace triqs::mc_tools {

  template <DoubleOrComplex MCSignType> class move_set {
    public:
    move_set(random_generator &rng) : rng_(rng) {}

    move_set(move_set const &) = delete;

    move_set &operator=(move_set const &) = delete;

    move_set(move_set &&) = default;

    move_set &operator=(move_set &&) = default;

    template <typename T>
      requires(MCMove<T, MCSignType>)
    void add(T &&m, std::string name, double weight) {
      if (weight < 0.0) throw std::runtime_error(fmt::format("Error in move_set: Weight of move {} is negative: {} < 0", name, weight));
      moves_.emplace_back(std::forward<T>(m));
      weights_.emplace_back(weight);
      names_.emplace_back(std::move(name));
      initialize();
    }

    double attempt();

    MCSignType accept() { return moves_[current_].accept() * attempt_sign_; }

    void reject() { moves_[current_].reject(); }

    void clear_statistics();

    void collect_statistics(mpi::communicator const &c);

    void calibrate(mpi::communicator const &c);

    [[nodiscard]] std::map<std::string, double> get_acceptance_rates() const;

    [[nodiscard]] std::string get_statistics(std::string const &prefix = "") const;

    [[nodiscard]] std::string get_timings(std::string const &prefix = "") const;

    [[nodiscard]] double total_duration() const {
      double total = 0.0;
      for (auto const &m : moves_) total += m.duration();
      return total;
    }

    [[nodiscard]] auto current() const { return current_; }

    [[nodiscard]] auto const &probabilities() const { return probs_; }

    [[nodiscard]] auto attempt_sign() const { return attempt_sign_; }

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

    friend void h5_write(h5::group g, std::string const &name, move_set const &ms) {
      auto gr = g.create_group(name);
      h5::write_hdf5_format(gr, ms);
      for (std::size_t i = 0; i < ms.moves_.size(); ++i) h5::write(gr, ms.names_[i], ms.moves_[i]);
    }

    friend void h5_read(h5::group g, std::string const &name, move_set &ms) {
      auto gr = g.open_group(name);
      h5::assert_hdf5_format(gr, ms);
      for (std::size_t i = 0; i < ms.moves_.size(); ++i) h5::read(gr, ms.names_[i], ms.moves_[i]);
    }
 
    private:
    // Normalize the weights and calculate the accumulated probabilities.
    void initialize();
 
    // Check the acceptance ratio, handle possible infinites, set the attemped sign and return the absolute value.
    double check_ratio(MCSignType ratio);
 
    private:
    std::vector<move<MCSignType>> moves_;
    std::vector<std::string> names_;
    std::size_t current_{0};
    std::reference_wrapper<random_generator> rng_;
    std::vector<double> weights_;
    std::vector<double> probs_;
    std::vector<double> acc_probs_;
    MCSignType attempt_sign_{1.0};
  };
 
  // Explicit template instantiation declarations.
  extern template class move_set<double>;
  extern template class move_set<std::complex<double>>;
 
} // namespace triqs::mc_tools