結果
| 問題 |
No.5007 Steiner Space Travel
|
| コンテスト | |
| ユーザー |
 Jiro_tech15
|
| 提出日時 | 2023-04-30 12:13:17 |
| 言語 | C++17(gcc12) (gcc 12.3.0 + boost 1.87.0) |
| 結果 |
TLE
|
| 実行時間 | - |
| コード長 | 21,618 bytes |
| コンパイル時間 | 4,669 ms |
| コンパイル使用メモリ | 217,176 KB |
| 実行使用メモリ | 153,080 KB |
| スコア | 0 |
| 最終ジャッジ日時 | 2023-04-30 12:14:04 |
| 合計ジャッジ時間 | 45,771 ms |
|
ジャッジサーバーID （参考情報） |
judge11 / judge12 |
（要ログイン）
| ファイルパターン | 結果 |
|---|---|
| other | TLE * 30 |
コンパイルメッセージ
main.cpp: 静的メンバ関数 ‘static void Pos::StaticInit(int)’ 内:
main.cpp:426:24: 警告: ‘dir’ may be used uninitialized [-Wmaybe-uninitialized]
426 | move_to[dir][p] = {adj_x, adj_y};
| ^
main.cpp:413:17: 備考: ‘dir’ はここで宣言されています
413 | Dir dir;
| ^~~
ソースコード
namespace atcoder {}
#ifdef LOCAL
#define dbg(x) cerr << __LINE__ << " : " << #x << " = " << (x) << endl;
#else
#define NDEBUG
#define dbg(x) true;
#pragma GCC target("avx2")
#pragma GCC optimize("O3")
#pragma GCC optimize("unroll-loops")
#endif
#ifdef GTEST
#include <gtest/gtest.h>
#endif
#include <math.h>
#include <algorithm>
#include <array>
#include <bitset>
#include <cassert>
#include <chrono>
#include <cstdlib>
#include <cstring>
#include <functional>
#include <iomanip>
#include <iostream>
#include <limits>
#include <list>
#include <map>
#include <numeric>
#include <queue>
#include <set>
#include <sstream>
#include <stack>
#include <string>
#include <tuple>
#include <unordered_map>
#include <unordered_set>
#include <vector>
#ifdef PERF
#include <gperftools/profiler.h>
#endif
using namespace std;
using namespace atcoder;
#define fast_io \
ios_base::sync_with_stdio(false); \
cin.tie(0); \
cout.tie(0);
#define ll long long int
#define rep(i, n) for (int i = 0; i < (int)(n); i++)
#define reps(i, n) for (int i = 1; i <= (int)(n); i++)
#define REP(i, n) for (int i = n - 1; i >= 0; i--)
#define REPS(i, n) for (int i = n; i > 0; i--)
#define MOD (long long int)(1e9 + 7)
#define INF (int)(1e9)
#define LINF (long long int)(1e18)
#define chmax(a, b) a = (((a) < (b)) ? (b) : (a))
#define chmin(a, b) a = (((a) > (b)) ? (b) : (a))
#define all(v) v.begin(), v.end()
typedef pair<int, int> Pii;
typedef pair<ll, ll> Pll;
constexpr double PI = acos(-1);
#ifdef NDEBUG
#define CHECK(v1, op, v2)
#else
#define CHECK(v1, op, v2) \
if (!((v1)op(v2))) { \
cerr << "ERROR:" << (v1) << " " << (v2) << endl; \
assert((v1)op(v2)); \
}
#endif
long double nCr(const int n, const int r) {
long double ret = 1;
rep(t, r) {
ret *= (n - t);
ret /= (r - t);
}
return ret;
}
template <typename T>
string to_string(const vector<T>& vec) {
string ret = "";
rep(i, vec.size()) {
ret += vec[i].to_string();
if (i + 1 != vec.size()) {
ret += ",";
}
}
return ret;
}
template <typename T>
ostream& operator<<(ostream& os, const vector<T>& vec) {
os << to_string(vec);
return os;
}
uint32_t xorshift() {
static uint32_t x = 12345789;
static uint32_t y = 362436069;
static uint32_t z = 521288629;
static uint32_t w = 88675123;
uint32_t t;
t = x ^ (x << 11);
x = y;
y = z;
z = w;
w ^= t ^ (t >> 8) ^ (w >> 19);
return w;
}
int rand(const uint32_t l, const uint32_t r) {
return xorshift() % (r - l) + l;
}
uint32_t rand_other_than(const uint32_t l, const uint32_t r,
const uint32_t other) {
const uint32_t num = rand(l, r - 1);
return num + (num >= other);
}
template <typename T>
const T& rand_vec(const vector<T>& vec) {
assert(vec.size() > 0);
return vec[rand(0, vec.size())];
}
template <typename T>
void shuffle(vector<T>& vec) {
rep(l, (int)vec.size() - 1) {
const int idx = rand(l, vec.size());
swap(vec[idx], vec[l]);
}
}
/* start */
vector<double> PARAMS = {};
class Timer {
chrono::system_clock::time_point _start, _end;
ll _sum = 0, _count = 0;
public:
void start() { _start = chrono::system_clock::now(); }
void stop() { _end = chrono::system_clock::now(); }
void add() {
const chrono::system_clock::time_point now = chrono::system_clock::now();
_sum += static_cast<double>(
chrono::duration_cast<chrono::nanoseconds>(now - _start).count());
_count++;
}
ll sum() const { return _sum / 1000; }
int count() const { return _count; }
string average() const {
if (_count == 0) {
return "NaN";
}
return to_string(_sum / 1000 / _count);
}
void reset() {
_start = chrono::system_clock::now();
_sum = 0;
_count = 0;
}
inline int ms() const {
const chrono::system_clock::time_point now = chrono::system_clock::now();
return static_cast<double>(
chrono::duration_cast<chrono::microseconds>(now - _start).count() /
1000);
}
inline int ns() const {
const chrono::system_clock::time_point now = chrono::system_clock::now();
return static_cast<double>(
chrono::duration_cast<chrono::microseconds>(now - _start).count());
}
};
#ifdef LOCAL
struct Timers : unordered_map<string, Timer> {
friend ostream& operator<<(ostream& os, const Timers& timers) {
for (const auto& pa : timers) {
os << pa.first << " time: " << pa.second.sum() / 1000
<< " count: " << pa.second.count() << endl;
}
return os;
}
};
#else
struct Timers {
struct Dummy {
void start() const {}
void add() const {}
};
Dummy dummy;
const Dummy& operator[](const std::string& str) { return dummy; }
friend ostream& operator<<(ostream& os, const Timers& timers) { return os; }
};
#endif
Timers global_timers;
// NBD = neighborhood
template <class Solution, class NBDGenerator>
struct SimulatedAnnealing {
public:
SimulatedAnnealing(double end_time, double start_temp, double end_temp)
: end_time_(end_time),
inv_time_(1.0 / end_time_),
cur_time_(0.0),
start_temp_(start_temp),
end_temp_(end_temp),
diff_temp_(end_temp - start_temp),
cur_temp_(start_temp) {
assert(start_temp >= end_temp);
timer_.start();
}
Solution run(const Solution& initial_sol, NBDGenerator& nbd_generator) {
Solution sol(initial_sol);
int acc = 0;
int g = 0;
for (;; ++g) {
if ((g & 0xf) == 0) {
#ifdef SA_MAX_G
if (g >= SA_MAX_G) {
break;
}
#endif
UpdateTime();
UpdateTemp();
if (cur_time_ > end_time_) {
break;
}
}
const auto nbd = nbd_generator.Generate(&sol, g);
constexpr int mod = 1'000'000'000;
const int r = rand(0, mod);
if (static_cast<double>(r) / mod < exp(-nbd->Diff() / cur_temp_)) {
acc++;
nbd->Update(&sol);
} else {
nbd->Restore(&sol);
}
}
cerr << "g,acc: " << g << " " << acc << endl;
return sol;
}
private:
double end_time_, inv_time_, cur_time_;
double start_temp_, end_temp_, diff_temp_, cur_temp_;
Timer timer_;
void UpdateTime() { cur_time_ = timer_.ms(); }
void UpdateTemp() {
cur_temp_ = max(0.0, start_temp_ - cur_time_ * diff_temp_ * inv_time_);
}
};
/* start */
template <typename Point, typename Dist_T>
class Graph_T {
public:
struct Edge {
Point to;
Dist_T dist;
Edge(const Point _to, const Dist_T _dist) : to(_to), dist(_dist) {}
};
Graph_T() {}
explicit Graph_T(const int n)
: n_(n), adjs_(n), dists_(n, vector<Dist_T>(n)) {}
int N() const { return n_; }
Dist_T Dist(const Point& from, const Point& to) const {
return dists_[from][to];
}
void AddEdge(const Point a, const Point b, const Dist_T dist) {
adjs_[a].emplace_back(b, dist);
dists_[a][b] = dist;
}
const vector<Edge>& Adj(const Point a) const { return adjs_[a]; }
void RemoveEdge(const Point a, const Point b) {
{
const auto itr = std::find_if(
all(adjs_[a]), [&](const Edge& edge) { return edge.to == b; });
assert(itr != adjs_[a].end());
adjs_[a].erase(itr);
dists_[a][b] = INF;
}
}
private:
int n_;
vector<vector<Edge>> adjs_;
vector<vector<Dist_T>> dists_;
};
enum Dir {
// y-1, y+1, x-1, x+1
kU,
kD,
kL,
kR
};
struct Pos {
int idx_;
Pos() {}
explicit Pos(const int _idx) : idx_(_idx) {}
Pos(int _x, int _y) : idx_(_y * N + _x) { assert(N > 0); }
int X() const { return pos_2_x[*this]; }
int Y() const { return pos_2_y[*this]; }
int Idx() const { return idx_; }
operator int() const { return Idx(); }
operator size_t() const { return Idx(); }
const vector<Pos>& Adj() const { return adj_poses[*this]; }
const vector<Dir>& AdjDirs() const { return adj_dirs[*this]; }
int Manhattan(const Pos& other) const {
return abs(X() - other.X()) + abs(Y() - other.Y());
}
int Euclid2(const Pos& other) const {
const int dx = X() - other.X();
const int dy = Y() - other.Y();
return dx * dx + dy * dy;
}
bool Move(const Dir dir) {
if (move_to[dir][*this].IsDummy()) {
return false;
} else {
*this = move_to[dir][*this];
return true;
}
}
bool operator<(const Pos& other) const { return this->Idx() < other.Idx(); }
bool operator==(const Pos& other) const { return this->Idx() == other.Idx(); }
bool operator!=(const Pos& other) const { return this->Idx() != other.Idx(); }
friend ostream& operator<<(ostream& os, const Pos& pos) {
os << pos.X() << " " << pos.Y();
return os;
}
bool IsDummy() const { return this->idx_ < 0; }
static Pos Dummy() {
Pos p;
p.idx_ = -1;
return p;
}
static Pos TryCreate(const int x, const int y, const int z) {
if (y < 0 || y >= N || x < 0 || x >= N) {
return Pos::Dummy();
} else {
return Pos(x, y);
}
}
static void StaticInit(const int n) {
N = n;
N2 = N * N;
N4 = N2 * N2;
adj_poses = vector<vector<Pos>>(N2);
adj_dirs = vector<vector<Dir>>(N2);
move_to = vector<vector<Pos>>(4, vector<Pos>(N2, Pos::Dummy()));
pos_2_y.resize(N2);
pos_2_x.resize(N2);
rep(y, N) {
rep(x, N) {
const Pos p{x, y};
pos_2_y[p] = y;
pos_2_x[p] = x;
for (int dy = -1; dy <= 1; ++dy) {
for (int dx = -1; dx <= 1; ++dx) {
if (abs(dy) + abs(dx) != 1) continue;
const int adj_y = y + dy;
const int adj_x = x + dx;
if (adj_y < 0 || adj_y >= N) continue;
if (adj_x < 0 || adj_x >= N) continue;
adj_poses[p].emplace_back(adj_x, adj_y);
Dir dir;
if (dy == -1) {
dir = kU;
} else if (dy == 1) {
dir = kD;
} else if (dx == -1) {
dir = kL;
} else if (dx == 1) {
dir = kR;
} else {
assert(false);
}
move_to[dir][p] = {adj_x, adj_y};
adj_dirs[p].emplace_back(dir);
}
}
}
}
}
static int N, N2, N4;
static vector<vector<Pos>> adj_poses;
static vector<vector<Dir>> adj_dirs;
static vector<vector<Pos>> move_to;
static vector<int> pos_2_x, pos_2_y;
};
vector<vector<Pos>> Pos::adj_poses;
vector<vector<Dir>> Pos::adj_dirs;
vector<vector<Pos>> Pos::move_to;
vector<int> Pos::pos_2_x, Pos::pos_2_y;
int Pos::N, Pos::N2, Pos::N4;
/* start */
constexpr int kPlanetType = 1;
constexpr int kStationType = 2;
struct Solution {
using Graph = Graph_T<int, int>;
Solution() {
orders.resize(N);
std::iota(all(orders), 0);
stations.resize(M);
rep(i, M) { stations[i] = Pos{rand(0, Pos::N + 1), rand(0, Pos::N + 1)}; }
}
static void StaticInit(const int _N, const int _M,
const vector<Pos>& _planets, const Graph& _graph) {
N = _N;
M = _M;
planets = _planets;
graph = _graph;
}
// [l, r)のreverse
int DiffSwap(const int l, const int r) const {
// l-1, l, r-1, r
// l-1, r-1, l, r
const int i0 = (l - 1 + N) % N;
const int i1 = l % N;
const int i2 = (r - 1 + N) % N;
const int i3 = r % N;
const int prev_dist =
graph.Dist(orders[i0], orders[i1]) + graph.Dist(orders[i2], orders[i3]);
const int next_dist =
graph.Dist(orders[i0], orders[i2]) + graph.Dist(orders[i1], orders[i3]);
return next_dist - prev_dist;
}
void Swap(int l, int r) {
if (l > r) {
swap(l, r);
}
reverse(orders.begin() + l, orders.begin() + r);
}
friend ostream& operator<<(ostream& os, const Solution& sol) {
int start_i = 0;
for (start_i = 0; start_i < N; start_i++) {
if (sol.orders[start_i] == 0) {
break;
}
}
#if 0
// その順序固定だとして、10^6通りの駅について、「その駅を通ったほうが嬉しい」辺のbitsetは何通りあるか？
using BS = bitset<100>;
unordered_set<BS> S;
rep(y, 1001) {
rep(x, 1001) {
const Pos station(x, y);
BS bs = 0;
rep(i, N) {
const int next_i = (i + 1) % N;
const auto p0 = planets[sol.orders[i]];
const auto p1 = planets[sol.orders[next_i]];
const int default_dist =
graph.Dist(sol.orders[i], sol.orders[next_i]);
const int station_dist =
p0.Euclid2(station) * alpha + station.Euclid2(p1) * alpha;
bs[i] = station_dist <= default_dist;
}
S.insert(bs);
}
}
cerr << endl;
cerr << "size: " << S.size() << endl;
#endif
vector<pair<int, int>> ans;
ans.emplace_back(kPlanetType, sol.orders[start_i]);
rep(add_i, N) {
const int i = (start_i + add_i) % N;
const int next_i = (i + 1) % N;
const Pos& from_planet = planets[sol.orders[i]];
const Pos& to_planet = planets[sol.orders[next_i]];
// stationを検討
int min_dist = Solution::graph.Dist(sol.orders[i], sol.orders[next_i]);
vector<int> best_stations;
// i,jの順に経由
rep(i, M) {
const Pos& station_i = sol.stations[i];
rep(j, M) {
if (i == j) continue;
const Pos& station_j = sol.stations[j];
const int dist = from_planet.Euclid2(station_i) * alpha +
station_i.Euclid2(station_j) +
station_j.Euclid2(to_planet) * alpha;
if (dist < min_dist) {
min_dist = dist;
best_stations = {i, j};
}
}
}
// iのみ経由
rep(i, M) {
const Pos& station = sol.stations[i];
const int dist = from_planet.Euclid2(station) * alpha +
station.Euclid2(to_planet) * alpha;
if (dist < min_dist) {
min_dist = dist;
best_stations = {i};
}
}
if (best_stations.size() == 0) {
// 惑星to惑星
int from = sol.orders[i];
int to = sol.orders[next_i];
vector<int> rev_planets;
while (from != to) {
rev_planets.emplace_back(to);
rep(i, N) {
if (i == to) continue;
if (graph.Dist(from, i) + graph.Dist(i, to) ==
graph.Dist(from, to)) {
to = i;
break;
}
}
}
REP(i, rev_planets.size()) {
ans.emplace_back(kPlanetType, rev_planets[i]);
}
} else {
for (const auto station : best_stations) {
ans.emplace_back(kStationType, station);
}
ans.emplace_back(kPlanetType, sol.orders[next_i]);
}
}
rep(i, M) { os << sol.stations[i] << endl; }
os << ans.size() << endl;
rep(i, ans.size()) {
const auto& [t, r] = ans[i];
os << t << " " << r + 1;
if (i + 1 != (int)ans.size()) {
os << endl;
}
}
int dist_sum = 0;
rep(i, (int)ans.size() - 1) {
const auto& [t0, r0] = ans[i];
const auto& [t1, r1] = ans[(i + 1) % ans.size()];
const Pos& from = (t0 == kPlanetType ? planets[r0] : sol.stations[r0]);
const Pos& to = (t1 == kPlanetType ? planets[r1] : sol.stations[r1]);
dist_sum += from.Euclid2(to) * (t0 == kPlanetType ? alpha : 1) *
(t1 == kPlanetType ? alpha : 1);
}
cerr << "score: " << std::round(1e9 / (1000 + sqrt(dist_sum))) << endl;
return os;
}
static vector<Pos> planets;
static int N, M;
static Graph graph;
static constexpr int alpha = 5;
static constexpr int alpha2 = 25;
vector<int> orders;
vector<Pos> stations;
};
int Solution::N;
int Solution::M;
vector<Pos> Solution::planets;
Graph_T<int, int> Solution::graph;
/* start */
struct NBD {
int Diff() const { return diff_; }
virtual void Update(Solution*) const = 0;
void Restore(Solution*) const {}
int diff_;
};
struct NBD_swap : public NBD {
NBD_swap() {}
NBD_swap(Solution* const sol, const int l, const int r) {
l_ = l;
r_ = r;
diff_ = sol->DiffSwap(l, r);
}
void Update(Solution* sol) const override { sol->Swap(l_, r_); }
int l_;
int r_;
};
struct NBDGenerator {
NBDGenerator() {}
NBD* Generate(Solution* const sol, const int g) {
// swap
const int l = rand(0, Solution::N);
const int r = rand_other_than(0, Solution::N, l);
nbd_swap0 = NBD_swap(sol, l, r);
return &nbd_swap0;
}
NBD_swap nbd_swap0;
};
class Solver {
public:
using Graph = Graph_T<int, int>;
Solver(istream& is) {
is >> N >> M;
Pos::StaticInit(1001);
planets.reserve(N);
rep(i, N) {
int a, b;
is >> a >> b;
planets.emplace_back(a, b);
}
graph = InitGraph();
Solution::StaticInit(N, M, planets, graph);
}
Graph InitGraph() const {
vector<vector<int>> dists(N, vector<int>(N, INF));
rep(i, N) {
rep(j, N) { dists[i][j] = alpha2 * planets[i].Euclid2(planets[j]); }
}
rep(k, N) {
rep(i, N) {
rep(j, N) { chmin(dists[i][j], dists[i][k] + dists[k][j]); }
}
}
Graph graph(N);
rep(i, N) {
rep(j, N) {
const int dist = dists[i][j];
graph.AddEdge(i, j, dist);
}
}
return graph;
}
using BS = bitset<100>;
unordered_map<BS, pair<Pos, ll>> CalcBestBSs(
const vector<int>& orders) const {
// その順序固定だとして、10^6通りの駅について、「その駅を通ったほうが嬉しい」辺のbitsetを求めて、
// 各bitsetごとに、距離の改善量が最も大きいものを保持する
// TODO: 駅to駅を考慮
// TODO: intにするかも
unordered_map<BS, pair<Pos, ll>> mp;
rep(y, 1001) {
rep(x, 1001) {
const Pos station(x, y);
BS bs = 0;
ll diff_sum = 0;
rep(i, N) {
const int next_i = (i + 1) % N;
const auto p0 = planets[orders[i]];
const auto p1 = planets[orders[next_i]];
const int default_dist = graph.Dist(orders[i], orders[next_i]);
const int station_dist =
p0.Euclid2(station) * alpha + station.Euclid2(p1) * alpha;
const ll diff = default_dist - station_dist;
if (diff > 0) {
bs[i] = true;
diff_sum += diff;
}
}
if (diff_sum == 0) continue;
if (mp[bs].second < diff_sum) {
mp[bs] = {Pos{x, y}, diff_sum};
}
}
}
cerr << "size: " << mp.size() << endl;
return mp;
}
Solution DPStations(const Solution& order_sol) const {
const auto orders = order_sol.orders;
const auto mp = CalcBestBSs(orders);
// dp[j][k] := i個目まで見て、j個駅を置いて、bit_countがk個
using P = pair<BS, ll>;
using IDX = tuple<int, int, int>;
vector dp(M + 1, vector(N + 1, P{BS{0}, (int)-INF}));
vector prev(mp.size() + 1, vector(M + 1, vector(101, IDX{0, 0, 0})));
dp[0][0].second = 0;
vector<Pos> pos_candidates;
pos_candidates.reserve(mp.size());
int i = -1;
for (const auto& [bs, pos_diff] : mp) {
++i;
const int bs_count = bs.count();
const auto& [pos, diff] = pos_diff;
pos_candidates.emplace_back(pos);
prev[i + 1] = prev[i];
REP(j, M) {
const auto next_j = j + 1;
REP(k, N + 1 - bs_count) {
if (dp[j][k].second < 0) continue;
if ((dp[j][k].first & bs) != 0) continue;
const auto next_bs = (dp[j][k].first | bs);
const auto next_k = next_bs.count();
const auto next_val = dp[j][k].second + diff;
if (dp[next_j][next_k].second < next_val) {
dp[next_j][next_k].first = next_bs;
dp[next_j][next_k].second = next_val;
prev[i + 1][next_j][next_k] = make_tuple(i, j, k);
}
}
}
}
IDX best_idx{0, 0, 0};
ll best_diff = 0;
rep(j, M + 1) {
rep(k, N + 1) {
const auto diff = dp[j][k].second;
if (best_diff < diff) {
best_diff = diff;
best_idx = make_tuple((int)mp.size(), j, k);
}
}
}
// 復元
vector<Pos> stations;
IDX idx = best_idx;
while (idx != make_tuple(0, 0, 0)) {
const auto& [i, j, k] = idx;
const auto& prev_idx = prev[i][j][k];
const auto& [prev_i, prev_j, prev_k] = prev_idx;
if (prev_j != j) {
stations.emplace_back(pos_candidates[prev_i]);
}
idx = prev_idx;
}
cerr << "station size: " << stations.size() << endl;
while (stations.size() < M) {
stations.emplace_back(Pos{500, 500});
}
auto sol = order_sol;
sol.stations = stations;
return sol;
}
Solution Solve(const int time_limit) const {
const Solution initial_sol;
SimulatedAnnealing<Solution, NBDGenerator> sa(time_limit - 500, 1, 0.01);
NBDGenerator nbd_generator;
const auto order_sol = sa.run(initial_sol, nbd_generator);
const auto final_sol = DPStations(order_sol);
return final_sol;
}
private:
int N, M;
vector<Pos> planets;
Graph graph;
static constexpr int alpha = 5;
static constexpr int alpha2 = 25;
};
int main(int argc, char* argv[]) {
fast_io;
if (argc >= 2) {
int idx = 0;
for (int i = 1; i < argc; ++i) {
PARAMS[idx++] = std::stod(argv[i]);
}
}
Timer timer;
timer.start();
Solver solver(cin);
auto sol = solver.Solve(900 - timer.ms());
cout << sol << endl;
return 0;
}