#1049648 (C++23) No.5007 Steiner Space Travel

• 問題: No.5007 Steiner Space Travel
• ユーザー: ぴぃいいいい
• 提出日時: 2025-03-06 03:34:44
• 言語: C++23 (gcc 13.3.0 + boost 1.87.0)
• 結果: AC
• 実行時間: 955 ms / 1,000 ms
• コード長: 8,179 bytes
• コンパイル時間: 7,657 ms
• コンパイル使用メモリ: 357,664 KB
• 実行使用メモリ: 8,608 KB
• スコア: 8,354,470
• 最終ジャッジ日時: 2025-03-06 03:35:23
• 合計ジャッジ時間: 39,077 ms
• ジャッジサーバーID （参考情報）: judge5 / judge2

ファイルパターン	結果
other	AC * 30

ソースコード

#include<bits/stdc++.h>
#include<atcoder/all>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/tree_policy.hpp>
using namespace std;
using namespace atcoder;
using namespace __gnu_pbds;
using mint = modint998244353;
#define Yes(n) cout << ((n) ? "Yes" : "No"  ) << endl
#define print(var) std::cout<<#var<<"="<<(var)<<std::endl
#define all(a) (a).begin(), (a).end()
#define vi vector<int>
#define vvi vector<vi>
#define vvvi vector<vvi>
#define ll long long
#define vll vector<ll>
#define vvll vector<vll>
#define vvvll vector<vvll>
#define vvvvll vector<vvvll>
#define vmi vector<mint>
#define vvmi vector<vmi>
#define vvvmi vector<vvmi>
#define vvvvmi vector<vvvmi>
#define vvvvvmi vector<vvvvmi>
#define vs vector<string>
#define pii pair<int,int>
#define vpii vector<pii>
#define vvpii vector<vpii>
#define bit(x,i)(((x)>>(i))&1)
#define inf (1<<30)
#define INF (1ll<<60)
#define X first
#define Y second
template<typename T> inline bool chmax(T &a, T b) { return ((a < b) ? (a = b, true) : (false)); }
template<typename T> inline bool chmin(T &a, T b) { return ((a > b) ? (a = b, true) : (false)); }
template<class T, class F> T nibutan(T ok, T ng, const F &f){while(abs(ok-ng)>1){T mid = (ok+ng)/2;(f(mid)?ok:ng) = mid;}return ok;}
template<class T> vector<T> digit(T x){vector<T> res; while(x>0){res.push_back(x%10); x/=10;} return res;}
ostream& operator<<(ostream& os, const mint& x){ os << x.val(); return os; }
template<class T> istream &operator>>(istream &is, vector<T> &vec){ for(auto &v:vec) is >> v; return is; }
template<class T> void coutvector(vector<T> x){ for(int i=0;i<(int)x.size();i++){if(i>0) cout<<" ";cout<<x[i];}cout<<endl;}
using Set = tree<int, null_type, less<int>, rb_tree_tag, tree_order_statistics_node_update>;
template<class S, class T> using Unordered_Map = gp_hash_table<S,T>;

#ifdef LOCAL
#include<dump.hpp>
CPP_DUMP_DEFINE_EXPORT_OBJECT(mint,val());
#else
#define cpp_dump
#endif

inline double get_time() {
    #ifdef ONLINE_JUDGE
    timeval tv; gettimeofday(&tv, NULL);
    return tv.tv_sec * 1e3 + tv.tv_usec * 1e-3;
    #else
    using namespace std::chrono;
    return duration_cast<milliseconds>(system_clock::now().time_since_epoch()).count();
    #endif
}
double start_time = get_time();
inline double elapsed(){ return get_time() - start_time; }

static uint32_t xor128(void){
    static uint32_t x=123456789,y=362436069,z=521288629,w=88675123;
    uint32_t t;
    t=(x^(x<<11));x=y;y=z;z=w; return( w=(w^(w>>19))^(t^(t>>8)) );
}
inline int randrange(int a){return (uint64_t(xor128()) * a >> 32); }
inline int randrange(int a, int b){ return randrange(b-a)+a; }
inline double randDouble(double a,double b){return a+(b-a)*xor128()/(double)ULONG_MAX;}

const int N = 100, M = 8;
vector<int> a(N),b(N);

void in(){
    int tmp; cin>>tmp>>tmp;
    for(int i=0;i<N;i++){
        cin>> a[i] >> b[i];
    }
}
int dist(int x, int y){
    return (a[x]-a[y])*(a[x]-a[y]) + (b[x]-b[y])*(b[x]-b[y]);
}
int dist(int x, pii p){
    return (a[x]-p.X)*(a[x]-p.X) + (b[x]-p.Y)*(b[x]-p.Y);
}
int dist(pii p1, pii p2){
    return (p1.X-p2.X)*(p1.X-p2.X) + (p1.Y-p2.Y)*(p1.Y-p2.Y);
}

vi init_sol(){
    vi ret(N); iota(all(ret),0);
    ret.push_back(0);
    return ret;
}

