引言
德洛内三角剖分（Delaunay Triangulation）是计算几何中的核心算法，它将平面上一组离散点连接成三角形网格，使得任意三角形的外接圆内不包含其他点。这一性质被称为”空圆性质”（Empty Circle Property），它保证了生成的三角形尽可能”规则”（接近等边），避免了细长条三角形的出现。德洛内三角剖分在有限元分析、地形建模、计算机图形学、地理信息系统等领域有着不可替代的作用。本文将讲解其数学原理，实现经典的Bowyer-Watson增量插入算法，并提供完整的Java代码。
德洛内三角剖分的数学原理
空圆性质
对于三角剖分中的任意三角形 $\triangle ABC$，其外接圆内（不含边界）不包含其他任何输入点。形式化地：
$$\forall P \in S \setminus {A, B, C}, \quad P \notin \text{circumcircle}(A, B, C)$$
最大化最小角性质
在所有可能的三角剖分中，德洛内三角剖分最大化最小内角。这意味着它尽可能避免了”扁平”的三角形，使得网格质量最优。
与Voronoi图的对偶关系
德洛内三角剖分与Voronoi图互为对偶图：
– 德洛内三角剖分的顶点 ↔ Voronoi图的生成点
– 德洛内三角剖分的边 ↔ Voronoi图的边（垂直平分线）
– 德洛内三角剖分的三角形中心 ↔ Voronoi图的顶点（外接圆圆心）
Bowyer-Watson增量插入算法
Bowyer-Watson算法是构造德洛内三角剖分的经典增量算法，其核心思想是逐个插入点，并在每次插入后维护德洛内性质。
算法步骤
1. 构造一个足够大的"超级三角形"，包含所有输入点
2. 初始化三角剖分：triangles = {超级三角形}
3. 对于每个输入点 p:
a. 找出所有外接圆包含 p 的三角形（"坏三角形"）
b. 找出坏三角形的边界（即不被两个坏三角形共享的边）
c. 删除所有坏三角形
d. 将 p 与边界上的每条边连接，形成新三角形
4. 删除与超级三角形顶点相关的所有三角形
5. 返回最终的三角剖分
点在外接圆内的判定
给定三角形顶点 A(x₁,y₁), B(x₂,y₂), C(x₃,y₃)，待测点 P(x,y)，使用行列式判定：
$$\begin{vmatrix} x_1 – x & y_1 – y & (x_1 – x)^2 + (y_1 – y)^2 \ x_2 – x & y_2 – y & (x_2 – x)^2 + (y_2 – y)^2 \ x_3 – x & y_3 – y & (x_3 – x)^2 + (y_3 – y)^2 \end{vmatrix} > 0$$
行列式大于0表示P在三角形ABC外接圆内部（假设ABC为逆时针顺序）。
Java实现
import javax.swing.*;
import java.awt.*;
import java.util.*;
import java.util.List;
/**
* 德洛内三角剖分 - Bowyer-Watson增量插入算法
*/
public class DelaunayTriangulation extends JPanel {
private List<Point> points; // 输入点集
private List<Triangle> triangles; // 三角剖分结果
private static final double EPS = 1e-9;
public DelaunayTriangulation(List<Point> points) {
this.points = new ArrayList<>(points);
this.triangles = computeDelaunay(points);
}
/**
* Bowyer-Watson增量插入算法主流程
*/
public List<Triangle> computeDelaunay(List<Point> inputPoints) {
List<Triangle> triangulation = new ArrayList<>();
// 步骤1：构造超级三角形，确保包含所有点
double minX = Double.MAX_VALUE, minY = Double.MAX_VALUE;
double maxX = -Double.MAX_VALUE, maxY = -Double.MAX_VALUE;
for (Point p : inputPoints) {
minX = Math.min(minX, p.x); minY = Math.min(minY, p.y);
maxX = Math.max(maxX, p.x); maxY = Math.max(maxY, p.y);
}
double dx = maxX - minX, dy = maxY - minY;
double deltaMax = Math.max(dx, dy);
double midX = (minX + maxX) / 2, midY = (minY + maxY) / 2;
// 超级三角形的三个顶点
Point p1 = new Point(midX - 20 * deltaMax, midY - deltaMax);
Point p2 = new Point(midX, midY + 20 * deltaMax);
Point p3 = new Point(midX + 20 * deltaMax, midY - deltaMax);
triangulation.add(new Triangle(p1, p2, p3));
// 步骤2-3：逐点插入
for (Point point : inputPoints) {
List<Triangle> badTriangles = new ArrayList<>();
// 找出所有外接圆包含当前点的三角形
for (Triangle t : triangulation) {
if (t.pointInCircumcircle(point)) {
badTriangles.add(t);
}
}
// 找出坏三角形的边界（多边形空洞）
List<Edge> polygon = new ArrayList<>();
for (Triangle t : badTriangles) {
for (Edge e : t.getEdges()) {
boolean shared = false;
for (Triangle other : badTriangles) {
if (other == t) continue;
if (other.hasEdge(e)) {
shared = true;
break;
}
}
if (!shared) {
polygon.add(e);
}
}
}
// 删除坏三角形
triangulation.removeAll(badTriangles);
// 将当前点与边界连接，形成新三角形
for (Edge e : polygon) {
triangulation.add(new Triangle(e.p1, e.p2, point));
}
}
// 步骤4：删除与超级三角形顶点相关的三角形
List<Triangle> result = new ArrayList<>();
for (Triangle t : triangulation) {
if (!t.containsVertex(p1) && !t.containsVertex(p2) && !t.containsVertex(p3)) {
result.add(t);
}
}
return result;
}
/**
* 三角形类
*/
static class Triangle {
Point a, b, c; // 三个顶点
Point circumcenter; // 外接圆圆心
double circumradius; // 外接圆半径
Triangle(Point a, Point b, Point c) {
this.a = a; this.b = b; this.c = c;
computeCircumcircle();
}
/**
* 计算外接圆（圆心和半径）
*/
private void computeCircumcircle() {
double d = 2 * (a.x * (b.y - c.y) + b.x * (c.y - a.y) + c.x * (a.y - b.y));
if (Math.abs(d) < EPS) {
circumcenter = new Point((a.x + b.x + c.x) / 3, (a.y + b.y + c.y) / 3);
circumradius = Double.MAX_VALUE;
return;
}
double ux = ((a.x * a.x + a.y * a.y) * (b.y - c.y)
+ (b.x * b.x + b.y * b.y) * (c.y - a.y)
+ (c.x * c.x + c.y * c.y) * (a.y - b.y)) / d;
double uy = ((a.x * a.x + a.y * a.y) * (c.x - b.x)
+ (b.x * b.x + b.y * b.y) * (a.x - c.x)
+ (c.x * c.x + c.y * c.y) * (b.x - a.x)) / d;
circumcenter = new Point(ux, uy);
circumradius = dist(circumcenter, a);
}
/**
* 判定点是否在外接圆内部（严格内部）
* 使用行列式法，假设三角形顶点为逆时针顺序
*/
boolean pointInCircumcircle(Point p) {
double ax = a.x - p.x, ay = a.y - p.y;
double bx = b.x - p.x, by = b.y - p.y;
double cx = c.x - p.x, cy = c.y - p.y;
double det = ax * (by * (cx * cx + cy * cy) - (bx * bx + by * by) * cy)
- ay * (bx * (cx * cx + cy * cy) - (bx * bx + by * by) * cx)
+ (ax * ax + ay * ay) * (bx * cy - by * cx);
return det > EPS;
}
List<Edge> getEdges() {
return Arrays.asList(new Edge(a, b), new Edge(b, c), new Edge(c, a));
}
boolean hasEdge(Edge e) {
return (samePoint(a, e.p1) && samePoint(b, e.p2)) ||
(samePoint(b, e.p1) && samePoint(c, e.p2)) ||
(samePoint(c, e.p1) && samePoint(a, e.p2)) ||
(samePoint(a, e.p2) && samePoint(b, e.p1)) ||
(samePoint(b, e.p2) && samePoint(c, e.p1)) ||
(samePoint(c, e.p2) && samePoint(a, e.p1));
}
boolean containsVertex(Point p) {
return samePoint(a, p) || samePoint(b, p) || samePoint(c, p);
}
private boolean samePoint(Point p1, Point p2) {
return Math.abs(p1.x - p2.x) < EPS && Math.abs(p1.y - p2.y) < EPS;
}
}
/**
* 边类（无向边）
*/
static class Edge {
Point p1, p2;
Edge(Point p1, Point p2) { this.p1 = p1; this.p2 = p2; }
}
/**
* 点类
*/
static class Point {
double x, y;
Point(double x, double y) { this.x = x; this.y = y; }
}
static double dist(Point p1, Point p2) {
return Math.sqrt((p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y));
}
// ================== 可视化 ==================
@Override
protected void paintComponent(Graphics g) {
super.paintComponent(g);
Graphics2D g2d = (Graphics2D) g;
g2d.setRenderingHint(RenderingHints.KEY_ANTIALIASING,
RenderingHints.VALUE_ANTIALIAS_ON);
// 绘制三角形
g2d.setColor(new Color(100, 150, 255, 80));
for (Triangle t : triangles) {
int[] xs = {(int) t.a.x, (int) t.b.x, (int) t.c.x};
int[] ys = {(int) t.a.y, (int) t.b.y, (int) t.c.y};
g2d.fillPolygon(xs, ys, 3);
}
g2d.setColor(new Color(50, 100, 200));
g2d.setStroke(new BasicStroke(1.5f));
for (Triangle t : triangles) {
g2d.drawLine((int) t.a.x, (int) t.a.y, (int) t.b.x, (int) t.b.y);
g2d.drawLine((int) t.b.x, (int) t.b.y, (int) t.c.x, (int) t.c.y);
g2d.drawLine((int) t.c.x, (int) t.c.y, (int) t.a.x, (int) t.a.y);
}
// 绘制输入点
g2d.setColor(Color.RED);
for (Point p : points) {
g2d.fillOval((int) p.x - 4, (int) p.y - 4, 8, 8);
}
g2d.setColor(Color.BLACK);
g2d.drawString("三角形数量: " + triangles.size(), 10, 20);
}
public static void main(String[] args) {
// 随机生成点集
List<Point> pts = new ArrayList<>();
Random rand = new Random(123);
for (int i = 0; i < 50; i++) {
pts.add(new Point(50 + rand.nextInt(500), 50 + rand.nextInt(400)));
}
JFrame frame = new JFrame("德洛内三角剖分 - Bowyer-Watson算法");
frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
frame.setSize(600, 520);
frame.add(new DelaunayTriangulation(pts));
frame.setVisible(true);
}
}
边翻转优化（Lawson算法）
在增量插入后，可以通过局部边翻转进一步优化网格质量。如果两个相邻三角形组成的四边形中，对角线交换后满足德洛内空圆性质，则执行翻转。
/**
* Lawson边翻转优化
* 对于每条内部边，检查翻转是否改善最小角
*/
public void lawsonFlip(List<Triangle> triangulation) {
boolean flipped;
do {
flipped = false;
// 遍历所有边，检查是否需要翻转
for (int i = 0; i < triangulation.size(); i++) {
Triangle t1 = triangulation.get(i);
for (int j = i + 1; j < triangulation.size(); j++) {
Triangle t2 = triangulation.get(j);
// 如果两个三角形共享一条边，检查空圆性质
// 若不满足，则翻转共享边
// 翻转操作：将四边形的另一条对角线作为共享边
// 此处为简化示意
}
}
} while (flipped);
}
复杂度分析
| 阶段 | 时间复杂度 | 说明 |
|---|---|---|
| 超级三角形构造 | O(1) | 固定3个顶点 |
| 逐点插入（最坏） | O(n²) | 每次需检查所有已有三角形 |
| 逐点插入（平均/随机） | O(n log n) | 点分布均匀时 |
| 边翻转优化 | O(n) ~ O(n²) | 通常少量翻转即可收敛 |
| 删除超级三角形 | O(n) | 遍历所有三角形 |
对于 n 个点的均匀分布，Bowyer-Watson算法的平均时间复杂度为 O(n log n)。若使用带空间索引的优化实现（如结合Delaunay树的增量更新），可达到稳定的最优复杂度。
应用场景
- 有限元分析（FEA）：将连续域离散化为三角形/四面体网格进行数值求解
- 地形建模：数字高程模型（DEM）的TIN（不规则三角网）表示
- 计算机图形学：三维模型的三角化、纹理映射参数化
- 地理信息系统：自然邻域插值、流域分析
- 路径规划：可通行区域的多边形近似与导航网格（NavMesh）生成
- 计算机视觉：三维重建中的点云表面重建
总结
德洛内三角剖分是计算几何中最优雅实用的算法之一。本文讲解了其空圆性质和最大化最小角特性，详细实现了Bowyer-Watson增量插入算法，并提供了完整的Java可视化代码。理解德洛内三角剖分与Voronoi图的对偶关系，有助于建立更完整的空间数据结构知识体系。对于大规模点集，可进一步研究分治法或结合空间索引的增量优化策略。
