动画原理：表面形变算法的思考与总结

前言：

之前我的文章 Mesh形变算法_mesh算法-CSDN博客就有大致的讨论过，介绍的也比较粗略！现在主要是想在Triangulated Surface Mesh Deformation方向上更深入的讨论一下！结合今年我对这一块的学习谈谈我的理解~

下面要介绍大致几种算法，不管是拉普拉斯形变算法或者它的形变还是ARAP的形变，本质上都是，通过在网格人为的定义一些控制点或控制边，然后根据控制点或控制边的移动来影响网格的形状。相应的优势很明显可以对网格的细节进行精细的控制，比如arap是一种刚性形变，它可以在形变工程中尽可能地保持细节！效果虽然好但是带来的也是需要考虑网格的拓扑结构和形变的质量，以及控制点或控制边等等的分布和数量的设计！在此之前推荐一本计算几何的经典书籍《Polygon Mesh Processing》：
里面的第七章就有很详细的讲述各种常见经典的形变算法，非常值得参考。

前置知识- 矩阵的迹：

矩阵的迹（Trace）是一个线性代数中的概念，它指的是矩阵主对角线上所有元素的总和。对于一个方形矩阵 A，其迹被表示为 Tr(A) 或者 trace(A)。
主对角线是指矩阵从左上角到右下角的一条线上的元素，这些元素的行索引和列索引相同。
例如，给定一个 3x3 的矩阵：

A = | a11 a12 a13 |
    | a21 a22 a23 |
    | a31 a32 a33 |

矩阵 A 的迹就是：Tr(A) = a11 + a22 + a33即，只需将主对角线上的元素加起来就可以了。
假设我们有如下的 2x2 矩阵：

B = | 4 9 |
    | 2 6 |

矩阵 B 的迹计算如下：

Tr(B) = 4 + 6 = 10

很简单对吧（~~）？直接将主对角线上的 4 和 6 相加得到迹为 10。
矩阵的迹性质：1、矩阵迹具有线性属性ext{tr}(A + B) = ext{tr}(A) + ext{tr}(B)；2、循环置换（Cyclic permutation） 对于多个矩阵的乘积，迹操作满足循环置换性，可以改变乘法中的项的顺序而不影响迹： ext{tr}(ABC) = ext{tr}(BCA) = ext{tr}(CAB)
特别注意并不意味着矩阵乘法是交换的，仅适用于迹。3、 特征值之和（Sum of eigenvalues）: 矩阵的迹等于其所有特征值之和。如果 A 的特征值是 λ1, λ2, …, λn，则： ext{tr}(A) = sum_{i=1}^{n} λi

计算机图形学、几何形变中矩阵的迹大致有以下作用：

• 物理模拟：在表征应力和应变等物理量时，矩阵的迹可以用来衡量表面或体积的变形程度。
• 网格形变：在网格形变（例如 ARAP 形变）中，局部旋转的最佳拟合通常通过 SVD 来寻找，而迹被用来评价变换的刚性程度，确保局部结构尽可能保持刚性。
• 几何处理：迹用来计算几何对象的局部描述符，例如用于估计曲面的均值曲率或高斯曲率。

前置知识-基于曲面的变形(Surface-Based Deformations)：

根据《Polygon Mesh Processing》第七章的形变定义，如下图对一个正方形的曲面S进行变形，固定曲面S浅蓝色的部分F，然后选取黄色部分H的顶点作为控制顶，将其向上拉动。可以看到经过变形后，没有被固定的部分(R，即深蓝色部分)的顶点的位置发生了相应的变换。 根本的问题就是如何选取合适的位移函数d，使得变形的结果符合需求。

sbd：位移函数d是从曲面S到三位空间的映射，即计算时在三角形网格上进行的。这类方法能够有很高的控制度，能够对每一个顶点都加以控制，但是其鲁棒性不好，运行效率往往取决网格的的复杂度和质量。
只需设置个别离散点的新位置来表达他所想要的形变，计算出剩余离散点应有的位置，改变一个模型的姿势或者说大体形态（全局变形），而不是模型的表面细节和组成部分（局部不变）。

前置知识- 奇异值SVD的分解：

SVD（Singular Value Decomposition）将任意一个复数或实数矩阵分解为三个特定的矩阵乘积的形式。假设我们有一个矩阵A m×n。SVD就是找到三个矩阵U、Σ和V^T，使得它们相乘的结果等于A。形式上来说就是：A = UΣV^T 这里的三个矩阵具有以下特点：

