Delaunay 삼각분할 — 빈 외접원 조건과 Bowyer–Watson 점 삽입

1934년, Boris Delaunay는 스승 Georgy Voronoi를 기리며 삼각분할 하나를 정의했다. 규칙은 단 한 줄이었다. 어떤 삼각형의 외접원 안에도 다른 점이 들어와선 안 된다. 이 조건 하나가 오늘날 거의 모든 비정렬 격자 생성기의 심장이 됐다.

빈 외접원 하나가 격자 품질을 정한다. Bowyer–Watson은 점을 하나씩 꽂아 그 성질을 유지한다. 부동소수점은 그 규칙을 조용히 배신한다.

슬라이더로 점 B를 위로 올리면 삼각형 A–B–C의 외접원이 커진다.

점 D가 원 안으로 들어오는 순간 대각선 A–C는 불법이 되고, Delaunay 해는 B–D로 뒤집힌다. 두 min-angle 숫자 중 큰 쪽이 선택된 대각선에 붙는지 보면 된다.

빈 외접원 하나로 정의되는 격자#

같은 점 집합을 삼각형으로 잇는 방법은 무수히 많다. Delaunay는 그중 단 하나를 고른다. 모든 삼각형의 외접원(세 꼭짓점을 지나는 원)이 비어 있는 분할이다.

삼각형 $T = \{a, b, c\}$의 외접원을 $C(T)$라 하면 조건은 이렇게 적힌다.

여기서 $\mathrm{int}\,C(T)$는 외접원의 내부, $p$는 격자의 다른 모든 점이다. 어떤 점도 다른 삼각형의 외접원 안쪽에 있으면 안 된다.

왜 하필 이 규칙인가. Delaunay 분할은 같은 점 집합의 모든 삼각분할 중 최소 내각을 최대화한다. 가늘고 긴 삼각형(sliver)을 최대한 피한다는 뜻이다. 가는 삼각형은 유한체적·유한요소에서 기울기 재구성과 행렬 조건수를 망가뜨린다. 그래서 격자 생성기는 Delaunay를 기본값으로 삼는다.

덤으로, Delaunay 삼각분할의 쌍대(dual)가 바로 Voronoi 다이어그램이다. 한 구조를 얻으면 다른 하나가 공짜로 따라온다.

InCircle 행렬식 — 원 안인지 한 번에 판정#

"점이 외접원 안에 있는가"를 외접원의 중심과 반지름을 구하지 않고 한 번에 답할 수 있다. 반시계 방향으로 정렬된 삼각형 $(a, b, c)$에 대해, 점 $d$가 외접원 내부에 있을 필요충분조건은 다음 행렬식이 양수인 것이다.

각 행은 한 꼭짓점을 $d$ 기준 상대좌표로 옮긴 뒤, 그 좌표와 거리 제곱을 나란히 둔 것이다. 이 값의 부호는 삼각형의 회전 방향(orientation)에 의존한다. 시계 방향이면 부호가 뒤집히므로, 실제 코드에서는 삼각형의 방향 부호를 곱해 정규화한다.

앞에서 B를 올렸을 때 D가 원 안으로 들어온 순간이, 이 행렬식의 부호가 바뀌는 자리다.

Bowyer–Watson — 점을 꽂고 구멍을 다시 메운다#

이제 격자를 만든다. Bowyer–Watson은 점을 하나씩 삽입하면서 매 단계 Delaunay 성질을 유지하는 증분(incremental) 알고리즘이다. 점을 덮는 다른 철학은 전진 프론트 쪽이다. 여기는 구멍을 지우고 다시 메운다. 흐름은 다섯 단계다.

1. 모든 점을 품을 만큼 큰 슈퍼 삼각형으로 시작한다.
2. 새 점 $p$에 대해, 외접원이 $p$를 품는 삼각형을 전부 찾는다 — 이들을 나쁜 삼각형이라 부른다.
3. 나쁜 삼각형들은 항상 별 모양(star-shaped) 구멍 하나를 이룬다. 이 구멍의 경계는 나쁜 삼각형 정확히 하나에만 속한 변들이다.
4. 나쁜 삼각형을 모두 지우고, $p$를 경계의 각 변과 이어 새 삼각형을 만든다.
5. 모든 점을 넣은 뒤, 슈퍼 삼각형의 꼭짓점을 공유하는 삼각형을 제거한다.

