Merge branch 'master' of https://github.com/libigl/libigl

Former-commit-id: 02f522d28b8e7772adccba1d25a32d5b07254f3f

include/igl/barycentric_coordinates.cpp

@@ -106,4 +106,5 @@ template void igl::barycentric_coordinates<Eigen::Matrix<double, -1, -1, 0, -1,
 template void igl::barycentric_coordinates<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3> >(Eigen::MatrixBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> >&);
 template void igl::barycentric_coordinates<Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Matrix<double, 1, 3, 1, 1, 3> >(Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> >&);
 template void igl::barycentric_coordinates<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false>, Eigen::Matrix<double, 1, 3, 1, 1, 3> >(Eigen::MatrixBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::MatrixBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1> const, 1, -1, false> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> >&);
+template void igl::barycentric_coordinates<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/copyleft/cgal/outer_element.cpp

@@ -1,9 +1,9 @@
 // This file is part of libigl, a simple c++ geometry processing library.
-// 
+//
 // Copyright (C) 2015 Qingnan Zhou <qnzhou@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 
+//
+// 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 "outer_element.h"
 #include <iostream>
@@ -23,7 +23,7 @@ IGL_INLINE void igl::copyleft::cgal::outer_vertex(
         IndexType & v_index,
         Eigen::PlainObjectBase<DerivedA> & A)
 {
-    // Algorithm: 
+    // Algorithm:
     //    Find an outer vertex (i.e. vertex reachable from infinity)
     //    Return the vertex with the largest X value.
     //    If there is a tie, pick the one with largest Y value.
@@ -112,7 +112,7 @@ IGL_INLINE void igl::copyleft::cgal::outer_edge(
     const ScalarArray3& outer_v = V.row(outer_vid);
     assert(candidate_faces.size() > 0);
 
-    auto get_vertex_index = [&](const IndexArray3& f, Index vid) -> Index 
+    auto get_vertex_index = [&](const IndexArray3& f, Index vid) -> Index
     {
         if (f[0] == vid) return 0;
         if (f[1] == vid) return 1;
@@ -204,80 +204,14 @@ IGL_INLINE void igl::copyleft::cgal::outer_edge(
     std::copy(incident_faces.begin(), incident_faces.end(), A.data());
 }
 
-template<
-    typename DerivedV,
-    typename DerivedF,
-    typename DerivedN,
-    typename DerivedI,
-    typename IndexType
-    >
-IGL_INLINE void igl::copyleft::cgal::outer_facet(
-        const Eigen::PlainObjectBase<DerivedV> & V,
-        const Eigen::PlainObjectBase<DerivedF> & F,
-        const Eigen::PlainObjectBase<DerivedN> & N,
-        const Eigen::PlainObjectBase<DerivedI> & I,
-        IndexType & f,
-        bool & flipped) {
-    // Algorithm:
-    //    Find an outer edge.
-    //    Find the incident facet with the largest absolute X normal component.
-    //    If there is a tie, keep the one with positive X component.
-    //    If there is still a tie, pick the face with the larger signed index
-    //    (flipped face has negative index).
-    typedef typename DerivedV::Scalar Scalar;
-    typedef typename DerivedV::Index Index;
-    const size_t INVALID = std::numeric_limits<size_t>::max();
-
-    Index v1,v2;
-    Eigen::Matrix<Index,Eigen::Dynamic,1> incident_faces;
-    outer_edge(V, F, I, v1, v2, incident_faces);
-    assert(incident_faces.size() > 0);
-
-    auto generic_fabs = [&](const Scalar& val) -> const Scalar {
-        if (val >= 0) return val;
-        else return -val;
-    };
-
-    Scalar max_nx = 0;
-    size_t outer_fid = INVALID;
-    const size_t num_incident_faces = incident_faces.size();
-    for (size_t i=0; i<num_incident_faces; i++) 
-    {
-        const auto& fid = incident_faces(i);
-        const Scalar nx = N(fid, 0);
-        if (outer_fid == INVALID) {
-            max_nx = nx;
-            outer_fid = fid;
-        } else {
-            if (generic_fabs(nx) > generic_fabs(max_nx)) {
-                max_nx = nx;
-                outer_fid = fid;
-            } else if (nx == -max_nx && nx > 0) {
-                max_nx = nx;
-                outer_fid = fid;
-            } else if (nx == max_nx) {
-                if ((max_nx >= 0 && outer_fid < fid) ||
-                    (max_nx <  0 && outer_fid > fid)) {
-                    max_nx = nx;
-                    outer_fid = fid;
-                }
-            }
-        }
-    }
-
-    assert(outer_fid != INVALID);
-    f = outer_fid;
-    flipped = max_nx < 0;
-}
 
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 #include <CGAL/Exact_predicates_exact_constructions_kernel.h>
-template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, int&, bool&);
-template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, unsigned long>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, unsigned long&, bool&);
-template void igl::copyleft::cgal::outer_facet<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::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int&, bool&);
-template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int&, bool&);
 template void igl::copyleft::cgal::outer_edge<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, long, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, long&, long&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
 template void igl::copyleft::cgal::outer_edge<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, long, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, long&, long&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
+template void igl::copyleft::cgal::outer_edge<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, long, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, long&, long&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
+template void igl::copyleft::cgal::outer_edge<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, long, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, long&, long&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
+template void igl::copyleft::cgal::outer_edge<Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, long, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, long&, long&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
 #endif

include/igl/copyleft/cgal/outer_element.h

@@ -73,38 +73,6 @@ namespace igl
               IndexType & v2,
               Eigen::PlainObjectBase<DerivedA> & A);
 
-
-      // Find a facet that is reachable from infinity without crossing any faces.
-      // Such facet is called "outer facet."
-      //
-      // Precondition: The input mesh must have all self-intersection resolved.
-      // I.e there is no duplicated vertices, no overlapping edge and no
-      // intersecting faces (the only exception is there could be topologically
-      // duplicated faces).  See cgal::remesh_self_intersections.h for how to
-      // obtain such input.
-      //
-      // Inputs:
-      //   V  #V by 3 list of vertex positions
-      //   F  #F by 3 list of triangle indices into V
-      //   N  #N by 3 list of face normals
-      //   I  #I list of facets to consider
-      // Outputs:
-      //   f  Index of the outer facet.
-      //   flipped  true iff the normal of f points inwards.
-      template<
-          typename DerivedV,
-          typename DerivedF,
-          typename DerivedN,
-          typename DerivedI,
-          typename IndexType
-          >
-      IGL_INLINE void outer_facet(
-              const Eigen::PlainObjectBase<DerivedV> & V,
-              const Eigen::PlainObjectBase<DerivedF> & F,
-              const Eigen::PlainObjectBase<DerivedN> & N,
-              const Eigen::PlainObjectBase<DerivedI> & I,
-              IndexType & f,
-              bool & flipped);
     }
   }
 }

include/igl/copyleft/cgal/outer_facet.cpp

@@ -79,10 +79,82 @@ IGL_INLINE void igl::copyleft::cgal::outer_facet(
     flipped = adj_faces[order[0]] > 0;
 }
 
+
+
+template<
+    typename DerivedV,
+    typename DerivedF,
+    typename DerivedN,
+    typename DerivedI,
+    typename IndexType
+    >
+IGL_INLINE void igl::copyleft::cgal::outer_facet(
+        const Eigen::PlainObjectBase<DerivedV> & V,
+        const Eigen::PlainObjectBase<DerivedF> & F,
+        const Eigen::PlainObjectBase<DerivedN> & N,
+        const Eigen::PlainObjectBase<DerivedI> & I,
+        IndexType & f,
+        bool & flipped) {
+    // Algorithm:
+    //    Find an outer edge.
+    //    Find the incident facet with the largest absolute X normal component.
+    //    If there is a tie, keep the one with positive X component.
+    //    If there is still a tie, pick the face with the larger signed index
+    //    (flipped face has negative index).
+    typedef typename DerivedV::Scalar Scalar;
+    typedef typename DerivedV::Index Index;
+    const size_t INVALID = std::numeric_limits<size_t>::max();
+
+    Index v1,v2;
+    Eigen::Matrix<Index,Eigen::Dynamic,1> incident_faces;
+    outer_edge(V, F, I, v1, v2, incident_faces);
+    assert(incident_faces.size() > 0);
+
+    auto generic_fabs = [&](const Scalar& val) -> const Scalar {
+        if (val >= 0) return val;
+        else return -val;
+    };
+
+    Scalar max_nx = 0;
+    size_t outer_fid = INVALID;
+    const size_t num_incident_faces = incident_faces.size();
+    for (size_t i=0; i<num_incident_faces; i++) 
+    {
+        const auto& fid = incident_faces(i);
+        const Scalar nx = N(fid, 0);
+        if (outer_fid == INVALID) {
+            max_nx = nx;
+            outer_fid = fid;
+        } else {
+            if (generic_fabs(nx) > generic_fabs(max_nx)) {
+                max_nx = nx;
+                outer_fid = fid;
+            } else if (nx == -max_nx && nx > 0) {
+                max_nx = nx;
+                outer_fid = fid;
+            } else if (nx == max_nx) {
+                if ((max_nx >= 0 && outer_fid < fid) ||
+                    (max_nx <  0 && outer_fid > fid)) {
+                    max_nx = nx;
+                    outer_fid = fid;
+                }
+            }
+        }
+    }
+
+    assert(outer_fid != INVALID);
+    f = outer_fid;
+    flipped = max_nx < 0;
+}
+
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 #include <CGAL/Exact_predicates_exact_constructions_kernel.h>
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, int&, bool&);
+template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, -1, 1, -1, -1>, unsigned long>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 1, -1, -1> > const&, unsigned long&, bool&);
+template void igl::copyleft::cgal::outer_facet<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::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int&, bool&);
+template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int&, bool&);
 template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, int&, bool&);
 template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> > const&, int&, bool&);
 template void igl::copyleft::cgal::outer_facet<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int&, bool&);

include/igl/copyleft/cgal/outer_facet.h

@@ -47,6 +47,39 @@ namespace igl
                 const Eigen::PlainObjectBase<DerivedI> & I,
                 IndexType & f,
                 bool & flipped);
