Nyaan's Library

This documentation is automatically generated by online-judge-tools/verification-helper

View on GitHub

:heavy_check_mark: Multipoint Simulated Annealing
(marathon/sa-manager.hpp)

焼きなまし法

焼きなまし法のライブラリを自分が使いやすいように改良したもの。

テンプレート

// 状態を持つ型
using state_t = vector<int>;
// スコアを持つ型
using score_t = double;
using pair_t = pair<state_t, score_t>;

// 焼きなましの初期解を入れる関数
pair_t get_init_state() {
  state_t ans;
  score_t score;
  // ans, score に初期解を入れる。
  // 多点焼きなましを行う場合は初期解に多様性が生まれるようにする。

  return {ans, score};
}

// 状態を更新する関数
void update_state(state_t& ans, score_t& score, double th) {
  // 遷移を作る
  // 基本的には以下のように値を複製して遷移先を作っていい
  // (差分更新したい場合は適宜頑張る)
  state_t nxt{ans};
  score_t nxt_score{score};

  // スコアの伸びを変数 diff に格納して値を元に遷移を採用するか決める。
  double diff = -1;
  if (diff > th) {
    // 遷移を採用して ans, score を更新 or 差分更新する。

  } else {
    // (必要があれば)元に戻す。
  }
}

pair_t sa(Timer& t) {
  SA_manager<state_t, score_t> sa(t, 1000, 10, 1900, 1);
  return sa.run<get_init_state,update_state>();
}

注意点

お気持ち

Depends on

Verified with

Code

#pragma once

#include <cassert>
#include <iostream>
#include <vector>
using namespace std;

#include "../misc/rng.hpp"
#include "../misc/timer.hpp"
#include "log_table.hpp"

template <typename state_t, typename score_t>
struct SA_manager {
 private:
  log_table rand_log;
  Timer timer;
  double start_time, end_time, inv_time, cur_time;
  double start_temp, diff_temp, cur_temp;
  int state_max, state_size;

  // 時間をセットする関数
  void set_time() { cur_time = timer.elapsed(); }

  // 時間を色々初期化する関数
  void init_time() {
    set_time();
    start_time = cur_time;
    inv_time = 1. / (end_time - start_time);
  }

  // 気温を設定して cur_temp に代入する関数
  void set_temp() {
    double passed_time = max(0., cur_time - start_time);
    cur_temp = max(0., start_temp - passed_time * inv_time * diff_temp);
  }

  // 焼きなましで持つ状態の個数の max を代入する関数
  void set_state_size() {
    double passed_time = max(0., cur_time - start_time);
    double n = 1.0 - max(0.0, passed_time * inv_time * 1.05 - 0.05);
    state_size = int(max(1., state_max * n) + 0.5);
  }

  // update すべきかどうかのしきい値を計算する
  double get_thresfold() { return cur_temp * rand_log(iter); }

 public:
  using pair_t = pair<state_t, score_t>;
  long long iter = 0, loop_count = 0;
  vector<pair_t> states;

#ifdef NyaanLocal
#define debug true
#else
#define debug false
#endif

  void dump() {
    // 重いので適宜コメントアウトする
    /**/
    if constexpr (debug) {
      if ((loop_count & 0x3FFF) == 0) {
        if (loop_count == 0) {
          cerr << "     iter | ";
          cerr << "    loops | ";
          cerr << "   time | ";
          cerr << "  cur temp | ";
          cerr << "size | ";
          cerr << "      cur score | ";
          cerr << endl;
        }
        cerr << setprecision(3);
        cerr << setw(9) << iter << " | ";
        cerr << setw(9) << loop_count << " | ";
        cerr << setw(7) << int(cur_time) << " | ";
        cerr << setw(10) << cur_temp << " | ";
        cerr << setw(4) << states.size() << " | ";
        auto comp_pair = [](const pair_t& lhs, const pair_t& rhs) {
          return lhs.second > rhs.second;
        };
        sort(begin(states), end(states), comp_pair);
        for (int i = 0; i < min<int>(states.size(), 5); i++) {
          cerr << " " << setw(14) << states[i].second << " | ";
        }
        cerr << endl;
      }
    }
    //*/
  }

