3. 微分几何

曲线

• 光滑平面曲线，可认为是可微的一维流形，参数化为：

$$\begin{gathered} \mathbf{x}:[a,b]\to {\mathbb{R}}^2, \\ \mathbf{x}(u)=(x(u),y(u){)}^T,u\in [a,b]\subset \mathbb{R} \end{gathered}$$

• 曲线的切向量${\mathbf{x}}^{\prime }(u)=(x^{\prime }(u),y^{\prime }(u){)}^T$，是一阶微分

• 假设参数化是正则的（regular），就是说${\mathbf{x}}^{\prime }(u)\ne 0$，那么法向量$\mathbf{n}(u)=x^{\prime }(u{)}^{\perp }/||x^{\prime }(u{)}^{\perp }||$，$\perp$是旋转90°的意思

• 同一条曲线可能有不同的参数化方式；但是对一个固定的$u$来说，参数化的结果（值$\mathbf{x}(u)$)是不一样的。通过重参数化可以实现不改变曲线形状下，更改参数化结果（值）

• 微分几何关心的是与参数化方式无关的曲线性质，如长度、曲率等

弧长

• 弧长参数化，是一种保留长度的映射，也是具有**等距性（isometry）**的重参数化

$$s=s(u)={\int }_a^u||{\mathbf{x}}^{\prime }(t)||dt$$

• 这种参数化与参数化方式无关，总是实现从参数域$[a,b]$到长度域$[0,L]$的映射，其中$L=\int a^b||x^{\prime }(u)||du$是曲线总长度
• 任意参数化方式，对曲线上的同一点，其到起点的弧长都是相等的
• 但是曲线不满足regular条件时，这个性质就不再成立了

曲率

• 设一个正则曲线用弧长参数化，可定义某点的曲率为：

$$\kappa (s):=||x^{\prime \prime }(s)||$$

曲率可以描述切向量的倒数与法向量的关系：
$${\mathbf{x}}^{\prime \prime }(s)=\kappa (s)\mathbf{n}(s)$$

处处曲率为0的曲线，是线段；

处处曲率相等的曲线，是圆

曲率还可以定义为内切圆半径的倒数

表面

弧长和曲率是曲线的欧式不变量，也就是对刚性运动无关。类似的，光滑表面也有度量

表面的参数化表示

• 例子：球面坐标系

$$\begin{gathered} (\theta ,\varphi )\to (x,y,z) \\ \theta \in [0,2\pi ],\varphi \in [-\pi ,\pi ] \\ \mathbf{x}(\theta ,\varphi )=(x(\theta ,\varphi ),y(\theta ,\varphi ),z(\theta ,\varphi ){)}^T=(R\cos (\theta )\cos (\varphi ),R\sin (\theta )\cos (\varphi ),R\sin (\varphi ){)}^T \end{gathered}$$

• 正式形式：$\mathbf{x}(u,v)=(x(u,v),y(u,v),z(u,v){)}^T$

• 以球面坐标系为例，表面的参数化具有以下定义：

  • 点$\mathbf{p}(x,y,z)$的球面坐标定义为$(\theta ,\varphi )$
  • 参数域中每条垂直线，也就是$\theta =c$，对应3D表面的一条曲线，iso-$\theta$曲线，就是连通两个极点的经线；
  • 参数域中每条水平先，也就是$\varphi =c$，对应3D表面的纬线，iso-$\varphi$曲线，$\varphi =0$对应的就是赤道

• 注意到参数域到表面出现的失真

度量性质

• 给出连续表面的参数化形式：

$$\mathbf{x}(u,v)=\left(\begin{array}{c}x(u,v)\\ y(u,v)\\ z(u,v)\end{array}\right),(u,v)\in \Omega \subset {\mathbb{R}}^2$$

• 对一点$\mathbf{x}(u_0,v_0)\in \mathcal{S}$，有两条***iso-parameter曲线***：

$$\begin{gathered} {\mathbf{C}}_{\mathbf{u}}(t)=\mathbf{x}(u_0+t,v_0) \\ {\mathbf{C}}_{\mathbf{v}}(t)=\mathbf{x}(u_0,v_0+t) \end{gathered}$$