+
+      // Find a facet that is reachable from infinity without crossing any faces.
+      // Such facet is called "outer facet."
+      //
+      // Precondition: The input mesh must have all self-intersection resolved.
+      // I.e there is no duplicated vertices, no overlapping edge and no
+      // intersecting faces (the only exception is there could be topologically
+      // duplicated faces).  See cgal::remesh_self_intersections.h for how to
+      // obtain such input.
+      //
+      // Inputs:
+      //   V  #V by 3 list of vertex positions
+      //   F  #F by 3 list of triangle indices into V
+      //   N  #N by 3 list of face normals
+      //   I  #I list of facets to consider
+      // Outputs:
+      //   f  Index of the outer facet.
+      //   flipped  true iff the normal of f points inwards.
+      template<
+          typename DerivedV,
+          typename DerivedF,
+          typename DerivedN,
+          typename DerivedI,
+          typename IndexType
+          >
+      IGL_INLINE void outer_facet(
+              const Eigen::PlainObjectBase<DerivedV> & V,
+              const Eigen::PlainObjectBase<DerivedF> & F,
+              const Eigen::PlainObjectBase<DerivedN> & N,
+              const Eigen::PlainObjectBase<DerivedI> & I,
+              IndexType & f,
+              bool & flipped);
+
     }
 
   }

include/igl/in_element.h

@@ -11,6 +11,7 @@
 #include "igl_inline.h"
 #include "AABB.h"
 #include <Eigen/Core>
+#include <Eigen/Sparse>
 
 namespace igl
 {

include/igl/parula.cpp

@@ -66,6 +66,7 @@ IGL_INLINE void igl::parula(
 // generated by autoexplicit.sh
 template void igl::parula<Eigen::Array<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Array<int, -1, 1, 0, -1, 1> > const&, bool, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::parula<double>(double, double*);
+template void igl::parula<double>(double, double&, double&, double&);
 template void igl::parula<Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, bool, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::parula<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, bool, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/parula.h

@@ -44,8 +44,8 @@ namespace igl
     const bool normalize,
     Eigen::PlainObjectBase<DerivedC> & C);
   // Inputs:
-  //   min_z  value at blue
-  //   max_z  value at red
+  //   min_Z  value at blue
+  //   max_Z  value at red
   template <typename DerivedZ, typename DerivedC>
   IGL_INLINE void parula(
     const Eigen::PlainObjectBase<DerivedZ> & Z,

include/igl/per_edge_normals.cpp

@@ -116,4 +116,7 @@ IGL_INLINE void igl::per_edge_normals(
 // Explicit template specialization
 template void igl::per_edge_normals<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 template void igl::per_edge_normals<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::PerEdgeNormalsWeightingType, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+template void igl::per_edge_normals<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::PerEdgeNormalsWeightingType, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::per_edge_normals<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::PerEdgeNormalsWeightingType, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::per_edge_normals<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/per_edge_normals.h

@@ -25,7 +25,8 @@ namespace igl
   // Inputs:
   //   V  #V by 3 eigen Matrix of mesh vertex 3D positions
   //   F  #F by 3 eigen Matrix of face (triangle) indices
-  //   weighting  weighting type
+  //   weight  weighting type
+  //   FN  #F by 3 matrix of 3D face normals per face
   // Output:
   //   N  #2 by 3 matrix of mesh edge 3D normals per row
   //   E  #E by 2 matrix of edge indices per row

include/igl/point_mesh_squared_distance.cpp

@@ -1,9 +1,9 @@
 // This file is part of libigl, a simple c++ geometry processing library.
-// 
+//
 // Copyright (C) 2014 Alec Jacobson <alecjacobson@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 
+//
+// 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 "point_mesh_squared_distance.h"
 #include "AABB.h"
@@ -17,7 +17,7 @@ template <
 IGL_INLINE void igl::point_mesh_squared_distance(
   const Eigen::PlainObjectBase<DerivedP> & P,
   const Eigen::PlainObjectBase<DerivedV> & V,
-  const Eigen::MatrixXi & Ele, 
+  const Eigen::MatrixXi & Ele,
   Eigen::PlainObjectBase<DerivedsqrD> & sqrD,
   Eigen::PlainObjectBase<DerivedI> & I,
   Eigen::PlainObjectBase<DerivedC> & C)
@@ -46,5 +46,6 @@ IGL_INLINE void igl::point_mesh_squared_distance(
 }
 
 #ifdef IGL_STATIC_LIBRARY
+template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/point_mesh_squared_distance.h

@@ -10,6 +10,7 @@
 #include "igl_inline.h"
 #include <Eigen/Core>
 #include <vector>
+
 namespace igl
 {
   // Compute distances from a set of points P to a triangle mesh (V,F)
@@ -17,7 +18,7 @@ namespace igl
   // Inputs:
   //   P  #P by 3 list of query point positions
   //   V  #V by 3 list of vertex positions
-  //   F  #F by (3|2|1) list of (triangle|edge|point) indices
+  //   Ele  #Ele by (3|2|1) list of (triangle|edge|point) indices
   // Outputs:
   //   sqrD  #P list of smallest squared distances
   //   I  #P list of primitive indices corresponding to smallest distances

include/igl/polar_svd.cpp

@@ -75,4 +75,6 @@ template void igl::polar_svd<Eigen::Matrix<float,3,3,0,3,3>,Eigen::Matrix<float,
 template void igl::polar_svd<Eigen::Matrix<double, 3, 3, 0, 3, 3>, Eigen::Matrix<double, 3, 3, 0, 3, 3>, Eigen::Matrix<double, 3, 3, 0, 3, 3>, Eigen::Matrix<double, 3, 3, 0, 3, 3>, Eigen::Matrix<double, 3, 1, 0, 3, 1>, Eigen::Matrix<double, 3, 3, 0, 3, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 3, 0, 3, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 3, 0, 3, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 3, 0, 3, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 3, 0, 3, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 1, 0, 3, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 3, 0, 3, 3> >&);
 template void igl::polar_svd<Eigen::Matrix<float, 2, 2, 0, 2, 2>, Eigen::Matrix<float, 2, 2, 0, 2, 2>, Eigen::Matrix<float, 2, 2, 0, 2, 2>, Eigen::Matrix<float, 2, 2, 0, 2, 2>, Eigen::Matrix<float, 2, 1, 0, 2, 1>, Eigen::Matrix<float, 2, 2, 0, 2, 2> >(Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 2, 0, 2, 2> > const&, Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 1, 0, 2, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<float, 2, 2, 0, 2, 2> >&);
 template void igl::polar_svd<Eigen::Matrix<double, 2, 2, 0, 2, 2>, Eigen::Matrix<double, 2, 2, 0, 2, 2>, Eigen::Matrix<double, 2, 2, 0, 2, 2>, Eigen::Matrix<double, 2, 2, 0, 2, 2>, Eigen::Matrix<double, 2, 1, 0, 2, 1>, Eigen::Matrix<double, 2, 2, 0, 2, 2> >(Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 2, 0, 2, 2> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 2, 0, 2, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 2, 0, 2, 2> >&);
+template void igl::polar_svd<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
+template void igl::polar_svd<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/readMESH.cpp

@@ -470,4 +470,6 @@ template bool igl::readMESH<Eigen::Matrix<float, -1, 3, 1, -1, 3>, Eigen::Matrix
 template bool igl::readMESH<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
 // generated by autoexplicit.sh
 template bool igl::readMESH<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template bool igl::readMESH<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template bool igl::readMESH<double, int>(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > >&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&);
 #endif

include/igl/slice_tets.cpp

@@ -226,4 +226,5 @@ IGL_INLINE void igl::slice_tets(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 template void igl::slice_tets<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, 4, 1, 0, 4, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 4, 1, 0, 4, 1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&);
+template void igl::slice_tets<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, double>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::SparseMatrix<double, 0, int>&);
 #endif

include/igl/unique.cpp

@@ -267,4 +267,5 @@ template void igl::unique<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<i
 template void igl::unique<double>(std::vector<double, std::allocator<double> > const&, std::vector<double, std::allocator<double> >&);
 template void igl::unique<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 template void igl::unique<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::unique_rows<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/unique_simplices.cpp

@@ -63,4 +63,5 @@ template void igl::unique_simplices<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen
 template void igl::unique_simplices<Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 template void igl::unique_simplices<Eigen::Matrix<double, -1, 2, 0, -1, 2>, Eigen::Matrix<double, -1, 2, 0, -1, 2>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, Eigen::Matrix<long, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&);
 template void igl::unique_simplices<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 2, 0, -1, 2>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+template void igl::unique_simplices<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/writeMESH.cpp

@@ -8,7 +8,7 @@
 #include "writeMESH.h"
 
 #include "verbose.h"
-#include "matrix_to_list.h"
+#include "list_to_matrix.h"
 #include <Eigen/Core>
 
 #include <iostream>
@@ -140,4 +140,7 @@ IGL_INLINE bool igl::writeMESH(
 template bool igl::writeMESH<Eigen::Matrix<double, 8, 3, 0, 8, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, 12, 3, 0, 12, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, 8, 3, 0, 8, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, 12, 3, 0, 12, 3> > const&);
 //template bool igl::writeMESH<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&);
 template bool igl::writeMESH<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&);
+
+template bool igl::writeMESH<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&);
+template bool igl::writeMESH<double, int>(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&);
 #endif

python/py_doc.cpp

@@ -150,6 +150,21 @@ const char *__doc_igl_jet = R"igl_Qu8mg5v7(// JET like MATLAB's jet
   //   r  red value
   //   g  green value
   //   b  blue value)igl_Qu8mg5v7";
+const char *__doc_igl_cat = R"igl_Qu8mg5v7(// Perform concatenation of a two matrices along a single dimension
+  // If dim == 1, then C = [A;B]. If dim == 2 then C = [A B]
+  // 
+  // Template:
+  //   Scalar  scalar data type for sparse matrices like double or int
+  //   Mat  matrix type for all matrices (e.g. MatrixXd, SparseMatrix)
+  //   MatC  matrix type for ouput matrix (e.g. MatrixXd) needs to support
+  //     resize
+  // Inputs:
+  //   A  first input matrix
+  //   B  second input matrix
+  //   dim  dimension along which to concatenate, 0 or 1
+  // Outputs:
+  //   C  output matrix
+  //   )igl_Qu8mg5v7";
 const char *__doc_igl_eigs = R"igl_Qu8mg5v7(See eigs for the documentation.)igl_Qu8mg5v7";
 const char *__doc_igl_per_corner_normals = R"igl_Qu8mg5v7(// Compute vertex normals via vertex position list, face list
   // Inputs:
