From c6aacee36e234a292884fd31e16351863cc9d071 Mon Sep 17 00:00:00 2001
From: "Junrui (John) Xu" <johnxu@cs.toronto.edu>
Date: Tue, 30 Oct 2018 12:47:01 -0400
Subject: [PATCH 1/2] Finished tutte

---
 src/tutte.cpp | 26 ++++++++++++++++++++++++--
 1 file changed, 24 insertions(+), 2 deletions(-)

diff --git a/src/tutte.cpp b/src/tutte.cpp
index 6f6109b..468bccf 100644
--- a/src/tutte.cpp
+++ b/src/tutte.cpp
@@ -1,11 +1,33 @@
 #include "tutte.h"
+#include "igl/boundary_loop.h"
+#include "igl/repdiag.h"
+#include "igl/map_vertices_to_circle.h"
+#include "igl/cotmatrix.h"
+#include "igl/min_quad_with_fixed.h"
 
 void tutte(
   const Eigen::MatrixXd & V,
   const Eigen::MatrixXi & F,
   Eigen::MatrixXd & U)
 {
-  // Replace with your code
-  U = V.leftCols(2);
+    Eigen::VectorXi b;
+    igl::boundary_loop(F,b);
+
+    Eigen::MatrixXd bc;
+    igl::map_vertices_to_circle(V,b,bc);
+
+    Eigen::SparseMatrix<double> L;
+    igl::cotmatrix(V,F,L);
+
+    Eigen::SparseMatrix<double> Q = -L;
+
+    igl::min_quad_with_fixed_data<double> data;
+    igl::min_quad_with_fixed_precompute(Q, b, Eigen::SparseMatrix<double>(), false, data);
+
+    const int dims = b.cols();
+    const int size = data.n;
+
+    igl::min_quad_with_fixed_solve(data, Eigen::MatrixXd::Zero(size, dims),
+                                     bc, Eigen::MatrixXd::Zero(0, dims), U);
 }
 

From 890659e96b97fa18eafead435460b1f64b29199e Mon Sep 17 00:00:00 2001
From: "Junrui (John) Xu" <johnxu@cs.toronto.edu>
Date: Wed, 31 Oct 2018 01:29:14 -0400
Subject: [PATCH 2/2] Finish lscm

---
 src/lscm.cpp               | 36 ++++++++++++++++++++++++++++++++++--
 src/vector_area_matrix.cpp | 28 +++++++++++++++++++++++++---
 2 files changed, 59 insertions(+), 5 deletions(-)

diff --git a/src/lscm.cpp b/src/lscm.cpp
index d6e1b6c..22aa13d 100644
--- a/src/lscm.cpp
+++ b/src/lscm.cpp
@@ -1,10 +1,42 @@
 #include "lscm.h"
+#include "vector_area_matrix.h"
+#include "igl/cotmatrix.h"
+#include "igl/massmatrix.h"
+#include "igl/repdiag.h"
+#include "igl/eigs.h"
 
 void lscm(
   const Eigen::MatrixXd & V,
   const Eigen::MatrixXi & F,
   Eigen::MatrixXd & U)
 {
-  // Replace with your code
-  U = V.leftCols(2);
+    Eigen::SparseMatrix<double> A;
+    vector_area_matrix(F, A);
+
+    Eigen::SparseMatrix<double> L;
+    igl::cotmatrix(V, F, L);
+
+    Eigen::SparseMatrix<double> L_diag;
+    igl::repdiag(L, 2, L_diag);
+
+    Eigen::SparseMatrix<double> Q = L_diag - 2. * A;
+
+    Eigen::SparseMatrix<double> M;
+    igl::massmatrix(V, F, igl::MASSMATRIX_TYPE_DEFAULT, M);
+
+    Eigen::SparseMatrix<double> B;
+    igl::repdiag(M, 2, B);
+
+    Eigen::MatrixXd v;
+    Eigen::VectorXd lambda;
+    igl::eigs(Q, B, 3, igl::EIGS_TYPE_SM, v, lambda);
+
+    const int size = V.rows();
+    U.resize(size, 2);
+    U.col(0) = v.col(2).head(size);
+    U.col(1) = v.col(2).tail(size);
+
+    Eigen::JacobiSVD<Eigen::MatrixXd> svd(U.transpose() * U, Eigen::ComputeThinU | Eigen::ComputeThinV);
+
+    U = U * svd.matrixU();
 }
diff --git a/src/vector_area_matrix.cpp b/src/vector_area_matrix.cpp
index 1d60ef8..3c95871 100644
--- a/src/vector_area_matrix.cpp
+++ b/src/vector_area_matrix.cpp
@@ -1,11 +1,33 @@
 #include "vector_area_matrix.h"
+#include "igl/boundary_loop.h"
 
 void vector_area_matrix(
   const Eigen::MatrixXi & F,
   Eigen::SparseMatrix<double>& A)
 {
-  // Replace with your code
-  int V_size = F.maxCoeff()+1;
-  A.resize(V_size*2,V_size*2);
+    int V_size = F.maxCoeff() + 1;
+    A.resize(V_size * 2, V_size * 2);
+
+    std::vector<std::vector<int>> b;
+    igl::boundary_loop(F, b);
+
+    std::vector<Eigen::Triplet<double>> tripletList;
+
+    for (int i = 0; i < b.size(); i++) {
+        std::vector<int> loop = b[i];
+
+        for (int j = 0; j < loop.size(); j++) {
+
+            int ui = loop[j];
+            int uj = loop[(j + 1) % loop.size()];
+
+            tripletList.push_back(Eigen::Triplet<double>(ui + V_size, uj, -0.25));
+            tripletList.push_back(Eigen::Triplet<double>(uj, ui + V_size, -0.25));
+            tripletList.push_back(Eigen::Triplet<double>(ui, uj + V_size, 0.25));
+            tripletList.push_back(Eigen::Triplet<double>(uj + V_size, ui, 0.25));
+        }
+    }
+
+    A.setFromTriplets(tripletList.begin(), tripletList.end());
 }