2024国际网格生成圆桌会议短课程摘要：网格生成的微分几何理论基础

停靠在巴尔的摩内湾的帆船 BAP Union号。

2024年3月5号这一周，美国工业与应用数学协会（SIAM）主办的国际网格圆桌会议（International Meshing Roundtable）在巴尔的摩召开。去年同样的会议在荷兰阿姆斯特丹举行。这一年来，国际工业软件行业发生了两次巨变。2024年一月中旬，芯片设计软件EDA龙头之一的Synopsis斥资350亿美元收购了仿真软件巨头Ansys；三月五号这一天，EDA另一龙头Cadence宣布斥资12.4亿美元收购了仿真领域的新秀BETA CAE。BETA CAE起源于希腊亚里士多德大学，由四个师生团队于1999年创立，其网格生成软件ANSA具有很多独特的技术特色。不同于其他有限元前处理软件，ANSA在处理几何模型和有限元网格方面具有很好的速度与质量，适应性和可定制性强，模型规模没有软件限制。这两次天价收购震动了国际工业软件行业，证明了这一领域的勃勃生机！参会者也都神采飞扬地热烈讨论。恰逢3月5号这一天，笔者在大会上给了两个小时的报告，题目是“微分几何应用于网格生成”。在报告中，笔者展示了与美国计算力学之父Tom Hughes合作，为BETA CAE完成的一个项目：世界上首次用等几何分析方法进行汽车整车白体（white body）的碰撞模拟实验。众所周知，汽车碰撞模拟是数值模拟之王，等几何分析方法用几何建模的样条函数空间作为物理场的函数空间。等几何分析与传统有限元法相比有很多优势，例如消除了几何近似误差；在同样精度下减少了面片的数目，从而提高了计算效率；在同样面片数目的情况下提高了计算精度，等等。但是传统的CAD模型包含大量的裁剪样条，从而缺乏水密性，无法直接用于等几何分析。因此，我们需要将原始CAD模型重新计算结构化四边形网格剖分，然后建立全局几何光滑、具有水密性的样条曲面。而曲面整体四边形网格理论，正是笔者演讲的重点之一。听众们表现出了巨大的热情，提出了很多有见地的问题。

图1. 汽车白体碰撞数值模拟实验。（与Tom Hughes, Kendrick Shepherd, BETA-CAE合作）。

CAE的有限元方法发源于1930年代，开始也是工程师根据直觉经验建立的经验性方法，后来职业数学家经过长期努力，建立了有限元方法的理论基础，例如纽约的克朗，北京的冯康等等。其理论来源于经典的偏微分方程理论，例如索伯列夫空间理论，里兹法，伽辽金法等等。几何建模方法起源于1970年代，主要是基于样条理论（开花理论Blossom），和逼近理论。工业软件涵盖了汽车、船舶、航空航天等广阔领域，发展迅猛。在历史发展中，很多基本问题并没有完全彻底解决，而是被广大工程师们用各种经验性实用技巧搪塞过去。例如，经典的NURBS样条曲面本质上是建立在仿射几何不变量上的，但是封闭曲面本身往往并不具备整体的共形结构，因此NURBS曲面不可避免地具有所谓的奇异点，在奇异点领域，保证曲面几何正则性的构造异常复杂。在CAD工业初期，人们很少考虑为拓扑复杂的曲面整体建模，因此并没有深入考虑仿射几何的整体存在性问题，这为后继发展带来巨大困难。一般曲面上都具有实或复射影结构，基于射影不变量构造的样条曲面将不会有奇异点的困扰。但是整个工业生产线都已经固定下来，更为科学的几何概念无法推翻现有的工业标准。再如，在计算流体力学应用中，或者大形变模拟问题中，曲面和体的整体结构化网格具有强烈的优点。但是整体生成结构化网格本身具有很强的挑战性，目前工业界的通用手法是从曲面边缘生成结构化四边形网格，然后逐次向内部推进，直至无法推进一步，最后总会遗漏一些单元，迫不得已用三角形面片填补，即生成所谓的混合网格。其实这些基本问题的理论在微分几何领域也是近期才建立起来，解决这些问题的算法也是近期才被发明出来。笔者团队希望用现代拓扑、微分几何、复几何理论，来重新审视古老的网格生成领域，夯实理论基础；应用近期发展出来的几何分析算法，来解决这些古老问题，以切实推进应用的发展。

