Triangulación de Delaunay — la regla del circuncírculo vacío e inserción Bowyer

El circuncírculo vacío se mantiene borrando solo los triángulos malos; un signo cerca de cero pasa a un predicado exacto.

En 1934, Boris Delaunay definió una triangulación en honor a su mentor Georgy Voronoi. La regla cabía en una sola línea: ningún otro punto puede quedar dentro del circuncírculo de un triángulo. Esa única condición se volvió el corazón de casi todo generador de mallas no estructuradas que se usa hoy.

Un solo circuncírculo vacío decide si la malla sirve. Bowyer–Watson inserta puntos uno a uno y conserva esa propiedad. El punto flotante traiciona la regla en silencio.

Al subir el vértice B, crece el circuncírculo del triángulo A–B–C.

vertex B height (drag B up toward the diagonal)

The circle is the circumcircle of triangle A–B–C. Slide B upward: when D crosses into the circle the diagonal A–C becomes illegal and the Delaunay solution flips to B–D. Notice the larger of the two min-angle numbers always belongs to the chosen diagonal — Delaunay maximizes the smallest angle.

En el instante en que el punto D entra en el círculo, la diagonal A–C se vuelve ilegal y la solución de Delaunay salta a B–D. El mayor de los dos números de min-angle queda en la diagonal elegida.

Una malla definida por un solo círculo vacío#

Hay infinitas maneras de conectar un conjunto de puntos en triángulos. Delaunay elige exactamente una: la triangulación en la que el circuncírculo de cada triángulo (el círculo que pasa por sus tres vértices) está vacío.

Sea $C(T)$ el circuncírculo del triángulo $T = \{a, b, c\}$ . La condición se escribe así:

\forall\, p \notin T : \quad p \notin \mathrm{int}\, C(T)

Aquí $\mathrm{int}\,C(T)$ es el interior del circuncírculo y $p$ es cualquier otro punto de la malla. Ningún punto puede quedar dentro del circuncírculo de otro triángulo.

¿Por qué esta regla? Entre todas las triangulaciones de un conjunto de puntos, la de Delaunay maximiza el ángulo mínimo. Es decir, evita en lo posible los triángulos largos y delgados (slivers). Los triángulos finos arruinan la reconstrucción de gradientes y el número de condición de la matriz en métodos de volúmenes y elementos finitos. Por eso los generadores de mallas toman Delaunay como opción por defecto.

Como bono, el dual de una triangulación de Delaunay es justamente el diagrama de Voronoi. Al obtener una estructura, la otra llega gratis.

El determinante InCircle — probar "dentro" de un solo golpe#

Se puede responder "¿está el punto dentro del circuncírculo?" sin calcular nunca el centro ni el radio del círculo. Para un triángulo $(a, b, c)$ ordenado en sentido antihorario, el punto $d$ está dentro del circuncírculo si y solo si este determinante es positivo:

\det \begin{pmatrix} a_x - d_x & a_y - d_y & (a_x - d_x)^2 + (a_y - d_y)^2 \\ b_x - d_x & b_y - d_y & (b_x - d_x)^2 + (b_y - d_y)^2 \\ c_x - d_x & c_y - d_y & (c_x - d_x)^2 + (c_y - d_y)^2 \end{pmatrix} > 0

Cada fila traslada un vértice a coordenadas relativas a $d$ y luego añade su distancia al cuadrado. El signo de este valor depende de la orientación del triángulo: se invierte con un giro horario, así que el código real lo normaliza multiplicando por el signo de orientación.

El instante en que D entró al círculo, al subir B, es el cambio de signo de este determinante.

Bowyer–Watson — suelta un punto y rellena el hueco#

Ahora se construye la malla. Otra forma de cubrir los mismos puntos es un frente de avance. Bowyer–Watson es un algoritmo incremental que inserta puntos uno a uno y mantiene la propiedad de Delaunay en cada paso. El flujo tiene cinco pasos.

Empezar con un supertriángulo lo bastante grande como para contener todos los puntos.
Para un punto nuevo $p$ , encontrar todos los triángulos cuyo circuncírculo contiene a $p$ — los llamamos triángulos malos.
Los triángulos malos siempre forman una sola cavidad en forma de estrella (star-shaped). Su frontera es el conjunto de aristas que pertenecen a exactamente un triángulo malo.
Borrar los triángulos malos y conectar $p$ con cada arista de la frontera, formando triángulos nuevos.
Tras insertar todos los puntos, quitar cualquier triángulo que comparta un vértice con el supertriángulo.

La clave es la propiedad del paso 3. Al insertar $p$ , los únicos triángulos que rompen la condición del circuncírculo vacío son "los triángulos cuyo circuncírculo contiene a $p$ ", y esos siempre forman una cavidad conexa. Basta con rellenar esa cavidad para que toda la malla vuelva a ser Delaunay.

El lienzo que se reorganiza con cada punto queda justo debajo de la salida de 60 puntos.

Volteo de Lawson — la elección local de una diagonal#

Donde Bowyer–Watson rellena toda una cavidad, el método de Lawson arregla una arista a la vez. Divide en tres el triángulo que contiene al punto nuevo y luego aplica la prueba in-circle a cada arista afectada. Si el punto opuesto cae dentro del circuncírculo, voltea la arista compartida (la diagonal). Repite hasta que no quede ninguna arista ilegal.

