import copy import random from time import perf_counter import argparse import sys import math class TimeKeeper: def __init__(self): self.start_time = perf_counter() def is_time_over(self, LIMIT): return (perf_counter() - self.start_time) >= LIMIT def time_now(self): return (perf_counter() - self.start_time) ############## def main(DEBUG): # ユークリッド距離 def euclidean_distance(point1, point2): return sum((a - b) ** 2 for a, b in zip(point1, point2)) ** 0.5 # K-meansクラスタリング def initialize_centroids(points, k): return random.sample(points, k) def assign_clusters(points, centroids): clusters = [] for point in points: distances = [euclidean_distance(point, centroid) for centroid in centroids] cluster = distances.index(min(distances)) clusters.append(cluster) return clusters def update_centroids(points, clusters, k): new_centroids = [] for i in range(k): cluster_points = [points[j] for j in range(len(points)) if clusters[j] == i] if cluster_points: new_centroid = [sum(dim) / len(cluster_points) for dim in zip(*cluster_points)] new_centroids.append(tuple(new_centroid)) return new_centroids def kmeans(points, k, max_iters=100): centroids = initialize_centroids(points, k) for _ in range(max_iters): clusters = assign_clusters(points, centroids) new_centroids = update_centroids(points, clusters, k) if centroids == new_centroids: break centroids = new_centroids return clusters, centroids def HeldKarp(dists): N = len(dists) dp = [[10**18] * N for _ in range(2**N)] for i in range(N): dp[1< 10**15: return [] now = start path = [] while now != goal: path.append(now) for i in range(N): if i != now and d[now][i] + d[i][goal] == d[now][goal]: now = i break path.append(goal) return path ######################## tk = TimeKeeper() LIMIT = 0.8 random.seed(0) ALPHA = 5 # 定数 N, M = map(int, input().split()) points = [] for i in range(N): x, y = map(int, input().split()) points.append((x, y)) id_to_xy = dict() for i in range(N): id_to_xy[i] = points[i] def solve(N, M, ALPHA, points): def calscore(stations, tour, ALPHA): s = 0 NN = len(tour) for i in range(NN-1): type0, id0 = tour[i] type1, id1 = tour[i+1] id0 -= 1 id1 -= 1 if type0 == 1: x0, y0 = id_to_xy[id0] else: x0, y0 = stations[id0] if type1 == 1: x1, y1 = id_to_xy[id1] else: x1, y1 = stations[id1] d2 = (x0 - x1)**2 + (y0 - y1)**2 if type0 + type1 == 2: coef = ALPHA**2 elif type0 + type1 == 3: coef = ALPHA elif type0 + type1 == 4: coef = 1 s += d2 * coef sc = int(10**9 / (1000 + s**0.5) + 0.5) return sc ############ # K-meansクラスタリングでM個のクラスに分類 while True: labels, centroids = kmeans(points, M) stations = [] for x, y in centroids: ix, iy = int(x+0.5), int(y+0.5) stations.append((ix, iy)) if len(stations) == M: break # スタート地点のグループが0になるように回転 stt = labels[0] labels = [(x-stt)%M for x in labels] for _ in range(stt): z = stations.pop(0) stations.append(z) # グループごとにp, stationの座標を # ps[i], stations[i] ps = [[] for _ in range(M)] id_to_xy = dict() for i, (x, y) in enumerate(points): ps[labels[i]].append((x, y, i)) id_to_xy[i] = (x, y) # print(ps) # 最短巡回ルート dists = [[0]*M for _ in range(M)] for i in range(M): for j in range(M): dists[i][j] = euclidean_distance(stations[i], stations[j]) tot_dist, tour_grp = HeldKarp(dists) # グループごとに解く tour = [] grp = 0 now = 0 # star_1 tour.append((1, now+1)) # type, id x, y = id_to_xy[now] ps[grp].remove((x, y, now)) # リストから削除 for grp in tour_grp: # star_1 -> station stx, sty = stations[grp] # stationのx,y tour.append((2, grp+1)) # type, id while ps[grp]: mndist = 10**18 for x, y, id in ps[grp]: d = euclidean_distance(stations[grp], (x, y)) if d < mndist: mndist = d nxt = id # station -> nearest star # move(grp, nxt) tour.append((1, nxt+1)) x, y = id_to_xy[nxt] ps[grp].remove((x, y, nxt)) # リストから削除 nowx, nowy = x, y now = nxt while True: # 最も近い地点とその距離を調べる mndpp = 10**18 for x, y, id in ps[grp]: d = euclidean_distance((nowx, nowy), (x, y)) if d < mndist: mndpp = d nxt = id nxtx, nxty = id_to_xy[nxt] d1 = euclidean_distance(stations[grp], (nowx, nowy)) d2 = euclidean_distance(stations[grp], (nxtx, nxty)) if ALPHA * mndpp < d1 + d2: # 次の地点へ直接行く tour.append((1, nxt+1)) ps[grp].remove((nxtx, nxty, nxt)) # リストから削除 nowx, nowy = nxtx, nxty now = nxt else: break # ステーションに戻る # move(now, grp) tour.append((2, grp+1)) # type, id tour.append((2, 1)) tour.append((1, 1)) sc = calscore(stations, tour, ALPHA) # print(sc, file=sys.stderr) return stations, tour, sc bestsc = 0 loop = 0 while True: # print(loop, file=sys.stderr) if tk.is_time_over(LIMIT): break ########### stations, tour, sc = solve(N, M, ALPHA, points) ########### if sc > bestsc: print('best', loop, sc, file=sys.stderr) bestsc = sc besttour = tour beststations = stations loop += 1 # スコア print('_SC', bestsc, file=sys.stderr) # 訪問箇所を追加してスコアが上がれば採用 # 全点対最短経路を求める C = [[10**18]*(N+M) for _ in range(N+M)] for i in range(N+M): C[i][i] = 0 for i0, (x0, y0) in enumerate(points): for i1, (x1, y1) in enumerate(points): C[i0][i1] = ALPHA**2 * ((x0 - x1)**2 + (y0 - y1)**2) C[i1][i0] = ALPHA**2 * ((x0 - x1)**2 + (y0 - y1)**2) for j in range(M): xs, ys = stations[j] C[i0][N+j] = ALPHA * ((x0 - xs)**2 + (y0 - ys)**2) C[N+j][i0] = ALPHA * ((x0 - xs)**2 + (y0 - ys)**2) for i in range(M): xs0, ys0 = stations[i] for j in range(j): xs1, ys1 = stations[j] C[i][j] = ((xs0 - xs1)**2 + (ys0 - ys1)**2) C = WF(C) anstour = [] anstour.append(besttour[0]) Nt = len(besttour) for i in range(Nt-1): type0, id0 = besttour[i] type1, id1 = besttour[i+1] id0 -= 1 id1 -= 1 # ★から★ if type0 == 1 and type1 == 1: nowc = C[id0][id1] path = get_path(C, id0, id1) if len(path) > 1: for j in range(1, len(path)): if path[j] < N: anstour.append((1, path[j]+1)) else: anstour.append((2, path[j]+1-N)) else: anstour.append((type1, id1+1)) else: anstour.append((type1, id1+1)) def calscore(stations, tour, ALPHA): s = 0 NN = len(tour) for i in range(NN-1): type0, id0 = tour[i] type1, id1 = tour[i+1] id0 -= 1 id1 -= 1 if type0 == 1: x0, y0 = id_to_xy[id0] else: x0, y0 = stations[id0] if type1 == 1: x1, y1 = id_to_xy[id1] else: x1, y1 = stations[id1] d2 = (x0 - x1)**2 + (y0 - y1)**2 if type0 + type1 == 2: coef = ALPHA**2 elif type0 + type1 == 3: coef = ALPHA elif type0 + type1 == 4: coef = 1 s += d2 * coef sc = int(10**9 / (1000 + s**0.5) + 0.5) return sc anssc = calscore(beststations, anstour, ALPHA) print('SC', anssc, file=sys.stderr) # 出力(ステーション) for x, y in beststations: print(x, y) # 出力(訪問箇所) print(len(anstour)) for type, id in anstour: print(type, id) if __name__ == '__main__': parser = argparse.ArgumentParser(description='Debug mode') parser.add_argument('--debug', action='store_true', help='Enable debug mode') args = parser.parse_args() main(args.debug)