1. U：左奇异向量矩阵矩阵大小为m×m，它是正交矩阵，即U的列向量互相正交，并且都是单位向量（长度为1）。
2. Σ（西格玛）：m×n的矩形对角矩阵。它的对角线上的元素是奇异值，并且奇异值按照从大到小的顺序排列。非对角线上的元素都是0。在实际应用中，奇异值越大表示那个方向上数据的变化越大。
3. V^T：右奇异向量矩阵的转置矩阵，大小为n×n的。原来的V是正交矩阵， 它就是其转置形式。和U类似，V中的列向量也是互相正交的，并且是单位向量。
   SVD分解在信号处理、统计学、计算机图形学以及机器学习中的降维压缩数据（如PCA主成分分析）等场合都会用到SVD。

举一个压缩ML数据的例子：假设A是一个人脸图片的矩阵，SVD能够识别哪些特征（比如脸部轮廓、眼睛位置等）是构成这些人脸最重要的。

这是因为较大的奇异值所对应的方向上，数据的变化量大，蕴含的信息丰富。通过保留这些主要特征，甚至可以将数据进行压缩，即用更少的信息来近似描述原始数据。

SVD性质：

1. 存在与唯一性：对于任何一个实数或复数矩阵A（不必是方阵），都存在一个奇异值分解，而且奇异值在Σ矩阵中是唯一确定的，即如果奇异值都是不相等的，那么这些奇异值是唯一的。然而，U和V矩阵的列向量可以有符号上的改变（即方向的镜像），因为奇异向量是根据A的A^TA 和 A^TA和AAT 计算得到的，并且仅确定到正负号。
2. 能量分布：A的Frobenius范数（即矩阵元素平方和的平方根）等于其奇异值的平方和如下：如果A的SVD为A = UΣV^T，并且Σ的对角线元素（即奇异值）为σ_1, σ_2, ..., σ_p，其中p是A的秩，则Frobenius范数满足
   ||A||_F = sqrt(σ_1^2 + σ_2^2 + ... + σ_p^2)。当你截取前k个最大奇异值及其对应的奇异向量时，这就给出了A的最佳k秩近似。
3. 数据压缩与去噪：由于SVD可以提供关于数据能量分布的信息，保留最大的几个奇异值及其对应的奇异向量，可以用于图像压缩和去噪，同时损失的信息相对较少。
4. 线性方程组的求解：SVD可用来计算矩阵的伪逆，进而找到线性方程组的最小二乘解。
5. 矩阵的秩：矩阵A的秩等于非零奇异值的个数，这个性质可以帮助我们判断矩阵是否满秩或找出其秩的大小。

SVD 3*3的例子：
奇异值分解（SVD）
奇异值分解(SVD)原理详解及推导

前置知识-形变能量

第一次接触形变的时候天天提到了能量一开始云里雾里，看的很懵懂很懵逼后来理解了所以我特意在这里强调！计算几何的形变能量是指数学上定义的一个函数，它衡量形变后的模型和原始模型之间的差异程度。能量通常用于指导网格形变的过程控制形变行为，并确保形变结果满足某些期望的性质，比如局部细节特征的保持、体积的不变性等。形变能量通常被用作最优化问题中的目标函数，使得通过找到最小化该能量函数的解（最小二乘法等）来获得最终的形变结果。
比如将要讲的能量：

1. 平滑能量：要求形变后的网格保持平滑，常用于避免过度扭曲或皱折。（Mesh平滑处理的几种算法比较）
2. ARAP 能量：As-Rigid-As-Possible 形变希望每个局部区域都能尽可能地保持原始的形状和大小，这就涉及了计算每个小区域的局部旋转，并试图最小化它们相对于全局变形的偏差。
3. 拉普拉斯形变能量：使用拉普拉斯坐标定义的能量，它试图保持原始网格和形变网格顶点的拉普拉斯坐标一致，以此保持网格的局部属性。

Laplacian Deformation形变及其改进：

