C++实现增量Delaunay三角剖分：数据结构、翻边与避坑指南-平芜编程栈

简介：本资源面向计算几何、图形学与C++算法学习者，提供一套用C++实现的增量Delaunay三角剖分完整工程，核心目标是最大化剖分中所有三角形的最小内角，从而获得形态更优的三角网格。项目采用顾贤峰教授提出的“半边”数据结构组织网格拓扑，主要实现集中在GaussCurvature_Sample模块的main.cpp中，适合作为算法课程设计、网格生成实验或图形学入门的参考代码。压缩包共89个文件，约4.17MB，包含13个h头文件、1个cpp源文件、1个sln解决方案与vcxproj工程文件，以及exe、pdb、obj等编译产物和tlog、log日志，另有gif演示、md说明与xml配置，便于直接编译运行与调试。目前已有1002人学习下载。读者可借此理解增量插入、边翻转与空外接圆准则的具体落地方式，掌握半边结构的遍历与维护技巧，并对照工程目录快速搭建自己的三角剖分实验环境。

1. 从一堆散点到一张三角网：增量 Delaunay 到底在解决什么

如果你做过点云重建、有限元网格生成、地形建模或者路径规划，大概率绕不开一个基础问题：给平面上几百到几百万个散点，怎么把它们连成一张“质量最好”的三角网。Delaunay 三角剖分就是这个问题最经典的答案，它的核心性质是空外接圆准则——任意一个三角形的外接圆内部不包含其他点。这条性质直接带来两个好处：三角形尽量接近等边，不会出现又细又长的“针状”三角形；同时它是对偶 Voronoi 图的基础，很多几何算法都建立在它之上。

而“增量”两个字，指的是构造方式。相比分治法、扫描线法，增量法更符合工程直觉：先拿三个点搭一个初始三角形，然后一个点一个点地插进去，每插一个点，就把它影响范围内的三角形“翻边”修正，直到整张网重新满足空外接圆准则。用 C++ 实现增量 Delaunay，难点不在数学推导，而在数据结构选型、点定位效率、退化情况处理这三件事上。这篇笔记就按“能跑起来、能调参数、能排错”的顺序，把这条路走一遍。

2. 增量 Delaunay 的算法骨架与数据结构选型

2.1 为什么选增量法而不是分治或扫描线

分治法把点集递归二分，理论复杂度 O(n log n)，但合并两个子三角网时的边拼接逻辑非常绕，调试成本高。扫描线法对点集排序依赖强，遇到重复点、共线点容易翻车。增量法单次插入最坏 O(n)，朴素实现整体 O(n²)，但配合点定位优化后平均能压到接近 O(n log n)，而且每一步状态都是合法的三角网，随时可以中断、可视化、增量更新。

我一般会这样判断：如果点集是静态的、一次性几百万点，分治更划算；如果点集是流式的、需要边插边查、或者要支持后续删点，增量法是首选。工程里绝大多数交互式建模场景都属于后者，所以增量法的实用价值更高。

2.2 半边结构与三角形邻接表怎么选

增量法最频繁的操作是“给定一个三角形，找它的三个邻居”和“翻转一条边”。用半边（half-edge）结构最自然：每条边拆成两条有向半边，每个半边记录起点、孪生半边、所属面、下一条半边。翻转操作就是改几个指针，O(1) 完成。

但半边结构实现门槛高，指针一多就容易出空指针。另一种做法是三角形数组加邻接索引：每个三角形存三个顶点索引和三个邻居三角形索引，翻转时同步更新两侧。这种写法内存连续、缓存友好，调试时打印也直观。下面是我常用的三角形结构定义：

