面向脏几何的曲面网格自动生成算法

孟祥桥，刘智伟，刘陶然，陈建军

(浙江大学 航空航天学院，杭州 310027)

0 引 言

由于建模误差、格式转换等一系列问题，模型往往存在法向错误、穿插、狭缝、贴合等一系列问题，我们将这类无法直接应用于数值计算的几何模型称为脏几何[1]。因此，在传统的网格生成流程中，网格生成前往往需要对几何进行清理操作，用户需要识别所有的几何错误，并依次修复[2]。这一修复过程往往需要大量人工交互，并借助复杂的图形界面，如CADfix[3]等专业的几何处理软件，且十分依赖于操作人员的经验。对于大量部件组成的复杂模型，人工修复的成本是难以忍受的。

Shephard 等提出了一种基于八叉树的全自动网格生成方法[4]，基于八叉树的网格生成方法可以一定程度上容忍几何错误，但仍存在一些固有的缺陷。模型的几何特征并非初始构造，该方法需要向基网格嵌入特征，使网格逼近真实几何[5]，基网格的质量严重影响几何特征嵌入的效果。在嵌入过程中，往往需要引入大量的单元来捕捉细小的几何特征，导致单元数目的大幅上升。

刘等人基于混合曲面表征[2,6-7]实现了装配体的“曲面嵌入”，可以解决多个部件之间可能存在的相交和重叠问题。其方法是在离散曲面上计算边界相交图，得到离散嵌入结果，之后根据连续-离散的映射关系完成连续曲面边界表征上的曲面自动嵌入[3]。该算法具有较高的时间效率，也保证了计算结果的几何精度，但其适用范围有限，鲁棒性不足，针对模型贯穿的问题，由于贯穿处缺少边界表征，曲面嵌入失效。

Hu 等通过将三角形-三角形求交升维至四面体-三角形求交，在初始三角形表面约束下生成纯四面体网格，并在有理数范围内实现了精确求交[8]。这种方法在一定程度上可以解决脏几何模型中出现的相交问题，但其局限性亦存在。一是相交计算十分耗时，二是此类方法仅能设定全局容差值，当缝隙间距大于其允许的容差值时，使用此方法不仅无法消除几何误差，反而可能急剧影响性能和网格质量，而容差值过大时，几何特征不能被完全保留。刘[9]也采用四面体化的方法来完成脏几何的网格生成，并创新性地引入了边界恢复算法[10]，从而大幅减少相交计算量。两者的四面体化方法中均要求输入边界法向正确，对于存在法向错误的几何，仍然需要大量的手工交互调整输入模型的法向。

本文针对脏几何中可能出现的法向错误、穿插、狭缝、贴合问题，通过四面体化方法自动修复脏几何，并生成高质量的曲面网格。其核心思想是一个有效的四面体网格的表面一定是一个有效的曲面网格。本文在输入模型的离散表征的基础上直接生成四面体网格，通过离散边界切割四面体解决了模型中存在的穿插问题，基于四面体优化解决狭缝和贴合问题，此时该四面体网格的表面便是有效的曲面网格，提取该表面并应用曲面网格重构提升表面网格质量，得到最终的曲面网格。

相较于八叉树方法，本方法具有更高的几何精度，可以生成贴体的曲面网格；相较于基于混合曲面表征的曲面嵌入方法，本方法具有更强的普适性和健壮性；相较于Hu 等的四面体化方法，本方法的主要优点有：

1）应用基于视图采样的方法完成法向修复，允许输入的脏几何中存在法向错误；

2）实现了基于容差场的局部修复与特征保持，避免四面体优化导致的特征丢失或修复不完全；

3）实现了基于多级Hausdorff 距离[11]的网格重构技术，可以在修复的曲面网格的基础上进行高精度的网格重构，避免网格重构导致的表面自相交。

1 算法流程

图1 和图2 给出了本文所提出的基于四面体化方法的面向脏几何曲面网格生成算法流程，该算法的输入为脏几何模型，输出为保形的高质量曲面网格。该算法主要包含以下步骤：

图1 算法流程简图Fig. 1 Diagram for the algorithm illustration