@@ -198,6 +213,29 @@ const char *__doc_igl_colon = R"igl_Qu8mg5v7(// Colon operator like matlab's col
   //     than hi, vice versa if hi<low
   // Output:
   //   I  list of values from low to hi with step size step)igl_Qu8mg5v7";
+const char *__doc_igl_fit_rotations = R"igl_Qu8mg5v7(// Known issues: This seems to be implemented in Eigen/Geometry:
+  // Eigen::umeyama
+  //
+  // FIT_ROTATIONS Given an input mesh and new positions find rotations for
+  // every covariance matrix in a stack of covariance matrices
+  // 
+  // Inputs:
+  //   S  nr*dim by dim stack of covariance matrices
+  //   single_precision  whether to use single precision (faster)
+  // Outputs:
+  //   R  dim by dim * nr list of rotations
+  //)igl_Qu8mg5v7";
+const char *__doc_igl_fit_rotations_planar = R"igl_Qu8mg5v7(// FIT_ROTATIONS Given an input mesh and new positions find 2D rotations for
+  // every vertex that best maps its one ring to the new one ring
+  // 
+  // Inputs:
+  //   S  nr*dim by dim stack of covariance matrices, third column and every
+  //   third row will be ignored
+  // Outputs:
+  //   R  dim by dim * nr list of rotations, third row and third column of each
+  //   rotation will just be identity
+  //)igl_Qu8mg5v7";
+const char *__doc_igl_fit_rotations_SSE = R"igl_Qu8mg5v7(See fit_rotations_SSE for the documentation.)igl_Qu8mg5v7";
 const char *__doc_igl_rotate_vectors = R"igl_Qu8mg5v7(// Rotate the vectors V by A radiants on the tangent plane spanned by B1 and
   // B2
   //
@@ -245,6 +283,17 @@ const char *__doc_igl_avg_edge_length = R"igl_Qu8mg5v7(// Compute the average ed
   //   l  average edge length
   //
   // See also: adjacency_matrix)igl_Qu8mg5v7";
+const char *__doc_igl_barycentric_coordinates = R"igl_Qu8mg5v7(// Compute barycentric coordinates in a tet
+  //
+  // Inputs:
+  //   P  #P by 3 Query points in 3d
+  //   A  #P by 3 Tet corners in 3d
+  //   B  #P by 3 Tet corners in 3d
+  //   C  #P by 3 Tet corners in 3d
+  //   D  #P by 3 Tet corners in 3d
+  // Outputs:
+  //   L  #P by 4 list of barycentric coordinates
+  //   )igl_Qu8mg5v7";
 const char *__doc_igl_lscm = R"igl_Qu8mg5v7(// Compute a Least-squares conformal map parametrization (equivalently
   // derived in "Intrinsic Parameterizations of Surface Meshes" [Desbrun et al.
   // 2002] and "Least Squares Conformal Maps for Automatic Texture Atlas
@@ -273,6 +322,41 @@ const char *__doc_igl_find_cross_field_singularities = R"igl_Qu8mg5v7(// Inputs:
   //   isSingularity    #V by 1 boolean eigen Vector indicating the presence of a singularity on a vertex
   //   singularityIndex #V by 1 integer eigen Vector containing the singularity indices
   //)igl_Qu8mg5v7";
+const char *__doc_igl_upsample = R"igl_Qu8mg5v7(// Subdivide a mesh without moving vertices: loop subdivision but odd
+  // vertices stay put and even vertices are just edge midpoints
+  // 
+  // Templates:
+  //   MatV  matrix for vertex positions, e.g. MatrixXd
+  //   MatF  matrix for vertex positions, e.g. MatrixXi
+  // Inputs:
+  //   V  #V by dim  mesh vertices
+  //   F  #F by 3  mesh triangles
+  // Outputs:
+  //   NV new vertex positions, V is guaranteed to be at top
+  //   NF new list of face indices
+  //
+  // NOTE: V should not be the same as NV,
+  // NOTE: F should not be the same as NF, use other proto
+  //
+  // Known issues:
+  //   - assumes (V,F) is edge-manifold.)igl_Qu8mg5v7";
+const char *__doc_igl_point_mesh_squared_distance = R"igl_Qu8mg5v7(// Compute distances from a set of points P to a triangle mesh (V,F)
+  //
+  // Inputs:
+  //   P  #P by 3 list of query point positions
+  //   V  #V by 3 list of vertex positions
+  //   Ele  #Ele by (3|2|1) list of (triangle|edge|point) indices
+  // Outputs:
+  //   sqrD  #P list of smallest squared distances
+  //   I  #P list of primitive indices corresponding to smallest distances
+  //   C  #P by 3 list of closest points
+  //
+  // Known bugs: This only computes distances to given primitivess. So
+  // unreferenced vertices are ignored. However, degenerate primitives are
+  // handled correctly: triangle [1 2 2] is treated as a segment [1 2], and
+  // triangle [1 1 1] is treated as a point. So one _could_ add extra
+  // combinatorially degenerate rows to Ele for all unreferenced vertices to
+  // also get distances to points.)igl_Qu8mg5v7";
 const char *__doc_igl_parula = R"igl_Qu8mg5v7(// PARULA like MATLAB's parula
   //
   // Inputs:
@@ -361,6 +445,27 @@ const char *__doc_igl_active_set = R"igl_Qu8mg5v7(// Known Bugs: rows of [Aeq;Ai
   // Benchmark: For a harmonic solve on a mesh with 325K facets, matlab 2.2
   // secs, igl/min_quad_with_fixed.h 7.1 secs
   //)igl_Qu8mg5v7";
+const char *__doc_igl_per_edge_normals = R"igl_Qu8mg5v7(// Compute face normals via vertex position list, face list
+  // Inputs:
+  //   V  #V by 3 eigen Matrix of mesh vertex 3D positions
+  //   F  #F by 3 eigen Matrix of face (triangle) indices
+  //   weight  weighting type
+  //   FN  #F by 3 matrix of 3D face normals per face
+  // Output:
+  //   N  #2 by 3 matrix of mesh edge 3D normals per row
+  //   E  #E by 2 matrix of edge indices per row
+  //   EMAP  #E by 1 matrix of indices from all edges to E
+  //)igl_Qu8mg5v7";
+const char *__doc_igl_covariance_scatter_matrix = R"igl_Qu8mg5v7(// Construct the covariance scatter matrix for a given arap energy
+  // Inputs:
+  //   V  #V by Vdim list of initial domain positions
+  //   F  #F by 3 list of triangle indices into V
+  //   energy  ARAPEnergyType enum value defining which energy is being used.
+  //     See ARAPEnergyType.h for valid options and explanations.
+  // Outputs:
+  //   CSM dim*#V/#F by dim*#V sparse matrix containing special laplacians along
+  //     the diagonal so that when multiplied by V gives covariance matrix
+  //     elements, can be used to speed up covariance matrix computation)igl_Qu8mg5v7";
 const char *__doc_igl_boundary_facets = R"igl_Qu8mg5v7(// BOUNDARY_FACETS Determine boundary faces (edges) of tetrahedra (triangles)
   // stored in T (analogous to qptoolbox's `outline` and `boundary_faces`).
   //
@@ -385,6 +490,27 @@ const char *__doc_igl_compute_frame_field_bisectors = R"igl_Qu8mg5v7(// Compute
   //   BIS1  #F by 3 eigen Matrix of the first per face frame field bisector
   //   BIS2  #F by 3 eigen Matrix of the second per face frame field bisector
   //)igl_Qu8mg5v7";
+const char *__doc_igl_edge_lengths = R"igl_Qu8mg5v7(// Constructs a list of lengths of edges opposite each index in a face
+  // (triangle/tet) list
+  //
+  // Templates:
+  //   DerivedV derived from vertex positions matrix type: i.e. MatrixXd
+  //   DerivedF derived from face indices matrix type: i.e. MatrixXi
+  //   DerivedL derived from edge lengths matrix type: i.e. MatrixXd
+  // Inputs:
+  //   V  eigen matrix #V by 3
+  //   F  #F by 2 list of mesh edges
+  //    or
+  //   F  #F by 3 list of mesh faces (must be triangles)
+  //    or
+  //   T  #T by 4 list of mesh elements (must be tets)
+  // Outputs:
+  //   L  #F by {1|3|6} list of edge lengths 
+  //     for edges, column of lengths
+  //     for triangles, columns correspond to edges [1,2],[2,0],[0,1]
+  //     for tets, columns correspond to edges
+  //     [3 0],[3 1],[3 2],[1 2],[2 0],[0 1]
+  //)igl_Qu8mg5v7";
 const char *__doc_igl_readOBJ = R"igl_Qu8mg5v7(// Read a mesh from an ascii obj file, filling in vertex positions, normals
   // and texture coordinates. Mesh may have faces of any number of degree
   //
@@ -487,6 +613,18 @@ const char *__doc_igl_min_quad_with_fixed_solve = R"igl_Qu8mg5v7(// Solves a sys
   //   sol  #unknowns+#lagrange by cols solution to linear system
   // Returns true on success, false on error)igl_Qu8mg5v7";
 const char *__doc_igl_min_quad_with_fixed = R"igl_Qu8mg5v7(See min_quad_with_fixed for the documentation.)igl_Qu8mg5v7";
+const char *__doc_igl_writeMESH = R"igl_Qu8mg5v7(// save a tetrahedral volume mesh to a .mesh file
+  //
+  // Templates:
+  //   Scalar  type for positions and vectors (will be cast as double)
+  //   Index  type for indices (will be cast to int)
+  // Input:
+  //   mesh_file_name  path of .mesh file
+  //   V  double matrix of vertex positions  #V by 3
+  //   T  #T list of tet indices into vertex positions
+  //   F  #F list of face indices into vertex positions
+  //
+  // Known bugs: Holes and regions are not supported)igl_Qu8mg5v7";
 const char *__doc_igl_unique = R"igl_Qu8mg5v7(// Act like matlab's [C,IA,IC] = unique(X)
   //
   // Templates:
@@ -550,6 +688,23 @@ const char *__doc_igl_slice_into = R"igl_Qu8mg5v7(// Act like the matlab Y(row_i
   //   Y  ym by yn lhs matrix
   // Output:
   //   Y  ym by yn lhs matrix, same as input but Y(R,C) = X)igl_Qu8mg5v7";
+const char *__doc_igl_slice_tets = R"igl_Qu8mg5v7(// SLICE_TETS Slice through a tet mesh (V,T) along a given plane (via its
+  // implicit equation).
+  //
+  // Inputs:
+  //   V  #V by 3 list of tet mesh vertices
+  //   T  #T by 4 list of tet indices into V 
+  //   plane  list of 4 coefficients in the plane equation: [x y z 1]'*plane = 0
+  //   Optional:
+  //     'Manifold' followed by whether to stitch together triangles into a
+  //       manifold mesh {true}: results in more compact U but slightly slower.
+  // Outputs:
+  //   U  #U by 3 list of triangle mesh vertices along slice
+  //   G  #G by 3 list of triangles indices into U
+  //   J  #G list of indices into T revealing from which tet each faces comes
+  //   BC  #U by #V list of barycentric coordinates (or more generally: linear
+  //     interpolation coordinates) so that U = BC*V
+  // )igl_Qu8mg5v7";
 const char *__doc_igl_n_polyvector = R"igl_Qu8mg5v7(// Inputs:
   //   v0, v1         the two #3 by 1 vectors
   //   normalized     boolean, if false, then the vectors are normalized prior to the calculation
@@ -577,6 +732,19 @@ const char *__doc_igl_boundary_loop = R"igl_Qu8mg5v7(// Compute list of ordered
   // Outputs:
   //   L  list of loops where L[i] = ordered list of boundary vertices in loop i
   //)igl_Qu8mg5v7";
+const char *__doc_igl_polar_svd = R"igl_Qu8mg5v7(// Computes the polar decomposition (R,T) of a matrix A using SVD singular
+  // value decomposition
+  //
+  // Inputs:
+  //   A  3 by 3 matrix to be decomposed
+  // Outputs:
+  //   R  3 by 3 rotation matrix part of decomposition (**always rotataion**)
+  //   T  3 by 3 stretch matrix part of decomposition
+  //   U  3 by 3 left-singular vectors
+  //   S  3 by 1 singular values
+  //   V  3 by 3 right-singular vectors
+  //
+  //)igl_Qu8mg5v7";
 const char *__doc_igl_comb_cross_field = R"igl_Qu8mg5v7(// Inputs:
   //   V          #V by 3 eigen Matrix of mesh vertex 3D positions
   //   F          #F by 4 eigen Matrix of face (quad) indices
@@ -592,6 +760,20 @@ const char *__doc_igl_invert_diag = R"igl_Qu8mg5v7(// Templates:
   //   X  an m by n sparse matrix
   // Outputs:
   //   Y  an m by n sparse matrix)igl_Qu8mg5v7";
+const char *__doc_igl_readMESH = R"igl_Qu8mg5v7(// load a tetrahedral volume mesh from a .mesh file
+  //
+  // Templates:
+  //   Scalar  type for positions and vectors (will be read as double and cast
+  //     to Scalar)
+  //   Index  type for indices (will be read as int and cast to Index)
+  // Input:
+  //   mesh_file_name  path of .mesh file
+  // Outputs:
+  //   V  double matrix of vertex positions  #V by 3
+  //   T  #T list of tet indices into vertex positions
+  //   F  #F list of face indices into vertex positions
+  //
+  // Known bugs: Holes and regions are not supported)igl_Qu8mg5v7";
 const char *__doc_igl_copyleft_comiso_miq = R"igl_Qu8mg5v7(// Inputs:
     //   V              #V by 3 list of mesh vertex 3D positions
     //   F              #F by 3 list of faces indices in V

python/py_doc.h

@@ -10,24 +10,34 @@ extern const char *__doc_igl_per_vertex_normals;
 extern const char *__doc_igl_sortrows;
 extern const char *__doc_igl_barycenter;
 extern const char *__doc_igl_jet;
+extern const char *__doc_igl_cat;
 extern const char *__doc_igl_eigs;
 extern const char *__doc_igl_per_corner_normals;
 extern const char *__doc_igl_massmatrix;
 extern const char *__doc_igl_colon;
+extern const char *__doc_igl_fit_rotations;
+extern const char *__doc_igl_fit_rotations_planar;
+extern const char *__doc_igl_fit_rotations_SSE;
 extern const char *__doc_igl_rotate_vectors;
 extern const char *__doc_igl_read_triangle_mesh;
 extern const char *__doc_igl_gaussian_curvature;
 extern const char *__doc_igl_avg_edge_length;
+extern const char *__doc_igl_barycentric_coordinates;
 extern const char *__doc_igl_lscm;
 extern const char *__doc_igl_find_cross_field_singularities;
+extern const char *__doc_igl_upsample;
+extern const char *__doc_igl_point_mesh_squared_distance;
 extern const char *__doc_igl_parula;
 extern const char *__doc_igl_setdiff;
 extern const char *__doc_igl_comb_frame_field;
 extern const char *__doc_igl_map_vertices_to_circle;
 extern const char *__doc_igl_writeOBJ;
 extern const char *__doc_igl_active_set;
+extern const char *__doc_igl_per_edge_normals;
+extern const char *__doc_igl_covariance_scatter_matrix;
 extern const char *__doc_igl_boundary_facets;
 extern const char *__doc_igl_compute_frame_field_bisectors;
+extern const char *__doc_igl_edge_lengths;
 extern const char *__doc_igl_readOBJ;
 extern const char *__doc_igl_cut_mesh_from_singularities;
 extern const char *__doc_igl_readDMAT;
@@ -37,6 +47,7 @@ extern const char *__doc_igl_doublearea_quad;
 extern const char *__doc_igl_min_quad_with_fixed_precompute;
 extern const char *__doc_igl_min_quad_with_fixed_solve;
 extern const char *__doc_igl_min_quad_with_fixed;
+extern const char *__doc_igl_writeMESH;
 extern const char *__doc_igl_unique;
 extern const char *__doc_igl_unique_rows;
 extern const char *__doc_igl_arap_precomputation;
@@ -44,10 +55,13 @@ extern const char *__doc_igl_arap_solve;
 extern const char *__doc_igl_cross_field_missmatch;
 extern const char *__doc_igl_grad;
 extern const char *__doc_igl_slice_into;
+extern const char *__doc_igl_slice_tets;
 extern const char *__doc_igl_n_polyvector;
 extern const char *__doc_igl_harmonic;
 extern const char *__doc_igl_boundary_loop;
+extern const char *__doc_igl_polar_svd;
 extern const char *__doc_igl_comb_cross_field;
 extern const char *__doc_igl_invert_diag;
+extern const char *__doc_igl_readMESH;
 extern const char *__doc_igl_copyleft_comiso_miq;
 extern const char *__doc_igl_copyleft_comiso_nrosy;

python/py_igl.cpp

@@ -49,6 +49,20 @@
 #include <igl/comb_frame_field.h>
 #include <igl/n_polyvector.h>
 
+#include <igl/point_mesh_squared_distance.h>
+#include <igl/AABB.h>
+#include <igl/readMESH.h>
+#include <igl/writeMESH.h>
+#include <igl/slice_tets.h>
+#include <igl/edge_lengths.h>
+#include <igl/upsample.h>
+#include <igl/cat.h>
+#include <igl/per_edge_normals.h>
+#include <igl/barycentric_coordinates.h>
+#include <igl/fit_rotations.h>
+#include <igl/polar_svd.h>
+#include <igl/covariance_scatter_matrix.h>
+
 
 void python_export_igl(py::module &m)
 {
@@ -99,4 +113,18 @@ void python_export_igl(py::module &m)
 #include "py_igl/py_comb_frame_field.cpp"
 #include "py_igl/py_n_polyvector.cpp"
 
+#include "py_igl/py_point_mesh_squared_distance.cpp"
+#include "py_igl/py_AABB.cpp"
+#include "py_igl/py_readMESH.cpp"
+#include "py_igl/py_writeMESH.cpp"
+#include "py_igl/py_slice_tets.cpp"
+#include "py_igl/py_edge_lengths.cpp"
+#include "py_igl/py_upsample.cpp"
+#include "py_igl/py_cat.cpp"
+#include "py_igl/py_per_edge_normals.cpp"
+#include "py_igl/py_barycentric_coordinates.cpp"
+#include "py_igl/py_fit_rotations.cpp"
+#include "py_igl/py_polar_svd.cpp"
+#include "py_igl/py_covariance_scatter_matrix.cpp"
+
 }

python/py_igl/py_AABB.cpp

@@ -0,0 +1,14 @@
+py::class_<igl::AABB<Eigen::MatrixXd,3> > AABB(m, "AABB");
+
+AABB
+.def(py::init<>())
+.def(py::init<const igl::AABB<Eigen::MatrixXd,3>& >())
+.def("init",[](igl::AABB<Eigen::MatrixXd,3>& tree, const Eigen::MatrixXd& V, const Eigen::MatrixXi& Ele)
+{
+    return tree.init(V, Ele, Eigen::Matrix<double, Eigen::Dynamic, 3>(), Eigen::Matrix<double, Eigen::Dynamic, 3>(), Eigen::MatrixXi(), 0); 
+})
+.def("squared_distance", [](const igl::AABB<Eigen::MatrixXd,3>& tree, const Eigen::MatrixXd& V, const Eigen::MatrixXi& Ele, const Eigen::MatrixXd& P, Eigen::MatrixXd& sqrD, Eigen::MatrixXi& I, Eigen::MatrixXd& C)
+{
+    return tree.squared_distance(V, Ele, P, sqrD, I, C);
+})
+;

python/py_igl/py_barycentric_coordinates.cpp