struct Triangle { int v[3]; // 三个顶点在点数组中的索引，逆时针排列 int n[3]; // n[i] 表示与边 (v[i], v[(i+1)%3]) 相对的邻居三角形索引 bool alive; // 标记是否已被删除，增量法里删除比物理移除更省事 }; struct Point { double x, y; int id; // 原始输入序号，便于回溯 };

v逆时针排列是约定，保证有向面积恒为正，后续 incircle 判断才不会符号混乱。n[i]的语义要记牢：它对应的是“边 v[i]→v[i+1] 对面的那个三角形”，不是共享这条边的三角形本身。这个索引关系搞反，翻边时就会连环错。

2.3 超级三角形与点定位策略

增量法需要一个“超级三角形”把所有点包住，避免插入边界点时处理无穷远。常见做法是取点集包围盒，向外扩一个足够大的边距，构造一个覆盖全场的三角形。注意这个三角形要足够大，否则后续插入外围点时，新三角形会超出超级三角形范围，导致翻边逻辑失效。

点定位决定整体效率。朴素做法是遍历所有存活三角形，判断点是否在内部，O(n) 一次插入。优化方向有两个：一是维护一个“最后插入的三角形”，利用空间局部性，从它出发沿邻接关系走向目标三角形（walk 算法）；二是建网格桶或跳表做粗定位。我一般先用 walk，代码量小、效果立竿见影。

// 从 start 三角形出发，沿邻接关系走向包含点 p 的三角形 int locate(const Point& p, int start, const std::vector<Triangle>& tris) { int cur = start; while (true) { const Triangle& t = tris[cur]; bool moved = false; for (int i = 0; i < 3; ++i) { int a = t.v[i], b = t.v[(i + 1) % 3]; // 若 p 在有向边 a->b 的右侧，说明目标在邻居方向 if (orient(points[a], points[b], p) < 0) { cur = t.n[i]; moved = true; break; } } if (!moved) return cur; // 三条边都在左侧，p 落在当前三角形内 } }

orient返回有向面积的两倍，正值表示 p 在 a→b 左侧。这个 walk 过程在点分布均匀时平均步数是个小常数，但点集极度聚集时会退化，所以生产代码里通常加一个步数上限，超限就回退到暴力遍历。

3. 插入一个点：翻边修正的完整实现

3.1 插入流程拆成四步

插入点 p 的完整流程是：定位到包含 p 的三角形；把该三角形一分为三；对三条新边分别检查是否满足空外接圆准则；不满足就翻转，翻转后继续检查受影响的边，直到稳定。这四步里，第三步和第四步是循环嵌套的，用一个栈来存“待检查边”最清晰。

void insertPoint(const Point& p, int startTri) { int t = locate(p, startTri, tris); // 第一步：一分为三，记录三条新边 int a = tris[t].v[0], b = tris[t].v[1], c = tris[t].v[2]; int t1 = newTriangle(a, b, p.id); int t2 = newTriangle(b, c, p.id); int t3 = newTriangle(c, a, p.id); // 建立三个新三角形之间的邻接，以及与原邻居的邻接 link(t1, t2, t3, t, a, b, c); tris[t].alive = false; // 第二步：三条新边入栈 std::stack<Edge> stk; stk.push({a, b, t1}); stk.push({b, c, t2}); stk.push({c, a, t3}); // 第三步：循环翻边 while (!stk.empty()) { Edge e = stk.top(); stk.pop(); if (!isDelaunay(e)) { flipAndPush(e, stk); } } }

newTriangle负责分配新三角形并初始化邻接为 -1，link负责把新三角形和原三角形的三个邻居接上。这里最容易出错的是邻接更新顺序：先接新三角形之间的边，再接原邻居，最后把原三角形标记删除。顺序反了会出现邻居指向已删除三角形的情况。

3.2 空外接圆判断与翻边操作

判断一条边是否 Delaunay，本质是看这条边对面的两个三角形，其中一个的外接圆是否包含另一个的第三个顶点。用 incircle 行列式判断，符号决定是否需要翻转。

// 返回正值表示 d 在 a,b,c 的外接圆内（假设 a,b,c 逆时针） double incircle(const Point& a, const Point& b, const Point& c, const Point& d) { double ax = a.x - d.x, ay = a.y - d.y; double bx = b.x - d.x, by = b.y - d.y; double cx = c.x - d.x, cy = c.y - d.y; return (ax * ax + ay * ay) * (bx * cy - by * cx) - (bx * bx + by * by) * (ax * cy - ay * cx) + (cx * cx + cy * cy) * (ax * by - ay * bx); }

翻边操作是把两个共享边的三角形组成的四边形，换成另一条对角线。设共享边为 (u,v)，两侧三角形第三个顶点分别是 p、q，翻转后新边是 (p,q)。翻边时要同步更新四个邻居的邻接索引，漏一个就会在后续 walk 时走进死循环。

void flipAndPush(const Edge& e, std::stack<Edge>& stk) { int t1 = e.tri, t2 = tris[t1].n[edgeIndex(t1, e.u, e.v)]; int p = thirdVertex(t1, e.u, e.v); int q = thirdVertex(t2, e.u, e.v); // 重连：t1 变为 (p, q, u)，t2 变为 (q, p, v) rebuildTriangle(t1, p, q, e.u); rebuildTriangle(t2, q, p, e.v); // 四条外围边重新入栈检查 stk.push({p, q, t1}); stk.push({q, u, t1}); stk.push({p, q, t2}); stk.push({q, v, t2}); }

edgeIndex根据顶点索引反查边在三角形中的位置，thirdVertex返回除 u、v 外的那个顶点。这两个辅助函数看着简单，但边界条件多，建议单独写单元测试。

3.3 用随机点集验证最小闭环

跑通的最小验证方式是：生成 100 个随机点，插入后遍历所有三角形，检查每个三角形的外接圆内是否还有别的点。如果全部通过，说明翻边逻辑正确。

bool validate() { for (const auto& t : tris) { if (!t.alive) continue; Point a = points[t.v[0]], b = points[t.v[1]], c = points[t.v[2]]; for (const auto& p : points) { if (p.id == t.v[0] || p.id == t.v[1] || p.id == t.v[2]) continue; if (incircle(a, b, c, p) > 1e-9) return false; } } return true; }

这个 O(n²) 的校验只适合小规模调试，但它是你写完翻边后最该先跑的东西。我见过太多人直接上大规模点集，结果三角形重叠、空洞，回头查半天，不如先用 100 个点把逻辑钉死。

4. 避坑与排查：增量 Delaunay 最容易翻车的五个地方

4.1 浮点误差导致翻边死循环

现象：插入某个点后程序卡死，栈里反复压入同一条边。原因：incircle 判断在临界情况下符号抖动，翻转后又判定需要翻回。解决：给 incircle 加一个容差阈值，绝对值小于 1e-9 时视为共圆，不翻转；同时给翻边循环加最大迭代次数，超限就打印当前边和四个顶点坐标，人工确认。

4.2 共线点导致初始三角形退化

现象：前三个点共线，初始三角形面积为 0，后续所有 orient 判断失效。原因：增量法默认初始三角形非退化。解决：插入前先扫描点集，找到第一个与已选点不共线的点再构造初始三角形；如果整个点集共线，直接返回空结果并提示。

4.3 超级三角形太小导致外围点插入失败

现象：插入靠近包围盒边缘的点时，新三角形跑到超级三角形外面，walk 找不到目标。原因：超级三角形边距不够。解决：边距取包围盒对角线长度的 10 到 20 倍，或者干脆用三个远端点构造。代价是数值范围变大，incircle 的浮点误差也会放大，所以边距不是越大越好，够用即可。

4.4 邻接索引更新遗漏导致 walk 走进已删除三角形

现象：locate 时访问到 alive 为 false 的三角形，或者邻居索引指向自己。原因：翻边时只更新了两个主三角形的邻接，忘了更新外围邻居的反向指针。解决：每次翻边后，对四个外围邻居做一次反向索引校验，确保它们的 n 数组里指向的是新三角形而不是旧索引。这个校验可以只在 debug 模式开，release 关掉。