CGAL笔记之网格生成——2D三角剖分和网格

最新推荐文章于 2024-01-23 10:24:02 发布

3333yyt

最新推荐文章于 2024-01-23 10:24:02 发布

阅读量1.7k

点赞数 2

分类专栏： CGAL 文章标签： c++ 算法图形渲染数据结构

本文链接：https://blog.csdn.net/younggung/article/details/129232094

版权

CGAL 专栏收录该内容

13 篇文章 13 订阅

订阅专栏

CGAL笔记之网格生成——2D三角剖分和网格

1 三角剖分和网格
2 网格

1 三角剖分和网格

1.1 定义

如果三角剖分的任何面的外接圆在其内部不包含顶点，则三角剖分是Delaunay 三角剖分。约束Delaunay 三角剖分是尽可能多的 Delaunay 的约束三角剖分。约束 Delaunay 三角剖分的任何小平面的外接圆在其内部不包含从该小平面可见的数据点。

如果一条边内切在一个空圆中（其内部不包含数据点），则称该边为Delaunay 边。如果该边的直径圆为空，则称该边为加布里埃尔边。

如果每条约束边都是 Delaunay 边，则称约束 Delaunay 三角剖分是一致的Delaunay 三角剖分。由于约束 Delaunay 三角剖分中的任何边要么是 Delaunay 边，要么是约束边，因此一致的 Delaunay 三角剖分实际上是 Delaunay 三角剖分。唯一的区别是一些边被标记为约束边。

如果每个约束边都是 Gabriel 边，则称约束 Delaunay 三角剖分是*一致的 Gabriel 三角剖分。*Gabriel 性质比 Delaunay 性质强，每条 Gabriel 边都是 Delaunay 边。因此，一致的 Gabriel 三角剖分也是一致的 Delaunay 三角剖分。

任何受约束的 Delaunay 三角剖分都可以细化为一致的 Delaunay 三角剖分或一致的 Gabriel 三角剖分，方法是在受约束的边上添加称为Steiner vertices 的顶点，直到它们被分解为小到足以成为 Delaunay 或 Gabriel 边的子约束。

1.2 建立一致的三角剖分

受约束的 Delaunay 三角剖分可以通过以下两个全局函数细化为一致的三角剖分：

template<class CDT>
void make_conforming_Delaunay_2 (CDT& t)

template<class CDT>
void make_conforming_Gabriel_2 (CDT& t)

在这两种情况下，模板参数CDT必须由受约束的 Delaunay 三角剖分类实例化（请参阅第2D 三角剖分）。

用于实例化参数的约束 Delaunay 三角剖分的几何特征CDT必须是概念的模型ConformingDelaunayTriangulationTraits_2。

受约束的 Delaunay 三角剖分t通过引用传递，并通过添加顶点细化为一致的 Delaunay 三角剖分或一致的 Gabriel 三角剖分。如果必须保留原始三角剖分以用于其他计算，建议用户复制输入三角剖分

make_conforming_Delaunay_2()和构建内部数据结构所使用的算法make_conforming_Gabriel_2()，如果在同一个三角剖分上连续调用这两个函数，则该内部数据结构将被计算两次。为了避免这些数据被构造两次，高级用户可以使用该类Triangulation_conformer_2<CDT>将受约束的 Delaunay 三角剖分细化为一致的 Delaunay 三角剖分，然后再细化为一致的 Gabriel 三角剖分。为了进一步控制细化算法，此类还提供了一次插入一个斯坦纳点的单独函数。

1.3 示例：进行符合 Delaunay 然后符合 Gabriel 的三角剖分

