1"""
2기하 알고리즘 (Computational Geometry)
3Computational Geometry Algorithms
4
5점, 선, 다각형 등의 기하학적 객체를 다루는 알고리즘입니다.
6"""
7
8from typing import List, Tuple, Optional
9from math import sqrt, atan2, pi, inf
10from functools import cmp_to_key
11
12
13# =============================================================================
14# 1. 기본 기하 연산
15# =============================================================================
16
17class Point:
18 """2D 점"""
19 def __init__(self, x: float, y: float):
20 self.x = x
21 self.y = y
22
23 def __add__(self, other: 'Point') -> 'Point':
24 return Point(self.x + other.x, self.y + other.y)
25
26 def __sub__(self, other: 'Point') -> 'Point':
27 return Point(self.x - other.x, self.y - other.y)
28
29 def __mul__(self, scalar: float) -> 'Point':
30 return Point(self.x * scalar, self.y * scalar)
31
32 def __eq__(self, other: 'Point') -> bool:
33 return abs(self.x - other.x) < 1e-9 and abs(self.y - other.y) < 1e-9
34
35 def __repr__(self) -> str:
36 return f"({self.x}, {self.y})"
37
38 def dot(self, other: 'Point') -> float:
39 """내적"""
40 return self.x * other.x + self.y * other.y
41
42 def cross(self, other: 'Point') -> float:
43 """외적 (z 성분)"""
44 return self.x * other.y - self.y * other.x
45
46 def norm(self) -> float:
47 """벡터 크기"""
48 return sqrt(self.x ** 2 + self.y ** 2)
49
50 def normalize(self) -> 'Point':
51 """단위 벡터"""
52 n = self.norm()
53 return Point(self.x / n, self.y / n) if n > 0 else Point(0, 0)
54
55
56def distance(p1: Point, p2: Point) -> float:
57 """두 점 사이 거리"""
58 return (p2 - p1).norm()
59
60
61# =============================================================================
62# 2. CCW (Counter-Clockwise)
63# =============================================================================
64
65def ccw(a: Point, b: Point, c: Point) -> int:
66 """
67 세 점의 방향 판별
68 반환: 1 (반시계), -1 (시계), 0 (일직선)
69 """
70 cross = (b - a).cross(c - a)
71 if cross > 1e-9:
72 return 1 # 반시계
73 elif cross < -1e-9:
74 return -1 # 시계
75 return 0 # 일직선
76
77
78def ccw_tuple(a: Tuple[float, float], b: Tuple[float, float], c: Tuple[float, float]) -> int:
79 """튜플 버전 CCW"""
80 cross = (b[0] - a[0]) * (c[1] - a[1]) - (b[1] - a[1]) * (c[0] - a[0])
81 if cross > 1e-9:
82 return 1
83 elif cross < -1e-9:
84 return -1
85 return 0
86
87
88# =============================================================================
89# 3. 선분 교차 판정
90# =============================================================================
91
92def segments_intersect(p1: Point, p2: Point, p3: Point, p4: Point) -> bool:
93 """
94 선분 p1-p2와 p3-p4의 교차 여부
95 """
96 d1 = ccw(p3, p4, p1)
97 d2 = ccw(p3, p4, p2)
98 d3 = ccw(p1, p2, p3)
99 d4 = ccw(p1, p2, p4)
100
101 # 일반적인 교차
102 if d1 * d2 < 0 and d3 * d4 < 0:
103 return True
104
105 # 경계 케이스 (점이 선분 위에 있는 경우)
106 def on_segment(p: Point, q: Point, r: Point) -> bool:
107 return (min(p.x, r.x) <= q.x <= max(p.x, r.x) and
108 min(p.y, r.y) <= q.y <= max(p.y, r.y))
109
110 if d1 == 0 and on_segment(p3, p1, p4):
111 return True
112 if d2 == 0 and on_segment(p3, p2, p4):
113 return True
114 if d3 == 0 and on_segment(p1, p3, p2):
115 return True
116 if d4 == 0 and on_segment(p1, p4, p2):
117 return True
118
119 return False
120
121
122def line_intersection(p1: Point, p2: Point, p3: Point, p4: Point) -> Optional[Point]:
123 """
124 두 직선의 교점 계산
125 직선 p1-p2와 직선 p3-p4
126 """
127 d1 = p2 - p1
128 d2 = p4 - p3
129 cross = d1.cross(d2)
130
131 if abs(cross) < 1e-9:
132 return None # 평행
133
134 t = (p3 - p1).cross(d2) / cross
135 return p1 + d1 * t
136
137
138# =============================================================================
139# 4. 볼록 껍질 (Convex Hull)
140# =============================================================================
141
142def convex_hull_graham(points: List[Point]) -> List[Point]:
143 """
144 Graham Scan 알고리즘
145 시간복잡도: O(n log n)
146 """
147 if len(points) < 3:
148 return points[:]
149
150 # 가장 아래, 왼쪽 점 찾기
151 start = min(points, key=lambda p: (p.y, p.x))
152
153 # 각도로 정렬
154 def polar_angle(p: Point) -> float:
155 return atan2(p.y - start.y, p.x - start.x)
156
157 def dist_sq(p: Point) -> float:
158 return (p.x - start.x) ** 2 + (p.y - start.y) ** 2
159
160 sorted_points = sorted(points, key=lambda p: (polar_angle(p), dist_sq(p)))
161
162 # 스택으로 볼록 껍질 구성
163 hull = []
164 for p in sorted_points:
165 while len(hull) >= 2 and ccw(hull[-2], hull[-1], p) <= 0:
166 hull.pop()
167 hull.append(p)
168
169 return hull
170
171
172def convex_hull_monotone_chain(points: List[Tuple[float, float]]) -> List[Tuple[float, float]]:
173 """
174 Monotone Chain 알고리즘
175 시간복잡도: O(n log n)
176 """
177 points = sorted(set(points))
178 if len(points) <= 1:
179 return points
180
181 # 하부 껍질
182 lower = []
183 for p in points:
184 while len(lower) >= 2 and ccw_tuple(lower[-2], lower[-1], p) <= 0:
185 lower.pop()
186 lower.append(p)
187
188 # 상부 껍질
189 upper = []
190 for p in reversed(points):
191 while len(upper) >= 2 and ccw_tuple(upper[-2], upper[-1], p) <= 0:
192 upper.pop()
193 upper.append(p)
194
195 return lower[:-1] + upper[:-1]
196
197
198# =============================================================================
199# 5. 다각형 연산
200# =============================================================================
201
202def polygon_area(vertices: List[Point]) -> float:
203 """
204 다각형 넓이 (신발끈 공식)
205 정점이 반시계 방향으로 정렬되어 있어야 함
206 """
207 n = len(vertices)
208 if n < 3:
209 return 0
210
211 area = 0
212 for i in range(n):
213 j = (i + 1) % n
214 area += vertices[i].cross(vertices[j])
215
216 return abs(area) / 2
217
218
219def polygon_perimeter(vertices: List[Point]) -> float:
220 """다각형 둘레"""
221 n = len(vertices)
222 perimeter = 0
223 for i in range(n):
224 j = (i + 1) % n
225 perimeter += distance(vertices[i], vertices[j])
226 return perimeter
227
228
229def point_in_polygon(point: Point, polygon: List[Point]) -> int:
230 """
231 점이 다각형 내부에 있는지 판정
232 반환: 1 (내부), 0 (경계), -1 (외부)
233 """
234 n = len(polygon)
235 winding = 0
236
237 for i in range(n):
238 j = (i + 1) % n
239 p1, p2 = polygon[i], polygon[j]
240
241 # 경계 위의 점 검사
242 if ccw(p1, p2, point) == 0:
243 if (min(p1.x, p2.x) <= point.x <= max(p1.x, p2.x) and
244 min(p1.y, p2.y) <= point.y <= max(p1.y, p2.y)):
245 return 0
246
247 # Winding number 계산
248 if p1.y <= point.y:
249 if p2.y > point.y and ccw(p1, p2, point) > 0:
250 winding += 1
251 else:
252 if p2.y <= point.y and ccw(p1, p2, point) < 0:
253 winding -= 1
254
255 return 1 if winding != 0 else -1
256
257
258def is_convex(polygon: List[Point]) -> bool:
259 """다각형이 볼록인지 검사"""
260 n = len(polygon)
261 if n < 3:
262 return False
263
264 sign = 0
265 for i in range(n):
266 d = ccw(polygon[i], polygon[(i + 1) % n], polygon[(i + 2) % n])
267 if d != 0:
268 if sign == 0:
269 sign = d
270 elif sign != d:
271 return False
272
273 return True
274
275
276# =============================================================================
277# 6. 점과 직선/선분 거리
278# =============================================================================
279
280def point_to_line_distance(point: Point, line_p1: Point, line_p2: Point) -> float:
281 """점에서 직선까지의 거리"""
282 v = line_p2 - line_p1
283 w = point - line_p1
284 return abs(v.cross(w)) / v.norm()
285
286
287def point_to_segment_distance(point: Point, seg_p1: Point, seg_p2: Point) -> float:
288 """점에서 선분까지의 거리"""
289 v = seg_p2 - seg_p1
290 w = point - seg_p1
291
292 c1 = w.dot(v)
293 if c1 <= 0:
294 return distance(point, seg_p1)
295
296 c2 = v.dot(v)
297 if c2 <= c1:
298 return distance(point, seg_p2)
299
300 t = c1 / c2
301 proj = seg_p1 + v * t
302 return distance(point, proj)
303
304
305# =============================================================================
306# 7. 가장 가까운 점 쌍 (Closest Pair)
307# =============================================================================
308
309def closest_pair(points: List[Point]) -> Tuple[Point, Point, float]:
310 """
311 가장 가까운 점 쌍 찾기
312 분할 정복: O(n log n)
313 """
314 def brute_force(pts: List[Point]) -> Tuple[Point, Point, float]:
315 min_dist = inf
316 p1, p2 = None, None
317 for i in range(len(pts)):
318 for j in range(i + 1, len(pts)):
319 d = distance(pts[i], pts[j])
320 if d < min_dist:
321 min_dist = d
322 p1, p2 = pts[i], pts[j]
323 return p1, p2, min_dist
324
325 def closest_split(pts_y: List[Point], mid_x: float, delta: float) -> Tuple[Point, Point, float]:
326 # 중간선 근처의 점들만
327 strip = [p for p in pts_y if abs(p.x - mid_x) < delta]
328
329 min_dist = delta
330 p1, p2 = None, None
331
332 for i in range(len(strip)):
333 j = i + 1
334 while j < len(strip) and strip[j].y - strip[i].y < min_dist:
335 d = distance(strip[i], strip[j])
336 if d < min_dist:
337 min_dist = d
338 p1, p2 = strip[i], strip[j]
339 j += 1
340
341 return p1, p2, min_dist
342
343 def divide_conquer(pts_x: List[Point], pts_y: List[Point]) -> Tuple[Point, Point, float]:
344 if len(pts_x) <= 3:
345 return brute_force(pts_x)
346
347 mid = len(pts_x) // 2
348 mid_point = pts_x[mid]
349
350 # 왼쪽/오른쪽 분할
351 left_x = pts_x[:mid]
352 right_x = pts_x[mid:]
353
354 left_set = set(id(p) for p in left_x)
355 left_y = [p for p in pts_y if id(p) in left_set]
356 right_y = [p for p in pts_y if id(p) not in left_set]
357
358 # 재귀 호출
359 l1, l2, left_dist = divide_conquer(left_x, left_y)
360 r1, r2, right_dist = divide_conquer(right_x, right_y)
361
362 if left_dist < right_dist:
363 best = (l1, l2, left_dist)
364 else:
365 best = (r1, r2, right_dist)
366
367 # 분할선 근처 검사
368 s1, s2, split_dist = closest_split(pts_y, mid_point.x, best[2])
369 if split_dist < best[2]:
370 return s1, s2, split_dist
371
372 return best
373
374 if len(points) < 2:
375 return None, None, inf
376
377 pts_x = sorted(points, key=lambda p: p.x)
378 pts_y = sorted(points, key=lambda p: p.y)
379
380 return divide_conquer(pts_x, pts_y)
381
382
383# =============================================================================
384# 8. 회전하는 캘리퍼스 (Rotating Calipers)
385# =============================================================================
386
387def rotating_calipers_diameter(hull: List[Point]) -> Tuple[Point, Point, float]:
388 """
389 볼록 껍질의 지름 (가장 먼 점 쌍)
390 시간복잡도: O(n)
391 """
392 n = len(hull)
393 if n < 2:
394 return None, None, 0
395 if n == 2:
396 return hull[0], hull[1], distance(hull[0], hull[1])
397
398 # 가장 먼 점 쌍 찾기
399 max_dist = 0
400 p1, p2 = None, None
401
402 j = 1
403 for i in range(n):
404 # hull[i]-hull[(i+1)%n] 에지에 대해 가장 먼 점 찾기
405 while True:
406 next_j = (j + 1) % n
407 # 삼각형 넓이 비교
408 area1 = abs((hull[(i + 1) % n] - hull[i]).cross(hull[j] - hull[i]))
409 area2 = abs((hull[(i + 1) % n] - hull[i]).cross(hull[next_j] - hull[i]))
410 if area2 > area1:
411 j = next_j
412 else:
413 break
414
415 d1 = distance(hull[i], hull[j])
416 d2 = distance(hull[(i + 1) % n], hull[j])
417
418 if d1 > max_dist:
419 max_dist = d1
420 p1, p2 = hull[i], hull[j]
421 if d2 > max_dist:
422 max_dist = d2
423 p1, p2 = hull[(i + 1) % n], hull[j]
424
425 return p1, p2, max_dist
426
427
428# =============================================================================
429# 9. 반평면 교집합 (Half-Plane Intersection)
430# =============================================================================
431
432class HalfPlane:
433 """반평면: ax + by + c >= 0"""
434 def __init__(self, a: float, b: float, c: float):
435 self.a = a
436 self.b = b
437 self.c = c
438 self.angle = atan2(b, a)
439
440 @classmethod
441 def from_points(cls, p1: Point, p2: Point) -> 'HalfPlane':
442 """p1 → p2 방향의 왼쪽 반평면"""
443 a = p2.y - p1.y
444 b = p1.x - p2.x
445 c = -a * p1.x - b * p1.y
446 return cls(a, b, c)
447
448 def side(self, p: Point) -> float:
449 """점이 어느 쪽에 있는지 (양수: 내부, 음수: 외부)"""
450 return self.a * p.x + self.b * p.y + self.c
451
452
453def half_plane_intersection_point(h1: HalfPlane, h2: HalfPlane) -> Optional[Point]:
454 """두 반평면 경계선의 교점"""
455 det = h1.a * h2.b - h2.a * h1.b
456 if abs(det) < 1e-9:
457 return None
458 x = (h1.b * h2.c - h2.b * h1.c) / det
459 y = (h2.a * h1.c - h1.a * h2.c) / det
460 return Point(x, y)
461
462
463# =============================================================================
464# 10. 실전 문제: 삼각형 넓이
465# =============================================================================
466
467def triangle_area(p1: Point, p2: Point, p3: Point) -> float:
468 """삼각형 넓이"""
469 return abs((p2 - p1).cross(p3 - p1)) / 2
470
471
472def triangle_circumcircle(p1: Point, p2: Point, p3: Point) -> Tuple[Point, float]:
473 """삼각형의 외접원 (중심, 반지름)"""
474 ax, ay = p1.x, p1.y
475 bx, by = p2.x, p2.y
476 cx, cy = p3.x, p3.y
477
478 d = 2 * (ax * (by - cy) + bx * (cy - ay) + cx * (ay - by))
479 if abs(d) < 1e-9:
480 return None, 0
481
482 ux = ((ax**2 + ay**2) * (by - cy) + (bx**2 + by**2) * (cy - ay) + (cx**2 + cy**2) * (ay - by)) / d
483 uy = ((ax**2 + ay**2) * (cx - bx) + (bx**2 + by**2) * (ax - cx) + (cx**2 + cy**2) * (bx - ax)) / d
484
485 center = Point(ux, uy)
486 radius = distance(center, p1)
487
488 return center, radius
489
490
491# =============================================================================
492# 테스트
493# =============================================================================
494
495def main():
496 print("=" * 60)
497 print("기하 알고리즘 (Computational Geometry) 예제")
498 print("=" * 60)
499
500 # 1. CCW
501 print("\n[1] CCW (방향 판별)")
502 a, b, c = Point(0, 0), Point(4, 0), Point(2, 2)
503 d = Point(2, -2)
504 print(f" A={a}, B={b}, C={c}")
505 print(f" CCW(A,B,C) = {ccw(a, b, c)} (반시계)")
506 print(f" CCW(A,B,D) = {ccw(a, b, d)} (시계, D={d})")
507
508 # 2. 선분 교차
509 print("\n[2] 선분 교차 판정")
510 p1, p2 = Point(0, 0), Point(4, 4)
511 p3, p4 = Point(0, 4), Point(4, 0)
512 print(f" 선분1: {p1}-{p2}")
513 print(f" 선분2: {p3}-{p4}")
514 print(f" 교차: {segments_intersect(p1, p2, p3, p4)}")
515
516 intersection = line_intersection(p1, p2, p3, p4)
517 print(f" 교점: {intersection}")
518
519 # 3. 볼록 껍질
520 print("\n[3] 볼록 껍질 (Convex Hull)")
521 points = [Point(0, 0), Point(1, 1), Point(2, 2), Point(0, 2),
522 Point(2, 0), Point(1, 0), Point(0, 1), Point(2, 1)]
523 hull = convex_hull_graham(points)
524 print(f" 점들: {points}")
525 print(f" 볼록 껍질: {hull}")
526
527 # 튜플 버전
528 pts_tuple = [(0, 0), (1, 1), (2, 2), (0, 2), (2, 0), (1, 0), (0, 1), (2, 1)]
529 hull_tuple = convex_hull_monotone_chain(pts_tuple)
530 print(f" Monotone Chain: {hull_tuple}")
531
532 # 4. 다각형 연산
533 print("\n[4] 다각형 연산")
534 polygon = [Point(0, 0), Point(4, 0), Point(4, 3), Point(0, 3)]
535 print(f" 다각형: {polygon}")
536 print(f" 넓이: {polygon_area(polygon)}")
537 print(f" 둘레: {polygon_perimeter(polygon)}")
538 print(f" 볼록: {is_convex(polygon)}")
539
540 test_point = Point(2, 1)
541 print(f" 점 {test_point} 위치: {point_in_polygon(test_point, polygon)} (1=내부)")
542
543 # 5. 점-선분 거리
544 print("\n[5] 점에서 선분까지 거리")
545 point = Point(2, 3)
546 seg_start = Point(0, 0)
547 seg_end = Point(4, 0)
548 dist = point_to_segment_distance(point, seg_start, seg_end)
549 print(f" 점: {point}")
550 print(f" 선분: {seg_start}-{seg_end}")
551 print(f" 거리: {dist}")
552
553 # 6. 가장 가까운 점 쌍
554 print("\n[6] 가장 가까운 점 쌍")
555 rand_points = [Point(2, 3), Point(12, 30), Point(40, 50),
556 Point(5, 1), Point(12, 10), Point(3, 4)]
557 p1, p2, dist = closest_pair(rand_points)
558 print(f" 점들: {rand_points}")
559 print(f" 가장 가까운 쌍: {p1}, {p2}")
560 print(f" 거리: {dist:.4f}")
561
562 # 7. 회전 캘리퍼스
563 print("\n[7] 볼록 껍질 지름 (회전 캘리퍼스)")
564 hull_for_diameter = [Point(0, 0), Point(4, 0), Point(5, 2), Point(3, 4), Point(1, 3)]
565 p1, p2, diam = rotating_calipers_diameter(hull_for_diameter)
566 print(f" 볼록 껍질: {hull_for_diameter}")
567 print(f" 가장 먼 쌍: {p1}, {p2}")
568 print(f" 지름: {diam:.4f}")
569
570 # 8. 삼각형
571 print("\n[8] 삼각형 연산")
572 t1, t2, t3 = Point(0, 0), Point(4, 0), Point(2, 3)
573 print(f" 삼각형: {t1}, {t2}, {t3}")
574 print(f" 넓이: {triangle_area(t1, t2, t3)}")
575 center, radius = triangle_circumcircle(t1, t2, t3)
576 print(f" 외접원 중심: {center}")
577 print(f" 외접원 반지름: {radius:.4f}")
578
579 # 9. 알고리즘 복잡도
580 print("\n[9] 알고리즘 복잡도")
581 print(" | 알고리즘 | 시간 복잡도 |")
582 print(" |-------------------|----------------|")
583 print(" | CCW | O(1) |")
584 print(" | 선분 교차 | O(1) |")
585 print(" | 볼록 껍질 | O(n log n) |")
586 print(" | 점 in 다각형 | O(n) |")
587 print(" | 가장 가까운 쌍 | O(n log n) |")
588 print(" | 회전 캘리퍼스 | O(n) |")
589
590 print("\n" + "=" * 60)
591
592
593if __name__ == "__main__":
594 main()