求偏导
$$\begin{gathered} {\mathbf{x}}_u(u_0,v_0):=\frac{\partial \mathbf{x}}{\partial u}(u_0,v_0) \\ {\mathbf{x}}_v(u_0,v_0):=\frac{\partial \mathbf{x}}{\partial v}(u_0,v_0) \end{gathered}$$
设一个正则参数化，也就是说${\mathbf{x}}_u\times {\mathbf{x}}_v\ne 0$，切平面可由这两个切向量张成；表面法向量同时垂直于这两个切向量：
$$\mathbf{n}=\frac{{\mathbf{x}}_u\times {\mathbf{x}}_v}{||{\mathbf{x}}_u\times {\mathbf{x}}_v||}$$

• 进一步可以给出方向导数。设参数空间的向量$\vec{\mathbf{w}}=(u_w,v_w{)}^T$，考虑一条用$(t,\vec{w})$参数化的直线$(u,v)=(u_0,v_0)+t\vec{\mathbf{w}}$，这条直线在表面上就是通过点$\mathbf{x}$的一条曲线

$${\mathbf{C}}_{\mathbf{w}}(t)=\mathbf{x}(u_0+t{u}_w,v_0+t{v}_w)$$

$\mathbf{w}$关于$\mathbf{x}$在$(u_0,v_0)$的方向导数定义为曲线${\mathbf{C}}_{\mathbf{w}}$在$t=0$时的切向，$\mathbf{w}=\partial {\mathbf{C}}_{\mathbf{w}}(t)/\partial t$

应用链式法则，$\mathbf{w}=\mathbf{J}\vec{\mathbf{w}}$，
$$\mathbf{J}=\left[\begin{array}{cc}\frac{\partial x}{\partial u}& \frac{\partial x}{\partial v}\\ \frac{\partial y}{\partial u}& \frac{\partial y}{\partial v}\\ \frac{\partial z}{\partial u}& \frac{\partial z}{\partial v}\end{array}\right]=[{\mathbf{x}}_u,{\mathbf{x}}_v]$$

第一基本形式

• 参数化函数$\mathbf{x}$的Jacobian矩阵对应的是参数空间向量$\vec{\mathbf{w}}$到表面切向量$\mathbf{w}$的线性映射
• Jacobian矩阵编码了表面的度量：允许度量角度、距离和面积怎样从参数域变换到表面的

设$\overrightarrow{{\mathbf{w}}_1}$和$\overrightarrow{{\mathbf{w}}_2}$是参数空间的两个单位方向向量，夹角余弦为${\overrightarrow{{\mathbf{w}}_1}}^T\overrightarrow{{\mathbf{w}}_2}$，那么对应的表面的两个切向量的夹角余弦为
$${\mathbf{w}}_1^T{\mathbf{w}}_2=({\mathbf{J}{\vec{\mathbf{w}}}_1}^T)(\mathbf{J}{\vec{\mathbf{w}}}_2)={\mathbf{w}}_1^T({\mathbf{J}}^T\mathbf{J})\overrightarrow{{\mathbf{w}}_2}$$
${\mathbf{J}}^T\mathbf{J}$就是第一基本形式：
$$\mathbf{I}={\mathbf{J}}^T\mathbf{J}=\left[\begin{array}{cc}E& F\\ F& G\end{array}\right]:=\left[\begin{array}{cc}{\mathbf{x}}_u^T{\mathbf{x}}_u& {\mathbf{x}}_u^T{\mathbf{x}}_v\\ {\mathbf{x}}_u^T{\mathbf{x}}_v& {\mathbf{x}}_v^T{\mathbf{x}}_v\end{array}\right]$$
第一基本形式除度量角度外，还可以度量切向量的长度的平方$||\mathbf{w}|{|}^2={\vec{\mathbf{w}}}^T\mathbf{I}\vec{\mathbf{w}}$

• 使用第一基本形式可以度量曲线$\mathbf{x}(t)=\mathbf{x}(\mathbf{u}(t))$的长度，其中$\mathbf{u}(t)=(u(t),v(t))$在参数域

