// 1.惑星だけでNearstNeighbor法 // 2.k-means++法でクラスタリングして宇宙ステーションを配置 // 2.calc_dist_planetの距離計算が間違ていたので修正 // 3.各辺の間に宇宙ステーションの挿入を試みる(cuの所属クラスタの宇宙ステーションの経由のみを考慮する簡易版) // 4.2-opt // 5.焼きなまし,calc_score,2-opt // 6.Or-opt // 7.挿入 // 8.Or-optの間違いを修正 // 9.宇宙ステーションを移動 // 10.宇宙ステーションの配置を最適化 // 11.宇宙ステーションの再接続(このバージョンでは起動をやめた),insert_stationで2種類の宇宙ステーションを考慮するよう改造 // 12.上のバージョンで宇宙ステーションの再接続をやめていなかったのでこのバージョンでオフ // 13.右によりバージョン10にもどした,10reconnect_stationをオフにした,insert_stationで2種類の宇宙ステーションを考慮するよう改造したものを戻した // 14.get_initial_solutionでまず限界まで2-opt,次に宇宙ステーションを挿入するよう変更 #include #include #include #include #include #include #include #include #include #include #define rep(i,a,b) for(int i=a;i #define vvi vector #define vvvi vector #define vvvvi vector #define vp vector