图2 算法流程图Fig. 2 Flow chat of the algorithm

1）构建混合曲面表征。建立模型的离散-连续表征数据结构[2]，采用B-rep 表示模型的连续表征，采用三角网格表示模型的离散表征，同时维护连续模型和离散模型之间的映射关系。得到离散表征后，我们采用基于视图采样的方法，完成法向调整。

2）初始四面体化。本文采用Bowyer-Watson 方法[6,12]完成四面体化，我们提取离散表征中的所有点，作为B-W 算法的输入，构造一个Delaunay 三角化。单次B-W 增量插点的过程为：（1）在当前Delaunay 三角化中找到包含待插入点的基单元；（2）利用单元的相邻关系在基单元附近快速搜索所有打破外接圆准则的单元；（3）删除所有（2）中的单元，形成空腔，连接待插入点与空腔边界顶点，形成新的Delaunay 三角化。

3）边界恢复与相交计算。初始四面体化中并不保证包含所有的原始曲面边界，不能包含的边界单元可分为两类。第一类是由于Delaunay 三角化只能保证点存在于三角化中，不能保证网格边和网格面也都存在于三角化中，这类边界可以通过边界恢复[9-10,13-14]的策略来找回原始边界，第二类是由于脏几何模型本身存在相交或重叠的情形，得到的三角化必然无法包含此类边界，如图1（c）所示蓝色区域，这类边界需要通过相交计算对四面体进行裁剪以找回边界。

4）四面体网格质量优化。相交计算仅能解决部件相交的问题，但模型中仍可能存在狭缝等不必要的特征。我们采用四面体网格质量优化的方法，在容差范围内抹除此类特征。

5）曲面网格重构。提取上述四面体网格的表面得到修复后的表面网格，修复后的表面网格解决了原始模型的一系列错误，但网格质量与网格密度均难以达到数值计算的要求，因此需要对4）中得到的表面网格进行重构。

2 算法实现

2.1 法向调整

本文算法的核心思想是生成有效的四面体网格，提取其表面得到最终有效的曲面网格。如图1 所示，我们在输入模型的外部加入包围盒，在整个包围盒区域内生成四面体网格，我们计算每个四面体单元相对输入表面的绕卷数[15]，筛选得到内部的四面体单元，从而方便地提取出有效的表面。绕卷数计算正确的前提是输入表面的法向是正确的，此外，在曲面网格重构中，也需要保持每个面上单元的法向一致。因此在完成混合曲面表征的构造后，我们需要对离散表征的面法向进行修复。

若离散表征中所有的三角单元都是相连的，即两个三角单元之间存在一条共享边，我们可以通过BFS 算法[16]完成面法向调整，即锁定某一个三角单元的法向，以该单元为起点，不断搜索其邻居，使与其相连的单元的法向调整为与该起点单元同向。但考虑最坏情况下，若几乎所有的面片均不存在与其他面片的拓扑连接，此时我们无法通过BFS 算法快速地完成法向调整。为了尽可能地解决脏几何中可能出现的法向问题，我们采用基于视图采样的法向调整策略。通过绘制该离散模型在各个方向上的视图，可以计算该模型中任意一个面片法向的可能性。理想情况下，我们在任意一个观测点v观察该模型，对于观察到的某一个三角形中的某一点p，视线方向 (v-p)与该三角形的法向n应 满足 (v-p)·n≤0。

定义三角面片ti的 法向正确的概率为P(ti)。我们对该模型建立一个包围球，从包围球上均匀选取N个观测点对该模型进行绘制，得到N张1 024×1 024的图像。定义O+(ti) 为三角面片ti在所有视图中法向与视线夹角大于90°的像素点数，O-(ti) 为三角面片ti在所有视图中法向与视线夹角小于90°的像素点数。最终可计算得到三角面片ti法向正确的概率为：

若P(ti)＞0.5，将保持该面片的法向，否则将该三角面片的法向翻转。在实际应用中，首先采用BFS 算法尽可能地将面片组合，然后采用上述视图采样法向调整策略完成法向的全局调整。对于单张视图，我们采用扫描线Z 缓冲器算法[17]，完成视图的快速绘制。图3 为法向调整前后的对比，其中绿色表示该三角面片的法向指向体的内部。

