Multipoint Simulated Annealing
(marathon/sa-manager.hpp)
- View this file on GitHub
- Last update: 2024-04-28 09:13:11+09:00
- Include:
#include "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>();
}
注意点
- 入力はグローバル領域に置いてよしなにする
-
state_t
が高次元のときは適宜std::array
を使うのが楽だと思う -
score_t
は値が大きければ大きいほど憂愁になるようにする - 焼きなましのループ回数 / 実行回数によってはエラー出力が重いのでコメントアウトする
-
end_temp
は 負数を指定できる。 温度は0.0
とのmax
を取るので負である間は山登り法と等価の挙動をする。- 温度は線形に変化するので、
start_temp
,0.0
間との外分比を考えれば山登り法の動作時間を制御できる。 -
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);
sort(begin(ret), end(ret));
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
*/