template bool chmin(T &a, const T& b){if(a>b){a=b; return true;}return false;} template bool chmax(T &a, const T& b){if(a; template using min_priority_queue = priority_queue, greater>; //const long long INF=1e18; const int INF=1000000000; const int64_t CYCLES_PER_SEC = 2300000000; const double TIMELIMIT = 0.95; constexpr int RANDOM_SEED = 42; // 乱数のシード値 int my_clamp(int x, int l, int r){ if(xr) x=r; return x; } struct Timer { int64_t start; Timer() { reset(); } void reset() { start = getCycle(); } void plus(double a) { start -= (a * CYCLES_PER_SEC); } inline double get() { return (double)(getCycle() - start) / CYCLES_PER_SEC; } inline int64_t getCycle() { uint32_t low, high; __asm__ volatile ("rdtsc" : "=a" (low), "=d" (high)); return ((int64_t)low) | ((int64_t)high << 32); } }; mt19937 mt(RANDOM_SEED); uniform_int_distribution<> rand8(0,7); uniform_int_distribution<> rand2050(20,50); uniform_int_distribution<> rand01(0,1); struct Input{ int N,M; vector

coord; Input(int N=0,int M=0, vector

coord={}): N(N), M(M), coord(coord) {} }; struct Solver { Timer timer; // グローバル変数============================== Input input; ll best_score; vector

best_ans; // 経由順 vector

best_station; // 宇宙ステーションの座標 vvi dist_planet; // 全惑星間の距離 // 便利関数============================== // 距離の2乗を求める template T dist2(T a,T b){return a*a+b*b;} // ans上の2点間の距離を求める int ans_to_dist2(Input &input, vector

&ans, vector

&station, int i, int j){ int t1,idx1,t2,idx2; tie(t1,idx1)=ans[i]; tie(t2,idx2)=ans[j]; // 両方惑星なら早期return if(t1==1 and t2==1) return dist_planet[idx1-1][idx2-1]; int x1,y1,x2,y2; if(t1==1) tie(x1,y1)=input.coord[idx1-1]; else tie(x1,y1)=station[idx1-1]; if(t2==1) tie(x2,y2)=input.coord[idx2-1]; else tie(x2,y2)=station[idx2-1]; int d=dist2(x2-x1,y2-y1); if(t1==1) d*=5; // 惑星は5倍 if(t2==1) d*=5; // 惑星は5倍 return d; } // 区間を逆回りにする void reverse_section(vector

&ans, int l, int r){ while(l &ans, vector

&station, int i, int j){ // 距離を求める int d1 = ans_to_dist2(input,ans,station,i,i+1) + ans_to_dist2(input,ans,station,j,j+1); int d2 = ans_to_dist2(input,ans,station,i,j) + ans_to_dist2(input,ans,station,i+1,j+1); // 改善できるなら入れ替える if(d2 &ans, vector

&station){ int cnt=0; // 交換回数をカウント rep(i,0,sz(ans)-3){ rep(j,i+2,sz(ans)-1){ // 始点と終点は同一なので-1 if(two_opt(input,ans,station,i,j)) cnt++; } } return cnt; } // k-means++法でM個のクラスタ中心を求める vi select_cluster_centers(Input &input){ vi centers; // クラスタ中心の集合 vi dist(input.N,INF); // 最寄のクラスタ中心までの距離 // データ点からランダムに1つ選びクラスタ中心とする uniform_int_distribution<> randN(0,input.N-1); int center=randN(mt); centers.emplace_back(center); // クラスタ中心を追加 // 残りのM-1個のクラスタ中心を求める rep(i,0,input.M-1){ int cx,cy; tie(cx,cy)=input.coord[center]; // 各点について最寄のクラスタ中心までの距離を求める rep(j,0,input.N){ int x,y; tie(x,y)=input.coord[j]; chmin(dist[j],dist2(x-cx,y-cy)); } // 最寄のクラスタ中心までの距離の内,最長のものを持つ点を求める int d=0; rep(j,0,input.N){ if(chmax(d,dist[j])){ center=j; } } centers.emplace_back(center); // クラスタ中心を追加 } return centers; } // k-means法でクラスタリング vi k_means(Input &input, vector> ¢er_pos){ vi cluster_num(input.N,-1); // 各点が属するクラスタ番号 while(1){ vi cluster_num_buf(input.N,-1); vector mn_dist(input.N,INF); // 各点の最寄クラスタまでの距離 // 各点を最寄のクラスタ中心に分類する rep(i,0,input.N){ // 点をまわす int x,y; tie(x,y)=input.coord[i]; rep(j,0,input.M){ // クラスタ中心をまわす double cx,cy; tie(cx,cy)=center_pos[j]; // 距離を更新 if(chmin(mn_dist[i],dist2(x-cx,y-cy))){ cluster_num_buf[i]=j; } } } // 終了判定,更新 if(cluster_num==cluster_num_buf) break; // 変化がなくなったら抜ける swap(cluster_num,cluster_num_buf); // 新しいクラスタ中心を求める rep(i,0,input.M) center_pos[i]=make_pair(0,0); // 初期化 vi point_cnt(input.M,0); // そのクラスタに所属する点の個数 // 座標を加算 rep(i,0,input.N){ int c_idx=cluster_num[i]; int x,y; tie(x,y)=input.coord[i]; center_pos[c_idx].first+=x; center_pos[c_idx].second+=y; point_cnt[c_idx]++; // クラスタに所属する点の個数をインクリメント } // 点の数で除算 rep(i,0,input.M){ center_pos[i].first=round((double)center_pos[i].first/point_cnt[i]); center_pos[i].second=round((double)center_pos[i].second/point_cnt[i]); } } return cluster_num; } // 宇宙ステーションを配置,各惑星のクラスター番号も返す vector

locate_station(Input &input,vi &cluster_num){ // k-means++ vi centers=select_cluster_centers(input); // M個のクラスタ中心を求める vector> center_pos(input.M); // 各クラスタ中心の座標 // 0初期化されるか rep(i,0,input.M){ center_pos[i]=input.coord[centers[i]]; } // k-means cluster_num=k_means(input,center_pos); // cluster_numは各惑星のクラスター番号,center_posも変わる // center_posをdoubleからintにする vector

station_pos(input.M); rep(i,0,input.M){ station_pos[i].first=round(center_pos[i].first); station_pos[i].second=round(center_pos[i].second); } return station_pos; } // 各辺の間に宇宙ステーションの挿入を試みる(mode=0が従来の片側の宇宙ステーションのみ考慮するバージョン) vector

insert_station(Input &input, vector

&ans, vector

station, vi cluster_num, int mode){ vector

res; res.emplace_back(ans[0]); // スタートはそのまま追加 // 各径路間にステーションを挿入するか判定 rep(i,0,sz(ans)-1){ int planet_idx_cu=ans[i].second-1; int planet_idx_nx=ans[i+1].second-1; // 距離が有利なら挿入する int station_num_cu=cluster_num[planet_idx_cu]; // 宇宙ステーションの番号 int station_num_nx=cluster_num[planet_idx_nx]; int px1,py1,px2,py2,sx1,sy1,sx2,sy2; tie(px1,py1)=input.coord[planet_idx_cu]; // 惑星の座標 tie(px2,py2)=input.coord[planet_idx_nx]; // 惑星の座標 tie(sx1,sy1)=station[station_num_cu];// 宇宙ステーションの座標 tie(sx2,sy2)=station[station_num_nx]; int cu_d=dist_planet[planet_idx_cu][planet_idx_nx]; // 最寄の宇宙ステーションが同一の場合 if(mode==0 or station_num_cu==station_num_nx){ int d=dist2(px1-sx1,py1-sy1)*5 + dist2(px2-sx1,py2-sy1)*5; // 宇宙ステーションを経由する場合のコスト if(cu_d>d){ res.emplace_back(P(2,station_num_cu+1)); } // 最寄の宇宙ステーションが異なる場合 }else{ // 宇宙ステーション1のみ経由 int nx_d1=dist2(px1-sx1,py1-sy1)*5 + dist2(px2-sx1,py2-sy1)*5; // 宇宙ステーション2のみ経由 int nx_d2=dist2(px1-sx2,py1-sy2)*5 + dist2(px2-sx2,py2-sy2)*5; // 宇宙ステーション1,2を経由 int nx_d12=dist2(px1-sx1,py1-sy1)*5 + dist2(sx1-sx2,sy1-sy2) + dist2(px2-sx2,py2-sy2)*5; // 最小コストのパターンを求める int mn_d=cu_d; int pat=0; if(chmin(mn_d,nx_d1)) pat=1; if(chmin(mn_d,nx_d2)) pat=2; if(chmin(mn_d,nx_d12)) pat=3; if(pat==1){ res.emplace_back(P(2,station_num_cu+1)); }else if(pat==2){ res.emplace_back(P(2,station_num_nx+1)); }else if(pat==3){ res.emplace_back(P(2,station_num_cu+1)); res.emplace_back(P(2,station_num_nx+1)); } } // 次の惑星に行く res.emplace_back(ans[i+1]); } return res; } // スコア計算 ll calc_score(Input &input, vector

&ans, vector

&station){ int sum=0; rep(i,0,sz(ans)-1){ sum+=ans_to_dist2(input,ans,station,i,i+1); } return round(1000000000/(1000+sqrt(sum))); } // 全惑星間の距離(a^2*D^2)を求める vvi calc_dist_planet(Input &input){ vvi res(input.N,vi(input.N,0)); rep(i,0,input.N-1){ int x1,y1; tie(x1,y1)=input.coord[i]; rep(j,i+1,input.N){ int x2,y2; tie(x2,y2)=input.coord[j]; int d=dist2(x2-x1,y2-y1)*25; res[i][j]=d; res[j][i]=d; } } return res; } // 初期解を生成する vector

get_initial_solution(Input &input, vector

&station, vvi &dist_planet){ vector

res; // 宇宙ステーションの位置を決める,各惑星のクラスタ番号を得る===== vi cluster_num; // 各惑星のクラスタ番号 station=locate_station(input,cluster_num); // 惑星だけでNearstNeighbor法===== vi visited(input.N,0); visited[0]=1; int cu=0; res.emplace_back(P(1,1)); // 最初の惑星からスタート rep(i,0,input.N-1){ // 最短距離の惑星を求める int nx=-1; int mn_d=INF; rep(i,0,input.N){ if(visited[i]) continue; if(chmin(mn_d,dist_planet[cu][i])){ nx=i; } } res.emplace_back(P(1,nx+1)); visited[nx]=1; swap(cu,nx); } res.emplace_back(P(1,1)); // 最初の惑星に戻る // 2-opt while(1){ if(two_opt_all(input,res,station)==0) break; } // 各辺の間に宇宙ステーションを挿入 res=insert_station(input,res,station,cluster_num,0); return res; } // 2-opt(make_new_ans用) void mna_two_opt(vector

&ans, vector

&station){ uniform_int_distribution<> rand_ans(0,sz(ans)-2); int a=rand_ans(mt); int b=rand_ans(mt); while(abs(a-b)<=1) b=rand_ans(mt); two_opt(input,ans,station,a,b); } // Or-opt(make_new_ans用) void mna_or_opt(Input &input, vector

&ans, vector

&station){ uniform_int_distribution<> rand_ans(0,sz(ans)-2); // 移動対象の区間を選択(区間長は1に固定) // 0は選ばない,sz(ans)-2は選ばない int tar_section=rand_ans(mt); while(tar_section==0 or tar_section==sz(ans)-2) tar_section=rand_ans(mt); // 移動先を選択 // tar_section-1,tar_section,tar_section+1は選ばない int saki=rand_ans(mt); while(tar_section-1<=saki and saki<=tar_section+1) saki=rand_ans(mt); int a=tar_section-1; int b=tar_section; int c=tar_section+1; int d=tar_section+2; int e=saki; int f=saki+1; int cu_d=ans_to_dist2(input,ans,station,a,b)+ans_to_dist2(input,ans,station,c,d)+ans_to_dist2(input,ans,station,e,f); int nx_d1=ans_to_dist2(input,ans,station,a,d)+ans_to_dist2(input,ans,station,e,b)+ans_to_dist2(input,ans,station,c,f); int nx_d2=ans_to_dist2(input,ans,station,a,d)+ans_to_dist2(input,ans,station,e,c)+ans_to_dist2(input,ans,station,b,f); // 反転あり bool rev=nx_d1>nx_d2; // 区間反転した方が良い場合 int mn_nx_d=min(nx_d1,nx_d2); // 更新可能な場合 if(mn_nx_d < cu_d){ P tar1,tar2; // 反転を伴う場合 if(rev){ tar1=ans[c]; tar2=ans[b]; // 反転を伴わない場合 }else{ tar1=ans[b]; tar2=ans[c]; } // 更新 vector

res; rep(i,0,sz(ans)){ if(i==b or i==c) continue; // 移動対象はスキップ res.emplace_back(ans[i]); if(i==e){ res.emplace_back(tar1); res.emplace_back(tar2); } } swap(ans,res); } } // 挿入 void mna_insert(Input &input, vector

&ans, vector

&station){ // 挿入対象の点を選ぶ // 始点,終点以外の点をランダムに指定 uniform_int_distribution<> rand_ans(0,sz(ans)-2); int tar=rand_ans(mt); while(tar==0) tar=rand_ans(mt); // 被挿入対象の点を選ぶ int insert_tar=rand_ans(mt); while(insert_tar==tar-1 or insert_tar==tar) insert_tar=rand_ans(mt); // 1.その点をなくして別の場所に入れる // 2.その点はそのままに別の場所に入れる int a=tar-1; int b=tar; int c=tar+1; int d=insert_tar; int e=insert_tar+1; int cu_d=ans_to_dist2(input,ans,station,a,b)+ans_to_dist2(input,ans,station,b,c)+ans_to_dist2(input,ans,station,d,e); int nx_d1=ans_to_dist2(input,ans,station,a,c)+ans_to_dist2(input,ans,station,d,b)+ans_to_dist2(input,ans,station,b,e); // 削除あり int nx_d2=ans_to_dist2(input,ans,station,a,b)+ans_to_dist2(input,ans,station,b,c)+ans_to_dist2(input,ans,station,d,b)+ans_to_dist2(input,ans,station,b,e); // 削除なし bool is_remove_b=false; // 1.何もしない if(cu_d<=nx_d1 and cu_d<=nx_d2){ return; // 2.点をなくして入れる }else if(nx_d1 res; rep(i,0,sz(ans)){ if(is_remove_b and i==b) continue; // 挿入対象はスキップ res.emplace_back(ans[i]); if(i==d){ res.emplace_back(ans[b]); } } swap(ans,res); } // 宇宙ステーションの移動 void mna_move_station(vector

&station){ // 宇宙ステーションのうち1つを指定 int tar=rand8(mt); // 適当に動かす int _dx=rand2050(mt); // 移動距離 int _dy=rand2050(mt); int pm_x=rand01(mt); // 方向 int pm_y=rand01(mt); if(pm_x==0) _dx*=-1; if(pm_y==0) _dy*=-1; int x,y; tie(x,y)=station[tar]; x=my_clamp(x+_dx,0,1000); y=my_clamp(y+_dy,0,1000); station[tar]=P(x,y); } // 交換 // λ-opt // pre_ansからnew_ansを作る void make_new_ans(Input &input, vector

&ans, vector

&station){ uniform_int_distribution<> rand_mna(0,3); int sel=rand_mna(mt); if(sel==0) mna_two_opt(ans,station); else if(sel==1) mna_or_opt(input,ans,station); else if(sel==2) mna_insert(input,ans,station); else if(sel==3) mna_move_station(station); } // 最後に宇宙ステーションの配置を最適化する void opt_station(Input &input, vector

&ans, vector

&station){ // 宇宙ステーションに接続された惑星の配置を加算 vector

station_sum(8); vi station_cnt(8); // 各宇宙ステーションにいくつの経路が加算されたか記録 rep(i,0,sz(ans)-1){ int t1,idx1,t2,idx2; tie(t1,idx1)=ans[i]; tie(t2,idx2)=ans[i+1]; if(t1==t2) continue; // 惑星と宇宙ステーションの間の経路以外はスキップ // 惑星の座標を加算する int planet_idx=idx1-1; int station_idx=idx2-1; if(t2==1) swap(planet_idx,station_idx); // 2番目の方が惑星なら入れ替える int x,y; tie(x,y)=input.coord[planet_idx]; station_sum[station_idx].first+=x; station_sum[station_idx].second+=y; station_cnt[station_idx]++; } // 座標の総和を経路数で割る rep(i,0,8){ if(station_cnt[i]==0) continue; // 0除算回避 station[i].first=round((double)station_sum[i].first/station_cnt[i]); station[i].second=round((double)station_sum[i].second/station_cnt[i]); } } // 宇宙ステーションを経由すべきか考え直す void reconnect_station(Input &input, vector

&ans, vector

&station){ // 宇宙ステーションなしの経路を求める vector

res; rep(i,0,sz(ans)){ if(ans[i].first==2) continue; // 宇宙ステーションはスキップ res.emplace_back(ans[i]); // 惑星のみ追加 } // 各惑星の最寄の宇宙ステーションの番号を求める vi nearest_station=get_nearest_station(input,station); // 各辺の間に宇宙ステーションを挿入 res=insert_station(input,res,station,nearest_station,1); swap(ans,res); } // 各惑星の最寄の宇宙ステーションの番号を求める(station位置はint型) vi get_nearest_station(Input &input, vector

&station){ vi cluster_num(input.N,-1); // 各点を最寄の宇宙ステーションに分類する rep(i,0,input.N){ // 点をまわす int x,y; tie(x,y)=input.coord[i]; int mn_dist=INF; // 最寄の宇宙ステーションまでの距離 rep(j,0,input.M){ // 宇宙ステーションをまわす int cx,cy; tie(cx,cy)=station[j]; // 距離を更新 if(chmin(mn_dist,dist2(x-cx,y-cy))){ cluster_num[i]=j; } } } return cluster_num; } // 初期解生成 void init(){ dist_planet=calc_dist_planet(input); // 全惑星間の距離を求める vector

station; vector

ans=get_initial_solution(input,station,dist_planet); best_ans=ans; best_station=station; } void SA(){ double ti = timer.get(); double start_time=ti; best_score = calc_score(input,best_ans,best_station); vector

pre_ans=best_ans; vector

pre_station=best_station; ll pre_score = best_score; //cerr << best_score << endl; double kankaku = TIMELIMIT / 10; double ne = ti; int cnt = 0; int cnt2 = 0; double start_temp = 10; double end_temp = 1; bool finded=false; while(1){ ti = timer.get(); if(TIMELIMIT < ti-start_time) break; double temp = start_temp + (end_temp - start_temp) * (ti-start_time) / TIMELIMIT; if (ti > ne) { ne += kankaku; cerr << temp << " " << pre_score << " " << best_score << endl; } // 新たな解の生成============================== vector

new_ans=pre_ans; vector

new_station=pre_station; make_new_ans(input,new_ans,new_station); // new_scoreを計算 ll new_score=calc_score(input,new_ans,new_station); // スコア計算 // 採用するか否かを決定 double prob = exp((new_score-pre_score)/temp); //prob/=1000000000; if (prob*10000 > 10000 * (rand() % 10000)) { pre_ans=new_ans; pre_station=new_station; pre_score=new_score; if(chmax(best_score,new_score)){ best_ans=new_ans; best_station=new_station; best_score=new_score; } cnt2++; } cnt++; } dump(cnt); dump(cnt2); //ti = timer.get(); //cerr << ti << endl; } void end_proc(){ //cerr << "before opt" << endl; //dump(best_score); opt_station(input,best_ans,best_station); best_score = calc_score(input,best_ans,best_station); //cerr << "after opt" << endl; //dump(best_score); //reconnect_station(input,best_ans,best_station); //best_score = calc_score(input,best_ans,best_station); //cerr << "after reconnect" << endl; //dump(best_score); } void input_func() { cin.tie(0); ios::sync_with_stdio(false); int N,M; cin >> N >> M; vector

coord(N); rep(i,0,N){ int a,b; cin >> a >> b; coord[i]=P(a,b); } input={N,M,coord}; } int solve() { init(); SA(); end_proc(); return 0; } void output(){ dump(best_score); rep(i,0,input.M){ int x,y; tie(x,y)=best_station[i]; cout << x << ' ' << y << endl; } fo(sz(best_ans)); rep(i,0,sz(best_ans)){ int t,r; tie(t,r)=best_ans[i]; cout << t << ' ' << r << endl; } } }; Solver solver; int main(int argc, char *argv[]) { solver.input_func(); solver.solve(); // 入力に合うグラフを作成 solver.output(); // クエリに答える return 0; }