@@ -0,0 +1,29 @@
+
+
+m.def("barycentric_coordinates", []
+(
+  const Eigen::MatrixXd& P,
+  const Eigen::MatrixXd& A,
+  const Eigen::MatrixXd& B,
+  const Eigen::MatrixXd& C,
+  const Eigen::MatrixXd& D,
+  Eigen::MatrixXd& L
+)
+{
+  return igl::barycentric_coordinates(P, A, B, C, D, L);
+}, __doc_igl_barycentric_coordinates,
+py::arg("P"), py::arg("A"), py::arg("B"), py::arg("C"), py::arg("D"), py::arg("L"));
+
+m.def("barycentric_coordinates", []
+(
+  const Eigen::MatrixXd& P,
+  const Eigen::MatrixXd& A,
+  const Eigen::MatrixXd& B,
+  const Eigen::MatrixXd& C,
+  Eigen::MatrixXd& L
+)
+{
+  return igl::barycentric_coordinates(P, A, B, C, L);
+}, __doc_igl_barycentric_coordinates,
+py::arg("P"), py::arg("A"), py::arg("B"), py::arg("C"), py::arg("L"));
+

python/py_igl/py_cat.cpp

@@ -0,0 +1,56 @@
+m.def("cat", []
+(
+  const int dim,
+  const Eigen::MatrixXd& A,
+  const Eigen::MatrixXd& B,
+  Eigen::MatrixXd& C
+)
+{
+  return igl::cat(dim, A, B, C);
+}, __doc_igl_cat,
+py::arg("dim"), py::arg("A"), py::arg("B"), py::arg("C"));
+
+m.def("cat", []
+(
+  const int dim,
+  Eigen::MatrixXd& A,
+  Eigen::MatrixXd& B
+)
+{
+  return igl::cat(dim, A, B);
+}, __doc_igl_cat,
+py::arg("dim"), py::arg("A"), py::arg("B"));
+
+m.def("cat", []
+(
+  const int dim,
+  Eigen::MatrixXi& A,
+  Eigen::MatrixXi& B
+)
+{
+  return igl::cat(dim, A, B);
+}, __doc_igl_cat,
+py::arg("dim"), py::arg("A"), py::arg("B"));
+
+//m.def("cat", []
+//(
+//  const std::vector<std::vector<Eigen::MatrixXd > > & A, 
+//  Eigen::MatrixXd & C
+//)
+//{
+//  return igl::cat(A, C);
+//}, __doc_igl_cat,
+//py::arg("A"), py::arg("C"));
+
+m.def("cat", []
+(
+  const int dim,
+  const Eigen::SparseMatrix<double>& A,
+  const Eigen::SparseMatrix<double>& B,
+  Eigen::SparseMatrix<double>& C
+)
+{
+  return igl::cat(dim, A, B, C);
+}, __doc_igl_cat,
+py::arg("dim"), py::arg("A"), py::arg("B"), py::arg("C"));
+

python/py_igl/py_covariance_scatter_matrix.cpp

@@ -0,0 +1,11 @@
+m.def("covariance_scatter_matrix", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const igl::ARAPEnergyType energy,
+  Eigen::SparseMatrix<double>& CSM
+)
+{
+  return igl::covariance_scatter_matrix(V,F,energy,CSM);
+}, __doc_igl_covariance_scatter_matrix,
+py::arg("V"), py::arg("F"), py::arg("energy"), py::arg("CSM"));

python/py_igl/py_edge_lengths.cpp

@@ -0,0 +1,11 @@
+m.def("edge_lengths", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  Eigen::MatrixXd& L
+)
+{
+  return igl::edge_lengths(V, F, L);
+}, __doc_igl_edge_lengths,
+py::arg("V"), py::arg("F"), py::arg("L"));
+

python/py_igl/py_fit_rotations.cpp

@@ -0,0 +1,22 @@
+m.def("fit_rotations", []
+(
+  const Eigen::MatrixXd& S,
+  const bool single_precision,
+  Eigen::MatrixXd& R
+)
+{
+  return igl::fit_rotations(S, single_precision, R);
+}, __doc_igl_fit_rotations,
+py::arg("S"), py::arg("single_precision"), py::arg("R"));
+
+
+m.def("fit_rotations_planar", []
+(
+  const Eigen::MatrixXd& S,
+  Eigen::MatrixXd& R
+)
+{
+  return igl::fit_rotations_planar(S, R);
+}, __doc_igl_fit_rotations_planar,
+py::arg("S"), py::arg("R"));
+

python/py_igl/py_parula.cpp

@@ -1,3 +1,25 @@
+//m.def("parula", []
+//(
+//  const double f,
+//  T * rgb
+//)
+//{
+//  return igl::parula(f, rgb);
+//}, __doc_igl_parula,
+//py::arg("f"), py::arg("rgb"));
+
+m.def("parula", []
+(
+  const double f,
+  double & r,
+  double & g,
+  double & b
+)
+{
+  return igl::parula(f, r, g, b);
+}, __doc_igl_parula,
+py::arg("f"), py::arg("r"), py::arg("g"), py::arg("b"));
+
 m.def("parula", []
 (
   const Eigen::MatrixXd& Z,
@@ -6,7 +28,7 @@ m.def("parula", []
 )
 {
   assert_is_VectorX("Z",Z);
-  return igl::parula(Z,normalize,C);
+  return igl::parula(Z, normalize, C);
 }, __doc_igl_parula,
 py::arg("Z"), py::arg("normalize"), py::arg("C"));
 
@@ -19,6 +41,7 @@ m.def("parula", []
 )
 {
   assert_is_VectorX("Z",Z);
-  return igl::parula(Z,min_Z,max_Z,C);
+  return igl::parula(Z, min_Z, max_Z, C);
 }, __doc_igl_parula,
 py::arg("Z"), py::arg("min_Z"), py::arg("max_Z"), py::arg("C"));
+

python/py_igl/py_per_edge_normals.cpp

@@ -0,0 +1,50 @@
+py::enum_<igl::PerEdgeNormalsWeightingType>(m, "PerEdgeNormalsWeightingType")
+    .value("PER_EDGE_NORMALS_WEIGHTING_TYPE_UNIFORM", igl::PER_EDGE_NORMALS_WEIGHTING_TYPE_UNIFORM)
+    .value("PER_EDGE_NORMALS_WEIGHTING_TYPE_AREA", igl::PER_EDGE_NORMALS_WEIGHTING_TYPE_AREA)
+    .value("PER_EDGE_NORMALS_WEIGHTING_TYPE_DEFAULT", igl::PER_EDGE_NORMALS_WEIGHTING_TYPE_DEFAULT)
+    .value("NUM_PER_EDGE_NORMALS_WEIGHTING_TYPE", igl::NUM_PER_EDGE_NORMALS_WEIGHTING_TYPE)
+    .export_values();
+
+
+m.def("per_edge_normals", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const igl::PerEdgeNormalsWeightingType weight,
+  const Eigen::MatrixXd& FN,
+  Eigen::MatrixXd& N,
+  Eigen::MatrixXi& E,
+  Eigen::MatrixXi& EMAP
+)
+{
+  return igl::per_edge_normals(V, F, weight, FN, N, E, EMAP);
+}, __doc_igl_per_edge_normals,
+py::arg("V"), py::arg("F"), py::arg("weight"), py::arg("FN"), py::arg("N"), py::arg("E"), py::arg("EMAP"));
+
+m.def("per_edge_normals", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const igl::PerEdgeNormalsWeightingType weight,
+  Eigen::MatrixXd& N,
+  Eigen::MatrixXi& E,
+  Eigen::MatrixXi& EMAP
+)
+{
+  return igl::per_edge_normals(V, F, weight, N, E, EMAP);
+}, __doc_igl_per_edge_normals,
+py::arg("V"), py::arg("F"), py::arg("weight"), py::arg("N"), py::arg("E"), py::arg("EMAP"));
+
+m.def("per_edge_normals", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  Eigen::MatrixXd& N,
+  Eigen::MatrixXi& E,
+  Eigen::MatrixXi& EMAP
+)
+{
+  return igl::per_edge_normals(V, F, N, E, EMAP);
+}, __doc_igl_per_edge_normals,
+py::arg("V"), py::arg("F"), py::arg("N"), py::arg("E"), py::arg("EMAP"));
+

python/py_igl/py_point_mesh_squared_distance.cpp

@@ -0,0 +1,15 @@
+m.def("point_mesh_squared_distance", []
+(
+  const Eigen::MatrixXd& P,
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& Ele,
+  Eigen::MatrixXd& sqrD,
+  Eigen::MatrixXi& I,
+  Eigen::MatrixXd& C
+)
+{
+  assert_is_VectorX("I",I);
+  return igl::point_mesh_squared_distance(P, V, Ele, sqrD, I, C);
+}, __doc_igl_point_mesh_squared_distance,
+py::arg("P"), py::arg("V"), py::arg("Ele"), py::arg("sqrD"), py::arg("I"), py::arg("C"));
+

python/py_igl/py_polar_svd.cpp

@@ -0,0 +1,26 @@
+m.def("polar_svd", []
+(
+  const Eigen::MatrixXd& A,
+  Eigen::MatrixXd& R,
+  Eigen::MatrixXd& T,
+  Eigen::MatrixXd& U,
+  Eigen::MatrixXd& S,
+  Eigen::MatrixXd& V
+)
+{
+  return igl::polar_svd(A, R, T, U, S, V);
+}, __doc_igl_polar_svd,
+py::arg("A"), py::arg("R"), py::arg("T"), py::arg("U"), py::arg("S"), py::arg("V"));
+
+
+m.def("polar_svd", []
+(
+  const Eigen::MatrixXd& A,
+  Eigen::MatrixXd& R,
+  Eigen::MatrixXd& T
+)
+{
+  return igl::polar_svd(A, R, T);
+}, __doc_igl_polar_svd,
+py::arg("A"), py::arg("R"), py::arg("T"));
+

python/py_igl/py_readMESH.cpp

@@ -0,0 +1,24 @@
+m.def("readMESH", []
+(
+  const std::string mesh_file_name,
+  std::vector<std::vector<double> > & V,
+  std::vector<std::vector<int> > & T,
+  std::vector<std::vector<int> > & F
+)
+{
+  return igl::readMESH(mesh_file_name, V, T, F);
+}, __doc_igl_readMESH,
+py::arg("mesh_file_name"), py::arg("V"), py::arg("T"), py::arg("F"));
+
+m.def("readMESH", []
+(
+  const std::string mesh_file_name,
+  Eigen::MatrixXd& V,
+  Eigen::MatrixXi& T,
+  Eigen::MatrixXi& F
+)
+{
+  return igl::readMESH(mesh_file_name, V, T, F);
+}, __doc_igl_readMESH,
+py::arg("mesh_file_name"), py::arg("V"), py::arg("T"), py::arg("F"));
+