图3 法向调整效果Fig. 3 Meshes before and after the normal adjustment

2.2 曲面求交

为了恢复模型的所有原始边界，我们使用原始边界的三角面片所在的平面来切割四面体。非退化情况下四面体单元被原始表面单元所在平面切割如图4 所示，被切割的四面体单元的侧面将会被截断为一个四边形和一个三角形，该四面体单元将形成一个四面体单元和一个五面体单元。首先对侧面的三个四边形分别进行Delaunay 三角化，再对该五面体单元进行分割，在其重心处插入新点，与侧面的新三角形顶点连接后形成分裂后的新四面体单元。贴合问题将产生零体积的四面体单元，本算法中的三角形求交计算全部采用基于GMP 的高精度有理数计算[18]，退化情况下舍弃零体积的新单元即可。基于精确的求交计算，我们可以保证输入模型中的所有自相交问题可以得到修复。

图4 四面体单元切割的非退化情形Fig. 4 Non-degenerate case of tetrahedral cutting

图5 展示了利用原始边界切割四面体解决相交问题，图5（a）中红色四面体与蓝色的原始边界出现相交，图5（b）为该四面体的底面。利用原始边界切割四面体，其结果如图5（c）和图5（d）所示。

图5 曲面求交过程Fig. 5 Surface intersection process

2.3 四面体网格优化

上述四面体-平面求交可以在有理数的范围内修复脏几何的自相交问题，但模型仍可能存在狭缝等不必要的细小特征。由于前述步骤已经完成了输入模型原始表面的恢复，那么在狭缝处必然存在狭窄甚至是退化的四面体单元，因此优化此处的四面体单元质量，即可抹除此类特征。图6 展示了通过四面体优化解决狭缝的详细过程。图示两圆柱之间存在狭缝需要被抹除，经过初始四面体化后，如图6（b）所示，狭缝处产生了图示红色的四面体单元，这类四面体单元的质量较差，经过四面体质量提升后，这类单元被抹去，得到的表面如图6（c）所示。

图6 四面体网格优化解决狭缝问题Fig. 6 Tetrahedral mesh optimization for slits

2.3.1 四面体质量提升

四面体单元质量的提升主要操作为边分裂、边折叠、面交换与点优化四种局部操作[19]，其示意图如图7 所示。本文通过循环执行这四种局部操作来逐步优化四面体单元。

图7 四面体单元质量提升Fig. 7 Quality improvement of tetrahedral elements

本文采用Rabinovich 等[8,20]提出的能量函数作为四面体单元的质量标准，该函数具有天然的各向同性，可以有效地甄别出各向异性的和负体积的四面体单元，公式如下：

其中，Jt是 将单元t转化为正单元（三维时为正四面体）的雅可比矩阵，D为网格维度，衡量四面体单元质量时取3。当某一单元的能量函数值大于临界值时，将被视为差单元，并通过上述质量提升操作提升单元质量。

四面体网格质量优化的目的是解决狭缝或重叠处产生的狭窄的或退化的四面体单元，对于空间中的四面体，其质量评判标准不同。为了加速优化效率，我们采用了多级能量阈值判别方案，即每个四面体单元优化的目标能量值不同，这里我们采用四面体单元至边界的距离来计算该单元的能量阈值。我们令原始表面附近的单元的目标能量值为Ets，该三角化中允许的最大能量值为kEts，则空间中某一四面体的目标优化质量Et可以通过线性插值得到，即：

其中，ds为该四面体与输入表面的最短距离。靠近表面的四面体单元，其目标能量值较小，而远离所有表面的四面体单元，其目标能量值较大，我们允许在距离表面远处产生较差的四面体单元。

2.3.2 包络检测与容差场

上述四面体网格质量优化方法中，边折叠与点优化涉及网格点的位移操作，若网格点的偏移量过大，可能导致网格与输入模型之间产生较大的偏差，我们将禁止这样的边折叠与点优化操作。为了确保本算法生成的表面网格与输入模型之间的误差在可控范围内，每次边折叠与点优化的操作执行后，需要进行包络检测以确保生成保形的表面网格，包络检测的方法如下。