曲线的切向量可表示为：
$$\frac{d\mathbf{x}(\mathbf{u}(t))}{dt}=\frac{\partial \mathbf{x}}{\partial u}\frac{\partial u}{\partial t}+\frac{\partial \mathbf{x}}{\partial v}\frac{\partial v}{\partial t}={\mathbf{x}}_u{u}_t+{\mathbf{x}}_v{u}_v$$
曲线的长度可表示为：
$$l(a,b)={\int }_a^b\sqrt{(u_t,v_t)\mathbf{I}(u_v,v_t{)}^T}dt={\int }_a^b\sqrt{E{u}_t^2+2F{u}_t{v}_t+G{v}_t^2}dt$$
表面面积可以记为
$$A={\iint }_U\sqrt{det(\mathbf{I})}dudv={\iint }_U\sqrt{EG-F^2}dudv$$

• 因为第一基本形式可以度量角度、距离和面积，因此可被认为是一种几何工具，也被叫做度量张量（metric tensor）

各向异性

除把参数空间的方向变换为表面的切向量外，还可以把参数空间的圆变换成小所愿，anisotropy ellipse

• 各向异性随缘的主轴：${\mathbf{e}}_1=\mathbf{J}\overrightarrow{{\mathbf{e}}_1},{\mathbf{e}}_2=\mathbf{J}\overrightarrow{{\mathbf{e}}_2}$
• 主轴长度${\sigma }_1=\sqrt{{\lambda }_1},{\sigma }_2=\sqrt{{\lambda }_2}$

表面曲率

对表面一点$\mathbf{p}$，设切向量是$\mathbf{t}=u_t{\mathbf{x}}_t+v_t{\mathbf{x}}_v$，对应参数空间的$\vec{\mathbf{t}}=(u_t,v_t{)}^T$

则法曲率（normal curvature）${\kappa }_n(\vec{\mathbf{t}})$是平面曲线（由切向量和表面法向量张成的平面与表面求交得到的曲线）的曲率

$${\kappa }_n(\vec{\mathbf{t}})=\frac{{\vec{\mathbf{t}}}^T\mathbf{I}\mathbf{I}\vec{\mathbf{t}}}{{\vec{\mathbf{t}}}^T\mathbf{I}\vec{\mathbf{t}}}=\frac{e{u}_t^2+2f{u}_t{v}_t+gv+t^2}{E{u}_t^2+2F{u}_t{v}_t+G{v}_t^2}$$

• $\mathbf{I}\mathbf{I}$是第二基本形式：

$$\mathbf{I}\mathbf{I}=\left[\begin{array}{cc}e& f\\ f& g\end{array}\right]=\left[\begin{array}{cc}{\mathbf{x}}_{uu}^T\mathbf{n}& {\mathbf{x}}_{uv}^T\mathbf{n}\\ {\mathbf{x}}_{uv}^T\mathbf{n}& {\mathbf{x}}_{vv}^T\mathbf{n}\end{array}\right]$$

$uu$是二阶偏导

• 假如绕着点$\mathbf{p}$的法向旋转切向量，然后得到不同的法曲率，这个法曲率有两个极值，称为主曲率${\kappa }_1,{\kappa }_2$（一个是最大曲率，一个是最小曲率）
• 如果${\kappa }_1\ne {\kappa }_2$，对应两个单位切向量${\mathbf{t}}_1,{\mathbf{t}}_2$，称为主方向
• 表面上两个曲率相等的点，称为umbilical或locally spherical，这些点的曲率性质是各向同的。平面或球表面的任一点都是umbilical的

Euler定理

欧拉关于法曲率和主曲率的关系：
$${\kappa }_n(\vec{\mathbf{t}})={\kappa }_1{\cos }^2\psi +{\kappa }_2{\sin }^2\psi$$
式中$\psi$是切向量$\mathbf{t}$和主方向${\mathbf{t}}_1$的夹角

• 这个关系说明，表面的曲率完全由两个主曲率决定；任意法曲率都是这两个曲率的凸集组合；
• 两个主方向永远是正交的

曲率张量