简单来说拉普拉斯形变是基于网格顶点的拉普拉斯坐标，它记录了顶点与其邻接顶点之间相对位置关系的系统。拉普拉斯坐标的定义通常基于网格的拉普拉斯算子，个描述网格局部几何特性的微分算子。对于三角网格中的每个顶点，其拉普拉斯坐标是该顶点与邻接顶点位置差的加权平均。
拉普拉斯形变以前大致讲过（Mesh形变算法） 或者看大佬博客图形处理（三）简单拉普拉斯网格变形-Siggraph 2004以及图形处理（五）基于旋转不变量的网格变形-Siggraph 2007
如下图拉普拉斯形变在旋转过程中的缺点：

1. 体积损失：在大范围的旋转下，传统的拉普拉斯形变可能导致模型体积的损失，产生塌陷的现象，特别是在模型的细长部分。这是因为拉普拉斯算子是一个线性算子，它不直接编码旋转信息。
2. 形状扭曲：由于传统拉普拉斯形变保持局部差分坐标，这些坐标并不完全与旋转不变性兼容。因此，当网格发生较大的旋转时，原先的局部结构可能会扭曲。
3. 细节损失：拉普拉斯形变是非刚性形变！在旋转之后，如果相同的拉普拉斯坐标被用于重建形状，那么模型的一些细节特征可能会丢失，导致形状的全局或局部特征改变。
   

基于旋转不变性的拉普拉斯形变是拉普拉斯形变的一个改进。拉普拉斯坐标被分解成一个方向向量（长度）和一个大小（权重）。然后这个方向向量会被调整，使其适应网格的旋转变化。这个改进的方法的关键是，它增加了对网格局部旋转的考虑。不仅保持了拉普拉斯坐标的大小不变，还尽量保持了方向向量在局部旋转下的不变性。这种方法可以使变形的结果看起来更加自然，特别是当网格经历复杂的旋转运动时。这些方法的共同点在于它们都利用了拉普拉斯坐标来保持网格的某些局部属性，基于旋转不变性的拉普拉斯形变由于考虑了旋转，往往在处理需要保留体积和细节的场景时更为优越。

ARAP形变：

As-Rigid-As-Possible跟拉普拉斯形变也是基于表面的变形，arap如何保证刚体：

希望尽量使所有的边都只经历刚体变换。 已知三角形具有稳定性，即三条边长度确定时形状确定。 尽量去保持每条边的长度也就是在尽量地保持局部的形状

。

ARAP形变框架中有两个关键的子问题需要解决：

1. 位置变量（p）: 这些是变形后顶点的位置。p 代表网格上每个顶点移动后的位置坐标。在优化过程中，这些位置会被更新以寻找变形的最终状态。
2. 旋转变量（R）: 这些是局部旋转。R 表示与网格顶点相关联的旋转矩阵，它们用于描述网格局部区域的刚体旋转。对每个顶点或网格单元，都会找到一个最佳的局部旋转矩阵 R，使得变形后的局部结构在旋转意义上保持与原始局部结构的一致性。
   求解子问题1与2的这个过程通常会迭代进行，直到满足一定的收敛条件为止。

所以如下图是ARAP整体流程：分别固定 R 和 P ，求解都很容易：固定 R ，计算P ，固定旋转变量 R，找到新的位置变量 P，以降低全局的变形能量求解问题2 ； 固定P ，计算R ，固定位置变量 P，解决旋转变量 R，使局部结构尽可能地保持刚性求解问题1；

上图公式左边是arap的约束，希望所有边发生变化后形变最小，就是能量公式。左边为0（接近0） 是被控制点不管约束然后尽可能的靠近控制点！min 是希望这个着整个等式最小， 因为都是范数的平方，所以最好情况下等于0。 上图公式右边是用户手动调整的一些点（控制点），添加一个参数alpha，用于调整拖动的力度。alpha越大假如接近无穷大即，这样形变的被控制点一定会到用户拖动的点（控制点）的指定位置。

子问题一、

子问题二、

详细的子问题1 与2求解这里就不展开说看了呢一改比我讲的好！能量函数的最小化问题需要在局部和全局做迭代求解：全局过程使用新的旋转矩阵更新线性系统右边的向量，从而求解并更新每个顶点的位置P；局部过程则使用新的顶点位置更新协方差矩阵，并对协方差矩阵做奇异值分解，求得旋转矩阵R。反复迭代直至满足终止条件，获得变形后的网格。ARAP方法有一个性质就是最小化目标函数会使网格的边长变化最小，所以前面才说ARAP是如何保证刚性变化的就是根据这个来的！
特别注意：操作过程不能保证所有的网格点都只有刚性变换，所以只能在最优解下（不能算最小二乘法但是其实是类似的，因为论文发表的时间比最小二乘法的时间要早），找到对网格所有顶点尽可能接近的旋转矩阵。