图2. 基于Ricci流的汽车表面非结构化网格生成。

在数值仿真的算法流程中，真正的瓶颈在于网格生成，其花费时间占整个流程的75%以上。在有限元方法中，最为常用的是线性有限元，即将曲面或者将实体进行三角剖分。实体三角剖分的关键在于曲面三角剖分。用伽辽金法求解椭圆型偏微分方程的一个弱解，弱解到真解的L2距离正比于三角剖分边长上界的平方，系数上界等于边长除以内切圆（内切球）的半径。如果三角剖分中有些内角趋于零，则系数上界趋于无穷。因此曲面三角网格生成算法需要保证最小角最大化，所有内角有严格的下界；同时网格顶点密度可以动态调节，以适应局部的几何或者物理场。但是，在传统网格生成领域，一直没有理论严密的算法，保证生成满足上面要求的曲面三角剖分。

在经典计算几何领域，欧式空间区域的三角剖分已经被充分研究，几乎所有的算法都是基于Delauay三角剖分及其对偶Voronoi图。给定平面样本点集，一个三角剖分是Delaunay的，如果每个三角形的外接圆内部不包含任意样本点。例如经典的Chew算法在平面凸区域中构建Delaunay三角剖分，如果三角形外接圆半径过大，则插入三角形的外心，更新三角剖分至Delaunay。这种Delaunay Refinement算法可以保证最小角不小于30度。Rupert Delaunay refinement算法可以处理更为复杂的平面区域，并且保持区域中的事先指定的直线段，最小角不小于20.7度。这些Delaunay Refinement 算法可以被推广到球面和双曲情形，并且保证内角下有界。

图3. 闭黎曼面单值化定理。

在20世纪初期，微分几何（黎曼面）理论中的中心定理被证明：所谓的单值化定理。这一定理断言，任何带有黎曼度量的曲面，无论拓扑或者几何多么复杂，都可以保角变换（共形变换）成常值高斯曲率的曲面，即三种标准曲面中的一种：单位圆球，平直轮胎或者双曲曲面。我们在标准常曲率曲面上计算三角剖分，应用Delaunay Refinement算法，可以保证最小角有界。然后在将三角剖分拉回到原来曲面上面。所谓保角映射，就是说局部上，映射是相似变换，但是整体上相似比处处不同。相似变换将无穷小圆映成无穷小圆，因此保持Delaunay三角剖分。保角变换保持曲线交角不变，因此最小角有界。这样，我们从理论上证明了高质量曲面三角剖分的存在性。

图4. 带边界黎曼面单值化定理。

那么，如何将抽象理论变成实用算法呢？这里的关键在于计算单值化定理。虽然单值化定理在20世纪初被证明，其实用算法却是一百年之后才被笔者团队发明出来。这是因为单值化定理的计算归结为求解Yamabe方程，而Yamabe方程强烈非线性，传统的有限元方法无能为力。给定曲面带有黎曼度量，共形变换其度量意味着找到一个标量函数，新的度量为，新的高斯曲率为

这个方程被称为是Yamabe方程，未知函数出现在指数上，因此传统数值PDE方法无法求解。在1980年代，Hamilton教授发明了Ricci流方法，用于证明庞加莱猜想，最后Perelman用Ricci流完成了这一壮举。曲面上的Ricci流方程是：

