結果
| 問題 | No.5025 Forest Fire |
| コンテスト | |
| ユーザー |
 prussian_coder
|
| 提出日時 | 2026-10-24 23:17:00 |
| 言語 | C++23(gcc16) (gcc 16.2.0 + boost 1.92.0 + ACL) |
| 結果 |
AC
不安定
|
| 実行時間 | 1,796 ms / 2,000 ms |
| + 663µs | |
| コード長 | 31,356 bytes |
| 記録 | |
| コンパイル時間 | 6,666 ms |
| コンパイル使用メモリ | 331,584 KB |
| 実行使用メモリ | 34,060 KB |
| スコア | 191,001,625 |
| 最終ジャッジ日時 | 2026-10-24 23:19:45 |
| 合計ジャッジ時間 | 156,039 ms |
|
ジャッジサーバーID （参考情報） |
judge5_0 / <nil> |
| 純コード判定しない問題か言語 |
（要ログイン）
| ファイルパターン | 結果 |
|---|---|
| other | AC * 50 |
ソースコード
#include <algorithm>
#include <array>
#include <chrono>
#include <cmath>
#include <cstdint>
#include <iostream>
#include <limits>
#include <random>
#include <string>
#include <utility>
#include <unordered_set>
#include <vector>
using namespace std;
using Connection = pair<int, int>;
constexpr int SAFE_ZONE_CRITERIA = 4;
constexpr int CAND_SIZES_LIM = 30000;
constexpr double TIME_LIMIT = 1.8;
constexpr double ROUTE_TIME_LIMIT = 1.0;
constexpr double MINCUT_TIME_LIMIT = 0.18;
constexpr double MINCUT_PROB = 0.5;
constexpr int REFINE_SAMPLES = 512;
constexpr double START_TEMP = 10000;
constexpr double END_TEMP = 100;
constexpr double TWO_ADD_ONE_REMOVE_PROB = 0.1;
constexpr double TWO_REMOVE_ONE_ADD_PROB = 0.1;
constexpr double SHUFFLE_PROB = 0.1;
constexpr double TWO_ADD_TWO_REMOVE_PROB = 0.1;
int N, K, N2;
vector<vector<int>> board;
vector<int> danger_zones;
int to_1d(int i, int j) { return i * N2 + j; }
pair<int, int> to_2d(int v) { return {v / N2, v % N2}; }
int connection_dist(int v1, int v2) {
auto [x1, y1] = to_2d(v1);
auto [x2, y2] = to_2d(v2);
return min(abs(x1 - x2), abs(y1 - y2));
}
int connection_dist(Connection edge) {
return connection_dist(edge.first, edge.second);
}
vector<int> create_route(int v1, int v2) {
auto [x1, y1] = to_2d(v1);
auto [x2, y2] = to_2d(v2);
if (x1 > x2) {
swap(x1, x2);
swap(y1, y2);
}
vector<int> route;
int x = x1, y = y1;
while (x != x2 || y != y2) {
if (x < x2) ++x;
if (y < y2) ++y;
else if (y > y2) --y;
if (x != x2 || y != y2) route.push_back(to_1d(x, y));
}
return route;
}
struct Dinic {
struct Edge { int to, next, capacity; };
vector<int> head, level, cursor, queue;
vector<Edge> edges;
explicit Dinic(int n) : head(n, -1), level(n), cursor(n), queue(n) {
edges.reserve(n * 6);
}
void add_edge(int from, int to, int capacity) {
edges.push_back({to, head[from], capacity});
head[from] = static_cast<int>(edges.size()) - 1;
edges.push_back({from, head[to], 0});
head[to] = static_cast<int>(edges.size()) - 1;
}
bool bfs(int source, int sink) {
fill(level.begin(), level.end(), -1);
int begin = 0, end = 0;
queue[end++] = source;
level[source] = 0;
while (begin < end) {
int v = queue[begin++];
for (int e = head[v]; e >= 0; e = edges[e].next) {
const auto& edge = edges[e];
if (edge.capacity > 0 && level[edge.to] < 0) {
level[edge.to] = level[v] + 1;
queue[end++] = edge.to;
}
}
}
return level[sink] >= 0;
}
int dfs(int v, int sink, int flow) {
if (v == sink) return flow;
for (int& e = cursor[v]; e >= 0; e = edges[e].next) {
auto& edge = edges[e];
if (edge.capacity == 0 || level[edge.to] != level[v] + 1) continue;
int sent = dfs(edge.to, sink, min(flow, edge.capacity));
if (sent) {
edge.capacity -= sent;
edges[e ^ 1].capacity += sent;
return sent;
}
}
return 0;
}
int max_flow(int source, int sink, int limit) {
int total = 0;
while (total < limit && bfs(source, sink)) {
cursor = head;
while (total < limit) {
int sent = dfs(source, sink, limit - total);
if (!sent) break;
total += sent;
}
}
return total;
}
};
// 4近傍の危険領域を頂点カットで分割する。各木を in -> out に分け、
// この辺に伐採コスト、隣接辺に十分大きい容量を置く。
// 射影の両端を端子領域にして、1本だけを孤立させるカットへの偏りを抑える。
vector<vector<int>> generate_mincuts(uint32_t seed) {
const auto start = chrono::steady_clock::now();
vector<vector<int>> components, result;
vector<int> label(N2 * N2, -1), local(N2 * N2, -1);
for (int v : danger_zones) if (label[v] < 0) {
int id = static_cast<int>(components.size());
components.push_back({v});
label[v] = id;
for (int j = 0; j < static_cast<int>(components.back().size()); ++j) {
int u = components.back()[j];
for (int d : {-N2, N2, -1, 1}) {
int w = u + d;
if (board[w / N2][w % N2] > SAFE_ZONE_CRITERIA && label[w] < 0) {
label[w] = id;
components.back().push_back(w);
}
}
}
}
sort(components.begin(), components.end(), [](const auto& a, const auto& b) {
return a.size() > b.size();
});
while (!components.empty() && components.back().size() < 40) components.pop_back();
if (components.empty()) return result;
mt19937 rng(seed ^ 0x9e3779b9U);
int mass = 0;
for (const auto& component : components) mass += static_cast<int>(component.size() * component.size());
constexpr array<pair<int, int>, 8> directions = {{{1, 0}, {0, 1}, {1, 1}, {1, -1},
{2, 1}, {2, -1}, {1, 2}, {1, -2}}};
constexpr int INF = 100000000;
for (int trial = 0; trial < 256; ++trial) {
if (chrono::duration<double>(chrono::steady_clock::now() - start).count() >= MINCUT_TIME_LIMIT) break;
int pick = static_cast<int>(rng() % mass), id = 0;
while (pick >= static_cast<int>(components[id].size() * components[id].size())) {
pick -= static_cast<int>(components[id].size() * components[id].size());
++id;
}
const auto& vertices = components[id];
int n = static_cast<int>(vertices.size());
auto [dx, dy] = directions[trial % directions.size()];
vector<int> projection(n), sorted_projection;
for (int j = 0; j < n; ++j) {
int v = vertices[j];
local[v] = j;
projection[j] = dx * (v / N2) + dy * (v % N2);
}
sorted_projection = projection;
sort(sorted_projection.begin(), sorted_projection.end());
int center = 20 + static_cast<int>(rng() % 61);
int half_width = 10 + static_cast<int>(rng() % 21);
int low = sorted_projection[n * max(5, center - half_width) / 100];
int high = sorted_projection[n * min(95, center + half_width) / 100];
if (low >= high) continue;
int source = 2 * n, sink = source + 1;
Dinic graph(2 * n + 2);
for (int j = 0; j < n; ++j) {
bool terminal = projection[j] <= low || projection[j] >= high;
// 小さな容量の揺らぎで、ほぼ同じ伐採本数の異なるカットを得る。
graph.add_edge(2 * j, 2 * j + 1, terminal ? INF : 1000 + static_cast<int>(rng() % 300));
if (projection[j] <= low) graph.add_edge(source, 2 * j, INF);
if (projection[j] >= high) graph.add_edge(2 * j + 1, sink, INF);
int v = vertices[j];
for (int d : {-N2, N2, -1, 1}) {
int u = v + d;
// 別の連結成分には辺がなく、境界の番兵マスは常に湿っている。
if (board[u / N2][u % N2] > SAFE_ZONE_CRITERIA) {
graph.add_edge(2 * j + 1, 2 * local[u], INF);
}
}
}
if (graph.max_flow(source, sink, 1300 * K + 1) >= 1300 * K + 1) continue;
vector<int> cut;
for (int j = 0; j < n; ++j) {
if (graph.level[2 * j] >= 0 && graph.level[2 * j + 1] < 0) cut.push_back(vertices[j]);
}
if (cut.empty() || static_cast<int>(cut.size()) > K) continue;
sort(cut.begin(), cut.end());
if (find(result.begin(), result.end(), cut) == result.end()) result.push_back(std::move(cut));
}
return result;
}
// 入力と候補辺だけで決まる情報をキャッシュし、反復用メモリを再利用する。
struct SearchWorkspace {
vector<int> routes, route_offsets;
vector<int> distances, danger_idx, trees, parents;
vector<pair<int, int>> adjacent_dangers;
vector<unsigned int> used;
vector<uint64_t> removal_mask, cached_removal_mask;
vector<int> removed_indices;
long long cached_score = 0;
bool score_cached = false;
unsigned int tree_stamp = 0;
void add_cut(const vector<int>& cut) {
routes.insert(routes.end(), cut.begin(), cut.end());
route_offsets.push_back(static_cast<int>(routes.size()));
distances.push_back(static_cast<int>(cut.size()));
}
explicit SearchWorkspace(const vector<Connection>& candidates)
: danger_idx(N2 * N2, -1), parents(danger_zones.size()),
used(N2 * N2, 0), removal_mask((danger_zones.size() + 63) / 64),
cached_removal_mask(removal_mask.size()) {
routes.reserve(candidates.size() * N2);
route_offsets.reserve(candidates.size() + 1);
route_offsets.push_back(0);
distances.reserve(candidates.size());
trees.reserve(N2 * N2);
removed_indices.reserve(N2 * N2);
for (auto [v1, v2] : candidates) {
auto [x1, y1] = to_2d(v1);
auto [x2, y2] = to_2d(v2);
if (x1 > x2) {
swap(x1, x2);
swap(y1, y2);
}
int x = x1, y = y1;
while (x != x2 || y != y2) {
if (x < x2) ++x;
if (y < y2) ++y;
else if (y > y2) --y;
if ((x != x2 || y != y2) && board[x][y] != 0) {
routes.push_back(to_1d(x, y));
}
}
route_offsets.push_back(static_cast<int>(routes.size()));
distances.push_back(connection_dist(v1, v2));
}
for (int i = 0; i < static_cast<int>(danger_zones.size()); ++i) {
danger_idx[danger_zones[i]] = i;
}
for (int i = 0; i < static_cast<int>(danger_zones.size()); ++i) {
int v = danger_zones[i];
for (int d : {1, N2}) {
int neighbor = v + d;
if (neighbor < N2 * N2 && danger_idx[neighbor] >= 0) {
adjacent_dangers.emplace_back(i, danger_idx[neighbor]);
}
}
}
}
const vector<int>& build(const vector<int>& edges) {
if (++tree_stamp == 0) {
fill(used.begin(), used.end(), 0);
tree_stamp = 1;
}
trees.clear();
for (auto it = edges.rbegin(); it != edges.rend(); ++it) {
for (int i = route_offsets[*it]; i < route_offsets[*it + 1]; ++i) {
int v = routes[i];
if (used[v] == tree_stamp) continue;
used[v] = tree_stamp;
trees.push_back(v);
if (static_cast<int>(trees.size()) == K) return trees;
}
}
return trees;
}
int find_root(int x) {
if (parents[x] < 0) return x;
return parents[x] = find_root(parents[x]);
}
long long evaluate(const vector<int>& selected_trees) {
constexpr int REMOVED = numeric_limits<int>::min();
fill(removal_mask.begin(), removal_mask.end(), 0);
removed_indices.clear();
for (int v : selected_trees) {
int idx = danger_idx[v];
if (idx < 0) continue;
uint64_t bit = uint64_t{1} << (idx % 64);
if ((removal_mask[idx / 64] & bit) == 0) {
removal_mask[idx / 64] |= bit;
removed_indices.push_back(idx);
}
}
// 削除対象が前回と完全一致する場合だけ、同じ評価値を再利用する。
if (score_cached && removal_mask == cached_removal_mask) return cached_score;
fill(parents.begin(), parents.end(), -1);
for (int idx : removed_indices) parents[idx] = REMOVED;
// 無向辺は一度だけ処理する。連結成分と評価値は元の処理と同じ。
for (auto [a, b] : adjacent_dangers) {
if (parents[a] == REMOVED || parents[b] == REMOVED) continue;
a = find_root(a);
b = find_root(b);
if (a == b) continue;
if (parents[a] > parents[b]) swap(a, b);
parents[a] += parents[b];
parents[b] = a;
}
long long score = 0;
for (int size : parents) {
if (size < 0 && size != REMOVED) score += 1LL * size * size;
}
cached_removal_mask = removal_mask;
cached_score = score;
score_cached = true;
return score;
}
};
// 共通の乾燥サンプルで比較する。各頂点を削除したときの改善量は
// DFS 木の low-link と部分木サイズから、全頂点分を線形時間で求める。
struct SampleEvaluator {
static constexpr int SAMPLES = 128;
int cells;
vector<int> vertices;
vector<array<int, 4>> neighbors;
vector<vector<unsigned char>> dry;
vector<unsigned char> cut;
vector<int> tin, low, subtree, component, component_sizes, order, separated;
vector<long long> separated_sq, gain, restore_cost;
int timer = 0;
explicit SampleEvaluator(int sample_count = SAMPLES) : cells(N2 * N2), neighbors(cells),
dry(sample_count, vector<unsigned char>(cells)), cut(cells),
tin(cells), low(cells), subtree(cells), component(cells),
separated(cells), separated_sq(cells), gain(cells), restore_cost(cells) {
uint64_t seed = 1469598103934665603ULL;
for (int x = 1; x <= N; ++x) for (int y = 1; y <= N; ++y) {
int v = to_1d(x, y);
neighbors[v] = {v - N2, v + N2, v - 1, v + 1};
seed = (seed ^ uint64_t(board[x][y])) * 1099511628211ULL;
if (board[x][y]) vertices.push_back(v);
}
mt19937_64 random(seed);
for (auto& sample : dry) for (int v : vertices) {
sample[v] = random() % 10 < uint64_t(board[v / N2][v % N2]);
}
order.reserve(cells);
component_sizes.reserve(cells);
}
void dfs(int v, int parent, int id, const vector<unsigned char>& active) {
tin[v] = low[v] = ++timer;
subtree[v] = 1;
separated[v] = 0;
separated_sq[v] = 0;
component[v] = id;
order.push_back(v);
for (int u : neighbors[v]) {
if (!active[u] || cut[u] || u == parent) continue;
if (tin[u]) {
low[v] = min(low[v], tin[u]);
} else {
dfs(u, v, id, active);
subtree[v] += subtree[u];
low[v] = min(low[v], low[u]);
if (low[u] >= tin[v]) {
separated[v] += subtree[u];
separated_sq[v] += 1LL * subtree[u] * subtree[u];
}
}
}
}
long long analyze() {
fill(gain.begin(), gain.end(), 0);
fill(restore_cost.begin(), restore_cost.end(), 0);
long long score = 0;
for (const auto& sample : dry) {
fill(tin.begin(), tin.end(), 0);
component_sizes.clear();
order.clear();
timer = 0;
for (int v : vertices) if (sample[v] && !cut[v] && !tin[v]) {
int id = static_cast<int>(component_sizes.size());
int begin = static_cast<int>(order.size());
dfs(v, -1, id, sample);
int size = subtree[v];
component_sizes.push_back(size);
score += 1LL * size * size;
for (int j = begin; j < static_cast<int>(order.size()); ++j) {
int u = order[j];
int rest = size - 1 - separated[u];
gain[u] += 1LL * size * size - separated_sq[u] - 1LL * rest * rest;
}
}
for (int v : vertices) if (cut[v] && sample[v]) {
array<int, 4> ids;
int count = 0, sum = 1;
long long square_sum = 0;
for (int u : neighbors[v]) if (sample[u] && !cut[u]) {
int id = component[u];
bool duplicate = false;
for (int j = 0; j < count; ++j) duplicate |= ids[j] == id;
if (duplicate) continue;
ids[count++] = id;
int size = component_sizes[id];
sum += size;
square_sum += 1LL * size * size;
}
restore_cost[v] += 1LL * sum * sum - square_sum;
}
}
return score;
}
};
// 変更した木を含む成分と、復元した木に隣接する成分だけを作り直す。
// サンプルごとの頂点削除効果を保持し、全サンプルの和を差分更新する。
struct DynamicEvaluator {
struct State {
vector<int> component, next, size, gain;
explicit State(int n) : component(n, -1), next(n, -1), size(n), gain(n) {}
};
const SampleEvaluator& samples;
vector<State> states;
vector<unsigned char> cut;
vector<long long> gain, restore_cost;
vector<int> tin, low, subtree, separated, separated_sq, order, affected, roots, inserted;
int timer = 0;
long long score = 0;
DynamicEvaluator(const SampleEvaluator& samples, const vector<int>& trees)
: samples(samples), cut(samples.cells), gain(samples.cells), restore_cost(samples.cells),
tin(samples.cells), low(samples.cells), subtree(samples.cells),
separated(samples.cells), separated_sq(samples.cells) {
for (int v : trees) cut[v] = 1;
states.reserve(samples.dry.size());
order.reserve(samples.cells);
affected.reserve(samples.cells);
for (int s = 0; s < static_cast<int>(samples.dry.size()); ++s) {
states.emplace_back(samples.cells);
fill(tin.begin(), tin.end(), 0);
timer = 0;
for (int v : samples.vertices) if (samples.dry[s][v] && !cut[v] && !tin[v]) {
build_component(s, v);
}
}
}
void dfs(int s, int v, int parent, int root) {
auto& state = states[s];
tin[v] = low[v] = ++timer;
subtree[v] = 1;
separated[v] = separated_sq[v] = 0;
state.component[v] = root;
order.push_back(v);
for (int u : samples.neighbors[v]) {
if (!samples.dry[s][u] || cut[u] || u == parent) continue;
if (tin[u]) {
low[v] = min(low[v], tin[u]);
} else {
dfs(s, u, v, root);
subtree[v] += subtree[u];
low[v] = min(low[v], low[u]);
if (low[u] >= tin[v]) {
separated[v] += subtree[u];
separated_sq[v] += subtree[u] * subtree[u];
}
}
}
}
void build_component(int s, int root) {
auto& state = states[s];
order.clear();
dfs(s, root, -1, root);
int size = static_cast<int>(order.size());
state.size[root] = size;
score += 1LL * size * size;
for (int i = 0; i < size; ++i) {
int v = order[i];
state.next[v] = i + 1 < size ? order[i + 1] : -1;
int rest = size - 1 - separated[v];
state.gain[v] = size * size - separated_sq[v] - rest * rest;
gain[v] += state.gain[v];
}
}
void invalidate(int s, int root) {
auto& state = states[s];
score -= 1LL * state.size[root] * state.size[root];
state.size[root] = 0;
for (int v = root; v >= 0; v = state.next[v]) {
gain[v] -= state.gain[v];
state.gain[v] = 0;
state.component[v] = -1;
tin[v] = 0;
affected.push_back(v);
}
}
// changes は相異なる木と、その変更後の伐採状態の組。
void apply(const vector<pair<int, bool>>& changes) {
for (auto [v, removed] : changes) cut[v] = removed;
for (int s = 0; s < static_cast<int>(states.size()); ++s) {
auto& state = states[s];
roots.clear();
inserted.clear();
for (auto [v, removed] : changes) {
if (!samples.dry[s][v]) continue;
if (state.component[v] >= 0) roots.push_back(state.component[v]);
else if (!removed) inserted.push_back(v);
if (!removed) for (int u : samples.neighbors[v]) {
if (state.component[u] >= 0) roots.push_back(state.component[u]);
}
}
if (roots.empty() && inserted.empty()) continue;
sort(roots.begin(), roots.end());
roots.erase(unique(roots.begin(), roots.end()), roots.end());
affected.clear();
for (int root : roots) invalidate(s, root);
for (int v : inserted) {
tin[v] = 0;
affected.push_back(v);
}
timer = 0;
for (int v : affected) if (!cut[v] && !tin[v]) build_component(s, v);
}
}
void restoration_costs(const vector<int>& trees) {
fill(restore_cost.begin(), restore_cost.end(), 0);
for (int s = 0; s < static_cast<int>(states.size()); ++s) {
const auto& state = states[s];
for (int v : trees) if (samples.dry[s][v]) {
array<int, 4> ids;
int count = 0, sum = 1;
long long square_sum = 0;
for (int u : samples.neighbors[v]) {
int root = state.component[u];
if (root < 0) continue;
bool duplicate = false;
for (int j = 0; j < count; ++j) duplicate |= ids[j] == root;
if (duplicate) continue;
ids[count++] = root;
int size = state.size[root];
sum += size;
square_sum += 1LL * size * size;
}
restore_cost[v] += 1LL * sum * sum - square_sum;
}
}
}
};
void refine_trees(vector<int>& trees, chrono::steady_clock::time_point start) {
SampleEvaluator samples(REFINE_SAMPLES);
DynamicEvaluator evaluator(samples, trees);
vector<unsigned char> tried(samples.cells);
auto expired = [&]() {
return chrono::duration<double>(chrono::steady_clock::now() - start).count() >= TIME_LIMIT;
};
int iterations = 0, improvements = 0;
bool costs_dirty = true;
while (!expired()) {
if (static_cast<int>(trees.size()) < K) {
int add = -1;
for (int v : samples.vertices) if (!evaluator.cut[v] &&
(add < 0 || evaluator.gain[v] > evaluator.gain[add])) add = v;
if (add < 0) break;
evaluator.apply({{add, true}});
trees.push_back(add);
costs_dirty = true;
} else {
// 棄却した交換では状態が元に戻るため、復元コストも再利用できる。
if (costs_dirty) evaluator.restoration_costs(trees);
costs_dirty = false;
int remove = -1;
for (int v : trees) if (!tried[v] &&
(remove < 0 || evaluator.restore_cost[v] < evaluator.restore_cost[remove])) remove = v;
if (remove < 0) break;
long long previous_score = evaluator.score;
evaluator.apply({{remove, false}});
int add = remove;
for (int v : samples.vertices) if (!evaluator.cut[v] && evaluator.gain[v] > evaluator.gain[add]) add = v;
if (evaluator.score - evaluator.gain[add] < previous_score) {
evaluator.apply({{add, true}});
*find(trees.begin(), trees.end(), remove) = add;
fill(tried.begin(), tried.end(), 0);
++improvements;
costs_dirty = true;
} else {
evaluator.apply({{remove, true}});
tried[remove] = 1;
}
}
++iterations;
}
cerr << "dynamic_iterations=" << iterations << " improvements=" << improvements
<< " sample_score=" << evaluator.score << '\n';
}
int main() {
const auto start_time = chrono::steady_clock::now();
ios::sync_with_stdio(false);
cin.tie(nullptr);
if (!(cin >> N >> K)) return 0;
N2 = N + 2;
board.assign(N2, vector<int>(N2, 0));
for (int i = 0; i < N; ++i) {
string row;
cin >> row;
for (int j = 0; j < N; ++j) board[i + 1][j + 1] = row[j] - '0';
}
uint32_t input_seed = 2166136261U;
for (const auto& row : board) for (int d : row) {
input_seed = (input_seed ^ static_cast<uint32_t>(d)) * 16777619U;
}
mt19937 rng(input_seed);
auto randint = [&](int lo, int hi) {
return uniform_int_distribution<int>(lo, hi)(rng);
};
auto random_real = [&]() {
return uniform_real_distribution<double>(0.0, 1.0)(rng);
};
vector<int> safe_zones;
for (int i = 0; i < N2; ++i) {
for (int j = 0; j < N2; ++j) {
if (board[i][j] <= SAFE_ZONE_CRITERIA) safe_zones.push_back(to_1d(i, j));
else danger_zones.push_back(to_1d(i, j));
}
}
unordered_set<uint64_t> candidate_set;
candidate_set.reserve(CAND_SIZES_LIM + 1);
vector<Connection> edge_candidates;
edge_candidates.reserve(CAND_SIZES_LIM + 1);
// 候補総数が上限より少ない入力でも終了するよう、試行回数を制限する。
for (int attempts = 0; candidate_set.size() <= CAND_SIZES_LIM && attempts < 1000000; ++attempts) {
int v1 = safe_zones[randint(0, static_cast<int>(safe_zones.size()) - 1)];
int v2 = safe_zones[randint(0, static_cast<int>(safe_zones.size()) - 1)];
int dist = connection_dist(v1, v2);
if (2 <= dist && dist <= 50) {
uint64_t key = (uint64_t(static_cast<uint32_t>(v1)) << 32) |
static_cast<uint32_t>(v2);
if (candidate_set.insert(key).second) edge_candidates.emplace_back(v1, v2);
}
}
// 旧版の set と同じ候補順に戻し、乱数に対する選択結果を維持する。
sort(edge_candidates.begin(), edge_candidates.end());
SearchWorkspace workspace(edge_candidates);
const int mincut_begin = static_cast<int>(edge_candidates.size());
for (const auto& cut : generate_mincuts(input_seed)) workspace.add_cut(cut);
const int candidate_count = static_cast<int>(workspace.distances.size());
const int candidate_max = candidate_count - 1;
auto random_candidate = [&]() {
// 最小カット候補は経路候補より少ないので、選択確率を別に確保する。
if (candidate_count > mincut_begin && random_real() < MINCUT_PROB) {
return randint(mincut_begin, candidate_max);
}
return randint(0, candidate_max);
};
cerr << "mincut_candidates=" << candidate_count - mincut_begin << '\n';
vector<int> edges, prev_edges, best_trees;
edges.reserve(128);
prev_edges.reserve(128);
int est_sizes = 0;
double current_score = numeric_limits<double>::infinity();
long long best_score = numeric_limits<long long>::max();
auto contains_edge = [&](int edge) {
return find(edges.begin(), edges.end(), edge) != edges.end();
};
while (true) {
double elapsed = chrono::duration<double>(chrono::steady_clock::now() - start_time).count();
if (elapsed > ROUTE_TIME_LIMIT) break;
if (candidate_max < 0) break;
double temperature = START_TEMP + (END_TEMP - START_TEMP) * elapsed / ROUTE_TIME_LIMIT;
int x = random_candidate();
double r = random_real();
if (contains_edge(x)) continue;
prev_edges = edges;
int prev_est_sizes = est_sizes;
double move_r = est_sizes >= 90 ? random_real() : 0.0;
if (est_sizes >= 90 && move_r >= 1.0 - SHUFFLE_PROB) {
shuffle(edges.begin(), edges.end(), rng);
} else if (est_sizes >= 90 &&
move_r >= 1.0 - SHUFFLE_PROB - TWO_ADD_TWO_REMOVE_PROB) {
if (edges.size() < 2 || candidate_count - static_cast<int>(edges.size()) < 2) continue;
int z = random_candidate();
while (z == x || contains_edge(z)) z = random_candidate();
int y1 = randint(0, static_cast<int>(edges.size()) - 1);
int y2 = randint(0, static_cast<int>(edges.size()) - 2);
if (y2 >= y1) ++y2;
if (y1 < y2) swap(y1, y2);
for (int y : {y1, y2}) {
est_sizes -= workspace.distances[edges[y]];
edges.erase(edges.begin() + y);
}
for (int new_edge : {x, z}) {
edges.push_back(new_edge);
est_sizes += workspace.distances[new_edge];
}
} else if (est_sizes <= 90) {
est_sizes += workspace.distances[x];
edges.push_back(x);
} else if (move_r < TWO_ADD_ONE_REMOVE_PROB && est_sizes < 105) {
if (candidate_count - static_cast<int>(edges.size()) < 2) continue;
int z = random_candidate();
while (z == x || contains_edge(z)) z = random_candidate();
int y = randint(0, static_cast<int>(edges.size()) - 1);
est_sizes -= workspace.distances[edges[y]];
edges.erase(edges.begin() + y);
for (int new_edge : {x, z}) {
edges.push_back(new_edge);
est_sizes += workspace.distances[new_edge];
}
} else if (TWO_ADD_ONE_REMOVE_PROB <= move_r &&
move_r < TWO_ADD_ONE_REMOVE_PROB + TWO_REMOVE_ONE_ADD_PROB) {
if (edges.size() < 2) continue;
int y1 = randint(0, static_cast<int>(edges.size()) - 1);
int y2 = randint(0, static_cast<int>(edges.size()) - 2);
if (y2 >= y1) ++y2;
if (y1 < y2) swap(y1, y2);
for (int y : {y1, y2}) {
est_sizes -= workspace.distances[edges[y]];
edges.erase(edges.begin() + y);
}
edges.push_back(x);
est_sizes += workspace.distances[x];
} else {
int y = randint(0, static_cast<int>(edges.size()) - 1);
est_sizes -= workspace.distances[edges[y]];
edges[y] = x;
est_sizes += workspace.distances[x];
}
const auto& trees = workspace.build(edges);
long long new_score = workspace.evaluate(trees);
if (new_score <= current_score || r < exp((current_score - new_score) / temperature)) {
current_score = static_cast<double>(new_score);
} else {
edges.swap(prev_edges);
est_sizes = prev_est_sizes;
}
if (new_score < best_score) {
best_trees = trees;
best_score = new_score;
}
}
refine_trees(best_trees, start_time);
cout << best_trees.size() << '\n';
for (int v : best_trees) {
auto [x, y] = to_2d(v);
cout << x - 1 << ' ' << y - 1 << '\n';
}
return 0;
}