初始化问题：分解成子问题1与2求解它们（类似于先有鸡还是先有蛋的问题），对此需要一开始固定哪一组变量。原作者的做法是如下：先把R初始化为3*3的单位三维旋转矩阵。原文如下

ARAP形变的性能消耗init和render阶段：

init：

• 预计算邻接关系：需要确定网格中每个顶点的邻域，这包括找到与每个顶点相连的边和顶点，建立数据结构来存储这些信息。
• 系统矩阵构建：建立线性方程组的稀疏矩阵(一共三个)，其复杂度与网格的顶点数和连接性有关!

render 每一帧：

• 局部步骤（Local
  Step）：对于每个顶点，需要找到一个最佳的旋转矩阵，这个旋转矩阵使得其邻域顶点的形变尽可能刚性。通过奇异值分解（SVD）来实现占实时消耗的60-70%左右。
• 全局步骤（Global
  Step）：更新顶点位置以最小化整体形变能量，求解一个线性方程组，该方程组的大小与网格的顶点数成正比。求解这样的方程组通常使用共轭梯度法占比不是很大。

每次迭代就上面提到的分别更新P和R，性能上主要是更新R，每个顶点都对应一个旋转矩阵，都需要进行svd分解计算，这一步占迭代时间的60-70%。由于 SVD 在局部步骤中的计算成本较高，确实，很多性能开销是发生在执行SVD分解上的。特别是对于大型网格，SVD的计算会导致显著的性能消耗。可以考虑在GPU端实现ARAP算法。

Laplacian Deformation形变与ARAP的比较：

基于旋转不变量的拉普拉斯形变： 利用了拉普拉斯坐标，它是描述顶点及其相邻顶点位置关系的相对坐标系统（主要是保持每个顶点与其邻接点的相对位置不变）。拉普拉斯坐标捕获了局部的微分几何特征，并倾向于保持物体形状的光滑性和连续性，但是不太考虑保持物体的刚性。
优点：

• 拉普拉斯形变能够比较好地保持网格的几何细节，因此适合微调模型。
• 相对于ARAP，它解决的线性系统通常更简单，求解速度可能更快，不需要SVD等的求解。

缺点：

• 拉普拉斯形变在处理大的变形时可能产生体积损失，也就是网格可能会塌陷。
• 因为它直接使用了原始的拉普拉斯坐标，所以没有显式地考虑刚体的运动。

ARAP形变：ARAP形变则着眼于局部的刚体旋转，尝试在变形过程中最小化顶点与其邻接点之间旋转导致的变形保持每个小片区域尽可能刚性，使局部网格在旋转后尽可能保持形状。局部形状应该像刚体一样旋转，而有限的弯曲只发生在连接这些局部刚体的接缝处。因此更好地应对旋转和扭曲操作，有效地保持物体的刚性。
优点：

• ARAP形变非常适合处理大的形变，它能够更好地保持网格的体积不变性，因此广泛用于动画和复杂变形（替代骨骼形变方案）。
• ARAP对用户交互友好，因为它允许直观的拖拽操作并且产生预期的结果。

缺点：

• ARAP需要在每次迭代时解决一个非线性最优化问题，比拉普拉斯要慢。
• 算法实现相对复杂，并且需要额外的步骤来计算和更新局部旋转。

算法对比总结：

• 当需要精细调整网格或执行小范围变形时，Laplacian Deformation是一个合适的选择，因为它能快速操作并且很好地保持细节。
• 对于广泛的、直观的大规模形变，ARAP提供了更加自然的结果，尤其是在动画和实时编辑应用中表现突出。

参考资料：

Mesh形变算法
用数学编辑3D模型（一）- Mesh Deformation with Laplacian Coordinates
用数学编辑3D模型（二）- As-Rigid-As-Possible
图形处理（三）简单拉普拉斯网格变形-Siggraph 2004
图形处理（五）基于旋转不变量的网格变形-Siggraph 2007
《As-Rigid-As-Possible Surface Modeling》
经典论文推导: As-Rigid-As-Possible(ARAP) Surface Modeling