vpii kmeans(const vpii &p) {

    // 初期化：とりあえず先頭のM個の点を初期重心にしとくぜ！
    vector<pair<double, double>> center;
    for (int i=0;i<M;i++) {
        center.push_back(p[i]);
    }

    vector<int> assignments(N, -1);
    bool changed = true;
    int iter = 0;
    const int max_iter = 100;

    while (changed && iter < max_iter) {
        changed = false;
        for (int i = 0; i < N; i++) {
            double minDist = numeric_limits<double>::max();
            int best = -1;
            for (int j = 0; j < M; j++) {
                double dx = p[i].X - center[j].X;
                double dy = p[i].Y - center[j].Y;
                double dist = dx * dx + dy * dy;
                if (dist < minDist) {
                    minDist = dist;
                    best = j;
                }
            }
            if (assignments[i] != best) {
                assignments[i] = best;
                changed = true;
            }
        }
        vector<pair<double, double>> nxt_center(M, {0.0, 0.0});
        vector<int> cnt(M, 0);
        for (int i = 0; i < N; i++) {
            int cluster = assignments[i];
            nxt_center[cluster].X += p[i].X;
            nxt_center[cluster].Y += p[i].Y;
            cnt[cluster]++;
        }
        for (int j = 0; j < M; j++) {
            if (cnt[j] > 0) {
                nxt_center[j].X /= cnt[j];
                nxt_center[j].Y /= cnt[j];
            } else {
                nxt_center[j] = center[j];
            }
        }
        center = nxt_center;
        iter++;
    }

    vpii result;
    for (int j = 0; j < M; j++) {
        int x = static_cast<int>(round(center[j].X));
        int y = static_cast<int>(round(center[j].Y));
        result.push_back({x, y});
    }
    return result;
}

vvi dist_table(const vpii &sta){
    vvi ret(N+M,vi(N+M,0));
    for(int i=0;i<N;i++)for(int j=0;j<N;j++){
        ret[i][j] = 25*dist(i,j);
    }
    for(int i=0;i<N;i++)for(int j=0;j<M;j++){
        ret[i][N+j] = ret[N+j][i] = 5*dist(i,sta[j]);
    }
    for(int i=0;i<M;i++)for(int j=0;j<M;j++){
        ret[N+i][N+j] = dist(sta[i],sta[j]);
    }
    
    for(int k=N;k<N+M;k++)for(int i=0;i<N+M;i++)for(int j=0;j<N+M;j++){
        chmin(ret[i][j],ret[i][k]+ret[k][j]);
    }

    return ret;
}

int solve(const vpii &sta, bool output=false){
    vvi d = dist_table(sta);

    vi tour = init_sol();
    int now_score = 0;
    for(int i=0;i<N;i++){
        now_score += d[tour[i]][tour[i+1]];
    }
    int best_score = now_score;

    int iter = 0;
    int last_update = 0;
    while(elapsed() < 990){
        int l,r;
        do{
            l = randrange(1,N);
            r = randrange(1,N);
        }while(l==r);
        if(l>r) swap(l,r);
        int diff = 0;
        diff -= d[tour[l-1]][tour[l]] + d[tour[r]][tour[r+1]];
        diff += d[tour[l-1]][tour[r]] + d[tour[l]][tour[r+1]];

        if(diff < 0){
            reverse(tour.begin()+l,tour.begin()+r+1);
            now_score += diff;
            if(chmin(best_score,now_score)){
                best_score = now_score;
                last_update = iter;
            }
        }
        iter++;
        if(iter - last_update > 10000){
            break;
        }
    }

    if(!output) return best_score;

    auto raw_d = [&](int x, int y){
        if(x < N && y < N) return 25*dist(x,y);
        if(x < N && y >= N) return 5*dist(x,sta[y-N]);
        if(x >= N && y < N) return 5*dist(y,sta[x-N]);
        return dist(sta[x-N],sta[y-N]);
    };

    vi sol{0};
    for(int i=1;i<(int)tour.size();i++){
        // tour[i-1] -> tour[i]までのルートを計算
        int now = tour[i-1];
        while(now!=tour[i]){
            int nxt = -1;
            for(int j=0;j<N+M;j++){
                if(j==now) continue;
                if(raw_d(now,j) + d[j][tour[i]] == d[now][tour[i]]){
                    nxt = j;
                    break;
                }
            }
            assert(nxt!=-1);
            sol.push_back(nxt);
            now = nxt;
        }
    }

    cerr << "best dist : " << best_score << endl;
    cerr << "best raw score : " << int(1e9/(1000+sqrt(best_score))) << endl;

    for(int i=0;i<M;i++) cout << sta[i].X << " " << sta[i].Y << endl;
    cout << sol.size() << endl;
    for(int i=0;i<sol.size();i++){
        if(sol[i] < N) cout << 1 << " " << sol[i]+1 << endl;
        else cout << 2 << " " << sol[i]-N+1 << endl;
    }

    return best_score;
}

int main(){
    cin.tie(nullptr);
    ios::sync_with_stdio(false);
    in();
    int best_score = inf;
    vpii best_sta;
    do{
        vpii p(N);
        for(int i=0;i<N;i++) p[i] = make_pair(a[i],b[i]);
        vpii sta = kmeans(p);
        int nxt_score = solve(sta);

        if(chmin(best_score,nxt_score)){
            best_sta = sta;
            cpp_dump(nxt_score);
        }
    }while(elapsed() < 950);

    cpp_dump(best_score);
    solve(best_sta,true);
    cpp_dump(elapsed());
}