input surface
harmonic map image
conformality
conformality

Surface Ricci Flow Tutorial with C++ Source

Smooth Case

Suppose $S$ is a surface embedded in the Euclidean space $R^3$, it has a natural induced Euclidean metric $\mathbf{g}$. Suppose $K$ is the Gaussian curvature induced by the metric $\mathbf{g}$. Choose a local coordinates, the metric tensor is $\mathbf{g}=(g_{ij})$. Then Hamilton's surface Ricci flow is defined as \[ \frac{ dg_{ij}(t) }{dt} = -2K(t)g_{ij}(t). \] Ricci flow deforms the Riemannian metric in a conformal way, \[ \mathbf{g}(t) = e^{2\lambda(t)}\mathbf{g}(0), \] the flow is given by \[ \frac{d\lambda(t)}{dt} = -K(t). \] The normalized Ricci flow is given by \[ \frac{ dg_{ij}(t) }{dt} = 2(\bar{K} -K(t) ) g_{ij}(t), \] where \[ \bar{K} = \frac{2\pi \chi(S)}{A(0)}, \] $\chi(S)$ is the Euler characteristic of the surface, $A(0)$ is the total area at time $0$. Normalized surface Ricci flow converges to the metric which induces constant Gaussian curvature $\bar{K}$.

Discrete Case

The surface is approximated by simplicial complex, namely triangle mesh $M$. A vertex is denoted as $v_i$; an oriented edge (halfedge) $[v_i,v_j]$, where $v_i$ is the source, $v_j$ is the target; an oriented face is $[v_i,v_j,v_k]$, where $v_i$, $v_j$ and $v_k$ are sorted counter-clock-wisely.
input surface
input surface
The discrete Riemannian metric is the edge length. Let $l_{ij}$ be the edge length for edge $[v_i,v_j]$, which satisfies the triangle inequality \[ l_{ij} + l_{jk} > l_{ki}. \] The discrete Gaussian curvature is the angle deficit. Suppose $v_i$ is an interior vertex of the mesh, \[ K_i = 2\pi - \sum_{jk}\theta_{i}^{jk}, \] where $\theta_i^{jk}$ is the corner angle at vertex $v_i$ in triangle $[v_i,v_j,v_k]$. If $v_i$ is on the boundary, then \[ K_i = \pi - \sum_{jk}\theta_{i}^{jk}. \] Surface Ricci flow conformally deforms the Riemannian metric, which transforms infinitesimal circles to infinitesimal circles. In discrete case, we replace infinitesimal circles by circles with finite radii. Each vertex $v_i$ is associated with a circle, $c(v_i,r_i)$. The edge length \[ l_{ij}^2 = \gamma_i^2 + \gamma_j^2 + 2I_{ij} \gamma_i\gamma_j, \] The discrete conformal metric deformation varies the radii $\gamma_i$'s, but preserves all the weights $I_{ij}$'s. Let \[ u_i = \log \gamma_i. \] Let the target curvature at vertex $v_i$ be $\bar{K}_i$. Denote $\mathbf{u}=(u_1,u_2,\cdots,u_n)$. Define an Energy \[ E(\mathbf{u}) = \int^\mathbf{u} \sum_{i=1}^n (\bar{K}_i - K_i) du_i. \] Then \[ \nabla E(\mathbf{u}) = (\bar{K}_1 - K_1, \bar{K}_2 - K_2, \cdots, \bar{K}_n - K_n )^T. \] The Hessian matrix of $E$ is positive definite on the sublinear space \[ \sum_{i=1}^n u_i = 0. \] Each triangle $[v_i,v_j,v_k]$ has a circle, which is orthogonal to three vertex circles $c(v_i,\gamma_i)$, $c(v_j,\gamma_j)$ and $c(v_k,\gamma_k)$, which is called the power circle. Suppose an edge $[v_i,v_j]$ attaches to two triangles $[v_i,v_j,v_k]$ and $[v_j,v_i,v_k]$, two power circle centers are $o_k$ and $o_l$. Define the edge weight $w_{ij}$ as \[ w_{ij} = \frac{|o_l - o_k| }{|v_i-v_j|}. \] Then we obtain \[ d\mathbf{K} = \Delta d\mathbf{u}, \] where $\Delta=(\Delta_{ij})$, \[ \Delta_{ij} = \left\{ \begin{array}{cl} -w_{ij} & [v_i,v_j] \in M\\ \sum_k w_{ik} & i= j \\ 0 & otherwise \end{array} \right. \] Namely, for each vertex $v_i$, \[ dK_i = \sum_{[v_i,v_j]\in M} w_{ij} ( du_i - du_j). \]

Algorithm

The computational algorithm is to use Newton's method to optimize the convex energy. The figure on top shows a conformal mapping from the human face to a rectangle. On the boundary, four corner vertices are selected, denoted as $\{v_0,v_1,v_2,v_3\}$. The target curvature on these 4 vertices equal to $\frac{\pi}{2}$, and zeros for any other vertex.
  1. Set the target curvature $\bar{K}_i$ for each vertex $v_i$, which satisfies the Gauss-Bonnet theorem \[ \sum_{v_i \in M} K_i = 2\pi \chi(M). \]
  2. For each vertex $v_i$, estimate the circle radius $\gamma_i$; for each edge, estimate the weight $I_{ij}$.
  3. Compute the edge length $l_{ij}$ from $\gamma_i$, $\gamma_j$ and $I_{ij}$, compute corner angles using cosine law.
  4. Compute current Gaussian curvature on each vertex $K_i$. If $\max_i |\bar{K}_i - K_i | < \delta$, return.
  5. For each face, compute the power circle; for each edge, compute the edge weight $w_{ij}$. Construct the Hessain matrix $\Delta$.
  6. Solve linear system, \[ \mathbf{\bar{K}} - \mathbf{K} = \Delta \delta\mathbf{u}, \]
  7. Update the radius vector \[ \mathbf{u} = \mathbf{u}+\delta \mathbf{u}, \gamma_i = e^{u_i}. \]
  8. Repeat from step 3 through step 7.

Source Code

The code is written in C++ using VisualStudio on Windows platform.
  1. Download the RicciFlowExtremalLength.zip file, unzip it.
  2. Enter RicciFlowExtremalLength/RicciFlowExtremalLength directory, double click on RiemannMapper.sln file. Compile the project using Visual Studio.
  3. The output binary is RicciFlowExtremalLength/bin/RicciFlow.exe
  4. Go to RicciFlowExtremalLength/demo/Alex directory, double click on RicciExtremalLength_test.bat
  5. The code will start running, and the results will be displayed simpleviewer.exe. The output file is Alex.remesh.uv.m.
The triangle mesh viewer can be found here.

Contact

If you encounter any problem during the compilation and running the code, please contact the author David Gu at $gu@cs.stonybrook.edu$.