笔者在2002年与丘成桐先生访问加州大学洛杉矶分校，师兄刘克峰教授在家中为丘先生举行生日宴会。笔者遇到很多师兄，交谈中提到了如何建立离散Ricci流理论的想法。笔者团队在2010年左右开始发表离散Ricci流算法，但是直到2018年才在微分几何杂志上发表离散Ricci流解的存在性、唯一性、收敛性的理论结果。迄今为止，离散曲面Ricci流算法是唯一的算法能够根据用户指定的高斯曲率构造黎曼度量。

由此可见，1910年代的单值化定理，1930年代的有限元法，1970年代的平面Delaunay三角剖分算法，1980年代的Ricci流理论，2010年代的离散Ricci流算法，最终从理论和算法层面完美地回答了曲面三角网格生成问题。

图5. 基于Ricci流的汽车空气动力学数值仿真。

当然，在工程实践层面，依然存在大量的挑战。例如大幅度面积畸变诱导的数值精度问题，非线性优化带来的计算复杂度问题，复杂整体曲面自动分片问题等等。这些工程问题需要硬件和软件的同步提高，和工程师的聪明才智才能逐步解决。

网格生成领域的另外一个问题是如何构造光滑标架场（或者cross field），即在曲面的每一个点，选取其切平面上的一个单位正交标架，当曲面上的点光滑变动的时候，标架也光滑变动。如果我们只考虑标架的x-轴，那么我们就得到一个光滑单位切向量场。固定曲面上一个点，所有的单位切向量构成一个圆周。曲面上所有的单位切向量构成一个三维流形：曲面的单位切丛。曲面光滑单位向量场（标架场）是单位切丛的全局截面。由纤维丛的示性类理论，单位切丛中全局光滑截面并不存在，其拓扑障碍为切向量场中的奇异点。围绕一个奇异点，单位切向量旋转的圈数定义为奇异点的指标。单位切向量场所有奇异点的指标和等于曲面的欧拉示性数，这些奇异点的位置并不重要，而总指标为常值，可以被视为曲面上的二维上同调类，即纤维丛的陈省身示性类。

图6. 曲面上的标架场。

反之，给定曲面上的离散点，和相应的指标，如果总指标等于欧拉示性数，那么由纤维丛拓扑障碍类理论，必然存在光滑标架场，以指定点为奇异点。但是如何从算法层面构造这样的标架场，却是非常具有挑战性的问题。直观上讲，我们可以用Ricci流构造一个新的度量，使得所有的曲率都集中到奇异点上。例如指标为+1的奇异点对应的曲率为90度，指标为-1的奇异点对应的曲率为-90度，等等。用Ricci流得到一个带有锥奇异点的平直度量，曲面去掉奇异点得到一个带刺孔的曲面。我们沿着任意一条环路在平直度量上平行移动一个标架，回到起点处，标架旋转的角度被称为是环路的holonomy。因为黎曼度量平直，因此holonomy只依赖于环路的通伦类。围绕每个刺孔，环路的holonomy是90度的整数倍，cross标架和自身重合。但是，考虑原来曲面的通伦群基底，它们的holonomy并不是90的整数倍。我们需要寻找某种补偿机制。

图6. 黎曼面上调和微分群的基底。

我们考虑曲面的de Rham上同调群，经典的Hodge理论断言，每个微分形式都有唯一的分解方式：

这里是恰当微分形式，是余恰当形式，是调和形式。所有调和形式构成的群与曲面de Rham上同调群同构。Hodge分解等价于在曲面上求解Poisson方程。在2002年左右，丘成桐先生指导笔者发明了Hodge分解算法，本质上是将有限元方法推广到流形上的外微分。如此，我们可以构造一个调和形式，其在每个通伦群基底上的积分等于这个基底的holonomy。这个调和微分形式给出了补偿机制。

我们固定曲面上一个基点，任选此点切平面上的一个标准正交标架F。对于曲面上任意一点p，选取一条路径连接基点和点p，沿着此路径平行移动标架F，将平移后的标架顺时针旋转，其转角等于调和形式在路径上的积分。如此所有环路的holonomy都等于90的整数倍，我们得到全局光滑的cross场。