  template <pair_t (*get_init_state)(),
            void (*update_state)(state_t&, score_t&, double)>
  pair_t run() {
    auto comp_pair = [](const pair_t& lhs, const pair_t& rhs) {
      return lhs.second > rhs.second;
    };

    // states の初期化
    set_time();
    if constexpr (debug) cerr << "init start : " << int(cur_time) << endl;
    states.clear();
    for (int i = 0; i < state_max; i++) states.push_back(get_init_state());
    sort(begin(states), end(states), comp_pair);

    // 時間・温度の初期化
    init_time(), set_temp();
    if constexpr (debug) cerr << "SA start   : " << int(cur_time) << endl;
    if (end_time <= cur_time) return states[0];

    iter = loop_count = 0;
    state_size = state_max;  // 現在の states のサイズ
    constexpr int L = 256;

    while (state_size != 1) {
      if ((loop_count & (L - 1)) == 0) {
        set_time();
        if (cur_time >= end_time) break;
        set_temp();
      }
      dump();
      loop_count++;
      for (auto& state : states) {
        iter++;
        double thresfold = get_thresfold();
        update_state(state.first, state.second, thresfold);
      }
      set_state_size();
      while (state_size < (int)states.size()) {
        sort(begin(states), end(states), comp_pair);
        states.pop_back();
      }
    }
    loop_count &= ~(L - 1);
    while (true) {
      set_time();
      if (cur_time >= end_time) break;
      set_temp();
      dump();
      for (int t = 0; t < L; t++) {
        double thresfold = get_thresfold();
        update_state(states[0].first, states[0].second, thresfold);
      }
      loop_count += L, iter += L;
    }
    sort(begin(states), end(states), comp_pair);
    if constexpr (debug) {
      cerr << "SA end     : " << int(cur_time) << endl;
      cerr << "iteration  : " << iter << endl;
      cerr << "loops      : " << loop_count << endl;
      cerr << "score      : " << states[0].second << endl;
      cerr << string(50, '*') << endl;
    }
    return states[0];
  }

  SA_manager(Timer& _timer, double _start_temp, double _end_temp,
             double _end_time, int _state_max = 1)
      : timer(_timer),
        start_time(0.),
        end_time(_end_time),
        inv_time(1. / _end_time),
        cur_time(0.),
        start_temp(_start_temp),
        diff_temp(_start_temp - _end_temp),
        cur_temp(_start_temp),
        state_max(_state_max),
        state_size(_state_max) {
    assert(_start_temp >= _end_temp);
    assert(_state_max > 0);
  }
};

#undef debug

/**
 * @brief Multipoint Simulated Annealing
 * @docs docs/marathon/sa-manager.md
 */
#line 2 "marathon/sa-manager.hpp"

#include <cassert>
#include <iostream>
#include <vector>
using namespace std;

#line 2 "misc/rng.hpp"

#line 2 "internal/internal-seed.hpp"

#include <chrono>
using namespace std;

namespace internal {
unsigned long long non_deterministic_seed() {
  unsigned long long m =
      chrono::duration_cast<chrono::nanoseconds>(
          chrono::high_resolution_clock::now().time_since_epoch())
          .count();
  m ^= 9845834732710364265uLL;
  m ^= m << 24, m ^= m >> 31, m ^= m << 35;
  return m;
}
unsigned long long deterministic_seed() { return 88172645463325252UL; }

// 64 bit の seed 値を生成 (手元では seed 固定)
// 連続で呼び出すと同じ値が何度も返ってくるので注意
// #define RANDOMIZED_SEED するとシードがランダムになる
unsigned long long seed() {
#if defined(NyaanLocal) && !defined(RANDOMIZED_SEED)
  return deterministic_seed();
#else
  return non_deterministic_seed();
#endif
}

}  // namespace internal
#line 4 "misc/rng.hpp"

