added cached version of slice

Former-commit-id: affe8d554e16c5550de99bf2edace7a2bb3e026b

include/igl/slice.cpp

@@ -350,6 +350,7 @@ template void igl::slice<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<in
 template Eigen::Matrix<double, -1, -1, 0, -1, -1> igl::slice<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&);
 template void igl::slice<Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int, Eigen::Matrix<double, -1, 3, 0, -1, 3>&);
 template void igl::slice<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, 3, 0, -1, 3> >(Eigen::DenseBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&);
+template void igl::slice<unsigned int, unsigned int>(Eigen::SparseMatrix<unsigned int, 0, int> const&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::SparseMatrix<unsigned int, 0, int>&);
 #ifdef WIN32
 template void igl::slice<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>>>(class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> const &, class Eigen::DenseBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> const &, int, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> &);
 template void igl::slice<class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>>>(class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>> const &, class Eigen::DenseBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> const &, int, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>> &);

include/igl/slice_cached.cpp

@@ -0,0 +1,55 @@
+// 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/.
+#include "slice_cached.h"
+
+#include <iostream>
+#include <vector>
+#include <utility>
+#include <igl/slice.h>
+
+template <typename TX, typename TY>
+IGL_INLINE void igl::slice_cached_precompute(
+  const Eigen::SparseMatrix<TX>& X,
+  const Eigen::Matrix<int,Eigen::Dynamic,1> & R,
+  const Eigen::Matrix<int,Eigen::Dynamic,1> & C,
+  Eigen::SparseMatrix<TY>& Y,
+  Eigen::VectorXi& data)
+{
+  // Create a sparse matrix whose entries are the ids
+  Eigen::SparseMatrix<unsigned> TS = X.template cast<unsigned>();
+
+  TS.makeCompressed();
+  for (unsigned i=0;i<TS.nonZeros();++i)
+    *(TS.valuePtr() + i) = i;
+
+  Eigen::SparseMatrix<unsigned> TS_sliced;
+  igl::slice(TS,R,C,TS_sliced);
+  Y = TS_sliced.cast<TY>();
+
+  data.resize(TS_sliced.nonZeros());
+  for (unsigned i=0;i<data.size();++i)
+  {
+    data[i] = *(TS_sliced.valuePtr() + i);
+    *(Y.valuePtr() + i) = *(X.valuePtr() + data[i]);
+  }
+}
+
+template <typename TX, typename TY>
+IGL_INLINE void igl::slice_cached(
+  const Eigen::SparseMatrix<TX>& X,
+  Eigen::SparseMatrix<TY>& Y,
+  const Eigen::VectorXi& data)
+{
+  for (unsigned i=0; i<data.size(); ++i)
+    *(Y.valuePtr() + i) = *(X.valuePtr() + data[i]);
+}
+
+#ifdef IGL_STATIC_LIBRARY
+template void igl::slice_cached_precompute<double, double>(Eigen::SparseMatrix<double, 0, int> const&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::SparseMatrix<double, 0, int>&, Eigen::Matrix<int, -1, 1, 0, -1, 1>&);
+template void igl::slice_cached<double, double>(Eigen::SparseMatrix<double, 0, int> const&, Eigen::SparseMatrix<double, 0, int>&, Eigen::Matrix<int, -1, 1, 0, -1, 1> const&);
+#endif

include/igl/slice_cached.h

@@ -0,0 +1,64 @@
+// 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_SLICE_CACHED_H
+#define IGL_SLICE_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
+{  
+
+  // Act like the matlab X(row_indices,col_indices) operator, where
+  // row_indices, col_indices are non-negative integer indices. This is a fast version
+  // of igl::slice that can analyze and store the sparsity structure. It is slower at the // first evaluation (slice_cached_precompute), but faster on the subsequent ones.
+  // 
+  // Inputs:
+  //   X  m by n matrix
+  //   R  list of row indices
+  //   C  list of column indices
+  //   
+  // Output:
+  //   Y  #R by #C matrix
+  //   data Temporary data used by slice_cached to repeat this operation
+  //
+  // Usage:
+  //
+  // // Construct and slice up Laplacian
+  // SparseMatrix<double> L,L_sliced;
+  // igl::cotmatrix(V,F,L);
+
+  // // Normal igl::slice call
+  // igl::slice(L,in,in,L_in_in);
+
+  // // Fast version
+  // static VectorXi data; // static or saved in a global state
+  // if (data.size() == 0)
+  //   igl::slice_cached_precompute(L,in,in,L_sliced,data);
+  // else
+  //   igl::slice_cached(L,L_sliced,temp);
+
+  template <typename TX, typename TY>
+  IGL_INLINE void slice_cached_precompute(
+    const Eigen::SparseMatrix<TX>& X,
+    const Eigen::Matrix<int,Eigen::Dynamic,1> & R,
+    const Eigen::Matrix<int,Eigen::Dynamic,1> & C,
+    Eigen::SparseMatrix<TY>& Y,
+    Eigen::VectorXi& data);
+
+  template <typename TX, typename TY>
+  IGL_INLINE void slice_cached(
+    const Eigen::SparseMatrix<TX>& X,
+    Eigen::SparseMatrix<TY>& Y,
+    const Eigen::VectorXi& data);
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "slice_cached.cpp"
+#endif
+}
+#endif