三角自适应网格 python

【Literature Notes】三角网格曲面的参数化

• 1.三角网格曲面参数化简介
• 2.边界点的映射
• 3.内部点的映射
• 4.算法结果验证
• 参考文献

---

1.三角网格曲面参数化简介

三角网格参数化就是对于给定的一个三角网格曲面和一个参数域，寻求一个从参数域上的点到三角网格上点的一一对应映射，在保持参数域上的拓扑信息与原始网格相同的同时谋求某种几何度量的变形最小化$^{[1]}$。因此需要做的工作一方面是寻找理论上的映射方法，另一方面是建立变形度量的方法以寻找最优映射。
三角网格参数化方法有以下几种分类：

1. 根据参数域的不同可以分为平面参数化和球面参数化等;
2. 根据网格的拓扑信息可以分为带边界网格参数化和封闭网格参数化，甚至是任意拓扑的三角网格参数化;
3. 考察不同的内在几何变量的变形,可以分为保面积参数化、保角参数化和等距参数化等;
4. 根据求解方式的不同可以分为线性方法和非线性方法;
5. 根据最后得到的参数化结果，可以分为局部参数化和整体参数化方法等$^{[1]}$。

网格曲面参数化的两大处理思路可以简单分为固定边界法与活动边界法$^{[4]}$。
典型的网格参数化算法有：

1. Floater保形参数化方法
2. Harmonic参数化方法
3. Intrinsic参数化方法

下文中将对Harmonic参数化方法（调和映射）进行展开。该方法将三角网格模型中的每一条边都看作为一个弹簧，先将三角网格的边界映射到预先定义好的多边形上，然后通过最小化整个三角网格模型的弹性势能来参数化内部空间点$^{[1]}$。调和映射属于固定边界参数化，相较于活动边界法（Free Boundary），其计算效率更高，但参数化效果不若后者。

先尝试采用固定边界法(单位正方形域或者圆域)，固定边界法（Fixed Boundary）：顾名思义，使用固定边界法对网格进行参数化的时候，首先将网格的边界顶点的参数值固定在一个凸多边形上（比如矩形或圆形），然后列出最小变形能量方程，对内部顶点的参数值进行求解。

2.边界点的映射

网格的边界点的参数值需要固定带参数域中预先定义好的多边形上，这个平面参数域多边形可取单位正方形域和单位圆域（为了便于生成网格曲面的等参数线刀轨）。
Step1：设网格曲面$M$总共有$n$个顶点，其中有$k$个边界顶点$v_i (i=1,…,k)$ ，找出网格的边界顶点。网格曲面边界顶点$v_i$ 的弦长参数化为：
$t_{i}=\frac{\sum_{j=1}^{i-1}\left|\boldsymbol{v}_{j+1}-\boldsymbol{v}_{j}\right|}{\sum_{j=1}^{n-1}\left|\boldsymbol{v}_{j+1}-\boldsymbol{v}_{j}\right|}\tag{1}$
Step2：将顶点的弦长参数化结果按比例映射到参数域中。
单位正方形域，边界顶点的参数化结果为：
$\left(x_{i}, y_{i}\right)= \begin{cases}\left(0.5-4 t_{i}, 0.5\right) & , t_{i} \in[0,0.25) \\ \left(-0.5,1.5-4 t_{i}\right) & , t_{i} \in[0.25,0.5) \\ \left(4 t_{i}-2.5,-0.5\right) & , t_{i} \in[0.5,0.75) \\ \left(0.5,4 t_{i}-3.5\right) & , t_{i} \in[0.75,1]\end{cases}\tag{2}$
单位圆域，则边界起始点和边界点$v_i$之间的圆弧所对应的圆心角为$2t_i\pi$ ，则边界顶点参数化的结果为：
$\left(x_{i}, y_{i}\right)=\left(\cos 2 t_{i} \pi, \sin 2 t_{i} \pi\right)\tag{3}$

3.内部点的映射

在调和映射的假设前提下，可以列出映射的变形能为(类似弹性势能)：
$E_{\text {harm }}(\mathrm{u})=\frac{1}{2} \sum_{(i, j) \in E} k_{i, j}\left\|\mathrm{u}_{i}-\mathrm{u}_{j}\right\|^{2}\tag{4}$
$k_{i,j}$为弹性系数,E为三角网格的边集。
其中$u_i,v_i$分别为边$(i,j)$的两个顶点映射到参数域中的位置，边界点的位置已经确定。
$k_{i, j}=\left(L_{i, k_{1}}^{2}+L_{j, k_{1}}^{2}-L_{i, j}^{2}\right) / A_{i, j, k_{1}}+\left(L_{i, k_{2}}^{2}+L_{j, k_{2}}^{2}-L_{i, j}^{2}\right) / A_{i, j, k_{2}}\tag{5}$
其中: $k_{1}, k_{2}$ 为边 $(i, j)$ 的两个相邻面的另外两个顶点, $L_{i, j}$ 表示边 $(i, j)$ 的 长度, $A_{i, j, k}$ 表示三角片 $f(i, j, k)$ 的面积。
采用拉格朗日乘子法解式(4)使得$Eharm(u)$最小，得到网格曲面中所有顶点映射到参数域中后的坐标位置。
以边界点为已知条件，内部点为未知量，将上式对未知量求导可得
$\left\{\begin{array}{l} \sum_{v_{j} \in \operatorname{adj}\left(v_{k+1}\right)} s_{k+1, j}\left(x_{k+1}-x_{j}\right)^{2}=0 \\ \sum_{v_{j} \in \operatorname{adj}\left(v_{k+1}\right)} s_{k+1, j}\left(y_{k+1}-y_{j}\right)^{2}=0 \\ \mathrm{~L} \mathrm{~L} \\ \sum_{v_{j} \in \operatorname{adj}\left(v_{i}\right)} s_{i, j}\left(x_{i}-x_{j}\right)^{2}=0 \\ \sum_{v_{j} \in \operatorname{adj}\left(v_{i}\right)} s_{i, j}\left(y_{i}-y_{j}\right)^{2}=0 \\ \mathrm{~L} \mathrm{~L} \\ \sum_{v_{j} \in \operatorname{adj}\left(v_{n}\right)} s_{n, j}\left(x_{n}-x_{j}\right)^{2}=0 \\ \sum_{v_{j} \in \operatorname{adj}\left(v_{n}\right)} s_{n, j}\left(y_{n}-y_{j}\right)^{2}=0 \end{array}\right.\tag{6}$
$adj(v_i)$为$v_i$相邻顶点，上式需要求解大型稀疏矩阵,运用超松弛迭代法可求出$x_i,y_i$，从而将网格的内部顶点全部映射到平面多边形的内部$^{[2]}$。

4.算法结果验证

网格模型：

单位正方形域参数映射：

单位圆域参数映射：

在求解式（6）中的大型稀疏矩阵时，为了矩阵建立方便，实际算法实践中采用如下等价化简进行求解$^{[1]}$。

参考文献

[1] 晏冬梅. 三角网格曲面参数化技术研究[D].南京航空航天大学,2005.
[2] 陈晓兵. 口腔修复体高效数控加工编程技术研究与实现[D].南京航空航天大学,2011.
[3] Sun Yuwen et al. Iso-parametric tool path generation from triangular meshes for free-form surface machining[J]. The International Journal of Advanced Manufacturing Technology, 2006, 28(7-8) : 721-726.