表面的局部性质可以用曲率张量$\mathbf{C}$表示，一个3×3的对称矩阵，特征值为${\kappa }_1,{\kappa }_2,0$，特征向量${\mathbf{t}}_1,{\mathbf{t}}_2,\mathbf{n}$，构造方法是$\mathbf{C}={\mathbf{P}\mathbf{D}\mathbf{P}}^{-1}$,$\mathbf{P}=[{\mathbf{t}}_1,{\mathbf{t}}_2,\mathbf{n}]$, $\mathbf{D}=diag({\kappa }_1,{\kappa }_2,0)$

其他两个曲率度量：

• 平均曲率$H=\frac{{\kappa }_1+{\kappa }_2}{2}$
• 高斯曲率$K={\kappa }_1{\kappa }_2$

高斯曲率可以把表面点分为三类：

• elliptical椭圆点（$K>0$）
• hyperbolic双曲线点（$K<0$）
• parabolic抛物线点（$K=0$）

双曲线点表面是局部鞍型；椭圆点是局部凸性的；抛物线点一般是分离前两者的

固有几何（Intrinsic Geometry）

• 仅依赖第一基本形式的性质称为固有的。固有几何只靠表面2D物体即可获知，不需要第三维信息
• 高斯曲率是固有几何，可以直接从第一基本形式推导出来；但是平均曲率不是
• 固有几何与参数化方法无关

Laplace算子

• Laplace算子定义为梯度的散度，对2-参数函数，可写做二阶偏导的和：

$$\Delta f=div\nabla f=div\left(\begin{array}{c}{f}_u\\ f_v\end{array}\right)=f_{uu}+f_{vv}$$

• Laplace-Beltrami算子拓展了Laplace算子到二维流形：

$${\Delta }_{\mathcal{S}}f=di{v}_{\mathcal{S}}{\nabla }_{\mathcal{S}}f$$

需要给出流形上梯度和散度的定义。给定坐标$\mathbf{x}$，一种计算方式是：
$${\Delta }_{\mathcal{S}}\mathbf{x}=-2H\mathbf{n}$$

尽管Laplace-Beltrami算子计算的时候使用了平均曲率，但是Laplace-Beltrami算子本身是intrinsic的

离散微分算子

对多边形网格来说，上述概念不能直接使用，需要定义离散微分算子。假设网格可以被看做是光滑表面的分段线性逼近，那么目标就是计算underlying表面的微分性质的approximation。经典的Laplace-Beltrami算子的离散化方法：

局部平均区域

• 总的思想是在一个局部邻域$\mathcal{N}(\mathbf{x})$，计算空间平均
• 通常，$\mathbf{x}$是网格顶点，n-ring领域或者其他方式用做平均区域
• 局部区域的大小影响离散算子的稳定性和准确率；淋浴越大，平均操作月光花，计算结果越稳定，但是不准确

几种平均区域的变种：

• Barycentric cell，连接三角形的重心和边的中点；
• Voronoi cell，连接三角形外心和边的中点
• Mixed Voronoi cell，修改一下，内角是钝角的三角形，用对面边中点取代外心

法向量

对独立的三角形$T=({\mathbf{x}}_i,{\mathbf{x}}_j,{\mathbf{x}}_k)$，法向量$\mathbf{n}(T)=\frac{({\mathbf{x}}_j-{\mathbf{x}}_i)\times ({\mathbf{x}}_k-{\mathbf{x}}_i)}{||({\mathbf{x}}_j-{\mathbf{x}}_i)\times ({\mathbf{x}}_k-{\mathbf{x}}_i)||}$

通过计算one-ring邻居的法向量的平均值，可以得到
$$\mathbf{n}(v)=\frac{{\sum }_{T\in {\mathcal{N}}_1(v)}{\alpha }_T\mathbf{n}(T)}{||{\sum }_{T\in {\mathcal{N}}_1(v)}{\alpha }_T\mathbf{n}(T)||}$$