图8 包络检测的采样点Fig. 8 Sampling points for the envelope detection

容差 ε的设置将会影响包络检测的时间效率与精确性，同时也会影响模型缺陷的修复效果。若要修复模型中的狭缝，那么此处的容差值应当大于狭缝的宽度，但其它区域应当尽可能地保证边界不发生改变。因此，我们定义了基于离散表征的容差场，任意一点处的容差分为两部分，全局容差值 ε0与局部容差值εp，该点处的容差值取两者中较小的一个。某一部件的 ε0为该部件的包围盒对角线长度的千分之一，局部容差值 εp需根据模型需求，在需要修复的狭缝区域设置较大的容差值，在必须严格保持的特征处设置较小的容差值。如图9 所示，我们通过设置楔形体脊边处的容差值为一个较小的容差值，实现特征保持的效果。

图9 局部容差实现特征保持Fig. 9 Feature preserving based on local tolerance

曲面求交与四面体质量提升两个环节中均可能产生新的边界点，新的边界点都是由边分裂形成的，故新边界点的容差值应取其边分裂前的两端点处容差值的最小值。

2.4 曲面网格重构

本文中采用经典的网格重构策略完成曲面网格重构，其核心算法与四面体网格重构相同，对需要优化的单元依次执行边分裂、边折叠、边交换与点优化操作，直至所有的表面网格均满足质量要求。

为了生成保形的表面网格，我们希望在重构过程中由于边分裂产生的新点、点优化后的点都能尽可能地贴近原始表面。网格重构中的误差估计通过计算Hausdorff 距离[11-21]来衡量，我们要求所有的网格点距离输入边界的Hausdorff 距离均小于hmax。由于表面网格优化中存在点优化的操作，在狭窄的区域有可能出现网格自相交的情况。为了避免这种情况，我们采用多级误差来避免由表面网格重构引起的表面自相交。在执行完毕一次曲面网格重构后，检测所有的表面网格是否存在自相交的情况，对于存在自相交的表面，令h′max=hmax/2，重新进行曲面网格重构，直至表面不存在自相交的区域。

3 数值验证

本算法采用C++语言实现，借助GMP 库完成基于有理数的精确计算。测试平台为Intel Core i7-6700，3.4 GHz，内存24 GB。

3.1 网格生成效果

本文选取了靶球、格栅鳍、F-16 战斗机、AIM-54 导弹、Ford 发动机和某压气机六个典型的脏几何模型作为实验对象，模型网格如图10 所示。

图10 典型脏几何模型Fig. 10 Typical models for dirty geometries

为验证本文算法的有效性，这里详细展示三个复杂脏几何的网格生成结果，分别为：靶球模型、格栅鳍模型和F-16 战斗机模型。

图11 展示了靶球模型的几何及其生成结果，该模型包含115 个体部件，存在15 177 处相交，图11（a）为曲面离散后的局部放大结果，可以明显地看出该处存在部件之间的相交问题，经过本文的四面体化方法修复的表面网格如图11（b）所示，部件之间的相交问题已经得到处理。

图11 靶球模型及其网格生成结果Fig. 11 Target ball model and its mesh generation

图12 展示了格栅鳍模型的几何及其网格生成结果，该模型的各片格栅鳍之间存在大量的贯穿问题，如图12（a）所示，对于此类模型，曲面嵌入的算法由于贯穿处缺少对应的边界表征，故无法支持部件之间互相贯穿的情形。而本文算法由于不依赖原始的边界表征，可以很好地解决这一问题，网格生成结果如图12（b）所示。

图12 格栅鳍模型及其网格生成结果Fig. 12 Grid fin model and its mesh generation