비결은 3단계의 성질에 있다. $p$ 삽입으로 빈 외접원 조건을 깨는 삼각형은 정확히 "외접원이 $p$를 품는 삼각형"뿐이고, 이들은 언제나 연결된 구멍을 만든다. 그 구멍만 다시 채우면 전체가 다시 Delaunay가 된다.

점을 꽂을 때마다 격자가 다시 짜이는 화면은, 60점 출력 바로 뒤에 둔다.

Lawson 뒤집기 — 대각선 하나의 국소 선택#

Bowyer–Watson이 구멍을 통째로 다시 메운다면, Lawson의 방법은 변 하나씩 고친다. 새 점을 품은 삼각형을 셋으로 쪼갠 뒤, 영향을 받은 각 변에 대해 in-circle 검사를 한다. 인접한 반대편 점이 외접원 안에 있으면 그 공유 변(대각선)을 뒤집는다. 더 이상 불법 변이 없을 때까지 반복한다.

이 뒤집기 한 번은 두 삼각형의 최소 내각을 반드시 키운다. 그래서 반복이 끝나면 자동으로 최소 내각이 최대화된 Delaunay 격자가 된다. 위 그림에서 본 A–C ↔ B–D 전환이 바로 이 한 번의 뒤집기다.

경계선을 강제로 살려야 할 때(constrained Delaunay)도 같은 도구가 쓰인다. 복원할 변 AB가 가로지르는 삼각형들의 띠(pipe)를 찾아, 그 안을 대각선 교환으로 다시 분할하면 AB가 격자 변으로 복구된다.

Python — 무작위 60점을 Delaunay로#

증분 삽입 전체를 numpy로 옮긴다. 외부 라이브러리 없이 돌아간다.

import numpy as np
 
def in_circle(a, b, c, d):
    """점 d가 삼각형 a,b,c의 외접원 안에 strictly 들어가면 True.
    삼각형 방향 부호로 정규화해 winding 무관하게 만든다."""
    ax, ay = a[0] - d[0], a[1] - d[1]
    bx, by = b[0] - d[0], b[1] - d[1]
    cx, cy = c[0] - d[0], c[1] - d[1]
    det = np.linalg.det(np.array([
        [ax, ay, ax * ax + ay * ay],
        [bx, by, bx * bx + by * by],
        [cx, cy, cx * cx + cy * cy],
    ]))
    orient = (b[0] - a[0]) * (c[1] - a[1]) - (b[1] - a[1]) * (c[0] - a[0])
    return det * np.sign(orient) > 1e-12
 
def super_triangle(points):
    pts = np.asarray(points)
    cx, cy = pts.mean(axis=0)
    r = np.max(np.linalg.norm(pts - [cx, cy], axis=1)) + 1.0
    return [(cx - 3 * r, cy - r), (cx + 3 * r, cy - r), (cx, cy + 3 * r)]
 
def bowyer_watson(points):
    st = super_triangle(points)
    tris = [(st[0], st[1], st[2])]
    for p in points:
        bad = [t for t in tris if in_circle(t[0], t[1], t[2], p)]
        edge_n = {}                              # 변마다 등장 횟수 집계
        for t in bad:
            for e in [(t[0], t[1]), (t[1], t[2]), (t[2], t[0])]:
                key = frozenset((e[0], e[1]))
                edge_n[key] = edge_n.get(key, 0) + 1
        boundary = [tuple(k) for k, n in edge_n.items() if n == 1]
        tris = [t for t in tris if t not in bad]
        for e in boundary:                       # 점 p로 구멍을 다시 메운다
            tris.append((e[0], e[1], p))
    s = set(st)
    return [t for t in tris if not (set(t) & s)]
 
rng = np.random.default_rng(7)
pts = [tuple(p) for p in rng.random((60, 2))]
mesh = bowyer_watson(pts)
print(f"points = {len(pts)}   triangles = {len(mesh)}")

출력은 대략 이렇다.

points = 60   triangles = 110

점을 추가할 때마다 격자가 스스로 다시 짜인다. 얇은 시안색 원들이 각 삼각형의 외접원이다. show circumcircles를 켜 두면, 어떤 노란 점도 다른 삼각형의 원 안에 들어가지 않는다는 빈 외접원 성질이 눈으로 보인다.

삼각형 수는 우연이 아니다. 평면 Delaunay에서 삼각형 개수는 $2N - 2 - h$로 정해진다 ($N$은 점 수, $h$는 볼록껍질 위의 점 수). $N = 60$이면 볼록껍질 점이 8개일 때 정확히 110개가 나온다. 격자가 채 그려지기 전에 개수만으로 검산할 수 있다.