python/py_igl/py_slice_tets.cpp

@@ -0,0 +1,23 @@
+m.def("slice_tets", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& T,
+  const Eigen::MatrixXd& plane,
+  Eigen::MatrixXd& U,
+  Eigen::MatrixXi& G,
+  Eigen::MatrixXi& J,
+  Eigen::SparseMatrix<double> & BC
+)
+{
+  assert_is_VectorX("plane", plane);
+  Eigen::VectorXd pl;
+  if (plane.size() != 0)
+    pl = plane;
+  assert_is_VectorX("J", J);
+  Eigen::VectorXi Jv;
+  if (J.size() != 0)
+    Jv = J;
+  return igl::slice_tets(V, T, pl, U, G, Jv, BC);
+}, __doc_igl_slice_tets,
+py::arg("V"), py::arg("T"), py::arg("plane"), py::arg("U"), py::arg("G"), py::arg("J"), py::arg("BC"));
+

python/py_igl/py_upsample.cpp

@@ -0,0 +1,22 @@
+m.def("upsample", []
+(
+  Eigen::MatrixXd& V,
+  Eigen::MatrixXi& F
+)
+{
+  return igl::upsample(V, F);
+}, __doc_igl_upsample,
+py::arg("V"), py::arg("F"));
+
+m.def("upsample", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  Eigen::MatrixXd& NV,
+  Eigen::MatrixXi& NF
+)
+{
+  return igl::upsample(V, F, NV, NF);
+}, __doc_igl_upsample,
+py::arg("V"), py::arg("F"), py::arg("NV"), py::arg("NF"));
+

python/py_igl/py_writeMESH.cpp

@@ -0,0 +1,24 @@
+m.def("writeMESH", []
+(
+  const std::string mesh_file_name,
+  const std::vector<std::vector<double> > & V,
+  const std::vector<std::vector<int> > & T,
+  const std::vector<std::vector<int> > & F
+)
+{
+  return igl::writeMESH(mesh_file_name, V, T, F);
+}, __doc_igl_writeMESH,
+py::arg("mesh_file_name"), py::arg("V"), py::arg("T"), py::arg("F"));
+
+m.def("writeMESH", []
+(
+  const std::string str,
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXd& T,
+  const Eigen::MatrixXi& F
+)
+{
+  return igl::writeMESH(str, V, T, F);
+}, __doc_igl_writeMESH,
+py::arg("str"), py::arg("V"), py::arg("T"), py::arg("F"));
+

python/py_vector.cpp

