## Assignment 4

Download the assignment archive DGP2019-A4.zip, and unzip the compressed file. You should have a directory named DGP2019-A4 containing a README.txt file and several subdirectories. Follow these instructions. You should upload your completed assignment through Canvas. As in previous assignments, your submission may have bugs, but it has to compile and the application has to run. Otherwise I will not grade your assignment.

In this assignment you will implement a number of smoothing algorithms according to the optimization approach described in Introduction to Geometry Processing Through Optimization

For this you have to edit the following two files

  core/Smoothing.hpp
core/Smoothing.cpp


More precisely, you have to complete the implementation of the following methods

  void  Smoothing::laplacianSmoothingVertexCoordinatesRun();
void  Smoothing::laplacianSmoothingFaceNormalsRun();
void  Smoothing::jacobiRun();
void  Smoothing::integrateNormalsRun();
void  Smoothing::coordsNormalsRun();


You should start by copying the files that you edited in previous assignments into the correspondiong src subdirectories.

#### 1) Unconstrained Laplacian Smoothing of Vertex Coordinates

  void  Smoothing::laplacianSmoothingVertexCoordinatesRun();


This process can be described as the application of the Jacobi descent method to the energy function $E_1(x) = \frac{1}{nE}\sum_{(i,j)} \| x_i-x_j \|^2$ where $$x_i$$ and $$x_j$$ denote the 3D coordinates of the $$i$$-th and $$j$$-th vertices, $$x$$ the concatenation of all the vertex coordinates, $$nE$$ is the total number of edges, and the sum is over all the edges $$(i,j)$$ of the primal graph of the mesh. Face normal vectors are recomputed at the end of the process.

STEPS

1.1) Compute vertex coordinate displacement vectors $\Delta x_i = \frac{1}{|i^\star|} \sum_{j\in i^\star} (x_j-x_i)$

To implement this step you should allocate a scalar accumulator $$w_i$$ and a $$3$$-dimensional accumulator $$\Delta x_i$$ for each mesh vertex, then for each primal edge $$e=(i,j)$$, select a positive edge weight $$w_e$$ (such as $$w_e=1$$), and increment the accumulators $\left\{\; \begin{matrix} w_i \leftarrow w_i+w_e \hfill \\ \Delta x_i \leftarrow \Delta x_i+w_e(x_j-x_i) \\ w_j \leftarrow w_j+w_e \hfill\\ \Delta x_j \leftarrow \Delta x_j+w_e(x_i-x_j) \end{matrix} \right.$

After all the edges have been considered, compute the displacement vectors as follows $\Delta x_i \leftarrow \Delta x_i / w_i$

Don't forget to test whether or not $$w_i$$ is positive before performing this operation.

1.2) Apply displacements $x_i \leftarrow x_i+\lambda\;\Delta x_i$

#### 2) Unconstrained Laplacian Smoothing of Face Normals

  void  Smoothing::laplacianSmoothingFaceNormalsRun();


This process can be described as the application of the Jacobi descent method to the energy function $E_2(n) = \frac{1}{nE^{\star}}\sum_{(f,g)} \| n_f-n_g \|^2$ where $$n_f$$ and $$n_g$$ denote the 3D normal vectors corresponding to the $$f$$-th and $$g$$-th faces, $$n$$ the concatenation of all the normal vectors, $$nE^{\star}$$ is the number of edges of the dual graph of the mesh, and the sum is over all the edges $$f,g$$ of the dual graph of the mesh (corresponding to the regular edges of the primal graph of the mesh). Vertex coordinates are not recomputed at the end of the process.

STEPS

2.1) Compute face normal displacement vectors $\Delta n_f = \frac{1}{|f^\star|} \sum_{g\in f^\star} (n_g-n_f)$

2.2) Apply displacements $n_f \leftarrow n_f+\lambda\;\Delta n_f$

2.3) Normalize face normals to unit length

#### 3) Jacobi Smoothing of Vertex Coordinates with Soft Constraints

  void  Smoothing::jacobiRun();


This process can be described as the application of the Jacobi descent method to the energy function $E_3(x) = \frac{\lambda_{D}}{nV} \sum_{i} \|x_i-x_i^0\| + \frac{\lambda_{S}}{nE} \sum_{(i,j)} \| x_i-x_j \|^2$

STEPS

3.1) Accumulate vertex weights and displacements $w_i = \frac{\lambda_{D}}{nV} + \frac{\lambda_{S}}{nE}\;|i^\star|$ $\Sigma x_i = \frac{\lambda_{D}}{nV} (x_i^0-x_i) + \frac{\lambda_{S}}{nE} \sum_{j\in i^\star} (x_j-x_i)$

3.2) Normalize displacements $\Delta x_i = \frac{1}{w_i}\; \Sigma x_i$

3.3) Apply displacements $x_i \leftarrow x_i+\lambda\;\Delta x_i$

