fixed types in biharmonic_coordinates.cpp

include/igl/biharmonic_coordinates.cpp

@@ -44,24 +44,28 @@ IGL_INLINE bool igl::biharmonic_coordinates(
 {
   using namespace Eigen;
   using namespace std;
+
+  typedef typename DerivedV::Scalar Scalar;
+  typedef typename DerivedT::Scalar Integer;
+
   // This is not the most efficient way to build A, but follows "Linear
   // Subspace Design for Real-Time Shape Deformation" [Wang et al. 2015].
-  SparseMatrix<double> A;
+  SparseMatrix<Scalar> A;
   {
-    DiagonalMatrix<double,Dynamic> Minv;
-    SparseMatrix<double> L,K;
+    DiagonalMatrix<Scalar, Dynamic> Minv;
+    SparseMatrix<Scalar> L, K;
     Array<bool,Dynamic,Dynamic> C;
     {
       Array<bool,Dynamic,1> I;
       on_boundary(T,I,C);
     }
-#ifdef false 
+#ifdef false
     // Version described in paper is "wrong"
     // http://www.cs.toronto.edu/~jacobson/images/error-in-linear-subspace-design-for-real-time-shape-deformation-2017-wang-et-al.pdf
-    SparseMatrix<double> N,Z,M;
+    SparseMatrix<Scalar> N, Z, M;
     normal_derivative(V,T,N);
     {
-      std::vector<Triplet<double> >ZIJV;
+      std::vector<Triplet<Scalar>> ZIJV;
       for(int t =0;t<T.rows();t++)
       {
         for(int f =0;f<T.cols();f++)
@@ -84,37 +88,37 @@ IGL_INLINE bool igl::biharmonic_coordinates(
     K = N+L;
     massmatrix(V,T,MASSMATRIX_TYPE_DEFAULT,M);
     // normalize
-    M /= ((VectorXd)M.diagonal()).array().abs().maxCoeff();
+    M /= ((Matrix<Scalar, Dynamic, 1>)M.diagonal()).array().abs().maxCoeff();
     Minv =
-      ((VectorXd)M.diagonal().array().inverse()).asDiagonal();
+      ((Matrix<Scalar, Dynamic, 1>)M.diagonal().array().inverse()).asDiagonal();
 #else
-    Eigen::SparseMatrix<double> M;
-    Eigen::MatrixXi E;
-    Eigen::VectorXi EMAP;
+    Eigen::SparseMatrix<Scalar> M;
+    Eigen::Matrix<Integer, Dynamic, Dynamic> E;
+    Eigen::Matrix<Integer, Dynamic, 1> EMAP;
     crouzeix_raviart_massmatrix(V,T,M,E,EMAP);
     crouzeix_raviart_cotmatrix(V,T,E,EMAP,L);
     // Ad  #E by #V facet-vertex incidence matrix
-    Eigen::SparseMatrix<double> Ad(E.rows(),V.rows());
+    Eigen::SparseMatrix<Scalar> Ad(E.rows(),V.rows());
     {
-      std::vector<Eigen::Triplet<double> > AIJV(E.size());
+      std::vector<Eigen::Triplet<Scalar>> AIJV(E.size());
       for(int e = 0;e<E.rows();e++)
       {
         for(int c = 0;c<E.cols();c++)
         {
-          AIJV[e+c*E.rows()] = Eigen::Triplet<double>(e,E(e,c),1);
+          AIJV[e + c * E.rows()] = Eigen::Triplet<Scalar>(e, E(e, c), 1);
         }
       }
       Ad.setFromTriplets(AIJV.begin(),AIJV.end());
     }
     // Degrees
-    Eigen::VectorXd De;
+    Eigen::Matrix<Scalar, Dynamic, 1> De;
     sum(Ad,2,De);
-    Eigen::DiagonalMatrix<double,Eigen::Dynamic> De_diag = 
+    Eigen::DiagonalMatrix<Scalar,Eigen::Dynamic> De_diag =
       De.array().inverse().matrix().asDiagonal();
     K = L*(De_diag*Ad);
     // normalize
-    M /= ((VectorXd)M.diagonal()).array().abs().maxCoeff();
-    Minv = ((VectorXd)M.diagonal().array().inverse()).asDiagonal();
+    M /= ((Matrix<Scalar, Dynamic, 1>)M.diagonal()).array().abs().maxCoeff();
+    Minv = ((Matrix<Scalar, Dynamic, 1>)M.diagonal().array().inverse()).asDiagonal();
     // kill boundary edges
     for(int f = 0;f<T.rows();f++)
     {
@@ -159,9 +163,9 @@ IGL_INLINE bool igl::biharmonic_coordinates(
   }
   const size_t dim = T.cols()-1;
   // Might as well be dense... I think...
-  MatrixXd J = MatrixXd::Zero(mp+mr,mp+r*(dim+1));
-  VectorXi b(mp+mr);
-  MatrixXd H(mp+r*(dim+1),dim);
+  Matrix<Scalar, Dynamic, Dynamic> J = Matrix<Scalar, Dynamic, Dynamic>::Zero(mp+mr,mp+r*(dim+1));
+  Matrix<Integer, Dynamic, 1> b(mp+mr);
+  Matrix<Scalar, Dynamic, Dynamic> H(mp+r*(dim+1),dim);
   {
     int v = 0;
     int c = 0;
@@ -194,7 +198,7 @@ IGL_INLINE bool igl::biharmonic_coordinates(
   // minimize    ½ W' A W'
   // subject to  W(b,:) = J
   return min_quad_with_fixed(
-    A,VectorXd::Zero(A.rows()).eval(),b,J,SparseMatrix<double>(),VectorXd(),true,W);
+    A,Matrix<Scalar, Dynamic, 1>::Zero(A.rows()).eval(),b,J,SparseMatrix<Scalar>(),Matrix<Scalar, Dynamic, 1>(),true,W);
 }
 
 #ifdef IGL_STATIC_LIBRARY