added templates and documentation

Former-commit-id: df3828f371ba561d0843121f5accb5fc5d3adaf3

include/igl/sparse_AtA_fast.cpp → include/igl/AtA_cached.cpp

@@ -5,18 +5,17 @@
 // This Source Code Form is subject to the terms of the Mozilla Public License
 // v. 2.0. If a copy of the MPL was not distributed with this file, You can
 // obtain one at http://mozilla.org/MPL/2.0/.
-#include "sparse_AtA_fast.h"
+#include "AtA_cached.h"
 
 #include <iostream>
 #include <vector>
-#include <unordered_map>
-#include <map>
 #include <utility>
 
-IGL_INLINE void igl::sparse_AtA_fast_precompute(
-    const Eigen::SparseMatrix<double>& A,
-    Eigen::SparseMatrix<double>& AtA,
-    igl::sparse_AtA_fast_data& data)
+template <typename Scalar>
+IGL_INLINE void igl::AtA_cached_precompute(
+    const Eigen::SparseMatrix<Scalar>& A,
+    Eigen::SparseMatrix<Scalar>& AtA,
+    igl::AtA_cached_data& data)
 {
   // 1 Compute At (this could be avoided, but performance-wise it will not make a difference)
   std::vector<std::vector<int> > Col_RowPtr;
@@ -47,7 +46,7 @@ IGL_INLINE void igl::sparse_AtA_fast_precompute(
     }
   }
 
-  Eigen::SparseMatrix<double> At = A.transpose();
+  Eigen::SparseMatrix<Scalar> At = A.transpose();
   At.makeCompressed();
   AtA = At * A;
   AtA.makeCompressed();
@@ -107,13 +106,14 @@ IGL_INLINE void igl::sparse_AtA_fast_precompute(
   }
   data.I_outer.push_back(data.I_row.size()); // makes it more efficient to iterate later on
 
-  igl::sparse_AtA_fast(A,AtA,data);
+  igl::AtA_cached(A,AtA,data);
 }
 