• ${\alpha }_T$有不同的选择：
  • 常数，不考虑长度、三角形面积、角度，计算结构可能反直觉
  • 基于三角形面积加权，${\alpha }_T=|T|$，也可能反直觉
  • 基于相邻三角形的内角，${\alpha }_T={\theta }_T$，计算开销较大，但效果较好

梯度

设一个分段线性函数$f$，方便起见，记$f(v_i)=f({\mathbf{x}}_i)=f({\mathbf{u}}_i)=f_i$，对三角形$({\mathbf{x}}_i,{\mathbf{x}}_j,{\mathbf{x}}_k)$的一点，可线性插值为：
$$f(\mathbf{u})=f_i{B}_i(\mathbf{u})+f_j{B}_j(\mathbf{u})+f_k{B}_k(\mathbf{u})$$
$\mathbf{u}=(u,v)$是2D共性参数化中的参数对

$B$是重心基函数：

重心基函数大小为1，且满足$B_i(\mathbf{u})+B_j(\mathbf{u})+B_k(\mathbf{u})=1$；并且$\nabla B_i(\mathbf{u})+\nabla B_j(\mathbf{u})+\nabla B_k(\mathbf{u})=0$

因此$f$的梯度：
$$\begin{gathered} \nabla f(\mathbf{u})=f_i\nabla B_i(\mathbf{u})+f_j\nabla B_j(\mathbf{u})+f_k\nabla B_k(\mathbf{u}) \\ \nabla f(\mathbf{u})=(f_j-f_i)\nabla B_j(\mathbf{u})+(f_k-f_i)\nabla B_k(\mathbf{u}) \end{gathered}$$

$$\nabla B_i(\mathbf{u})=\frac{({\mathbf{x}}_k-x_j{)}^{\perp }}{2A_T}$$

沿着顶点$v_i$对边垂直方向，下降速度最快

最终，梯度表示为：
$$\nabla f(\mathbf{u})=(f_j-f_i)\frac{({\mathbf{x}}_i-x_k{)}^{\perp }}{2A_T}+(f_k-f_i)\frac{({\mathbf{x}}_j-x_i{)}^{\perp }}{2A_T}$$

离散Laplace-Meltrami算子

两种形式：Uniform graph Laplacian和余切形式

Uniform Laplacian

