# 定数 NUM_PLANETS = 100i64 NUM_STATIONS = 8i64 SUM_POINTS = NUM_PLANETS + NUM_STATIONS # 総地点数 PLANET_MOVE_COST = 5i64 # 問題文ではα COORD_MIN = 0i64 # 座標の下限値 COORD_MAX = 1000i64 # 座標の上限値 TYPE_PLANET = 1i64 # 経由地が惑星 TYPE_STATION = 2i64 # 経由地が宇宙ステーション INF = 10i64 ** 18 # 距離の初期値など、あり得ない大きな数に使う TIME_LIMIT = 0.95 # 時間制限 # 入力受け取り # 惑星の数と宇宙ステーションの数、問題文のNとMに対応 start_time = Time.utc n, m = read_line.split.map(&.to_i64) # 各惑星の座標 coord_planets = Hash(Int64, Tuple(Int64, Int64)).new n.times do |i| a, b = read_line.split.map(&.to_i64) coord_planets[i] = {a, b} end cluster_of_planet, centers, clusters = make_cluster(coord_planets) coord_stations = Hash(Int64, Tuple(Int64, Int64)).new m.times do |i| coord_stations[NUM_PLANETS + i] = centers[i] end STDERR.puts coord_stations dist_table = make_dist_table(coord_planets, coord_stations) # 解法を実行し出力 ans_route = solve(coord_planets, coord_stations, cluster_of_planet, clusters, dist_table, start_time) coord_stations.each_value do |xy| x, y = xy puts "#{x} #{y}" end puts ans_route.size ans_route.each do |point| type = get_point_type(point) point += 1 point -= NUM_PLANETS if type == TYPE_STATION puts "#{type} #{point}" end # ソルバ def solve(coord_planets, coord_stations, cluster_of_planet, clusters, dist_table, start_time) route = cluster_route(coord_planets, coord_stations, cluster_of_planet, clusters, dist_table) # route = greedy_route(coord_planets, coord_stations, dist_table) score, opt2_route = opt2(coord_planets, coord_stations, dist_table, route) # ans_route = [] of Int64 # route_size = opt2_route.size.to_i64 // 8 # station_number = NUM_PLANETS # station_change_cnt = 0i64 # opt2_route.each_with_index do |point, i| # ans_route << point # if i % route_size == 0 && station_change_cnt < 8 # coord_stations[station_number] = coord_planets[point] # ans_route << station_number # station_number += 1 # station_change_cnt += 1 # end # end # ans_route opt2_route end # 2点間の移動に要するエネルギーを計算 def calc_cost(pos1 : Tuple(Int64, Int64), pos2 : Tuple(Int64, Int64), type1 : Int64, type2 : Int64) : Int64 x1, y1 = pos1 x2, y2 = pos2 euc_dist = ((Math.sqrt((x1 - x2) ** 2 + (y1 - y2) ** 2)) ** 2).to_i64 if type1 == TYPE_PLANET && type2 == TYPE_PLANET return (PLANET_MOVE_COST ** 2) * euc_dist elsif type1 == TYPE_STATION && type2 == TYPE_STATION return euc_dist else return PLANET_MOVE_COST * euc_dist end end # 指定した地点が惑星か宇宙ステーションかを判別し、該当する判別番号を返す # ただし地点番号は問題文と異なり0-indexedであることに注意せよ def get_point_type(point : Int64) point < NUM_PLANETS ? TYPE_PLANET : TYPE_STATION end # 各経由地間の距離を計算し、二次元配列にまとめて返す def make_dist_table(coord_planets : Hash(Int64, Tuple(Int64, Int64)), coord_stations : Hash(Int64, Tuple(Int64, Int64))) : Array(Array(Int64)) table_size = NUM_PLANETS + NUM_STATIONS dist_table = Array.new(table_size) { [0i64] * table_size } table_size.times do |i| table_size.times do |j| pos1 = i < NUM_PLANETS ? coord_planets[i] : coord_stations[i] pos2 = j < NUM_PLANETS ? coord_planets[j] : coord_stations[j] type1 = get_point_type(i) type2 = get_point_type(j) dist_table[i][j] = calc_cost(pos1, pos2, type1, type2) end end # ワーシャルフロイド法 SUM_POINTS.times do |k| SUM_POINTS.times do |i| SUM_POINTS.times do |j| d = dist_table[i][k] + dist_table[k][j] dist_table[i][j] = dist_table[i][j] < d ? dist_table[i][j] : d end end end dist_table end # ダイクストラ法によって経由地iから経由地jへの最短経路を復元 def dijkstra(i : Int64, j : Int64, coord_planets : Hash(Int64, Tuple(Int64, Int64)), coord_stations : Hash(Int64, Tuple(Int64, Int64))) dists = [INF] * SUM_POINTS prev_points = [-1i64] * SUM_POINTS q = AtCoder::PriorityQueue(Tuple(Int64, Int64)).new { |a, b| a >= b } dists[i] = 0i64 q << {0i64, i} until q.empty? dist, point = q.pop.not_nil! next if dist > dists[point] SUM_POINTS.times do |next_point| type1 = get_point_type(point) type2 = get_point_type(next_point) pos1 = type1 == TYPE_PLANET ? coord_planets[point] : coord_stations[point] pos2 = type2 == TYPE_PLANET ? coord_planets[next_point] : coord_stations[next_point] next_dist = dist + calc_cost(pos1, pos2, type1, type2) if next_dist < dists[next_point] prev_points[next_point] = point dists[next_point] = next_dist q << {next_dist, next_point} end end end point = j path = [] of Int64 while point != i path << point point = prev_points[point] end path.reverse! end # 距離の近い惑星同士のクラスタを生成 def make_cluster(coord_planets) # k-means法もどき cluster_of_planet = Hash(Int64, Int64).new NUM_PLANETS.times do |planet| cluster_of_planet[planet] = Random.rand(NUM_STATIONS) end clusters = Array.new(NUM_STATIONS) { [] of Int64 } cluster_of_planet.each do |k, v| clusters[v] << k end centers = [] of Tuple(Int64, Int64) NUM_STATIONS.times do centers << coord_planets.sample[1] end finish = false until finish prev_centers = centers new_clusters = Array.new(NUM_STATIONS) { [] of Int64 } old_cluster_of_planet = cluster_of_planet.clone NUM_PLANETS.times do |planet| planet_pos = coord_planets[planet] dists = centers.map_with_index { |center, cluster_number| x2, y2 = center {calc_cost(planet_pos, center, TYPE_PLANET, TYPE_STATION), cluster_number} }.sort! dists.each do |dist_and_cluster_number| dist, cluster_number = dist_and_cluster_number new_clusters[cluster_number] << planet cluster_of_planet[planet] = cluster_number break end end if new_clusters.includes?([] of Int64) return {old_cluster_of_planet, centers, clusters} end clusters = new_clusters centers = [] of Tuple(Int64, Int64) clusters.each do |cluster| centers << calc_center_of_cluster(cluster, coord_planets) end finish = true if centers == prev_centers end STDERR.puts clusters STDERR.puts cluster_of_planet {cluster_of_planet, centers, clusters} end # クラスタの中心座標を求める def calc_center_of_cluster(cluster, coord_planets) center = [0.0, 0.0] cluster_size = cluster.size STDERR.puts cluster_size cluster.each do |point| center[0] += coord_planets[point][0] center[1] += coord_planets[point][1] end center[0] /= cluster_size unless center[0] == 0 center[1] /= cluster_size unless center[0] == 0 {center[0].to_i64, center[1].to_i64} end # 各宇宙ステーションを一周しながら周辺の惑星を訪問し初期解を作る def cluster_route(coord_planets, coord_stations, cluster_of_planet, clusters, dist_table) current_station = cluster_of_planet[0] visited = [false] * NUM_STATIONS visited[current_station] = true route = [0i64] (NUM_STATIONS).times do nearest_dist = INF nearest_station = -1i64 NUM_STATIONS.times do |next_station| next if visited[next_station] if dist_table[current_station + NUM_PLANETS][next_station + NUM_PLANETS] < nearest_dist nearest_dist = dist_table[current_station + NUM_PLANETS][next_station + NUM_PLANETS] nearest_station = next_station end end route << current_station + NUM_PLANETS clusters[current_station].each do |planet| route << planet if planet != 0 route << current_station + NUM_PLANETS end route << current_station + NUM_PLANETS current_station = nearest_station visited[current_station] = true end route << (cluster_of_planet[0] + NUM_PLANETS) route << 0i64 STDERR.puts route route end # 単純な貪欲法で初期解を作る def greedy_route(coord_planets : Hash(Int64, Tuple(Int64, Int64)), coord_stations : Hash(Int64, Tuple(Int64, Int64)), dist_table) point = 0i64 visited = [false] * SUM_POINTS visited[0] = true route = [0i64] (SUM_POINTS - 1).times do |i| nearest_dist = INF nearest_point = -1i64 SUM_POINTS.times do |next_point| next if visited[next_point] if dist_table[point][next_point] < nearest_dist nearest_dist = dist_table[point][next_point] nearest_point = next_point end end path = dijkstra(point, nearest_point, coord_planets, coord_stations) route += path point = nearest_point visited[point] = true end route += dijkstra(point, 0i64, coord_planets, coord_stations) STDERR.puts route route end # 2-opt法で経路を改善 def opt2(coord_planets : Hash(Int64, Tuple(Int64, Int64)), coord_stations : Hash(Int64, Tuple(Int64, Int64)), dist_table, route_arg) total_dist = 0i64 route = route_arg.clone loop do cnt = 0i64 (SUM_POINTS - 2).times do |i| i1 = i + 1 (i + 2).upto(SUM_POINTS - 1) do |j| j1 = j == SUM_POINTS - 1 ? 0i64 : j + 1 if i != 0 || j1 != 0 l1 = dist_table[route[i]][route[i1]] l2 = dist_table[route[j]][route[j1]] l3 = dist_table[route[i]][route[j]] l4 = dist_table[route[i1]][route[j1]] if l1 + l2 > l3 + l4 route[i1..j] = route[i1..j].reverse! cnt += 1 end end end end total_dist += cnt break if cnt == 0 end {total_dist, route} end # ac-library.cr by hakatashi https://github.com/google/ac-library.cr # # Copyright 2022 Google LLC # # Licensed under the Apache License, Version 2.0 (the "License"); # you may not use this file except in compliance with the License. # You may obtain a copy of the License at # # https://www.apache.org/licenses/LICENSE-2.0 # # Unless required by applicable law or agreed to in writing, software # distributed under the License is distributed on an "AS IS" BASIS, # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. # See the License for the specific language governing permissions and # limitations under the License. module AtCoder # Implements standard priority queue like [std::priority_queue](https://en.cppreference.com/w/cpp/container/priority_queue). # # ``` # q = AtCoder::PriorityQueue(Int64).new # q << 1_i64 # q << 3_i64 # q << 2_i64 # q.pop # => 3 # q.pop # => 2 # q.pop # => 1 # ``` class PriorityQueue(T) getter heap : Array(T) def initialize initialize { |a, b| a <= b } end # Initializes queue with the custom comperator. # # If the second argument `b` should be popped earlier than # the first argument `a`, return `true`. Else, return `false`. # # ``` # q = AtCoder::PriorityQueue(Int64).new { |a, b| a >= b } # q << 1_i64 # q << 3_i64 # q << 2_i64 # q.pop # => 1 # q.pop # => 2 # q.pop # => 3 # ``` def initialize(&block : T, T -> Bool) @heap = Array(T).new @compare_proc = block end # Pushes value into the queue. def push(v : T) @heap << v index = @heap.size - 1 while index != 0 parent = (index - 1) // 2 if @compare_proc.call(@heap[index], @heap[parent]) break end @heap[parent], @heap[index] = @heap[index], @heap[parent] index = parent end end # Alias of `push` def <<(v : T) push(v) end # Pops value from the queue. def pop if @heap.size == 0 return nil end if @heap.size == 1 return @heap.pop end ret = @heap.first @heap[0] = @heap.pop index = 0 while index * 2 + 1 < @heap.size child = if index * 2 + 2 < @heap.size && !@compare_proc.call(@heap[index * 2 + 2], @heap[index * 2 + 1]) index * 2 + 2 else index * 2 + 1 end if @compare_proc.call(@heap[child], @heap[index]) break end @heap[child], @heap[index] = @heap[index], @heap[child] index = child end ret end # Returns `true` if the queue is empty. delegate :empty?, to: @heap # Returns size of the queue. delegate :size, to: @heap end end