Un solo volteo siempre aumenta el ángulo mínimo de los dos triángulos implicados. Así que cuando el bucle termina, queda automáticamente una malla de Delaunay con el ángulo mínimo maximizado. El cambio A–C ↔ B–D de la figura de arriba es exactamente uno de esos volteos.

La misma herramienta resuelve la Delaunay con restricciones (constrained Delaunay), cuando hay que conservar una arista de frontera. Encuentra la tubería (pipe) de triángulos que cruza una arista de recuperación AB y vuelve a triangular esa franja intercambiando diagonales hasta que AB reaparezca como arista de la malla.

Python — 60 puntos aleatorios en una malla de Delaunay#

Toda la inserción incremental llevada a numpy. Corre sin ninguna biblioteca externa.

import numpy as np
 
def in_circle(a, b, c, d):
    """True si el punto d está estrictamente dentro del circuncírculo de a,b,c.
    Normalizado por el signo de orientación, así el winding no importa."""
    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 = {}                              # cuenta cuántas veces aparece cada arista
        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:                       # rellena la cavidad con el punto 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)}")

La salida es, aproximadamente:

points = 60   triangles = 110

show circumcircles

Each thin cyan ring is the circumcircle of one triangle. In a valid Delaunay mesh no vertex (yellow dot) ever sits inside another triangle's ring — that is the empty-circumcircle property. Add points and watch the mesh repair itself.

Cada vez que se agrega un punto, la malla se reorganiza sola. Los anillos cian delgados son los circuncírculos. Con show circumcircles activo, ningún punto amarillo cae dentro del anillo de otro triángulo.

El número de triángulos no es casualidad. En una malla de Delaunay plana la cantidad de triángulos queda fija en $2N - 2 - h$ , donde $N$ es el número de puntos y $h$ el número de puntos sobre la envolvente convexa. Con $N = 60$ y 8 puntos en la envolvente, salen exactamente 110. Se puede verificar el conteo incluso antes de dibujar la malla.

El 110 se cuenta a mano así. Dos por 60 es 120. Al restar 2 quedan 118. Al restar los 8 puntos de la envolvente quedan 110. Si se olvida el término de la envolvente, aparece 118, y una salida de 110 parece un fallo. Los 8 que faltan son h. La semilla 7, sesenta puntos en el cuadrado unitario, es lo que produjo ese 8.

La prueba in-circle cuenta un punto como interior solo cuando el determinante por el signo de orientación supera 1e-12. Cuatro puntos en un mismo círculo dan determinante 0, así que la prueba no se activa. Ninguna diagonal cuenta como interior, y se conserva la que se tendió primero. Un positivo minúsculo entre 0 y 1e-12 también se trata como exterior. Las cuatro esquinas de cada cuadrado caen en esa banda.

En el supertriángulo, r es la distancia máxima al centroide más 1. Los vértices están en (cx−3r, cy−r), (cx+3r, cy−r) y (cx, cy+3r). La base mide 6r de ancho y el tramo vertical mide 4r. Los puntos de muestra caben en un disco de radio r−1, y los vértices súper quedan a unas 3r del centro. Los 60 puntos ya están dentro antes de la inserción. Al borrar los triángulos que tocan un vértice súper, quedan 110.

Cuando los puntos casi caen sobre un mismo círculo#

La regla elegante tiene una trampa. Cuando cuatro puntos caen exactamente sobre un mismo círculo (cocirculares), el determinante in-circle vale cero. Cualquiera de las dos diagonales es Delaunay y la triangulación no es única. El peor caso es una malla cuadrada: las cuatro esquinas de cada cuadrado son cocirculares, así que la entrada más regular resulta ser la más ambigua.

El problema de verdad surge cerca de ese cero, por culpa del punto flotante. Cuando el determinante alcanza el límite de precisión de máquina, su signo se invierte al azar, la topología cae en contradicción y el código se queda en un bucle infinito o se cae. Por eso bibliotecas de producción como Triangle y CGAL usan predicados exactos adaptativos (adaptive exact predicates, Shewchuk) en lugar de un simple double. Calculan solo los dígitos necesarios para acertar el signo siempre.

En tres dimensiones aparece una trampa más. La condición de circunesfera vacía no prohíbe tetraedros casi planos (slivers). Por eso el mallado 3D agrega refinamiento de Delaunay (Ruppert, Chew), insertando puntos de Steiner para subir la calidad.

La última trampa práctica es la frontera. Una malla de Delaunay cruda no respeta la frontera geométrica dibujada. Para forzar que una arista de frontera sobreviva hace falta el paso de Delaunay con restricciones de arriba. "¿Por qué la malla atraviesa la pared?" es la señal de que ese paso se omitió.

Si el signo in-circle se cree tal cual en una malla cuadrada#

Cuatro esquinas sobre un mismo círculo dejan el determinante en 0. El umbral 1e-12 no llama a eso interior. Cualquiera de las dos diagonales es Delaunay, pero el código conserva la que apareció primero. Una malla cuadrada entera es ese empate, así que la entrada más regular es la que sacude la topología.

El conteo de triángulos es el filtro barato. En N=60, un resultado de 118 significa que se omitió el término de envolvente h. Un resultado de 110 coincide con h=8 de esta semilla. Si el conteo cuadra y la pared igual se corta, lo sospechoso es la recuperación de frontera que falta, no el circuncírculo vacío.

Triangulación de Delaunay — la regla del circuncírculo vacío e inserción Bowyer–Watson