그 110을 손으로 다시 세면 이렇다. 2 곱하기 60은 120이다. 120에서 2를 빼면 118, 118에서 껍질 점 8을 빼면 110이다. 껍질 항을 빼먹으면 118이 나오고, 코드가 110을 찍으면 버그처럼 보인다. 빠진 8이 $h$다. 반대로 삼각형이 110이면 118에서 110을 뺀 8이 껍질 점 수다. 점을 세어 보기 전에 검산이 끝난다. 시드 7의 단위 정사각형 60점이 그 8을 줬다.

in-circle은 행렬식에 방향 부호를 곱한 값이 1e-12보다 클 때만 내부로 친다. 네 점이 한 원이면 행렬식은 0이라 이 검사를 통과하지 못한다. 어느 대각선도 내부가 아니므로, 알고리즘은 먼저 깔린 대각선을 그대로 둔다. 0과 1e-12 사이의 아주 작은 양수도 바깥으로 처리된다. 정사각 격자는 네 꼭짓점이 늘 이 구간에 걸친다. 가장 규칙적인 입력이 부호를 가장 먼저 잃는다.

슈퍼 삼각형의 r은 무게중심에서 가장 먼 점까지의 거리에 1을 더한 값이다. 꼭짓점은 (cx−3r, cy−r), (cx+3r, cy−r), (cx, cy+3r)이다. 밑변 폭은 6r이고, 아래끝에서 위끝까지는 4r이다. 샘플 점은 반지름 r−1인 원반 안에 있고, 슈퍼 꼭짓점은 중심에서 대략 3r 밖에 있다. 60점을 넣기 전에 이미 그 안에 들어 있다. 마지막에 슈퍼 꼭짓점을 공유하는 삼각형을 지우면 남은 개수가 110이다.

좌표가 거의 한 원 위에 있을 때#

우아한 규칙에는 함정이 있다. 네 점이 정확히 한 원 위에 있으면(cocircular) in-circle 행렬식은 0이 된다. 어느 대각선을 골라도 Delaunay이고, 분할이 유일하지 않다. 최악은 정사각 격자다 — 모든 정사각형의 네 꼭짓점이 cocircular라, 가장 규칙적인 입력이 가장 모호한 입력이 된다.

진짜 문제는 그 0 근처에서 부동소수점이 일으킨다. 행렬식이 머신 정밀도 한계에 닿으면 부호가 무작위로 뒤집히고, 위상(topology)이 모순에 빠져 코드가 무한 루프에 걸리거나 죽는다. 그래서 Triangle, CGAL 같은 실전 라이브러리는 단순 double 대신 adaptive exact predicates(Shewchuk)를 쓴다. 필요한 자릿수만 정확히 계산해 부호만큼은 절대 틀리지 않게 한다.

3차원으로 가면 함정이 하나 더 늘어난다. 빈 외접구 조건은 거의 납작한 사면체(sliver)를 막지 못한다. 그래서 3D에서는 Delaunay refinement(Ruppert·Chew)로 Steiner 점을 추가해 품질을 끌어올린다. 경계 변을 살리는 쪽의 같은 도구는 제약 Delaunay에 따로 적혀 있다.

마지막 실무 함정은 경계다. 날것의 Delaunay는 당신이 그린 형상 경계를 존중하지 않는다. 경계 변을 강제로 살리려면 위에서 본 constrained Delaunay가 따로 필요하다. "왜 격자가 벽을 뚫고 지나가지?"는 이걸 빼먹었다는 신호다.

정사각 격자에서 in-circle 부호를 그대로 쓰면#

네 꼭짓점이 한 원 위에 있으면 행렬식은 0이다. 1e-12 문턱은 그 0을 내부로 인정하지 않는다. 어느 대각선을 골라도 Delaunay인데, 코드는 먼저 생긴 쪽만 남긴다. 정사각 격자 전체가 이 동률이라, 가장 규칙적인 입력이 위상부터 흔들린다.

삼각형 개수로 먼저 거르는 편이 싸다. N=60에서 118이 나오면 껍질 항 h를 빼먹은 것이고, 110이면 이 시드의 h=8과 맞다. 개수가 맞는데도 벽이 뚫리면, 빈 외접원이 아니라 경계 복원을 뺀 쪽을 의심한다.