namespace my_rand {
using i64 = long long;
using u64 = unsigned long long;

// [0, 2^64 - 1)
u64 rng() {
  static u64 _x = internal::seed();
  return _x ^= _x << 7, _x ^= _x >> 9;
}

// [l, r]
i64 rng(i64 l, i64 r) {
  assert(l <= r);
  return l + rng() % u64(r - l + 1);
}

// [l, r)
i64 randint(i64 l, i64 r) {
  assert(l < r);
  return l + rng() % u64(r - l);
}

// choose n numbers from [l, r) without overlapping
vector<i64> randset(i64 l, i64 r, i64 n) {
  assert(l <= r && n <= r - l);
  unordered_set<i64> s;
  for (i64 i = n; i; --i) {
    i64 m = randint(l, r + 1 - i);
    if (s.find(m) != s.end()) m = r - i;
    s.insert(m);
  }
  vector<i64> ret;
  for (auto& x : s) ret.push_back(x);
  return ret;
}

// [0.0, 1.0)
double rnd() { return rng() * 5.42101086242752217004e-20; }
// [l, r)
double rnd(double l, double r) {
  assert(l < r);
  return l + rnd() * (r - l);
}

template <typename T>
void randshf(vector<T>& v) {
  int n = v.size();
  for (int i = 1; i < n; i++) swap(v[i], v[randint(0, i + 1)]);
}

}  // namespace my_rand

using my_rand::randint;
using my_rand::randset;
using my_rand::randshf;
using my_rand::rnd;
using my_rand::rng;
#line 2 "misc/timer.hpp"

#line 4 "misc/timer.hpp"
using namespace std;

struct Timer {
  chrono::high_resolution_clock::time_point st;

  Timer() { reset(); }
  void reset() { st = chrono::high_resolution_clock::now(); }

  long long elapsed() {
    auto ed = chrono::high_resolution_clock::now();
    return chrono::duration_cast<chrono::milliseconds>(ed - st).count();
  }
  long long operator()() { return elapsed(); }
};
#line 2 "marathon/log_table.hpp"

struct log_table {
  static constexpr int M = 65536;
  static constexpr int mask = M - 1;
  double l[M];
  log_table() : l() {
    unsigned long long x = 88172645463325252ULL;
    double log_u64max = log(2) * 64;
    for (int i = 0; i < M; i++) {
      x = x ^ (x << 7);
      x = x ^ (x >> 9);
      l[i] = log(double(x)) - log_u64max;
    }
  }
  double operator()(int i) const { return l[i & mask]; }
};
#line 11 "marathon/sa-manager.hpp"

template <typename state_t, typename score_t>
struct SA_manager {
 private:
  log_table rand_log;
  Timer timer;
  double start_time, end_time, inv_time, cur_time;
  double start_temp, diff_temp, cur_temp;
  int state_max, state_size;

  // 時間をセットする関数
  void set_time() { cur_time = timer.elapsed(); }

  // 時間を色々初期化する関数
  void init_time() {
    set_time();
    start_time = cur_time;
    inv_time = 1. / (end_time - start_time);
  }

  // 気温を設定して cur_temp に代入する関数
  void set_temp() {
    double passed_time = max(0., cur_time - start_time);
    cur_temp = max(0., start_temp - passed_time * inv_time * diff_temp);
  }

  // 焼きなましで持つ状態の個数の max を代入する関数
  void set_state_size() {
    double passed_time = max(0., cur_time - start_time);
    double n = 1.0 - max(0.0, passed_time * inv_time * 1.05 - 0.05);
    state_size = int(max(1., state_max * n) + 0.5);
  }

  // update すべきかどうかのしきい値を計算する
  double get_thresfold() { return cur_temp * rand_log(iter); }