此示例将几个线段插入到受约束的 Delaunay 三角剖分中，使其符合 Delaunay，然后符合 Gabriel。在每一步，打印三角剖分的顶点数。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Triangulation_conformer_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Constrained_Delaunay_triangulation_2<K> CDT;
typedef CDT::Point Point;
typedef CDT::Vertex_handle Vertex_handle;
int main()
{
  CDT cdt;
  // construct a constrained triangulation
  Vertex_handle
    va = cdt.insert(Point( 5., 5.)),
    vb = cdt.insert(Point(-5., 5.)),
    vc = cdt.insert(Point( 4., 3.)),
    vd = cdt.insert(Point( 5.,-5.)),
    ve = cdt.insert(Point( 6., 6.)),
    vf = cdt.insert(Point(-6., 6.)),
    vg = cdt.insert(Point(-6.,-6.)),
    vh = cdt.insert(Point( 6.,-6.));
  cdt.insert_constraint(va,vb);
  cdt.insert_constraint(vb,vc);
  cdt.insert_constraint(vc,vd);
  cdt.insert_constraint(vd,va);
  cdt.insert_constraint(ve,vf);
  cdt.insert_constraint(vf,vg);
  cdt.insert_constraint(vg,vh);
  cdt.insert_constraint(vh,ve);
  std::cout << "Number of vertices before: "
            << cdt.number_of_vertices() << std::endl;
  // make it conforming Delaunay
  CGAL::make_conforming_Delaunay_2(cdt);
  std::cout << "Number of vertices after make_conforming_Delaunay_2: "
            << cdt.number_of_vertices() << std::endl;
  // then make it conforming Gabriel
  CGAL::make_conforming_Gabriel_2(cdt);
  std::cout << "Number of vertices after make_conforming_Gabriel_2: "
            << cdt.number_of_vertices() << std::endl;
}

在这里插入图片描述

从左到右：初始 Delaunay 三角剖分、相应的一致 Delaunay 三角剖分和相应的 Gabriel 三角剖分。

2 网格

2.1 定义

网格是将给定区域划分为形状和大小满足多个条件的单纯形。

域是用户想要划分网格的区域。它必须是平面的有界区域。该域由平面直线图定义，简称Pslg ，它是一组线段，该组中的两个线段要么不相交，要么共享一个端点。Pslg的段是将由网格中的边并集表示的约束。Pslg还可以包含将显示为网格顶点的孤立点。

Pslg的段是边界或内部约束的段。Pslg的段必须覆盖域的边界。

Pslg将平面划分为几个连接的组件。默认情况下，域是有界连通分量的并集。

下图未使用种子点定义的域示例及其可能的网格。
在这里插入图片描述

用户可以通过提供一组种子点来覆盖此默认值。种子点要么标记要划分网格的组件，要么标记不划分网格的（孔）组件。

下图使用相同Pslg定义的另一个域和用于定义孔的两个种子点。在相应的网格中，这两个孔被三角化但没有网格化。

在这里插入图片描述

2.2 形状和尺寸标准

三角形的形状标准是下界乙关于外接半径与最短边长的比值。这样的界限意味着下界反正弦1个2乙关于三角形的最小角和上界π− 2 ∗反正弦1个2乙在最大角度上。不幸的是，只有在以下情况下才能保证算法的终止≥ _2个–√，对应于20.7角度上的度数。

大小标准可以是倾向于选择小三角形的任何标准。例如，尺寸准则可以是三角形最长边长度的上限，或者外接圆半径的上限。大小界限可以在整个域中变化。例如，大小标准可以对与给定线相交的三角形施加较小的大小。

criteria两种类型的标准都在作为网格函数参数传递的对象中定义。

2.3 网格算法

网格划分问题的输入是一个Pslg和一组描述要划分网格的域的种子，以及一组大小和形状标准。在此包中实现的算法从输入Pslg的约束 Delaunay 三角剖分开始，并使用 Delaunay 细化方法生成网格。此方法将新顶点插入到三角剖分中，尽可能远离其他顶点，并在满足条件时停止。

如果输入Pslg的入射段之间的所有角度都大于60度数，如果外接半径/边缘比的界限大于2个–√，算法保证以满足尺寸和形状标准的网格终止。

