# 解法概要 # まず惑星をk-means法を用いて宇宙ステーションの個数分のクラスタに分ける # 次に惑星1を起点にして宇宙ステーションを一周して惑星1に戻る初期解を作る # 宇宙ステーションに寄る度に同クラスタに属する惑星も巡っておく # 2-opt法を用いて初期解の経路を改善する # k-means法の初期値を乱択で変更しながら上記を0.3秒繰り返し、より良い初期解を得る(山登り法) # 仕上げに焼きなまし法を用いて宇宙ステーションの位置を改善する(結局山登りと結果変わらないけど...) # 定数 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 ** 9 # 距離の初期値など、あり得ない大きな数に使う TIME_LIMIT = 0.95 # 時間制限 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 # 解法を実行し出力 ans_route, coord_stations = solve(coord_planets, coord_stations, cluster_of_planet, clusters, 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, start_time) loop_times = 0i64 change_times = 0i64 # より良い初期解を求めて山登り法 route = cluster_route(coord_planets, coord_stations, cluster_of_planet, clusters) cost = calc_route_total_cost(route, coord_planets, coord_stations) opt2_route = route opt2_route = opt2(coord_planets, coord_stations, route) opt2_cost = calc_route_total_cost(opt2_route, coord_planets, coord_stations) opt2_time_limit = 0.3 while (Time.utc - start_time).total_seconds <= opt2_time_limit loop_times += 1 new_cluster_of_planet, new_centers, new_clusters = make_cluster(coord_planets) new_coord_stations = Hash(Int64, Tuple(Int64, Int64)).new NUM_STATIONS.times do |station_number| new_coord_stations[station_number + NUM_PLANETS] = new_centers[station_number] end new_route = cluster_route(coord_planets, new_coord_stations, new_cluster_of_planet, new_clusters) new_cost = calc_route_total_cost(new_route, coord_planets, coord_stations) new_opt2_route = opt2(coord_planets, new_coord_stations, new_route) new_opt2_cost = calc_route_total_cost(new_opt2_route, coord_planets, new_coord_stations) if new_opt2_cost < opt2_cost coord_stations = new_coord_stations cluster_of_planet = new_cluster_of_planet clusters = new_clusters centers = new_centers route = new_route cost = new_cost opt2_cost = new_opt2_cost opt2_route = new_opt2_route change_times += 1 end end # 宇宙ステーションの位置を改善するための焼きなまし法 temp_start = 100i64 temp_end = 1i64 sa_time_limit = TIME_LIMIT score = calc_score(opt2_cost) change_limit = 10i64 best_score = score best_coord_stations = coord_stations while (Time.utc - start_time).total_seconds <= TIME_LIMIT loop_times += 1 target_station = Random.rand(NUM_STATIONS) x1, y1 = coord_stations[target_station + NUM_PLANETS] amount_of_change_x = Random.rand(-(change_limit)..change_limit) amount_of_change_y = Random.rand(-(change_limit)..change_limit) x2 = x1 + amount_of_change_x y2 = y1 + amount_of_change_y next unless 0i64 <= x2 <= COORD_MAX && 0i64 <= y2 <= COORD_MAX new_coord_stations = coord_stations.clone new_coord_stations[target_station + NUM_PLANETS] = {x2, y2} new_opt2_cost = calc_route_total_cost(opt2_route, coord_planets, new_coord_stations) new_score = calc_score(new_opt2_cost) temp = temp_start + (temp_end - temp_start) * (Time.utc - start_time).total_seconds / sa_time_limit prob = Math.exp((new_score - score) / temp) rand_mod = ((Random.rand(Int32::MIN..Int32::MAX) % INF) / INF) force_next = prob > rand_mod if new_score > score || force_next coord_stations = new_coord_stations opt2_cost = new_opt2_cost score = new_score change_times += 1 if score > best_score best_coord_stations = coord_stations best_score = score end end end STDERR.puts "loop_times: #{loop_times}" STDERR.puts "change_times: #{change_times}" {opt2_route, best_coord_stations} end # スコア計算(ざっくり) def calc_score(cost) ((10i64 ** 9) / (1000i64 + Math.sqrt(cost))).round 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 calc_route_total_cost(route, coord_planets, coord_stations) cost = 0i64 route.each_cons_pair do |a, b| type1 = get_point_type(a) type2 = get_point_type(b) pos1 = type1 == TYPE_PLANET ? coord_planets[a] : coord_stations[a] pos2 = type2 == TYPE_PLANET ? coord_planets[b] : coord_stations[b] cost += calc_cost(pos1, pos2, type1, type2) end cost 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) range = 30 NUM_STATIONS.times do x1, y1 = coord_planets.sample[1] x2 = Random.rand((x1 - range)..(x1 + range)) y2 = Random.rand((y1 - range)..(y1 + range)) x2 = x1 unless 0i64 <= x2 <= COORD_MAX y2 = y1 unless 0i64 <= y2 <= COORD_MAX centers << {x2, y2} 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 {cluster_of_planet, centers, clusters} end # クラスタの中心座標を求める def calc_center_of_cluster(cluster, coord_planets) center = [0.0, 0.0] cluster_size = 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 calc_cluster_total_cost(clusters, coord_planets, coord_stations) total_cost = 0i64 NUM_STATIONS.times do |cluster_number| pos1 = coord_stations[cluster_number + NUM_PLANETS] total_cost += clusters[cluster_number].sum { |planet| calc_cost(pos1, coord_planets[planet], TYPE_STATION, TYPE_PLANET) } end total_cost end # 各宇宙ステーションを一周しながら周辺の惑星を訪問し初期解を作る def cluster_route(coord_planets, coord_stations, cluster_of_planet, clusters) 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] dist = calc_cost(coord_stations[current_station + NUM_PLANETS], coord_stations[next_station + NUM_PLANETS], TYPE_STATION, TYPE_STATION) if dist < nearest_dist nearest_dist = dist 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 route end # 2-opt法で経路を改善 def opt2(coord_planets : Hash(Int64, Tuple(Int64, Int64)), coord_stations : Hash(Int64, Tuple(Int64, Int64)), route_arg : Array(Int64)) : Array(Int64) total_dist = 0i64 route = route_arg.clone route_size = route.size loop do cnt = 0i64 (route_size - 2).times do |i| i1 = i + 1 (i + 2).upto(route_size - 1) do |j| j1 = j == route_size - 1 ? 0i64 : j + 1 if i != 0 && j1 != 0 point1 = route[i] point2 = route[j] point3 = route[i1] point4 = route[j1] type1 = get_point_type(point1) type2 = get_point_type(point2) type3 = get_point_type(point3) type4 = get_point_type(point4) pos1 = type1 == TYPE_PLANET ? coord_planets[point1] : coord_stations[point1] pos2 = type2 == TYPE_PLANET ? coord_planets[point2] : coord_stations[point2] pos3 = type3 == TYPE_PLANET ? coord_planets[point3] : coord_stations[point3] pos4 = type4 == TYPE_PLANET ? coord_planets[point4] : coord_stations[point4] l1 = calc_cost(pos1, pos3, type1, type3) # dist_table[route[i]][route[i1]] l2 = calc_cost(pos2, pos4, type2, type4) # dist_table[route[j]][route[j1]] l3 = calc_cost(pos1, pos2, type1, type2) # dist_table[route[i]][route[j]] l4 = calc_cost(pos3, pos4, type3, type4) # dist_table[route[i1]][route[j1]] if l1 + l2 > l3 + l4 section_size = j - i1 + 1 reverse_times = section_size // 2 route[i1..j] = route[i1..j].reverse! cnt += 1 end end end end total_dist += cnt break if cnt == 0 end route end # 経路の部分破壊と再構築 def change_route(route, coord_planets, coord_stations) route_size = route.size section_size = Random.rand(2i64..4i64) change_start_i = Random.rand(route_size) section_limit = change_start_i + section_size < route_size - 1 ? change_start_i + section_size : route_size - 1 change_end_i = Random.rand(change_start_i..section_limit) section = route[change_start_i..change_end_i] cost = calc_route_total_cost(section, coord_planets, coord_stations) ans_section = section changed_route = section section.permutations(section.size).to_a.each do |new_seciton| new_cost = calc_route_total_cost(new_seciton, coord_planets, coord_stations) if new_cost < cost ans_section = new_seciton cost = new_cost end end new_route = route.clone new_route[change_start_i..change_end_i] = ans_section new_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