ひとりで勝手にはじめた蟻本読書会 平面走査 計算幾何 蟻本読書会 の続きです。
凸包とグラハムスキャン
二次元平面に N 箇所の牧場があります。i 番目の牧場は(x_i, y_i)にあります。同じ座標に複数の座標は存在しません。もっとも遠いところにある牧場の距離を求めよ。ただしNの最大値は50000とする。
このような問題を解く場合、N箇所の牧場の距離 N *(N – 1)/ 2 とおりすべてを調べていては制限時間内に間に合わせることはできません。省略できる計算を省略するにはどうすればよいかを考えます。
凸包上の点だけ調べる
ある3点を結んでできる三角形の内部にある点は最遠点対をなすことはありません。そこでこのような点を取り除く処理を繰り返していくと、与えられた点集合のうちもっとも外側にある点だけが残ります。このようなもっとも外側の点集合は、元の点集合を含む最小の凸多角形の頂点となり、この多角形を元の点集合の凸包といいます。凸包がわかればあとは凸包上のすべてのペアについてだけ距離を計算すればよいので、制限時間内に解を得ることができます。
凸包を求める方法
では凸包を得るにはどうすればよいのでしょうか?ここではグラハムスキャンを用いて凸包を求めます。まず点をX軸でソートします。ソート後の最初の点と最後の点はかならず凸包の頂点となります。あとは間の点をうえがわとしたがわにわけてそれぞれ求めます。
ソートされた点を繋げて凸包をつくっていきますが、点を追加することで凸ではない形になる場合があります。この場合は凹の部分を取り除いて処理を進めていきます。
上側であれば最初のふたつの点はリストに追加して、3つ目以降は最後と最後から二番目に追加された点で構成される線分の傾きと最後に追加された点とこれから追加されようとする点で構成される線分の傾きを比較します。後者のほうが大きいと凹の部分ができてしまうので、最後に追加された点を取り除き、もう一度同じような比較の処理をおこないます。これで凸包の上側を完成させることができます。
下側も同じような処理をおこない、両者を合わせれば凸包が完成します。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 |
using System; using System.Collections.Generic; using System.Linq; class Position { public Position(double x, double y) { X = x; Y = y; } public double X = 0; public double Y = 0; } class Program { static List<Position> GrahamScan(List<Position> positions) { List<Position> sorted = positions.OrderBy(_ => _.X).ThenBy(_ => _.Y).ToList(); // 凸包の上半分 List<Position> rets = new List<Position>(); foreach (Position pos in sorted) AddList(rets, pos, true); // 下半分 List<Position> rets2 = new List<Position>(); foreach (Position pos in sorted) AddList(rets2, pos, false); // 最初と最後には同じ点が格納されているので取り除く if (rets2[0] == rets[0]) rets2.RemoveAt(0); if (rets2.Last() == rets.Last()) rets2.RemoveAt(rets2.Count - 1); // 時計回りになるように両者を結合する(下半分を反転して上半分に追加) rets2.Reverse(); rets.AddRange(rets2); return rets; } // 凸包の上半分(または下半分)を構成する点のリストに点を追加する static void AddList(List<Position> list, Position pos, bool upper) { if (list.Count < 2) // 最初の2つはそのまま追加 list.Add(pos); else // 3つ目以降 { int lastIndex = list.Count - 1; double a = GetSlope(list[lastIndex - 1], list[lastIndex]) - GetSlope(list[lastIndex], pos); if ((upper && a > 0) || (!upper && a < 0)) list.Add(pos); else { // 凹になる場合はひとつ前に追加したものを削除してもう一度追加処理をやり直す list.RemoveAt(lastIndex); AddList(list, pos, upper); } } } // 点 p1 と p2 で構成される線分の傾きを求める static double GetSlope(Position p1, Position p2) { if (p2.Y == p1.Y) return 0; if (p2.X - p1.X == 0 && p2.Y > p1.Y) // p2.X - p1.X == 0 のとき除算違反にならないように対策 return double.MaxValue; if (p2.X - p1.X == 0 && p2.Y < p1.Y) return double.MinValue; return (p2.Y - p1.Y) / (p2.X - p1.X); } } |
応用問題
D – Big Bang
天文学者である高橋君はその宇宙の膨張の速度を計測することにしました。
高橋君はある 2 つの日について、同じ N 個の星の位置を観測しました。星の位置は座標平面上の点として記録されます。つまり各日の観測結果は座標平面上の N 個の点からなる点集合になります。
2 回の観測の結果を見比べてみると、1 回目の観測結果である点集合に対して以下の 3 つの操作を順に実行すると 2 回目の観測結果である点集合に一致することがわかりました。
同じ向きに同じ距離だけ平行移動する。
原点を中心に同じ角度だけ回転する。
原点を中心に P 倍 (1 ≦ P) に相似拡大する。つまり点 (a, b) を点 (a×P, b×P) に移すという操作をすべての点に実行する。P の値がわかれば膨張速度を求めることができそうです。1 回目と 2 回目の観測結果が与えられるので P を求めてください。
入力されるデータ
N
Ax_1 Ay_1
Ax_2 Ay_2
:
Ax_N Ay_N
Bx_1 By_1
Bx_2 By_2
:
Bx_N By_N(ただし、
2 ≦ N ≦ 10^5
-15,000 ≦ Ax_i, Ay_i ≦ 15,000
-10^9 ≦ Bx_i, By_i ≦ 10^91 回目も 2 回目も、同じ点に複数の星が観測されることはない。
1 回目に観測された i 番目の星と 2 回目に観測された i 番目の星が同一の星とは限らない )
「1 回目に観測された i 番目の星と 2 回目に観測された i 番目の星が同一の星とは限らない」ので、AとBをそのまま比較することはできません。しかし「原点を中心に P 倍に相似拡大する」というのであれば1回目の観測でもっとも離れていた星のペアと2回目の観測でもっとも離れている星のペアは同じ星のペアであり、ここを比較すれば P の値を求めることができます。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 |
using System; using System.Collections.Generic; using System.Linq; // Positionクラスは上記と同じなので省略 class Program { // GrahamScanメソッド等は上記と同じなので省略 static void Main() { int N = int.Parse(Console.ReadLine()); List<Position> A = new List<Position>(); for (int i = 0; i < N; i++) { double[] vs = Console.ReadLine().Split().Select(_ => double.Parse(_)).ToArray(); double x = vs[0]; double y = vs[1]; A.Add(new Position(x, y)); } List<Position> B = new List<Position>(); for (int i = 0; i < N; i++) { double[] vs = Console.ReadLine().Split().Select(_ => double.Parse(_)).ToArray(); double x = vs[0]; double y = vs[1]; B.Add(new Position(x, y)); } A = GrahamScan(A); B = GrahamScan(B); // A の最遠点対の距離(2乗)を求める double maxBefore2 = double.MinValue; for(int i = 1; i < A.Count; i++) maxBefore2 = Math.Max(maxBefore2, Math.Pow(A[i - 1].X - A[i].X, 2) + Math.Pow(A[i - 1].Y - A[i].Y, 2)); maxBefore2 = Math.Max(maxBefore2, Math.Pow(A[0].X - A[A.Count - 1].X, 2) + Math.Pow(A[0].Y - A[A.Count - 1].Y, 2)); // B の最遠点対の距離(2乗)を求める double maxAfter2 = double.MinValue; for (int i = 1; i < B.Count; i++) maxAfter2 = Math.Max(maxAfter2, Math.Pow(B[i - 1].X - B[i].X, 2) + Math.Pow(B[i - 1].Y - B[i].Y, 2)); maxAfter2 = Math.Max(maxAfter2, Math.Pow(B[0].X - B[B.Count - 1].X, 2) + Math.Pow(B[0].Y - B[B.Count - 1].Y, 2)); // 用いているのは距離の2乗なので除算したあと平方根を取る Console.WriteLine(Math.Sqrt(maxAfter2/ maxBefore2)); } } |
動的凸包
C – 泥棒
あなたは、川と緑に囲まれた小さいながらも豊かな王国の泥棒です。
新しい王が治安向上のために監視所を設置し始めたので監視に引っかからないように気をつける必要があります。
監視所では他の監視所との間をまっすぐ結んだ線の上に誰がいるかを調べています。
あなたはこの線を横切ったり、線上の家に盗みに入ることはできません。
もちろん、監視所を通過することも出来ません
盗みに入る家と監視所の位置が与えられるので監視に引っかからずにアジトから盗みに入る家までたどり着けるかを判定し、監視に引っかからないなら SAFE 、引っかかるなら DANGER と出力せよ。ただし、盗みに入る段階でまだ設置されていない監視所について監視に引っかかることはありません。
入力されるデータ
N
S_1 X_1 Y_1
S_2 X_2 Y_2
…
S_N X_N Y_Nただし、N は盗みに入る家と監視所の数の合計を表す整数 (2 ≦ N ≦ 100,000)である。
S_i は TARGET もしくは MONITOR であり、それぞれ盗みに入る家と監視所を表している。
盗みに入る家・監視所は X_i, Y_i(1 ≦ X_i, Y_i ≦ 100,000,000)であり、すべて異なる位置にある。
家に盗みに入る日時・監視所が設置される日時は入力で与えられる順の通りである。
すなわち、ある家に盗みに入る時に既に設置されている監視所とはそれより前の入力で与えられた監視所のことである。
アジトは原点(0, 0)にある。
アジトが原点にあり、盗みに入る家と監視所の座標の最小値が 1 であることから、盗みに入る家が監視所によって囲まれていたら DANGER ということになります。この問題も凸包を使えば解けそうです。
前問では最初にすべての座標が与えられていましたが、この問題では凸包を動的に構成していくことになります。上半分と下半分にわけて凸包を構築します。X座標をみてそれぞれのどの位置に追加するのかを決めます。N ≦ 100,000 と大きいので、追加する位置は二分探索法で求めます。
すでに構築されている凸包の内部には点を追加しません。また点を追加するまえに凹んだ部分ができるかどうかを確認し、凹みができる場合はそのような点を取り除いてから追加の処理をおこないます。あとは同じです。
なお、凸包の内部に相当するか?点を追加して凹んだ部分ができるかどうかは符号付き面積を利用して判定しています。
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 |
using System; using System.Collections.Generic; using System.Linq; class Position { public Position(double x, double y) { X = x; Y = y; } public double X = 0; public double Y = 0; } class Program { // 点を追加することで凹んだ部分ができた場合は取り除きながら凸包を構成する点を追加していく static void InsertList(Position pos, List<Position> list, bool isUpper) { if (list.Count < 2) { // 最初の点はそのまま追加する if (list.Count == 0) list.Add(pos); else if (list.Count == 1) { // ふたつめの点はX座標が後ろの場合、同じ場合はY座標が大きい場合、後ろに追加する // それ以外のときは最初の要素とする if (list[0].X < pos.X || (list[0].X == pos.X && list[0].Y < pos.Y)) list.Add(pos); else list.Insert(0, pos); } return; } // 二分探索法で挿入すべき位置を探す int idx = LowerBound(list, pos.X); if (idx == 0) Insert(0); else if (idx >= list.Count) Insert(list.Count); else { // 上側の凸包を構築しているときはそれより下(内側)に点を追加してはならない if (isUpper && IsLower(list[idx - 1], list[idx], pos)) return; // 下側の凸包を構築しているときはそれより上(内側)に点を追加してはならない if (!isUpper && IsUpper(list[idx - 1], list[idx], pos)) return; // それ以外は追加してOK Insert(idx); } // 再帰呼び出し。点を追加することで凹んだ部分ができた場合は取り除きながら凸包を構成する点を追加していく void Insert(int index) { // 点を追加することで凹んだ部分が左にできないかチェック if (index - 2 >= 0) { if (isUpper && IsLower(list[index - 2], pos, list[index - 1])) { list.RemoveAt(index - 1); Insert(index - 1); return; } if (!isUpper && IsUpper(list[index - 2], pos, list[index - 1])) { list.RemoveAt(index - 1); Insert(index - 1); return; } } // 点を追加することで凹んだ部分が右にできないかチェック if (index + 1 < list.Count) { if (isUpper && IsLower(pos, list[index + 1], list[index])) { list.RemoveAt(index); Insert(index); return; } if (!isUpper && IsUpper(pos, list[index + 1], list[index])) { list.RemoveAt(index); Insert(index); return; } } // 凹んだ部分ができない場合はリストに点を追加する list.Insert(index, pos); } } // 第一引数から第二引数に向かう線分より第三引数の点は上にあるか?(線上も含む) static bool IsUpper(Position start, Position end, Position point) { double x1 = end.X - start.X; double y1 = end.Y - start.Y; double x2 = point.X - start.X; double y2 = point.Y - start.Y; return x1 * y2 - x2 * y1 >= 0; } // 第一引数から第二引数に向かう線分より第三引数の点は下にあるか?(線上も含む) static bool IsLower(Position start, Position end, Position point) { double x1 = end.X - start.X; double y1 = end.Y - start.Y; double x2 = point.X - start.X; double y2 = point.Y - start.Y; return x1 * y2 - x2 * y1 <= 0; } // 二分探索法を点を追加する位置を探す static int LowerBound(List<Position> positions, double target) { if (positions[0].X >= target) return 0; int left = 0; int right = positions.Count; while (right - left > 1) { int mid = (left + right) / 2; if (target <= positions[mid].X) right = mid; else left = mid; } return right; } static void Main() { int N = int.Parse(Console.ReadLine()); List<Position> upper = new List<Position>(); List<Position> lower = new List<Position>(); for (int i = 0; i < N; i++) { string[] vs = Console.ReadLine().Split(); double x = int.Parse(vs[1]); double y = int.Parse(vs[2]); Position pos = new Position(x, y); if (vs[0] == "MONITOR") { InsertList(pos, upper, true); InsertList(pos, lower, false); } if (vs[0] == "TARGET") { bool inside = false; if (upper.Count > 0) { int index = LowerBound(upper, x); if (index - 1 >= 0 && index < upper.Count) inside = IsLower(upper[index - 1], upper[index], pos); } if (inside && lower.Count > 0) { int index = LowerBound(lower, x); if (index - 1 >= 0 && index < lower.Count) inside = IsUpper(lower[index - 1], lower[index], pos); } // 最初の2点が垂直に位置してその間に点がある場合だけ別に考える(index が 0 となって上の判定がうまくできない) if (upper.Count >= 2 && upper[0].X == upper[1].X && upper[0].X == pos.X) { if(upper[0].Y <= pos.Y && pos.Y <= upper[1].Y) inside = true; } if (lower.Count >= 2 && lower[0].X == lower[1].X && lower[0].X == pos.X) { if (lower[1].Y <= pos.Y && pos.Y <= lower[0].Y) inside = true; } if (inside) Console.WriteLine("DANGER"); else Console.WriteLine("SAFE"); } } } } |