图13 展示了F-16 战斗机模型的几何及其网格生成结果，该模型由8 个部件组装而成，部件之间存在相交与狭缝等问题。如图13（a）展示了该模型存在的部件之间的相交问题，粉色尾鳍与淡蓝色机身之间存在部件相交问题，图13（b）为该处网格的生成结果。图13（c）展示了该模型部件之间存在的狭缝问题，深蓝色为该模型机翼，淡蓝色区域为该模型机身，机翼与机身之间存在狭缝，图13（d）为该处网格生成结果，可以看出狭缝已经被抹除。

图13 F-16 模型及其网格生成结果Fig. 13 F-16 aircraft model and its mesh generation

F-16 战斗机模型中存在一些特殊的几何特征在网格生成中需要得到保留，如图14 所示机翼后缘，我们希望在此处能够保留直线特征，后缘处应用2.3.2 节所述的局部容差可以实现该处的特征保留，在机翼后缘处设置较小的容差值，实现了该处特征保持的效果，TetWild[8]在此处生成的网格如图14（a）所示，无法保留此类特征。

图14 F-16 战斗机局部特征保持对比Fig. 14 Comparison of local feature preservation for the F-16 fighter

此外，由于本文算法的网格重构过程中引入了尺寸场[21]，可以根据曲率特征实现网格的局部加密，即在曲率较大的区域自动生成较密的网格，图14（b）中的机翼前后缘区域的网格受曲率特征的控制，网格尺寸较小，这一特征也使得该网格更能满足后续数值分析的需要。

AIM-54 导弹、Ford 发动机和某压气机的部分模型错误与修复结果如图15 所示。图15（a）中AIM-54 的机翼与机身之间，图15（c）中Ford 发动机红圈处存在贴合错误，其修复结果如图15（b）和15（d）所示。图15（e）中压气机存在狭缝，修复结果如图15（f）所示。

图15 AIM-54、Ford 发动机、压气机模型与网格Fig. 15 Models and meshes of AIM-54, a Ford engine and a compressor

3.2 对比与数据统计

常用的网格生成算法均不具备表面修复的能力，Hu 等完成的TetWild[8]旨在任意离散表面上生成四面体网格，具备一定的表面修复能力，3.1 节中的结果表明，本文算法具有更好的修复效果，并且由于本文加入了曲面网格重构的策略，可以在保形的基础上尽可能地生成高质量的表面网格。

我们将前述6 个典型的脏几何模型作为实验对象，与TetWild[8]进行网格质量与时间效率对比。由于TetWild 只能接受离散输入，且无法处理法向错误的情况，我们将混合曲面表征中已经完成法向调整的离散表征结果作为TetWild 的输入。

我们选用平均最小角度和平均最大边长比作为网格质量对比项，平均最小角度为该网格所有单元中每个单元的最小角度的平均值，平均最大边长比为该网格所有单元中每个单元的最大边长比的平均值。前述6 个模型分别采用TetWild 与本文算法处理所得网格的质量统计与对比结果如表1 所示。对比表明，本文算法生成的网格的平均最小角度与平均最大边长比均优于TetWild。

表1 网格质量对比Table 1 Comparison of mesh quality

前述6 个模型分别采用TetWild 与本文算法处理所耗时间的数据统计与对比结果如表2 所示。对比表明，TetWild 所耗时间是本文算法的3 倍及以上，本文算法在时间效率上领先TetWild。

表2 时间效率对比Table 2 Comparison of time efficiency

4 结 论

本文针对脏几何存在的法向错误、穿插、狭缝、贴合问题，应用基于四面体化的曲面网格修复方法，开发了适用于脏几何的自动网格生成算法。本算法实现了以往四面体化方法无法完成的法向修复；应用边界恢复算法，可以大幅减少参与求交的三角形数量，提高求交性能；应用分级的误差检测，可以大幅提高包络检测精度，从而保证生成保形的曲面网格；引入容差场，实现了特征保持与特征消除；应用多级误差估计，避免曲面网格重构产生自相交，最终生成了高质量的表面网格。上述对比实验也展示了本算法在修复效果、时间效率与网格质量上的优势。

本文算法需要手动设置容差，从而达到好的修复效果，未来将持续探索更加简便与自动化的容差设置方式。此外，本文算法针对存在破洞的几何尚无法完成修复，需要进一步研究补面的算法。