#### 4) Integration of Face Normals with Soft Vertex Constraints

  void  Smoothing::integrateNormalsRun();


This process can be described as the application of the Jacobi descent method to the energy function $E_4(x) = \frac{\lambda_{D}}{nV} \sum_{i} \|x_i-x_i^0\| + \frac{\lambda_{L}}{nHE} \sum_{(f,i,j)} (n_f^t(x_i-x_j))^2 + \frac{\lambda_{S}}{nE} \sum_{(i,j)} \| x_i-x_j \|^2$

STEPS

4.1) Accumulate vertex weight matrices and displacements $W_i= \frac{\lambda_{D}}{nV}\;I + \frac{\lambda_{L}}{nHE} \sum_{(f,i,j)\in i^\star} (n_f n_f^t)+ \frac{\lambda_{S}}{nE}\;|i^\star|\;I$ $\Sigma x_i = \frac{\lambda_{D}}{nV} (x_i^0-x_i) + \frac{\lambda_{L}}{nHE} \sum_{(f,i,j)\in i^\star} (n_f n_f^t)(x_j-x_i) + \frac{\lambda_{S}}{nE} \sum_{j\in i^\star} (x_j-x_i)$

4.2) Solve $W_i\;\Delta x_i = \Sigma x_i$ Remember that since $$W_i$$ is a $$3\times 3$$ symmetric and positive definite matrix, the Cholesky decomposition theorem startes that there exists a unique lower triangular $$3\times 3$$ matrix $$L_i$$, so that $$W_i=L_i\,L_i^t$$. As a result, solving this system of $$3$$ linear equations in $$3$$ variables can be reduced to the solution of two triangular systems of equations $L_i\;u_i = \Sigma x_i$ $L_i^t\;\Delta x_i = u_i$ where $$u_i$$ is an auxiliary variable.

Since a $$3\times 3$$ symmetric and positive definite matrix $$W$$ only has $$6$$ degrees of freedom. it can be stored in compact form as a $$6$$-dimensional vector $$[W0,W1,W2,W3,W4,W5]$$ so that $W = \left[ \begin{matrix} W0 & W1 & W3 \\ W1 & W2 & W4 \\ W3 & W4 & W5 \end{matrix} \right]$

Since a lower triangular $$3\times 3$$ matrix $$L$$ also has only $$6$$ degrees of freedom. it can be stored in compact form as a $$6$$-dimensional vector $$[L0,L1,L2,L3,L4,L5]$$ so that $L = \left[ \begin{matrix} L0 & 0 & 0 \\ L1 & L2 & 0 \\ L3 & L4 & L5 \end{matrix} \right]$ Rather than using a library method to compute the Cholesky decomposition of the weight matrices, and since we are only dealing with $$3\times 3$$ matrices here, it is more efficient to hard-code the Cholesky decomposition algorithm as part of the smoothing inner loop. The equations to compute the coefficients of the matrix $$L$$ from the coefficients of the matrix $$W$$ can be derived by expanding the matrix product $\left[ \begin{matrix} W0 & W1 & W3 \\ W1 & W2 & W4 \\ W3 & W4 & W5 \end{matrix} \right] = \left[ \begin{matrix} L0 & 0 & 0 \\ L1 & L2 & 0 \\ L3 & L4 & L5 \end{matrix} \right] \left[ \begin{matrix} L0 & L1 & L3 \\ 0 & L2 & L4 \\ 0 & 0 & L5 \end{matrix} \right]$

Just to be safe, make sure that you test the sign of the argument before you compute a square root, since otherwise your program will probably crash if the matrix $$W$$ is not positive definite.

4.3) Apply displacements $x_i \leftarrow x_i+\lambda\;\Delta x_i$

#### 5) Jacobi Smoothing of Vertex Coordinates and Face Normals with Soft Vertex Constraints

  void  Smoothing::coordsNormalsRun();


This process can be described as the application of the Jacobi descent method to the energy function, where both the vertex coordinates $$x$$ and the face normal vectors $$n$$ are regarded as variables. $E_5(x,n) = \frac{\lambda_{D}}{nV} \sum_{i} \|x_i-x_i^0\| + \frac{\lambda_{L}}{nHE} \sum_{(f,i,j)} (n_f^t(x_i-x_j))^2 + \frac{\lambda_{S}}{nE^{\star}} \sum_{(f,g)} \| n_f-n_g \|^2$

STEPS

You should be able to figure out what the steps are in this case. But keep in mind that both the vertex coordinates and the face normal vectors are regarded as variables, and must be updated simultaneously. For this you need to accumulate weigths and displacement vectors for both vertex coordinates and face normals.

#### Test Files

A number of simple test files have been added to the assets directory. Since the polygon meshes described in these files are small, they can be used to debug your implementation. Some of these files contain scene graphs with multiple IndexedFaceSet nodes, some contain regular polygon meshes, and some contain singular polygon meshes. Some meshes have boundaries and others do not.