如果一些输入角度小于60度，该算法最终将得到一个网格，其中一些三角形违反了小输入角度附近的标准。这是不可避免的，因为无法抑制由输入段形成的小角度。此外，一些具有小输入角度的域不能用甚至小于小输入角度的角度进行网格划分。请注意，如果域是多边形区域，则生成的网格将满足尺寸和形状标准，但输入角度较小。此外，该算法可能会成功生成角度下界大于20.7度，但没有这样的保证。

2.4 构建网格

通过调用全局函数从受约束的 Delaunay 三角剖分中获得网格：

template<class CDT, class Criteria>
void refine_Delaunay_mesh_2 (CDT &t, const Criteria& criteria)

模板参数CDT必须由受约束的 Delaunay 三角剖分类实例化。为了覆盖域，此函数的一个版本有两个定义种子点序列的参数。

的几何特征类CDT必须是概念的模型DelaunayMeshTraits_2。ConformingDelaunayTriangulationTraits_2这个概念改进了添加几何谓词和构造函数的概念。模板参数Criteria必须是的模型MeshingCriteria_2。这个概念定义了三角形必须满足的标准。CGAL 为这个概念提供了两个模型：

Delaunay_mesh_criteria_2<CDT>，它定义了一个形状标准，它限制了三角形的最小角度，
Delaunay_mesh_size_criteria_2<CDT>，这为先前的标准添加了最大边长的界限。

refine_Delaunay_mesh_2()如果在同一个三角剖分上使用不同的标准多次调用该函数，则该算法会在每次调用时重建用于网格划分的内部数据结构。为了避免每次调用都重建数据结构，高级用户可以使用类Delaunay_mesher_2<CDT>。此类还提供分步功能。这些函数一次插入一个顶点。

任何类型的对象Delaunay_mesher_2<CDT>都是从对 a 的引用构造的CDT，并且有几个成员函数来定义要划分网格的域和划分CDT. 有关详细信息，请参见下面给出的示例和参考手册。请注意，CDT不应在对象的生命周期内对其进行外部修改Delaunay_mesher_2<CDT>。

构建网格后，可以使用is_in_domain()面类型的成员函数确定三角剖分的哪些面在网格域中。

2.5 使用 Lloyd 优化网格

该包还提供了一个全局函数，可以在 Delaunay 细化生成的网格上运行 Lloyd 优化迭代。此网格优化的目标是改善网格内部的角度，并使它们尽可能接近 60 度。

template< class CDT >
Mesh_optimization_return_code lloyd_optimize_mesh_2(CDT& cdt);

此优化过程交替将顶点重新定位到其 Voronoi 单元的质心，并更新三角剖分的 Delaunay 连通性。质心是根据大小函数计算的，该函数旨在保留 Delaunay 细化生成的网格中点的局部密度。

在这里插入图片描述

图中refine_Delaunay_mesh_2() （左）为统一尺寸标准生成的网格。（右）显示了经过 100 次 Lloyd 优化迭代后的相同网格。

下图显示了这些网格内部角度refine_Delaunay_mesh_2()、lloyd_optimize_mesh_2()的直方图。

在这里插入图片描述
Delaunay 细化后以及 Lloyd 优化 10 次和 100 次迭代后网格内角度的直方图。经过 Delaunay 细化后，角度在区间 [28.5; 121.9]度。经过 10 次 Lloyd 优化迭代后，它们在 [29.1; 110.8]。100 次迭代将它们带到 [29.3; 109.9]。

2.6 示例

2.6.1 使用全局函数的示例