【Taubin95】
$$\Delta f(v_i)=\frac{1}{|{\mathcal{N}}_1(v_i|}\sum _{v_j\in {\mathcal{N}}_1(v_i)}(f_j-f_i)$$
应用到坐标函数$\mathbf{x}$时，uniform Laplacian相当于一个从中心顶点${\mathbf{x}}_i$指向one-ring平均顶点的向量

这种算式，当所有顶点共面时，计算结果不为0；然而我们想让结果为0。因此，uniform Laplacian对非均匀网格不适用。为什么呢？因为算式只考虑了邻接关系，没有考虑顶点的空间分布。

从另一方面，这种不变形又可以用来改进局部顶点分布

余切形式

【Meyer03】

在一个局部平均区域$A_i=A(v_i)$对梯度的散度的积分；使用散度定理，对向量函数$\mathbf{F}$：
$${\int }_{A_i}div\mathbf{F}(\mathbf{u})dA={\int }_{\partial A_i}\mathbf{F}(\mathbf{u})\cdot \mathbf{n}(\mathbf{u})ds$$
这个方程把平均区域$A_i$的积分和边缘$\partial A_i$的积分关联起来

应用于Laplacian，有
$${\int }_{A_i}\Delta f(\mathbf{u})dA={\int }_{A_i}div\nabla f(\mathbf{u})dA={\int }_{\partial A_i}\nabla f(\mathbf{u})\cdot \mathbf{n}(\mathbf{u})ds$$
$\Delta f(\mathbf{u})$是是散度，$div$是偏导，平均区域的梯度的散度可以用梯度与法向量的点乘对边缘的积分来替代

我们一个一个三角形看过去，比如对一个三角形$T$的区域，$\partial A_i\cap T$，求一下梯度的散度的积分；又由于同一个三角形内，梯度处处相等，所以有：
$${\int }_{\partial A_i\cap T}\nabla f(\mathbf{u})\cdot \mathbf{n}(\mathbf{u})ds=\nabla f(\mathbf{u}){\int }_{\partial A_i\cap T}\mathbf{n}(\mathbf{u})ds$$
又，这个边缘是由两个边的中点加上顶点组成的，所以有：
$$上式=\nabla f(\mathbf{u})\cdot (\mathbf{a}-\mathbf{b}{)}^{\perp }=\frac{1}{2}\nabla f(\mathbf{u})\cdot ({\mathbf{x}}_j-{\mathbf{x}}_k{)}^{\perp }$$
把前面梯度公式带进来，有
$$上式=(f_j-f_i)\frac{({\mathbf{x}}_i-{\mathbf{x}}_k{)}^{\perp }\cdot ({\mathbf{x}}_j-{\mathbf{x}}_k{)}^{\perp }}{4A_T}+(f_k-f_i)\frac{({\mathbf{x}}_j-{\mathbf{x}}_i{)}^{\perp }\cdot ({\mathbf{x}}_j-{\mathbf{x}}_k{)}^{\perp }}{4A_T}$$
正好让${\gamma }_j,{\gamma }_k$表示另外两个顶点$v_j,v_k$的夹角，根据正弦定理和余弦定理，上面的式子可以进一步简化为：
$$上式=\frac{1}{2}(\cot {\gamma }_k(f_j-f_i)+\cot {\gamma }_j(f_k-f_i))$$
最后对整个平均区域积分，得到
$${\int }_{A_i}\Delta f(\mathbf{u})dA=\frac{1}{2}\sum _{v_j\in {\mathcal{N}}_1(v_i)}(\cot {\alpha }_{i,j}+\cot {\beta }_{i,j})(f_j-f_i)$$
则Laplace-Beltrami算子的余切形式为：
$$\Delta f(v_i):=\frac{1}{2A_i}\sum _{v_j\in {\mathcal{N}}_1(v_i)}(\cot {\alpha }_{i,j}+\cot {\beta }_{i,j})(f_j-f_i)$$

• 余切形式也有不足
  • 如果${\alpha }_{i,j}+{\beta }_{i,j}>\pi$时，余切权重变成负值，可能引起三角形翻转；
  • 余切形式不是intrinsic的
• 类似的，也可以定义一个散度算子（Tong2003）

考虑一个向量域$\mathbf{w}:\mathcal{S}\to {\mathbb{R}}^3$，每个三角形定义一个常向量${\mathbf{w}}_T$,离散散度就是每个顶点$v_i$计算一个标量值$div\mathbf{w}(v_i)$,
$$div\mathbf{w}(v_i)=\frac{1}{A_i}\sum _{T\in {\mathcal{N}}_1(v_i)}\nabla B_i{|}_T\cdot {\mathbf{w}}_T{A}_T$$
这里$\nabla B_i{|}_T$是顶点$v_i$在三角形$T$中的基函数的梯度向量（常向量）

离散曲率

Laplace-Beltrami算子提供了平均曲法向的离散逼近，于是定义离散平均曲率的绝对值为：
$$H(v_i)=\frac{1}{2}||\Delta {\mathbf{x}}_i||$$
高斯曲率的离散算子：
$$K(v_i)=\frac{1}{A_i}(2\pi -\sum _{v_j\in {\mathcal{N}}_1(v_i)}{\theta }_j)$$
由此，可以得到主曲率：
$${\kappa }_{1,2}(v_i)=H(v_i)\pm \sqrt{H(v_i{)}^2-K(v_i)}$$

离散曲率张量

怎么定义离散曲率张量？

对每个边，沿着边指定最小曲率为零；根据二面角指定最大曲率。对一个局部区域，计算和它相较的所有边的和，然后平均
$$\mathbf{C}(v)=\frac{1}{A(v)}\sum _{\mathbf{e}\in A(v)}\beta (\mathbf{e})||\mathbf{e}\cap A(v)||\bar{\mathbf{e}}{\bar{\mathbf{e}}}^T$$
$\beta (\mathbf{e})$是有向二面角；$\mathbf{e}$是边；$\bar{\mathbf{e}}=\mathbf{e}/||\mathbf{e}||$