由此可见，整体标架场构造的理论基础是纤维丛的示性类理论和Hodge理论，算法依赖于离散Ricci流和Hodge分解，而Hodge分解算法本质上是有限元方法在流形上外微分的直接推广。

在CAD/CAE领域，曲面结构化四边形网格生成具有根本的重要性。但是长期以来，全局结构化网格生成一直基于各种经验性方法，例如局部构造，全局推广的方法，即波前推进法（advancing wavefront），这种方法无法保证全局结构化，最后总会用一些三角形面片填补空隙。人们也发明了各种基于混合网格的数值模拟算法。但是自从等几何分析方法横空出世，人们迫切需要全局结构化四边形网格的算法。目前大多数方法需要用户手工输入奇异点，或者手动将复杂曲面分割成规则面片，再生成全局四边形网格。如何自动生成高质量整体结构化网格，成为是工业软件的核心难点。观察上图，我们看到四边形网格的奇异点构型并非随机，而是具有强烈的内在规律，如何寻找结构化网格奇异点构型的控制方程（govering equation）成为整个领域发展的瓶颈。

图7. 曲面上的结构化网格。

让我们回到最简单情形：考察轮胎的四边形网格剖分。我们可以轻易构造一个四边形网格，没有任何奇异点（即拓扑度不等于4的顶点）；但是我们无法构造一个四边形网格，只有两个奇异点，一个拓扑度为3，另一个拓扑度为5。假设这样的网格存在，那么它符合欧拉示性数的拓扑要求，也符号高斯-博内定理的黎曼几何要求。那么，为什么轮胎曲面不存在如此3-5四边形网格剖分？

我们需要考察其他几何结构，即复几何结构。笔者团队近期证明：四边形网格诱导曲面的一个共形结构和一个亚纯四次微分。具有共形结构的曲面是黎曼面，黎曼面的余切丛具有复结构，被称为全纯余切丛。全纯余切丛的四次幂是这一问题的背景几何。任何一个亚纯四次微分都是这一全纯线从的整体亚纯截面，其零点极点构型就是这一复纤维丛的示性类，满足特定的Abel方程，即亚纯四次微分诱导的除子的Jacobi像是Jacobi簇中一个特定点:

这里是一个亚纯四次微分，是一个全纯微分，是Abel-Jacobi映射。

图8. 黎曼面上的全纯微分群基底。

例如上图，我们计算亏格为2黎曼曲面上的全纯微分，其零点分布并非任意而是严格满足Abel方程。

由此可见，结构化曲面网格生成理论来自于代数几何与复几何。从算法角度而言，我们需要Ricci流，Hodge分解，和求解Abel-Jacobi方程。

CAD/CAE是一个巨大的工业，也是一个资本雄厚、利润丰厚的工业。数十年的发展和数百万工程师的努力，CAD/CAE已经具有成熟的工业软件体系。但是很多基本的问题被历史所遗留，一直没有得到根本解决，而是通过各种经验性方法和人工方法加以妥协变通。现代微分几何和拓扑有助于这些基本问题的解决。虽然这一领域历史悠久，资源丰厚，但是它并不傲慢，依然拥抱各种新思想、新理论，因此它能够一直在加速发展。

目前AI的发展颠覆了很多传统领域，ChatGPT终结了自然语言处理领域（NLP），深度学习的发展颠覆了计算机视觉、图像处理等领域，Sora的横空出世将会为计算机图形学带来巨变。但是CAD/CAE领域需要精确的数值计算结果，目前AI的概率统计方法依然无法取得精确的数值结果。CAD/CAE领域需要很多几何拓扑、物理偏微分方程的背景知识，近期无法被AI大模型所取代。我们期待现代几何的理论和方法更加深入地进入这一领域，真正为网格生成奠定理论基础，促进下一次产业革命！