-IGL_INLINE void igl::sparse_AtA_fast(
-    const Eigen::SparseMatrix<double>& A,
-    Eigen::SparseMatrix<double>& AtA,
-    const igl::sparse_AtA_fast_data& data)
+template <typename Scalar>
+IGL_INLINE void igl::AtA_cached(
+    const Eigen::SparseMatrix<Scalar>& A,
+    Eigen::SparseMatrix<Scalar>& AtA,
+    const igl::AtA_cached_data& data)
 {
   for (unsigned i=0; i<data.I_outer.size()-1; ++i)
   {
@@ -125,4 +125,6 @@ IGL_INLINE void igl::sparse_AtA_fast(
 
 
 #ifdef IGL_STATIC_LIBRARY
+template void igl::AtA_cached<double>(Eigen::SparseMatrix<double, 0, int> const&, Eigen::SparseMatrix<double, 0, int>&, igl::AtA_cached_data const&);
+template void igl::AtA_cached_precompute<double>(Eigen::SparseMatrix<double, 0, int> const&, Eigen::SparseMatrix<double, 0, int>&, igl::AtA_cached_data&);
 #endif

include/igl/AtA_cached.h

@@ -0,0 +1,68 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2017 Daniele Panozzo <daniele.panozzo@gmail.com>
+// 
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_ATA_CACHED_H
+#define IGL_ATA_CACHED_H
+#include "igl_inline.h"
+#define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
+#include <Eigen/Dense>
+#include <Eigen/Sparse>
+namespace igl
+{  
+  struct AtA_cached_data
+  {
+    // Weights
+    Eigen::VectorXd W;
+
+    // Flatten composition rules
+    std::vector<int> I_row;
+    std::vector<int> I_col;
+    std::vector<int> I_w;
+
+    // For each entry of AtA, points to the beginning
+    // of the composition rules
+    std::vector<int> I_outer;
+  };
+
+  // Computes At * W * A, where A is sparse and W is diagonal. Divides the 
+  // construction in two phases, one
+  // for fixing the sparsity pattern, and one to populate it with values. Compared to
+  // evaluating it directly, this version is slower for the first time (since it requires a
+  // precomputation), but faster to the subsequent evaluations.
+  //
+  // Input:
+  //   A m x n sparse matrix
+  //   data stores the precomputed sparsity pattern, data.W contains the optional diagonal weights (stored as a dense vector). If W is not provided, it is replaced by the identity.
+  // Outputs:
+  //   AtA  m by m matrix computed as AtA * W * A
+  //
+  // Example:
+  // AtA_data = igl::AtA_cached_data();
+  // AtA_data.W = W;
+  // if (s.AtA.rows() == 0)
+  //   igl::AtA_cached_precompute(s.A,s.AtA,s.AtA_data);
+  // else
+  //   igl::AtA_cached(s.A,s.AtA,s.AtA_data);
+  template <typename Scalar>
+  IGL_INLINE void AtA_cached_precompute(
+    const Eigen::SparseMatrix<Scalar>& A,
+    Eigen::SparseMatrix<Scalar>& AtA,
+    AtA_cached_data& data);
+
+  template <typename Scalar>
+  IGL_INLINE void AtA_cached(
+    const Eigen::SparseMatrix<Scalar>& A,
+    Eigen::SparseMatrix<Scalar>& AtA,
+    const AtA_cached_data& data);
+  
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "AtA_cached.cpp"
+#endif
+
+#endif

include/igl/slim.cpp

@@ -35,8 +35,8 @@
 #include <Eigen/IterativeLinearSolvers>
 
 #include <igl/Timer.h>
-#include <igl/sparse_fast.h>
-#include <igl/sparse_AtA_fast.h>
+#include <igl/sparse_cached.h>
+#include <igl/AtA_cached.h>
 
 #ifdef CHOLMOD
 #include <Eigen/CholmodSupport>
@@ -505,10 +505,10 @@ namespace igl
       if (s.A.rows() == 0)
       {
         s.A = Eigen::SparseMatrix<double>(s.dim * s.dim * s.f_n, s.dim * s.v_n);
-        igl::sparse_fast_precompute(IJV,s.A,s.A_data);
+        igl::sparse_cached_precompute(IJV,s.A,s.A_data);
       }
       else
-        igl::sparse_fast(IJV,s.A,s.A_data);
+        igl::sparse_cached(IJV,s.A,s.A_data);
       #else
       Eigen::SparseMatrix<double> A(s.dim * s.dim * s.f_n, s.dim * s.v_n);
       buildA(s,IJV);
@@ -534,9 +534,9 @@ namespace igl
       #ifdef SLIM_CACHED
       s.AtA_data.W = s.WGL_M;
       if (s.AtA.rows() == 0)
-        igl::sparse_AtA_fast_precompute(s.A,s.AtA,s.AtA_data);
+        igl::AtA_cached_precompute(s.A,s.AtA,s.AtA_data);
       else
-        igl::sparse_AtA_fast(s.A,s.AtA,s.AtA_data);
+        igl::AtA_cached(s.A,s.AtA,s.AtA_data);
 
       L = s.AtA + s.proximal_p * id_m; //add also a proximal 
       L.makeCompressed();

include/igl/slim.h

@@ -17,7 +17,7 @@
 #define SLIM_CACHED 
 
 #ifdef SLIM_CACHED
-#include <igl/sparse_AtA_fast.h>
+#include <igl/AtA_cached.h>
 #endif
 
 namespace igl
@@ -77,7 +77,7 @@ struct SLIMData
   Eigen::SparseMatrix<double> A;
   Eigen::VectorXi A_data;
   Eigen::SparseMatrix<double> AtA;
-  igl::sparse_AtA_fast_data AtA_data;
+  igl::AtA_cached_data AtA_data;
   #endif
 };
 

include/igl/sparse_AtA_fast.h

@@ -1,47 +0,0 @@
-// This file is part of libigl, a simple c++ geometry processing library.
-// 
-// Copyright (C) 2017 Daniele Panozzo <daniele.panozzo@gmail.com>
-// 
-// This Source Code Form is subject to the terms of the Mozilla Public License 
-// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
-// obtain one at http://mozilla.org/MPL/2.0/.
-#ifndef IGL_SPARSE_ATA_FAST_H
-#define IGL_SPARSE_ATA_FAST_H
-#include "igl_inline.h"
-#define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
-#include <Eigen/Dense>
-#include <Eigen/Sparse>
-namespace igl
-{  
-  struct sparse_AtA_fast_data
-  {
-    // Weights
-    Eigen::VectorXd W;
-
-    // Flatten composition rules
-    std::vector<int> I_row;
-    std::vector<int> I_col;
-    std::vector<int> I_w;
-
-    // For each entry of AtA, points to the beginning
-    // of the composition rules
-    std::vector<int> I_outer;
-  };
-
-  IGL_INLINE void sparse_AtA_fast_precompute(
-    const Eigen::SparseMatrix<double>& A,
-    Eigen::SparseMatrix<double>& AtA,
-    sparse_AtA_fast_data& data);
-
-  IGL_INLINE void sparse_AtA_fast(
-    const Eigen::SparseMatrix<double>& A,
-    Eigen::SparseMatrix<double>& AtA,
-    const sparse_AtA_fast_data& data);
-  
-}
-
-#ifndef IGL_STATIC_LIBRARY
-#  include "sparse_AtA_fast.cpp"
-#endif
-
-#endif

include/igl/sparse_fast.cpp → include/igl/sparse_cached.cpp

@@ -5,7 +5,7 @@
 // This Source Code Form is subject to the terms of the Mozilla Public License
 // v. 2.0. If a copy of the MPL was not distributed with this file, You can
 // obtain one at http://mozilla.org/MPL/2.0/.
-#include "sparse_fast.h"
+#include "sparse_cached.h"
 
 #include <iostream>
 #include <vector>
@@ -14,10 +14,11 @@
 #include <map>
 #include <utility>
 
-IGL_INLINE void igl::sparse_fast_precompute(
-  const Eigen::VectorXi & I,
-  const Eigen::VectorXi & J,
-  Eigen::SparseMatrix<double>& X,
+template <typename DerivedI, typename Scalar>
+IGL_INLINE void igl::sparse_cached_precompute(
+  const Eigen::MatrixBase<DerivedI> & I,
+  const Eigen::MatrixBase<DerivedI> & J,
+  Eigen::SparseMatrix<Scalar>& X,
   Eigen::VectorXi& data)
 {
   // Generates the triplets
@@ -26,12 +27,13 @@ IGL_INLINE void igl::sparse_fast_precompute(
     t[i] = Eigen::Triplet<double>(I[i],J[i],1);
 
   // Call the triplets version
-  sparse_fast_precompute(t,X,data);
+  sparse_cached_precompute(t,X,data);
 }
 
-IGL_INLINE void igl::sparse_fast_precompute(
-  const std::vector<Eigen::Triplet<double> >& triplets,
-  Eigen::SparseMatrix<double>& X,
+template <typename Scalar>
+IGL_INLINE void igl::sparse_cached_precompute(
+  const std::vector<Eigen::Triplet<Scalar> >& triplets,
+  Eigen::SparseMatrix<Scalar>& X,
   Eigen::VectorXi& data)
 {
     // Construct an empty sparse matrix
@@ -79,9 +81,10 @@ IGL_INLINE void igl::sparse_fast_precompute(
 
 }
   
-IGL_INLINE void igl::sparse_fast(
-  const std::vector<Eigen::Triplet<double> >& triplets,
-  Eigen::SparseMatrix<double>& X,
+template <typename Scalar>
+IGL_INLINE void igl::sparse_cached(
+  const std::vector<Eigen::Triplet<Scalar> >& triplets,
+  Eigen::SparseMatrix<Scalar>& X,
   const Eigen::VectorXi& data)
 {
   assert(triplets.size() == data.size());
@@ -95,9 +98,10 @@ IGL_INLINE void igl::sparse_fast(
     *(X.valuePtr() + data[i]) += triplets[i].value();
 }
 
-IGL_INLINE void igl::sparse_fast(
-  const Eigen::VectorXd & V,
-  Eigen::SparseMatrix<double>& X,
+template <typename DerivedV, typename Scalar>
+IGL_INLINE void igl::sparse_cached(
+  const Eigen::MatrixBase<DerivedV>& V,
+  Eigen::SparseMatrix<Scalar>& X,
   const Eigen::VectorXi& data)
 {
   assert(V.size() == data.size());
@@ -113,4 +117,6 @@ IGL_INLINE void igl::sparse_fast(
 
 
 #ifdef IGL_STATIC_LIBRARY
+template void igl::sparse_cached_precompute<double>(std::vector<Eigen::Triplet<double, Eigen::SparseMatrix<double, 0, int>::Index>, std::allocator<Eigen::Triplet<double, Eigen::SparseMatrix<double, 0, int>::Index> > > const&, Eigen::SparseMatrix<double, 0, int>&, Eigen::Matrix<int, -1, 1, 0, -1, 1>&);
+template void igl::sparse_cached<double>(std::vector<Eigen::Triplet<double, Eigen::SparseMatrix<double, 0, int>::Index>, std::allocator<Eigen::Triplet<double, Eigen::SparseMatrix<double, 0, int>::Index> > > const&, Eigen::SparseMatrix<double, 0, int>&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&);
 #endif

include/igl/sparse_fast.h → include/igl/sparse_cached.h

@@ -5,8 +5,8 @@
 // This Source Code Form is subject to the terms of the Mozilla Public License 
 // v. 2.0. If a copy of the MPL was not distributed with this file, You can 
 // obtain one at http://mozilla.org/MPL/2.0/.
-#ifndef IGL_SPARSE_FAST_H
-#define IGL_SPARSE_FAST_H
+#ifndef IGL_SPARSE_CACHED_H
+#define IGL_SPARSE_CACHED_H
 #include "igl_inline.h"
 #define EIGEN_YES_I_KNOW_SPARSE_MODULE_IS_NOT_STABLE_YET
 #include <Eigen/Dense>
@@ -36,31 +36,47 @@ namespace igl
   // Outputs:
   //   X  m by n matrix of type T whose entries are to be found 
   //
-  IGL_INLINE void sparse_fast_precompute(
-    const Eigen::VectorXi & I,
-    const Eigen::VectorXi & J,
-    Eigen::SparseMatrix<double>& X,
+  // Example:
+  //   Eigen::SparseMatrix<double> A;
+  //   std::vector<Eigen::Triplet<double> > IJV;
+  //   buildA(IJV);
+  //   if (A.rows() == 0)
+  //   {
+  //     A = Eigen::SparseMatrix<double>(rows,cols);
+  //     igl::sparse_cached_precompute(IJV,A,A_data);
+  //   }
+  //   else
+  //     igl::sparse_cached(IJV,s.A,s.A_data);
+
+  template <typename DerivedI, typename Scalar>
+  IGL_INLINE void sparse_cached_precompute(
+    const Eigen::MatrixBase<DerivedI> & I,
+    const Eigen::MatrixBase<DerivedI> & J,
+    Eigen::SparseMatrix<Scalar>& X,
     Eigen::VectorXi& data);
   
-  IGL_INLINE void sparse_fast_precompute(
-    const std::vector<Eigen::Triplet<double> >& triplets,
-    Eigen::SparseMatrix<double>& X,
+  template <typename Scalar>
+  IGL_INLINE void sparse_cached_precompute(
+    const std::vector<Eigen::Triplet<Scalar> >& triplets,
+    Eigen::SparseMatrix<Scalar>& X,
     Eigen::VectorXi& data);
 
-  IGL_INLINE void sparse_fast(
-    const std::vector<Eigen::Triplet<double> >& triplets,
-    Eigen::SparseMatrix<double>& X,
+  template <typename Scalar>
+  IGL_INLINE void sparse_cached(
+    const std::vector<Eigen::Triplet<Scalar> >& triplets,
+    Eigen::SparseMatrix<Scalar>& X,
     const Eigen::VectorXi& data);
 
-  IGL_INLINE void sparse_fast(
-    const Eigen::VectorXd & V,
-    Eigen::SparseMatrix<double>& X,
+  template <typename DerivedV, typename Scalar>
+  IGL_INLINE void sparse_cached(
+    const Eigen::MatrixBase<DerivedV>& V,
+    Eigen::SparseMatrix<Scalar>& X,
     const Eigen::VectorXi& data);
   
 }
 
 #ifndef IGL_STATIC_LIBRARY
-#  include "sparse_fast.cpp"
+#  include "sparse_cached.cpp"
 #endif
 
 #endif

tutorial/710_SLIM/main.cpp

@@ -152,7 +152,7 @@ void deform_3d_demo_iter(igl::viewer::Viewer& viewer) {
     double soft_const_p = 1e5;
     sData.exp_factor = 5.0;
     slim_precompute(V,F,V_0,sData,igl::SLIMData::EXP_CONFORMAL,b,bc,soft_const_p);
-    cout << "precomputed" << endl;
+    //cout << "precomputed" << endl;
 
     first_iter = false;
     display_3d_mesh(viewer);