 public:
  using pair_t = pair<state_t, score_t>;
  long long iter = 0, loop_count = 0;
  vector<pair_t> states;

#ifdef NyaanLocal
#define debug true
#else
#define debug false
#endif

  void dump() {
    // 重いので適宜コメントアウトする
    /**/
    if constexpr (debug) {
      if ((loop_count & 0x3FFF) == 0) {
        if (loop_count == 0) {
          cerr << "     iter | ";
          cerr << "    loops | ";
          cerr << "   time | ";
          cerr << "  cur temp | ";
          cerr << "size | ";
          cerr << "      cur score | ";
          cerr << endl;
        }
        cerr << setprecision(3);
        cerr << setw(9) << iter << " | ";
        cerr << setw(9) << loop_count << " | ";
        cerr << setw(7) << int(cur_time) << " | ";
        cerr << setw(10) << cur_temp << " | ";
        cerr << setw(4) << states.size() << " | ";
        auto comp_pair = [](const pair_t& lhs, const pair_t& rhs) {
          return lhs.second > rhs.second;
        };
        sort(begin(states), end(states), comp_pair);
        for (int i = 0; i < min<int>(states.size(), 5); i++) {
          cerr << " " << setw(14) << states[i].second << " | ";
        }
        cerr << endl;
      }
    }
    //*/
  }

  template <pair_t (*get_init_state)(),
            void (*update_state)(state_t&, score_t&, double)>
  pair_t run() {
    auto comp_pair = [](const pair_t& lhs, const pair_t& rhs) {
      return lhs.second > rhs.second;
    };

    // states の初期化
    set_time();
    if constexpr (debug) cerr << "init start : " << int(cur_time) << endl;
    states.clear();
    for (int i = 0; i < state_max; i++) states.push_back(get_init_state());
    sort(begin(states), end(states), comp_pair);

    // 時間・温度の初期化
    init_time(), set_temp();
    if constexpr (debug) cerr << "SA start   : " << int(cur_time) << endl;
    if (end_time <= cur_time) return states[0];

    iter = loop_count = 0;
    state_size = state_max;  // 現在の states のサイズ
    constexpr int L = 256;

    while (state_size != 1) {
      if ((loop_count & (L - 1)) == 0) {
        set_time();
        if (cur_time >= end_time) break;
        set_temp();
      }
      dump();
      loop_count++;
      for (auto& state : states) {
        iter++;
        double thresfold = get_thresfold();
        update_state(state.first, state.second, thresfold);
      }
      set_state_size();
      while (state_size < (int)states.size()) {
        sort(begin(states), end(states), comp_pair);
        states.pop_back();
      }
    }
    loop_count &= ~(L - 1);
    while (true) {
      set_time();
      if (cur_time >= end_time) break;
      set_temp();
      dump();
      for (int t = 0; t < L; t++) {
        double thresfold = get_thresfold();
        update_state(states[0].first, states[0].second, thresfold);
      }
      loop_count += L, iter += L;
    }
    sort(begin(states), end(states), comp_pair);
    if constexpr (debug) {
      cerr << "SA end     : " << int(cur_time) << endl;
      cerr << "iteration  : " << iter << endl;
      cerr << "loops      : " << loop_count << endl;
      cerr << "score      : " << states[0].second << endl;
      cerr << string(50, '*') << endl;
    }
    return states[0];
  }

  SA_manager(Timer& _timer, double _start_temp, double _end_temp,
             double _end_time, int _state_max = 1)
      : timer(_timer),
        start_time(0.),
        end_time(_end_time),
        inv_time(1. / _end_time),
        cur_time(0.),
        start_temp(_start_temp),
        diff_temp(_start_temp - _end_temp),
        cur_temp(_start_temp),
        state_max(_state_max),
        state_size(_state_max) {
    assert(_start_temp >= _end_temp);
    assert(_state_max > 0);
  }
};

#undef debug

/**
 * @brief Multipoint Simulated Annealing
 * @docs docs/marathon/sa-manager.md
 */
Back to top page