# 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.
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.
- 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).
\]
- For each vertex $v_i$, estimate the circle radius $\gamma_i$; for each edge, estimate the weight $I_{ij}$.
- Compute the edge length $l_{ij}$ from $\gamma_i$, $\gamma_j$ and $I_{ij}$, compute corner angles using cosine law.
- Compute current Gaussian curvature on each vertex $K_i$. If $\max_i |\bar{K}_i - K_i | < \delta$, return.
- For each face, compute the power circle; for each edge, compute the edge weight $w_{ij}$. Construct the Hessain matrix $\Delta$.
- Solve linear system,
\[
\mathbf{\bar{K}} - \mathbf{K} = \Delta \delta\mathbf{u},
\]
- Update the radius vector
\[
\mathbf{u} = \mathbf{u}+\delta \mathbf{u}, \gamma_i = e^{u_i}.
\]
- Repeat from step 3 through step 7.

### Source Code

The code is written in C++ using VisualStudio on Windows platform.
- Download the RicciFlowExtremalLength.zip file, unzip it.
- Enter RicciFlowExtremalLength/RicciFlowExtremalLength directory, double click on RiemannMapper.sln file. Compile the project using Visual Studio.
- The output binary is RicciFlowExtremalLength/bin/RicciFlow.exe
- Go to RicciFlowExtremalLength/demo/Alex directory, double click on RicciExtremalLength_test.bat
- 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$.