CGAL 2D约束Delaunay三角化remove操作抛出异常原因排查

这是《CGAL - 2D Delaunay triangulation - remove is not removing》的续篇。为解决点无法删除的问题，我尝试通过cdt.locate()定位目标点所在面，匹配到对应顶点后，插入该顶点的精确坐标，再调用cdt.remove()删除刚插入的点，但程序抛出异常终止。请问异常产生的原因是什么？

相关代码 #

//---TYPEDEFS

typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
typedef CGAL::Constrained_triangulation_face_base_2<K> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb,Fb> TDS;
typedef CGAL::Exact_predicates_tag Itag;
typedef CGAL::Constrained_Delaunay_triangulation_2<K,TDS,Itag> CDT;

typedef CDT::Point Point;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Vertex_circulator Vertex_circulator;
typedef CDT::Face_handle Face_handle;

...

int RemovePontosDesabilitados(CDT& cdt, tysetPonto3D& SetPonDesabilitados)
{
  for (auto& RPonto : SetPonDesabilitados)
  {
    Point PonRemover(RPonto.x, RPonto.y);

    auto Face(cdt.locate(PonRemover));

    if (Face != NULL)
    {
       int C(0);
       bool Achou(false);

       while (C < 3 && !Achou)
       {
         if((fabs(Face->vertex(C)->point().x()) - RPonto.x) < 1e-5 &&
            (fabs(Face->vertex(C)->point().y()) - RPonto.y) < 1e-5) 
         {
           Achou = true;
         }
         else
         {
           C++;
         }
       }

       if (Achou)
       {
         Point PonRemover(Face->vertex(C)->point().x(), Face->vertex(C)->point().y());

         auto Tam = cdt.number_of_vertices();
         auto PonIns = cdt.insert(PonRemover);
         Tam = cdt.number_of_vertices();
         cdt.remove(PonIns);
         Tam = cdt.number_of_vertices();
       }
    }
  }
  return 1;
}

异常原因分析 #

• 重复插入触发约束顶点删除限制
  CGAL的Constrained_Delaunay_triangulation_2中，插入已存在的点时，insert()不会创建新顶点，而是返回指向已有顶点的句柄。此时调用remove()删除该句柄，本质是尝试删除原本就存在的顶点——如果这个顶点是约束边的端点，CGAL不允许直接删除这类顶点，强行操作就会抛出异常。你的额外插入操作完全冗余，反而触发了这个限制。

• 逻辑绕路导致错误
  你已经通过locate()找到目标顶点的句柄（Face->vertex(C)），完全不需要再插入一次。直接使用这个已有句柄调用remove()才是正确路径。

• 浮点数匹配的潜在隐患
  手动用fabs判断点匹配容易受精度误差影响，建议直接使用CGAL内核提供的精确比较：Face->vertex(C)->point() == PonRemover，避免漏判或误判。

修正后的代码示例 #

int RemovePontosDesabilitados(CDT& cdt, tysetPonto3D& SetPonDesabilitados)
{
  for (auto& RPonto : SetPonDesabilitados)
  {
    Point PonRemover(RPonto.x, RPonto.y);

    auto Face(cdt.locate(PonRemover));

    if (Face != nullptr)
    {
       int C(0);
       bool Achou(false);

       while (C < 3 && !Achou)
       {
         if(Face->vertex(C)->point() == PonRemover) 
         {
           Achou = true;
         }
         else
         {
           C++;
         }
       }

       if (Achou)
       {
         Vertex_handle target_vertex = Face->vertex(C);
         // 先检查顶点是否关联约束边
         bool has_constraint = false;
         Vertex_circulator vc = cdt.incident_vertices(target_vertex);
         if (vc != nullptr) {
           do {
             if (cdt.is_constrained(target_vertex, vc)) {
               has_constraint = true;
               break;
             }
           } while (++vc != cdt.incident_vertices(target_vertex));
         }

         if (!has_constraint) {
           cdt.remove(target_vertex);
         } else {
           // 先删除关联的约束边，再删除顶点
           std::vector<std::pair<Point, Point>> constraints_to_remove;
           vc = cdt.incident_vertices(target_vertex);
           do {
             if (cdt.is_constrained(target_vertex, vc)) {
               constraints_to_remove.emplace_back(target_vertex->point(), vc->point());
             }
           } while (++vc != cdt.incident_vertices(target_vertex));

           for (auto& cons : constraints_to_remove) {
             cdt.remove_constraint(cons.first, cons.second);
           }
           cdt.remove(target_vertex);
         }
       }
    }
  }
  return 1;
}

内容的提问来源于stack exchange，提问作者Rabelo