以下示例将几个线段插入到约束三角剖分中，然后使用全局函数对其进行网格剖分refine_Delaunay_mesh_2()。大小和形状标准是标准类提供的默认标准Delaunay_mesh_criteria_2<K>。没有给出种子，这意味着网格域覆盖了除无界分量之外的整个平面。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Delaunay_mesher_2.h>
#include <CGAL/Delaunay_mesh_face_base_2.h>
#include <CGAL/Delaunay_mesh_size_criteria_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
typedef CGAL::Delaunay_mesh_face_base_2<K> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> Tds;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, Tds> CDT;
typedef CGAL::Delaunay_mesh_size_criteria_2<CDT> Criteria;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Point Point;
int main()
{
  CDT cdt;
  Vertex_handle va = cdt.insert(Point(-4,0));
  Vertex_handle vb = cdt.insert(Point(0,-1));
  Vertex_handle vc = cdt.insert(Point(4,0));
  Vertex_handle vd = cdt.insert(Point(0,1));
  cdt.insert(Point(2, 0.6));
  cdt.insert_constraint(va, vb);
  cdt.insert_constraint(vb, vc);
  cdt.insert_constraint(vc, vd);
  cdt.insert_constraint(vd, va);
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Meshing the triangulation..." << std::endl;
  CGAL::refine_Delaunay_mesh_2(cdt, Criteria(0.125, 0.5));
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
}

2.6.2 使用类 Delaunay_mesher_2<CDT> 的示例

此示例使用该类Delaunay_mesher_2<CDT>并两次调用refine_mesh()成员函数，更改其间的大小和形状标准。在这种情况下，使用两倍的全局函数refine_Delaunay_mesh_2()会降低效率，因为算法所需的一些内部结构将被构建两次。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Delaunay_mesher_2.h>
#include <CGAL/Delaunay_mesh_face_base_2.h>
#include <CGAL/Delaunay_mesh_size_criteria_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
typedef CGAL::Delaunay_mesh_face_base_2<K> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> Tds;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, Tds> CDT;
typedef CGAL::Delaunay_mesh_size_criteria_2<CDT> Criteria;
typedef CGAL::Delaunay_mesher_2<CDT, Criteria> Mesher;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Point Point;
int main()
{
  CDT cdt;
  Vertex_handle va = cdt.insert(Point(-4,0));
  Vertex_handle vb = cdt.insert(Point(0,-1));
  Vertex_handle vc = cdt.insert(Point(4,0));
  Vertex_handle vd = cdt.insert(Point(0,1));
  cdt.insert(Point(2, 0.6));
  cdt.insert_constraint(va, vb);
  cdt.insert_constraint(vb, vc);
  cdt.insert_constraint(vc, vd);
  cdt.insert_constraint(vd, va);
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Meshing the triangulation with default criterias..."
            << std::endl;
  Mesher mesher(cdt);
  mesher.refine_mesh();
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Meshing with new criterias..." << std::endl;
  // 0.125 is the default shape bound. It corresponds to abound 20.6 degree.
  // 0.5 is the upper bound on the length of the longuest edge.
  // See reference manual for Delaunay_mesh_size_traits_2<K>.
  mesher.set_criteria(Criteria(0.125, 0.5));
  mesher.refine_mesh();
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
}

2.6.3 使用种子的例子

此示例使用全局函数refine_Delaunay_mesh_2()，但使用一个种子定义域。大小和形状标准是标准类提供的默认标准Delaunay_mesh_criteria_2<K>。

网格构建完成后，is_in_domain()使用面的成员函数对其进行计数。