@@ -91,6 +91,7 @@ py::class_<Type> bind_eigen_2(py::module &m, const char *name,
         .def("size", [](const Type &m) { return m.size(); })
         .def("cols", [](const Type &m) { return m.cols(); })
         .def("rows", [](const Type &m) { return m.rows(); })
+        .def("shape", [](const Type &m) { return std::tuple<int,int>(m.rows(), m.cols()); })
 
         /* Extract rows and colums */
         .def("col", [](const Type &m, int i) {
@@ -356,6 +357,8 @@ py::class_<Type> bind_eigen_sparse_2(py::module &m, const char *name,
         .def("size", [](const Type &m) { return m.size(); })
         .def("cols", [](const Type &m) { return m.cols(); })
         .def("rows", [](const Type &m) { return m.rows(); })
+        .def("shape", [](const Type &m) { return std::tuple<int,int>(m.rows(), m.cols()); })
+
 
         /* Initialization */
         .def("setZero", [](Type &m) { m.setZero(); })
@@ -515,6 +518,7 @@ py::class_<Type> bind_eigen_diagonal_2(py::module &m, const char *name,
         .def("size", [](const Type &m) { return m.size(); })
         .def("cols", [](const Type &m) { return m.cols(); })
         .def("rows", [](const Type &m) { return m.rows(); })
+        .def("shape", [](const Type &m) { return std::tuple<int,int>(m.rows(), m.cols()); })
 
         /* Initialization */
         .def("setZero", [](Type &m) { m.setZero(); })
@@ -648,6 +652,8 @@ void python_export_vector(py::module &m) {
        else
           return "Numerical Issue";
     })
+    .def("analyzePattern",[](Eigen::SimplicialLLT<Eigen::SparseMatrix<double > >& s, const Eigen::SparseMatrix<double>& a) { return s.analyzePattern(a); })
+    .def("factorize",[](Eigen::SimplicialLLT<Eigen::SparseMatrix<double > >& s, const Eigen::SparseMatrix<double>& a) { return s.factorize(a); })
     .def("solve",[](const Eigen::SimplicialLLT<Eigen::SparseMatrix<double > >& s, const Eigen::MatrixXd& rhs) { return Eigen::MatrixXd(s.solve(rhs)); })
     ;
 

python/scripts/basic_function.mako

@@ -0,0 +1,46 @@
+% for enum in enums:
+py::enum_<\
+% for n in enum['namespaces']:
+${n}::\
+% endfor
+${enum['name']}>(m, "${enum['name']}")
+% for c in enum['constants']:
+    .value("${c}", \
+% for n in enum['namespaces']:
+${n}::\
+% endfor
+${c})
+% endfor
+    .export_values();
+% endfor
+
+
+% for func in functions:
+m.def("${func['name']}", []
+(
+  % for p in func['parameters'][:-1]:
+  ${p['type']} ${p['name']},
+  % endfor
+  ${func['parameters'][-1]['type']} ${func['parameters'][-1]['name']}
+)
+{
+  return \
+% for n in func['namespaces']:
+${n}::\
+% endfor
+${func['name']}(\
+% for p in func['parameters'][:-1]:
+${p['name']}, \
+% endfor
+${func['parameters'][-1]['name']});
+}, __doc_\
+% for n in func['namespaces']:
+${n}_\
+% endfor
+${func['name']},
+% for p in func['parameters'][:-1]:
+py::arg("${p['name']}"), \
+% endfor
+py::arg("${func['parameters'][-1]['name']}"));
+
+% endfor

python/scripts/generate_bindings.py

@@ -0,0 +1,298 @@
+#!/usr/bin/env python3
+#
+#  Syntax: generate_docstrings.py <path_to_c++_header_files> <path_to_python_files>
+#
+#  Extract documentation from C++ header files to use it in libiglPython bindings
+#
+
+import os, sys, glob
+import pickle
+
+import shutil
+from joblib import Parallel, delayed
+from multiprocessing import cpu_count
+from mako.template import Template
+from parser import parse
+
+
+# http://stackoverflow.com/questions/3207219/how-to-list-all-files-of-a-directory-in-python
+def get_filepaths(directory):
+    """
+    This function will generate the file names in a directory
+    tree by walking the tree either top-down or bottom-up. For each
+    directory in the tree rooted at directory top (including top itself),
+    it yields a 3-tuple (dirpath, dirnames, filenames).
+    """
+    file_paths = []  # List which will store all of the full filepaths.
+
+    # Walk the tree.
+    for root, directories, files in os.walk(directory):
+        for filename in files:
+            # Join the two strings in order to form the full filepath.
+            filepath = os.path.join(root, filename)
+            file_paths.append(filepath)  # Add it to the list.
+
+    return file_paths  # Self-explanatory.
+
+
+def get_name_from_path(path, basepath, prefix, postfix):
+    f_clean = path[len(basepath):]
+    f_clean = f_clean.replace(basepath, "")
+    f_clean = f_clean.replace(postfix, "")
+    f_clean = f_clean.replace(prefix, "")
+    f_clean = f_clean.replace("/", "_")
+    f_clean = f_clean.replace("\\", "_")
+    f_clean = f_clean.replace(" ", "_")
+    f_clean = f_clean.replace(".", "_")
+    return f_clean
+
+
+def map_parameter_types(name, cpp_type, parsed_types, errors, enum_types):
+    # TODO Replace with proper regex matching and derive types from templates, comment parsing, names in cpp files
+    # CAUTION: This is work in progress mapping code to get a grip of the problem
+    # Types to map
+    #    const int dim -> const int& dim ?
+    result = []
+
+    if cpp_type.startswith("const"):
+        result.append("const ")
+        cpp_type = cpp_type[6:]  # Strip const part
+
+    # Handle special types
+    skip_parsing = False
+    if cpp_type.startswith("MatY"):
+        result.append("Eigen::SparseMatrix<double>&")
+        skip_parsing = True
+    if cpp_type == "std::vector<std::vector<Scalar> > &":
+        result.append("std::vector<std::vector<double> > &")
+        skip_parsing = True
+    if cpp_type == "std::vector<std::vector<Index> > &":
+        result.append("std::vector<std::vector<int> > &")
+        skip_parsing = True
+    for constant in enum_types:
+        if cpp_type.endswith(constant):
+            result.append(cpp_type)
+            skip_parsing = True
+
+    if len(parsed_types) == 0:
+        errors.append("Empty typechain: %s" % cpp_type)
+        if cpp_type == "int" or cpp_type == "bool":
+            return cpp_type, True
+        else:
+            return cpp_type, False
+
+    # print(parsed_types, cpp_type)
+    if not skip_parsing:
+        for i, t in enumerate(parsed_types):
+
+            if t == "Eigen":
+                result.append("Eigen::")
+                continue
+            if t == "std":
+                result.append("std::")
+                continue
+
+            if t == "PlainObjectBase" or t == "MatrixBase":
+                if name == "F":
+                    result.append("MatrixXi&")
+                elif name == "V":
+                    result.append("MatrixXd&")
+                else:
+                    result.append("MatrixXd&")
+                break
+            if t == "MatrixXi":
+                result.append("MatrixXi&")
+                break
+            if t == "MatrixXd" or t == "VectorXd":
+                result.append("MatrixXd&")
+                break
+            if t == "SparseMatrix" and len(parsed_types) >= i + 2 and (
+                    parsed_types[i + 1] == "Scalar" or parsed_types[i + 1] == "T"):
+                result.append("SparseMatrix<double>&")
+                break
+            if t == "SparseVector" and len(parsed_types) >= i + 2 and (parsed_types[i + 1] == "Scalar" or parsed_types[
+                    i + 1] == "T"):
+                result.append("SparseMatrix<double>&")
+                break
+
+            if t == "bool" or t == "int" or t == "double" or t == "unsigned" or t == "string":
+                if cpp_type.endswith("&"):
+                    result.append(t + " &")
+                else:
+                    result.append(t)
+                break
+
+            else:
+                errors.append("Unknown typechain: %s" % cpp_type)
+                return cpp_type, False
+
+
+    return "".join(result), True
+
+
+if __name__ == '__main__':
+
+    if len(sys.argv) != 2:
+        print('Syntax: %s <path_to_c++_files>' % sys.argv[0])
+        exit(-1)
+
+    errors = {"missing": [], "empty": [], "others": [], "incorrect": [], "render": [], "various": []}
+    files = {"complete": [], "partial": [], "errors": [], "others": [], "empty": []}
+
+    # List all files in the given folder and subfolders
+    cpp_base_path = sys.argv[1]
+    cpp_file_paths = get_filepaths(cpp_base_path)
+
+    # Add all the .h filepaths to a dict
+    print("Collecting cpp files for parsing...")
+    mapping = {}
+    cppmapping = {}
+    for f in cpp_file_paths:
+        if f.endswith(".h"):
+            name = get_name_from_path(f, cpp_base_path, "", ".h")
+            mapping[name] = f
+
+        if f.endswith(".cpp"):
+            name = get_name_from_path(f, cpp_base_path, "", ".cpp")
+            cppmapping[name] = f
+
+    # Add all python binding files to a list
+    implemented_names = list(mapping.keys())  # ["point_mesh_squared_distance"]
+    implemented_names.sort()
+    single_postfix = ""
+    single_prefix = ""
+
+    # Create a list of all cpp header files
+    files_to_parse = []
+    cppfiles_to_parse = []
+    for n in implemented_names:
+        files_to_parse.append(mapping[n])
+
+        if n not in cppmapping:
+            errors["missing"].append("No cpp source file for function %s found." % n)
+        else:
+            cppfiles_to_parse.append(cppmapping[n])
+
+    # Parse c++ header files
+    print("Parsing header files...")
+    load_headers = False
+    if load_headers:
+        with open("headers.dat", 'rb') as fs:
+            dicts = pickle.load(fs)
+    else:
+        job_count = cpu_count()
+        dicts = Parallel(n_jobs=job_count)(delayed(parse)(path) for path in files_to_parse)
+
+    if not load_headers:
+        print("Saving parsed header files...")
+        with open("headers.dat", 'wb') as fs:
+            pickle.dump(dicts, fs)
+
+    # Not yet needed, as explicit template parsing does not seem to be supported in clang
+    # Parse c++ source files
+    # cppdicts = Parallel(n_jobs=job_count)(delayed(parse)(path) for path in cppfiles_to_parse)
+
+    # Change directory to become independent of execution directory
+    print("Generating directory tree for binding files...")
+    path = os.path.dirname(__file__)
+    if path != "":
+        os.chdir(path)
+    shutil.rmtree("generated")
+    os.makedirs("generated/complete")
+    os.mkdir("generated/partial")
+
+    print("Generating and writing binding files...")
+    for idx, n in enumerate(implemented_names):
+        d = dicts[idx]
+        contained_elements = sum(map(lambda x: len(x), d.values()))
+
+        # Skip files that don't contain functions/enums/classes
+        if contained_elements == 0:
+            errors["empty"].append("Function %s contains no parseable content in cpp header. Something might be wrong." % n)
+            files["empty"].append(n)
+            continue
+
+        # Add functions with classes to others
+        if len(d["classes"]) != 0 or len(d["structs"]) != 0:
+            errors["others"].append("Function %s contains classes/structs in cpp header. Skipping" % n)
+            files["others"].append(n)
+            continue
+
+        # Work on files that contain only functions/enums and namespaces
+        if len(d["functions"]) + len(d["namespaces"]) + len(d["enums"]) == contained_elements:
+            correct_functions = []
+            incorrect_functions = []
+
+            # Collect enums to generate binding files
+            enums = []
+            enum_types = []
+            for e in d["enums"]:
+                enums.append({"name": e.name, "namespaces": d["namespaces"], "constants": e.constants})
+                enum_types.append(e.name)
+
+            # Collect functions to generate binding files
+            for f in d["functions"]:
+                parameters = []
+                correct_function = True
+                f_errors = []
+                for p in f.parameters:
+                    typ, correct = map_parameter_types(p[0], p[1], p[2], f_errors, enum_types)
+                    correct_function &= correct
+                    parameters.append({"name": p[0], "type": typ})
+
+                if correct_function and len(parameters) > 0: #TODO add constants like EPS
+                    correct_functions.append({"parameters": parameters, "namespaces": d["namespaces"], "name": f.name})
+                elif len(parameters) > 0:
+                    incorrect_functions.append({"parameters": parameters, "namespaces": d["namespaces"], "name": f.name})
+                    errors["incorrect"].append("Incorrect function in %s: %s, %s\n" % (n, f.name, ",".join(f_errors)))
+                else:
+                    errors["various"].append("Function without pars in %s: %s, %s\n" % (n, f.name, ","
+                                                                                                     "".join(f_errors)))
+
+            # Write binding files
+            try:
+                tpl = Template(filename='basic_function.mako')
+                rendered = tpl.render(functions=correct_functions, enums=enums)
+                tpl1 = Template(filename='basic_function.mako')
+                rendered1 = tpl.render(functions=incorrect_functions, enums=enums)
+                path = "generated/"
+                if len(incorrect_functions) == 0 and (len(correct_functions) != 0 or len(enums) != 0):
+                    path += "complete/"
+                    with open(path + single_prefix + "py_" + n + ".cpp", 'w') as fs:
+                        fs.write(rendered)
+                    files["complete"].append(n)
+                else:
+                    path += "partial/"
+                    with open(path + single_prefix + "py_" + n + ".cpp", 'w') as fs:
+                        fs.write("// COMPLETE BINDINGS ========================\n")
+                        fs.write(rendered)
+                        fs.write("\n\n\n\n// INCOMPLETE BINDINGS ========================\n")
+                        fs.write(rendered1)
+
+                    if len(correct_functions) != 0:
+                        files["partial"].append(n)
+                    else:
+                        files["errors"].append(n)
+
+            except Exception as e:
+                files["errors"].append(n)
+                errors["render"].append("Template rendering failed:" + n + " " + str(correct_functions) + ", incorrect "
+                                                                                                "functions are " + str(
+                    incorrect_functions) + str(e) + "\n")
+
+    print("Writing error and overview files...")
+    with open("errors.txt" + single_postfix, 'w') as fs:
+        l = list(errors.keys())
+        l.sort()
+        for k in l:
+            fs.write("%s: %i \n" %(k, len(errors[k])))
+            fs.writelines("\n".join(errors[k]))
+            fs.write("\n\n\n")
+
+    with open("files.txt" + single_postfix, 'w') as fs:
+        l = list(files.keys())
+        l.sort()
+        for k in l:
+            fs.write("%s: %i \n" %(k, len(files[k])))
+            fs.writelines("\n".join(files[k]))
+            fs.write("\n\n\n")

python/scripts/generate_docstrings.py

@@ -87,16 +87,16 @@ if __name__ == '__main__':
     cpplines = []
 
     for idx, n in enumerate(implemented_names):
-        dict = dicts[idx]
-        contained_elements = sum(map(lambda x: len(x), dict.values()))
+        d = dicts[idx]
+        contained_elements = sum(map(lambda x: len(x), d.values()))
         # Check for files that don't contain functions/enums/classes
         if contained_elements == 0:
             print("Function %s contains no parseable content in cpp header. Something might be wrong." % n)
             continue
         else:
             names = []
-            namespaces = "_".join(dict["namespaces"])  # Assumption that all entities lie in deepest namespace
-            for f in dict["functions"]:
+            namespaces = "_".join(d["namespaces"])  # Assumption that all entities lie in deepest namespace
+            for f in d["functions"]:
                 h_string = "extern const char *__doc_" + namespaces + "_" + f.name + ";\n"
                 docu_string = "See " + f.name + " for the documentation."
                 if f.documentation != "":

python/scripts/parser.py

@@ -3,10 +3,16 @@ import os
 from threading import Thread
 
 import clang.cindex
+import ccsyspath
 import itertools
 from mako.template import Template
 
 
+
+
+
+
+
 def get_annotations(node):
     return [c.displayname for c in node.get_children()
             if c.kind == clang.cindex.CursorKind.ANNOTATE_ATTR]
@@ -17,33 +23,33 @@ class Function(object):
         self.name = cursor.spelling
         self.annotations = get_annotations(cursor)
         self.access = cursor.access_specifier
-        #        template_pars = [c.extent for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.TEMPLATE_TYPE_PARAMETER]
-        #        parameter_dec = [c for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.PARM_DECL]
-        #        print(parameter_dec, template_pars)
-        #        print(cursor.get_num_template_arguments(), cursor.get_template_argument_type(0), cursor.get_template_argument_value(0), template_pars, parameter_dec)
-        self.parameters = []
-        self.parnames = []
+#        template_pars = [c.extent for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.TEMPLATE_TYPE_PARAMETER]
+        parameter_dec = [c for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.PARM_DECL]
+        
+        parameters = []
+        for p in parameter_dec:
+            children = []
+            for c in p.get_children():
+#                print(c.spelling)
+                children.append(c.spelling)
+            parameters.append((p.spelling, p.type.spelling, children))
+
+        self.parameters = parameters
         self.documentation = cursor.raw_comment
 
 
 class Enum(object):
     def __init__(self, cursor):
         self.name = cursor.spelling
-        self.annotations = get_annotations(cursor)
-        self.access = cursor.access_specifier
-        #        template_pars = [c.extent for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.TEMPLATE_TYPE_PARAMETER]
-        #        parameter_dec = [c for c in cursor.get_children() if c.kind == clang.cindex.CursorKind.PARM_DECL]
-        #        print(parameter_dec, template_pars)
-        #        print(cursor.get_num_template_arguments(), cursor.get_template_argument_type(0), cursor.get_template_argument_value(0), template_pars, parameter_dec)
-        self.parameters = []
-        self.parnames = []
+        self.constants = [c.spelling for c in cursor.get_children() if c.kind ==
+                       clang.cindex.CursorKind.ENUM_CONSTANT_DECL]
         self.documentation = cursor.raw_comment
 
-# class Class(object):
-#    def __init__(self, cursor):
-#        self.name = cursor.spelling
+class Class(object):
+    def __init__(self, cursor):
+        self.name = cursor.spelling
 #        self.functions = []
-#        self.annotations = get_annotations(cursor)
+        self.annotations = get_annotations(cursor)
 
 #        for c in cursor.get_children():
 #            if (c.kind == clang.cindex.CursorKind.CXX_METHOD and
@@ -51,19 +57,16 @@ class Enum(object):
 #                f = Function(c)
 #                self.functions.append(f)
 
-def find_namespace_node(c):
-    if (c.kind == clang.cindex.CursorKind.NAMESPACE and c.spelling == "igl"):
-        return c
-    else:
-        for child_node in c.get_children():
-            return find_namespace_node(child_node)
 
 
 def traverse(c, path, objects):
     if c.location.file and not c.location.file.name.endswith(path):
         return
 
-    # print(c.kind, c.spelling)
+    if c.spelling == "PARULA_COLOR_MAP": # Fix to prevent python stack overflow from infinite recursion
+        return
+
+#    print(c.kind, c.spelling)
 
 
     if c.kind == clang.cindex.CursorKind.TRANSLATION_UNIT or c.kind == clang.cindex.CursorKind.UNEXPOSED_DECL:
@@ -76,7 +79,7 @@ def traverse(c, path, objects):
         pass
 
     elif c.kind == clang.cindex.CursorKind.FUNCTION_TEMPLATE:
-        # print("Function Template", c.spelling, c.raw_comment)
+#        print("Function Template", c.spelling, c.raw_comment)
         objects["functions"].append(Function(c))
         return
 
@@ -90,6 +93,18 @@ def traverse(c, path, objects):
         objects["enums"].append(Enum(c))
         return
 
+    elif c.kind == clang.cindex.CursorKind.CLASS_DECL:
+        objects["classes"].append(Class(c))
+        return
+
+    elif c.kind == clang.cindex.CursorKind.CLASS_TEMPLATE:
+        objects["classes"].append(Class(c))
+        return
+
+    elif c.kind == clang.cindex.CursorKind.STRUCT_DECL:
+        objects["structs"].append(Class(c))
+        return
+
     else:
         # print("Unknown", c.kind, c.spelling)
         pass
@@ -100,20 +115,16 @@ def traverse(c, path, objects):
 
 def parse(path):
     index = clang.cindex.Index.create()
-    tu = index.parse(path, ['-x', 'c++', '-std=c++11', '-fparse-all-comments', '-DIGL_STATIC_LIBRARY'])
-    # Clang can't parse files with missing definitions, add static library definition
-    objects = {"functions": [], "enums": [], "namespaces": [], "classes": []}
+    # Clang can't parse files with missing definitions, add static library definition or not?
+    args = ['-x', 'c++', '-std=c++11', '-fparse-all-comments', '-DIGL_STATIC_LIBRARY']
+    args.append('-I/usr/include/eigen3/') # TODO Properly add all needed includes
+    syspath = ccsyspath.system_include_paths('clang++-3.7') # Add the system libraries
+    incargs = [(b'-I' + inc).decode("utf-8") for inc in syspath]
+    args.extend(incargs)
+
+    tu = index.parse(path, args)
+    objects = {"functions": [], "enums": [], "namespaces": [], "classes": [], "structs": []}
     traverse(tu.cursor, path, objects)
-
-    #    tpl = Template(filename='bind.mako')
-    #    rendered = tpl.render(functions=functions)
-
-    #    OUTPUT_DIR = 'generated'
-
-    #    if not os.path.isdir(OUTPUT_DIR): os.mkdir(OUTPUT_DIR)
-
-    #    with open("generated/{}.bind.cc".format(sys.argv[1]), "w") as f:
-    #        f.write(rendered)
     return objects
 
 if __name__ == '__main__':

python/tutorial/704_SignedDistance.py

@@ -0,0 +1,140 @@
+from __future__ import print_function
+
+# Add the igl library to the modules search path
+import sys, os
+sys.path.insert(0, os.getcwd() + "/../")
+
+import pyigl as igl
+from iglhelpers import e2p
+import math
+TUTORIAL_SHARED_PATH = "../../tutorial/shared/"
+
+global V, F, T, tree, FN, VN, EN, E, EMAP, max_distance, slice_z, overlay
+
+V = igl.eigen.MatrixXd()
+F = igl.eigen.MatrixXi()
+T = igl.eigen.MatrixXi()
+tree = igl.AABB()
+FN = igl.eigen.MatrixXd()
+VN = igl.eigen.MatrixXd()
+EN = igl.eigen.MatrixXd()
+E = igl.eigen.MatrixXi()
+EMAP = igl.eigen.MatrixXi()
+
+max_distance = 1
+slice_z = 0.5
+overlay = False
+
+viewer = igl.viewer.Viewer()
+
+def update_visualization(viewer):
+    global V, F, T, tree, FN, VN, EN, E, EMAP, max_distance, slice_z, overlay
+    plane = igl.eigen.MatrixXd([0.0, 0.0, 1.0, -((1-slice_z) * V.col(2).minCoeff() + slice_z * V.col(2).maxCoeff())])
+    V_vis = igl.eigen.MatrixXd()
+    F_vis = igl.eigen.MatrixXi()
+
+    # Extract triangle mesh slice through volume mesh and subdivide nasty triangles
+    J = igl.eigen.MatrixXi()
+    bary = igl.eigen.SparseMatrixd()
+    igl.slice_tets(V, T, plane, V_vis, F_vis, J, bary)
+    max_l = 0.03
+#    while True:
+#        l = igl.eigen.MatrixXd()
+#        igl.edge_lengths(V_vis, F_vis, l)
+#        l /= (V_vis.colwise().maxCoeff() - V_vis.colwise().minCoeff()).norm()
+#        
+#        if l.maxCoeff() < max_l:
+#            break
+#        
+#        bad = e2p(l.rowwiseMaxCoeff())
+#        bad = bad > max_l
+#        F_vis_bad = igl.eigen.MatrixXi()
+#        F_vis_good = igl.eigen.MatrixXi()
+#        igl::slice_mask(F_vis, bad, 1, F_vis_bad);
+#        igl::slice_mask(F_vis, (bad!=true).eval(), 1, F_vis_good);
+#        igl.upsample(V_vis, F_vis_bad)
+#        F_vis = igl.cat(1, F_vis_bad, F_vis_good)
+
+
+#    #Compute signed distance
+#    S_vis = igl.eigen.MatrixXd()
+#    I = igl.eigen.MatrixXi()
+#    N = igl.eigen.MatrixXd()
+#    C = igl.eigen.MatrixXd()
+
+#    # Bunny is a watertight mesh so use pseudonormal for signing
+#    igl.signed_distance_pseudonormal(V_vis, V, F, tree, FN, VN, EN, EMAP, S_vis, I, C, N)
+
+#    # push to [0,1] range
+#    S_vis.array() = 0.5*(S_vis.array()/max_distance)+0.5;
+#    C_vis = igl.eigen.MatrixXi()
+#    # color without normalizing
+#    igl.parula(S_vis, False, C_vis)
+
+
+#    const auto & append_mesh = [&C_vis,&F_vis,&V_vis](const Eigen::MatrixXd & V, const Eigen::MatrixXi & F, const RowVector3d & color)
+
+#    F_vis.conservativeResize(F_vis.rows() + F.rows(), 3)
+#    F_vis.bottomRows(F.rows()) = F.array() + V_vis.rows()
+#    V_vis.conservativeResize(V_vis.rows() + V.rows(), 3)
+#    V_vis.bottomRows(V.rows()) = V
+#    C_vis.conservativeResize(C_vis.rows() + V.rows(), 3)
+#    C_vis.bottomRows(V.rows()).rowwise() = color
+
+#    if overlay:
+#        append_mesh(V, F, RowVector3d(0.8,0.8,0.8))
+
+    viewer.data.clear()
+    viewer.data.set_mesh(V_vis, F_vis)
+#    viewer.data.set_colors(C_vis)
+    viewer.core.lighting_factor = overlay
+
+
+
+
+def key_down(viewer, key, modifier):
+    global slice_z, overlay
+
+    if key == ord(' '):
+        overlay = not overlay
+    elif key == ord('.'):
+        slice_z = min(slice_z + 0.01, 0.99)
+    elif key == ord(','):
+        slice_z = max(slice_z - 0.01, 0.01)
+    else:
+        return False
+
+    update_visualization(viewer)
+    return True
+
+
+print("Press [space] to toggle showing surface.")
+print("Press '.'/',' to push back/pull forward slicing plane.")
+
+#Load mesh: (V,T) tet-mesh of convex hull, F contains original surface triangles
+igl.readMESH(TUTORIAL_SHARED_PATH + "bunny.mesh", V, T, F);
+
+#Call to point_mesh_squared_distance to determine bounds
+sqrD = igl.eigen.MatrixXd()
+I = igl.eigen.MatrixXi()
+C = igl.eigen.MatrixXd()
+igl.point_mesh_squared_distance(V, V, F, sqrD, I, C)
+max_distance = math.sqrt(sqrD.maxCoeff())
+
+#Precompute signed distance AABB tree
+tree.init(V, F)
+
+#Precompute vertex, edge and face normals
+igl.per_face_normals(V, F, FN)
+igl.per_vertex_normals(V, F, igl.PER_VERTEX_NORMALS_WEIGHTING_TYPE_ANGLE, FN, VN)
+igl.per_edge_normals(V, F, igl.PER_EDGE_NORMALS_WEIGHTING_TYPE_UNIFORM, FN, EN, E, EMAP)
+
+#Plot the generated mesh
+update_visualization(viewer);
+viewer.callback_key_down = key_down
+viewer.core.show_lines = False
+viewer.launch()
+
+
+
+