#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Delaunay_mesher_2.h>
#include <CGAL/Delaunay_mesh_face_base_2.h>
#include <CGAL/Delaunay_mesh_size_criteria_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Triangulation_vertex_base_2<K> Vb;
typedef CGAL::Delaunay_mesh_face_base_2<K> Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb> Tds;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, Tds> CDT;
typedef CGAL::Delaunay_mesh_size_criteria_2<CDT> Criteria;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Point Point;
int main()
{
  CDT cdt;
  Vertex_handle va = cdt.insert(Point(2,0));
  Vertex_handle vb = cdt.insert(Point(0,2));
  Vertex_handle vc = cdt.insert(Point(-2,0));
  Vertex_handle vd = cdt.insert(Point(0,-2));
  cdt.insert_constraint(va, vb);
  cdt.insert_constraint(vb, vc);
  cdt.insert_constraint(vc, vd);
  cdt.insert_constraint(vd, va);
  va = cdt.insert(Point(3,3));
  vb = cdt.insert(Point(-3,3));
  vc = cdt.insert(Point(-3,-3));
  vd = cdt.insert(Point(3,0-3));
  cdt.insert_constraint(va, vb);
  cdt.insert_constraint(vb, vc);
  cdt.insert_constraint(vc, vd);
  cdt.insert_constraint(vd, va);
  std::list<Point> list_of_seeds;
  list_of_seeds.push_back(Point(0, 0));
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Meshing the domain..." << std::endl;
  CGAL::refine_Delaunay_mesh_2(cdt, list_of_seeds.begin(), list_of_seeds.end(),
                               Criteria());
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Number of finite faces: " << cdt.number_of_faces() << std::endl;
  int mesh_faces_counter = 0;
  for(CDT::Finite_faces_iterator fit = cdt.finite_faces_begin();
      fit != cdt.finite_faces_end(); ++fit)
  {
    if(fit->is_in_domain()) ++mesh_faces_counter;
  }
  std::cout << "Number of faces in the mesh domain: " << mesh_faces_counter << std::endl;
}

2.6.4 使用 Lloyd 优化器的示例

此示例使用全局函数lloyd_optimize_mesh_2()。网格是使用的函数生成refine_Delaunay_mesh_2()的CGAL::Delaunay_mesher_2，然后使用进行优化lloyd_optimize_mesh_2()。优化将在max_iteration_number交替顶点重定位和 Delaunay 连接更新的 10 次（由设置）迭代后停止。可以使用更多终止条件，并在参考手册中进行了详细说明。

#define CGAL_MESH_2_OPTIMIZER_VERBOSE
//#define CGAL_MESH_2_OPTIMIZERS_DEBUG
//#define CGAL_MESH_2_SIZING_FIELD_USE_BARYCENTRIC_COORDINATES
#include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
#include <CGAL/Constrained_Delaunay_triangulation_2.h>
#include <CGAL/Delaunay_mesher_2.h>
#include <CGAL/Delaunay_mesh_face_base_2.h>
#include <CGAL/Delaunay_mesh_vertex_base_2.h>
#include <CGAL/Delaunay_mesh_size_criteria_2.h>
#include <CGAL/lloyd_optimize_mesh_2.h>
#include <iostream>
typedef CGAL::Exact_predicates_inexact_constructions_kernel K;
typedef CGAL::Delaunay_mesh_vertex_base_2<K>                Vb;
typedef CGAL::Delaunay_mesh_face_base_2<K>                  Fb;
typedef CGAL::Triangulation_data_structure_2<Vb, Fb>        Tds;
typedef CGAL::Constrained_Delaunay_triangulation_2<K, Tds>  CDT;
typedef CGAL::Delaunay_mesh_size_criteria_2<CDT>            Criteria;
typedef CGAL::Delaunay_mesher_2<CDT, Criteria>              Mesher;
typedef CDT::Vertex_handle Vertex_handle;
typedef CDT::Point Point;
int main()
{
  CDT cdt;
  Vertex_handle va = cdt.insert(Point(-2,0));
  Vertex_handle vb = cdt.insert(Point(0,-2));
  Vertex_handle vc = cdt.insert(Point(2,0));
  Vertex_handle vd = cdt.insert(Point(0,1));
  cdt.insert(Point(2, 0.6));
  cdt.insert_constraint(va, vb);
  cdt.insert_constraint(vb, vc);
  cdt.insert_constraint(vc, vd);
  cdt.insert_constraint(vd, va);
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Meshing..." << std::endl;
  Mesher mesher(cdt);
  mesher.set_criteria(Criteria(0.125, 0.05));
  mesher.refine_mesh();
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
  std::cout << "Run Lloyd optimization...";
  CGAL::lloyd_optimize_mesh_2(cdt,
    CGAL::parameters::max_iteration_number = 10);
  std::cout << " done." << std::endl;
  std::cout << "Number of vertices: " << cdt.number_of_vertices() << std::endl;
}