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

Former-commit-id: 37d21138d8c747e6d9310fac91b8009a7b8dbfbb

.gitignore

@@ -93,3 +93,5 @@ tests/bin
 python/build3
 *.pyc
 python/build4
+python/scripts/generated
+python/builddebug

.travis.yml

@@ -11,13 +11,13 @@ matrix:
         - cd python
         - mkdir build
         - cd build
-        - cmake -DCMAKE_CXX_COMPILER=g++-4.8 -DCMAKE_C_COMPILER=gcc-4.8 -DLIBIGL_WITH_EMBREE=OFF ../
+        - cmake -DCMAKE_CXX_COMPILER=g++-4.8 -DCMAKE_C_COMPILER=gcc-4.8 -DLIBIGL_WITH_EMBREE=OFF -DLIBIGL_USE_STATIC_LIBRARY=ON ../
         - make -j 2
         - cd ../../
         - cd tutorial
         - mkdir build
         - cd build
-        - cmake -DLIBIGL_USE_STATIC_LIBRARY=ON -DCMAKE_CXX_COMPILER=g++-4.8 -DCMAKE_C_COMPILER=gcc-4.8 -DLIBIGL_WITH_EMBREE=OFF ../
+        - cmake -DLIBIGL_USE_STATIC_LIBRARY=ON -DCMAKE_CXX_COMPILER=g++-4.8 -DCMAKE_C_COMPILER=gcc-4.8 ../
         - make -j 2
       addons:
         apt:

README.md

@@ -23,8 +23,8 @@ It is **a header-only library**. You do not need to compile anything to use,
 just include igl headers (e.g. `#include <igl/cotmatrix.h>`) and run.  Each
 header file contains a single function (e.g. `igl/cotmatrix.h` contains
 `igl::cotmatrix()`). Most are tailored to operate on a generic triangle mesh
-stored in an n-by-3 matrix of vertex positions V and an m-by-3 matrix of
-triangle indices F.
+stored in an n-by-3 matrix of vertex positions `V` and an m-by-3 matrix of
+triangle indices `F`.
 
 _Optionally_ the library may also be [pre-compiled](optional/) into a statically
 linked library, for faster compile times with your projects. This only effects
@@ -244,7 +244,7 @@ page](https://github.com/libigl/libigl/issues).
 
 ## Copyright
 2016 Alec Jacobson, Daniele Panozzo, Christian Schüller, Olga Diamanti, Qingnan
-Zhou, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
+Zhou, Sebastian Koch, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
 Giorgis, Luigi Rocca, Leonardo Sacht, Kevin Walliman, Olga Sorkine-Hornung, and others.
 
 Please see individual files for appropriate copyright notices.

include/igl/AABB.cpp

@@ -28,7 +28,7 @@
 
 template <typename DerivedV, int DIM>
   template <typename Derivedbb_mins, typename Derivedbb_maxs>
-inline void igl::AABB<DerivedV,DIM>::init(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::init(
     const Eigen::PlainObjectBase<DerivedV> & V,
     const Eigen::MatrixXi & Ele, 
     const Eigen::PlainObjectBase<Derivedbb_mins> & bb_mins,
@@ -92,7 +92,7 @@ inline void igl::AABB<DerivedV,DIM>::init(
 }
 
   template <typename DerivedV, int DIM>
-inline void igl::AABB<DerivedV,DIM>::init(
+void igl::AABB<DerivedV,DIM>::init(
     const Eigen::PlainObjectBase<DerivedV> & V,
     const Eigen::MatrixXi & Ele)
 {
@@ -102,7 +102,7 @@ inline void igl::AABB<DerivedV,DIM>::init(
 }
 
   template <typename DerivedV, int DIM>
-inline void igl::AABB<DerivedV,DIM>::init(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::init(
     const Eigen::PlainObjectBase<DerivedV> & V,
     const Eigen::MatrixXi & Ele, 
     const Eigen::MatrixXi & SI,
@@ -201,14 +201,14 @@ inline void igl::AABB<DerivedV,DIM>::init(
 }
 
 template <typename DerivedV, int DIM>
-inline bool igl::AABB<DerivedV,DIM>::is_leaf() const
+IGL_INLINE bool igl::AABB<DerivedV,DIM>::is_leaf() const
 {
   return m_primitive != -1;
 }
 
 template <typename DerivedV, int DIM>
 template <typename Derivedq>
-inline std::vector<int> igl::AABB<DerivedV,DIM>::find(
+IGL_INLINE std::vector<int> igl::AABB<DerivedV,DIM>::find(
     const Eigen::PlainObjectBase<DerivedV> & V,
     const Eigen::MatrixXi & Ele, 
     const Eigen::PlainObjectBase<Derivedq> & q,
@@ -298,7 +298,7 @@ inline std::vector<int> igl::AABB<DerivedV,DIM>::find(
 }
 
 template <typename DerivedV, int DIM>
-inline int igl::AABB<DerivedV,DIM>::subtree_size() const
+IGL_INLINE int igl::AABB<DerivedV,DIM>::subtree_size() const
 {
   // 1 for self
   int n = 1;
@@ -318,7 +318,7 @@ inline int igl::AABB<DerivedV,DIM>::subtree_size() const
 
 template <typename DerivedV, int DIM>
 template <typename Derivedbb_mins, typename Derivedbb_maxs>
-inline void igl::AABB<DerivedV,DIM>::serialize(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::serialize(
     Eigen::PlainObjectBase<Derivedbb_mins> & bb_mins,
     Eigen::PlainObjectBase<Derivedbb_maxs> & bb_maxs,
     Eigen::VectorXi & elements,
@@ -350,7 +350,7 @@ inline void igl::AABB<DerivedV,DIM>::serialize(
 }
 
 template <typename DerivedV, int DIM>
-inline typename igl::AABB<DerivedV,DIM>::Scalar 
+IGL_INLINE typename igl::AABB<DerivedV,DIM>::Scalar 
 igl::AABB<DerivedV,DIM>::squared_distance(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
@@ -363,7 +363,7 @@ igl::AABB<DerivedV,DIM>::squared_distance(
 
 
 template <typename DerivedV, int DIM>
-inline typename igl::AABB<DerivedV,DIM>::Scalar 
+IGL_INLINE typename igl::AABB<DerivedV,DIM>::Scalar 
 igl::AABB<DerivedV,DIM>::squared_distance(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
@@ -447,7 +447,7 @@ template <
   typename DerivedsqrD, 
   typename DerivedI, 
   typename DerivedC>
-inline void igl::AABB<DerivedV,DIM>::squared_distance(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::squared_distance(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
   const Eigen::PlainObjectBase<DerivedP> & P,
@@ -475,7 +475,7 @@ template <
   typename DerivedsqrD, 
   typename DerivedI, 
   typename DerivedC>
-inline void igl::AABB<DerivedV,DIM>::squared_distance(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::squared_distance(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
   const AABB<Derivedother_V,DIM> & other,
@@ -508,7 +508,7 @@ template <
   typename DerivedsqrD, 
   typename DerivedI, 
   typename DerivedC>
-inline typename igl::AABB<DerivedV,DIM>::Scalar igl::AABB<DerivedV,DIM>::squared_distance_helper(
+IGL_INLINE typename igl::AABB<DerivedV,DIM>::Scalar igl::AABB<DerivedV,DIM>::squared_distance_helper(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
   const AABB<Derivedother_V,DIM> * other,
@@ -733,7 +733,7 @@ inline typename igl::AABB<DerivedV,DIM>::Scalar igl::AABB<DerivedV,DIM>::squared
 }
 
 template <typename DerivedV, int DIM>
-inline void igl::AABB<DerivedV,DIM>::leaf_squared_distance(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::leaf_squared_distance(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
   const RowVectorDIMS & p,
@@ -752,7 +752,7 @@ inline void igl::AABB<DerivedV,DIM>::leaf_squared_distance(
 
 
 template <typename DerivedV, int DIM>
-inline void igl::AABB<DerivedV,DIM>::set_min(
+IGL_INLINE void igl::AABB<DerivedV,DIM>::set_min(
   const RowVectorDIMS & 
 #ifndef NDEBUG
   p
@@ -781,7 +781,7 @@ inline void igl::AABB<DerivedV,DIM>::set_min(
 
 
 template <typename DerivedV, int DIM>
-inline bool 
+IGL_INLINE bool 
 igl::AABB<DerivedV,DIM>::intersect_ray(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
@@ -816,7 +816,7 @@ igl::AABB<DerivedV,DIM>::intersect_ray(
 }
 
 template <typename DerivedV, int DIM>
-inline bool 
+IGL_INLINE bool 
 igl::AABB<DerivedV,DIM>::intersect_ray(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
@@ -871,7 +871,7 @@ igl::AABB<DerivedV,DIM>::intersect_ray(
 }
 
 template <typename DerivedV, int DIM>
-inline bool 
+IGL_INLINE bool 
 igl::AABB<DerivedV,DIM>::intersect_ray(
   const Eigen::PlainObjectBase<DerivedV> & V,
   const Eigen::MatrixXi & Ele, 
@@ -938,6 +938,7 @@ igl::AABB<DerivedV,DIM>::intersect_ray(
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::init(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&);
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 2>::init(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&);
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_distance<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::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<double, -1, -1, 0, -1, -1> >&) const;
@@ -946,8 +947,4 @@ template double igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_
 template double igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 2>::squared_distance(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::Matrix<double, 1, 2, 1, 1, 2> const&, int&, Eigen::Matrix<double, 1, 2, 1, 1, 2>&) const;
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 2>::squared_distance<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::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<double, -1, -1, 0, -1, -1> >&) const;
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_distance<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::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<double, -1, -1, 0, -1, -1> >&) const;
-
-
 #endif
-
-

include/igl/AABB.h

@@ -9,6 +9,7 @@
 #define IGL_AABB_H
 
 #include "Hit.h"
+#include "igl_inline.h"
 #include <Eigen/Core>
 #include <Eigen/Geometry>
 #include <vector>
@@ -81,7 +82,7 @@ public:
       }
       // Seems like there should have been an elegant solution to this using
       // the copy-swap idiom above:
-      inline void deinit()
+      IGL_INLINE void deinit()
       {
         m_primitive = -1;
         m_box = Eigen::AlignedBox<Scalar,DIM>();
@@ -105,7 +106,7 @@ public:
       //   elements  max_tree list of element or (not leaf id) indices into Ele
       //   i  recursive call index {0}
       template <typename Derivedbb_mins, typename Derivedbb_maxs>
-        inline void init(
+        IGL_INLINE void init(
             const Eigen::PlainObjectBase<DerivedV> & V,
             const Eigen::MatrixXi & Ele, 
             const Eigen::PlainObjectBase<Derivedbb_mins> & bb_mins,
@@ -113,7 +114,7 @@ public:
             const Eigen::VectorXi & elements,
             const int i = 0);
       // Wrapper for root with empty serialization
-      inline void init(
+      IGL_INLINE void init(
           const Eigen::PlainObjectBase<DerivedV> & V,
           const Eigen::MatrixXi & Ele);
       // Build an Axis-Aligned Bounding Box tree for a given mesh.
@@ -128,13 +129,13 @@ public:
       //   I  #I list of indices into Ele of elements to include (for recursive
       //     calls)
       // 
-      inline void init(
+      IGL_INLINE void init(
           const Eigen::PlainObjectBase<DerivedV> & V,
           const Eigen::MatrixXi & Ele, 
           const Eigen::MatrixXi & SI,
           const Eigen::VectorXi & I);
       // Return whether at leaf node
-      inline bool is_leaf() const;
+      IGL_INLINE bool is_leaf() const;
       // Find the indices of elements containing given point: this makes sense
       // when Ele is a co-dimension 0 simplex (tets in 3D, triangles in 2D).
       //
@@ -148,7 +149,7 @@ public:
       // Returns:
       //   list of indices of elements containing q
       template <typename Derivedq>
-      inline std::vector<int> find(
+      IGL_INLINE std::vector<int> find(
           const Eigen::PlainObjectBase<DerivedV> & V,
           const Eigen::MatrixXi & Ele, 
           const Eigen::PlainObjectBase<Derivedq> & q,
@@ -156,7 +157,7 @@ public:
 
       // If number of elements m then total tree size should be 2*h where h is
       // the deepest depth 2^ceil(log(#Ele*2-1))
-      inline int subtree_size() const;
+      IGL_INLINE int subtree_size() const;
 
       // Serialize this class into 3 arrays (so we can pass it pack to matlab)
       //
@@ -166,7 +167,7 @@ public:
       //   elements  max_tree list of element or (not leaf id) indices into Ele
       //   i  recursive call index into these arrays {0}
       template <typename Derivedbb_mins, typename Derivedbb_maxs>
-        inline void serialize(
+        IGL_INLINE void serialize(
             Eigen::PlainObjectBase<Derivedbb_mins> & bb_mins,
             Eigen::PlainObjectBase<Derivedbb_maxs> & bb_maxs,
             Eigen::VectorXi & elements,
@@ -186,14 +187,14 @@ public:
       //
       // Known bugs: currently assumes Elements are triangles regardless of
       // dimension.
-      inline Scalar squared_distance(
+      IGL_INLINE Scalar squared_distance(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & p,
         int & i,
         RowVectorDIMS & c) const;
 //private:
-      inline Scalar squared_distance(
+      IGL_INLINE Scalar squared_distance(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & p,
@@ -201,21 +202,21 @@ public:
         int & i,
         RowVectorDIMS & c) const;
       // All hits
-      inline bool intersect_ray(
+      IGL_INLINE bool intersect_ray(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & origin,
         const RowVectorDIMS & dir,
         std::vector<igl::Hit> & hits) const;
       // First hit
-      inline bool intersect_ray(
+      IGL_INLINE bool intersect_ray(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & origin,
         const RowVectorDIMS & dir,
         igl::Hit & hit) const;
 //private:
-      inline bool intersect_ray(
+      IGL_INLINE bool intersect_ray(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & origin,
@@ -230,7 +231,7 @@ public:
         typename DerivedsqrD, 
         typename DerivedI, 
         typename DerivedC>
-      inline void squared_distance(
+      IGL_INLINE void squared_distance(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const Eigen::PlainObjectBase<DerivedP> & P,
@@ -243,7 +244,7 @@ public:
         typename DerivedsqrD, 
         typename DerivedI, 
         typename DerivedC>
-      inline void squared_distance(
+      IGL_INLINE void squared_distance(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const AABB<Derivedother_V,DIM> & other,
@@ -258,7 +259,7 @@ private:
         typename DerivedsqrD, 
         typename DerivedI, 
         typename DerivedC>
-      inline Scalar squared_distance_helper(
+      IGL_INLINE Scalar squared_distance_helper(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const AABB<Derivedother_V,DIM> * other,
@@ -269,14 +270,14 @@ private:
         Eigen::PlainObjectBase<DerivedI> & I,
         Eigen::PlainObjectBase<DerivedC> & C) const;
       // Helper function for leaves: works in-place on sqr_d
-      inline void leaf_squared_distance(
+      IGL_INLINE void leaf_squared_distance(
         const Eigen::PlainObjectBase<DerivedV> & V,
         const Eigen::MatrixXi & Ele, 
         const RowVectorDIMS & p,
         Scalar & sqr_d,
         int & i,
         RowVectorDIMS & c) const;
-      inline void set_min(
+      IGL_INLINE void set_min(
         const RowVectorDIMS & p,
         const Scalar sqr_d_candidate,
         const int i_candidate,

include/igl/ambient_occlusion.cpp

@@ -10,7 +10,7 @@
 #include "ray_mesh_intersect.h"
 #include "EPS.h"
 #include "Hit.h"
-#include <thread>
+#include "parallel_for.h"
 #include <functional>
 #include <vector>
 #include <algorithm>
@@ -38,40 +38,27 @@ IGL_INLINE void igl::ambient_occlusion(
   // Embree seems to be parallel when constructing but not when tracing rays
   const MatrixXf D = random_dir_stratified(num_samples).cast<float>();
 
-  const size_t nthreads = n<1000?1:std::thread::hardware_concurrency();
+  const auto & inner = [&P,&N,&num_samples,&D,&S,&shoot_ray](const int p)
   {
-    std::vector<std::thread> threads(nthreads);
-    for(int t = 0;t<nthreads;t++)
+    const Vector3f origin = P.row(p).template cast<float>();
+    const Vector3f normal = N.row(p).template cast<float>();
+    int num_hits = 0;
+    for(int s = 0;s<num_samples;s++)
     {
-      threads[t] = std::thread(std::bind(
-        [&P,&N,&shoot_ray,&S,&num_samples,&D](const int bi, const int ei, const int t)
-        {
-          // loop over mesh vertices in this chunk
-          for(int p = bi;p<ei;p++)
-          {
-            const Vector3f origin = P.row(p).template cast<float>();
-            const Vector3f normal = N.row(p).template cast<float>();
-            int num_hits = 0;
-            for(int s = 0;s<num_samples;s++)
-            {
-              Vector3f d = D.row(s);
-              if(d.dot(normal) < 0)
-              {
-                // reverse ray
-                d *= -1;
-              }
-              if(shoot_ray(origin,d))
-              {
-                num_hits++;
-              }
-            }
-            S(p) = (double)num_hits/(double)num_samples;
-          }
-        },t*n/nthreads,(t+1)==nthreads?n:(t+1)*n/nthreads,t));
+      Vector3f d = D.row(s);
+      if(d.dot(normal) < 0)
+      {
+        // reverse ray
+        d *= -1;
+      }
+      if(shoot_ray(origin,d))
+      {
+        num_hits++;
+      }
     }
-    std::for_each(threads.begin(),threads.end(),[](std::thread& x){x.join();});
-  }
-
+    S(p) = (double)num_hits/(double)num_samples;
+  };
+  parallel_for(n,inner,1000);
 }
 
 template <
@@ -147,4 +134,5 @@ template void igl::ambient_occlusion<Eigen::Matrix<double, -1, -1, 0, -1, -1>, E
 template void igl::ambient_occlusion<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -1, 1, 0, -1, 1> >(std::function<bool (Eigen::Matrix<float, 3, 1, 0, 3, 1> const&, Eigen::Matrix<float, 3, 1, 0, 3, 1> const&)> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
 // generated by autoexplicit.sh
 template void igl::ambient_occlusion<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, -1, 1, 0, -1, 1> >(std::function<bool (Eigen::Matrix<float, 3, 1, 0, 3, 1> const&, Eigen::Matrix<float, 3, 1, 0, 3, 1> const&)> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
+template void igl::ambient_occlusion<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(std::function<bool (Eigen::Matrix<float, 3, 1, 0, 3, 1> const&, Eigen::Matrix<float, 3, 1, 0, 3, 1> const&)> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/cat.cpp

@@ -6,12 +6,14 @@
 // 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 "cat.h"
+
 #include <cstdio>
 
 // Bug in unsupported/Eigen/SparseExtra needs iostream first
 #include <iostream>
 #include <unsupported/Eigen/SparseExtra>
 
+
 // Sparse matrices need to be handled carefully. Because C++ does not 
 // Template:
 //   Scalar  sparse matrix scalar type, e.g. double
@@ -22,6 +24,7 @@ IGL_INLINE void igl::cat(
     const Eigen::SparseMatrix<Scalar> & B, 
     Eigen::SparseMatrix<Scalar> & C)
 {
+
   assert(dim == 1 || dim == 2);
   using namespace Eigen;
   // Special case if B or A is empty
@@ -36,6 +39,7 @@ IGL_INLINE void igl::cat(
     return;
   }
 
+#if false
   // This **must** be DynamicSparseMatrix, otherwise this implementation is
   // insanely slow
   DynamicSparseMatrix<Scalar, RowMajor> dyn_C;
@@ -77,6 +81,109 @@ IGL_INLINE void igl::cat(
   }
 
   C = SparseMatrix<Scalar>(dyn_C);
+#elif false
+  std::vector<Triplet<Scalar> > CIJV;
+  CIJV.reserve(A.nonZeros() + B.nonZeros());
+  {
+    // Iterate over outside of A
+    for(int k=0; k<A.outerSize(); ++k)
+    {
+      // Iterate over inside
+      for(typename SparseMatrix<Scalar>::InnerIterator it (A,k); it; ++it)
+      {
+        CIJV.emplace_back(it.row(),it.col(),it.value());
+      }
+    }
+    // Iterate over outside of B
+    for(int k=0; k<B.outerSize(); ++k)
+    {
+      // Iterate over inside
+      for(typename SparseMatrix<Scalar>::InnerIterator it (B,k); it; ++it)
+      {
+        int r = (dim == 1 ? A.rows()+it.row() : it.row());
+        int c = (dim == 2 ? A.cols()+it.col() : it.col());
+        CIJV.emplace_back(r,c,it.value());
+      }
+    }
+
+  }
+
+  C = SparseMatrix<Scalar>( 
+      dim == 1 ? A.rows()+B.rows() : A.rows(),
+      dim == 1 ? A.cols()          : A.cols()+B.cols());
+  C.reserve(A.nonZeros() + B.nonZeros());
+  C.setFromTriplets(CIJV.begin(),CIJV.end());
+#else
+  C = SparseMatrix<Scalar>( 
+      dim == 1 ? A.rows()+B.rows() : A.rows(),
+      dim == 1 ? A.cols()          : A.cols()+B.cols());
+  Eigen::VectorXi per_col = Eigen::VectorXi::Zero(C.cols());
+  if(dim == 1)
+  {
+    assert(A.outerSize() == B.outerSize());
+    for(int k = 0;k<A.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (A,k); it; ++it)
+      {
+        per_col(k)++;
+      }
+      for(typename SparseMatrix<Scalar>::InnerIterator it (B,k); it; ++it)
+      {
+        per_col(k)++;
+      }
+    }
+  }else
+  {
+    for(int k = 0;k<A.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (A,k); it; ++it)
+      {
+        per_col(k)++;
+      }
+    }
+    for(int k = 0;k<B.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (B,k); it; ++it)
+      {
+        per_col(A.cols() + k)++;
+      }
+    }
+  }
+  C.reserve(per_col);
+  if(dim == 1)
+  {
+    for(int k = 0;k<A.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (A,k); it; ++it)
+      {
+        C.insert(it.row(),k) = it.value();
+      }
+      for(typename SparseMatrix<Scalar>::InnerIterator it (B,k); it; ++it)
+      {
+        C.insert(A.rows()+it.row(),k) = it.value();
+      }
+    }
+  }else
+  {
+    for(int k = 0;k<A.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (A,k); it; ++it)
+      {
+        C.insert(it.row(),k) = it.value();
+      }
+    }
+    for(int k = 0;k<B.outerSize();++k)
+    {
+      for(typename SparseMatrix<Scalar>::InnerIterator it (B,k); it; ++it)
+      {
+        C.insert(it.row(),A.cols()+k) = it.value();
+      }
+    }
+  }
+  C.makeCompressed();
+
+#endif
+
 }
 
 template <typename Derived, class MatC>

include/igl/copyleft/cgal/half_space_box.cpp

@@ -14,7 +14,10 @@ IGL_INLINE void igl::copyleft::cgal::half_space_box(
   typedef CGAL::Point_3<CGAL::Epeck> Point;
   typedef CGAL::Vector_3<CGAL::Epeck> Vector;
   typedef CGAL::Epeck::FT EScalar;
-  Eigen::RowVector3d avg = V.colwise().mean();
+  Eigen::Matrix<typename DerivedV::Scalar,1,3> avg(0,0,0);
+  for(int v = 0;v<V.rows();v++) for(int c = 0;c<V.cols();c++) avg(c) += V(v,c);
+  avg /= V.rows();
+
   Point o3(avg(0),avg(1),avg(2));
   Point o2 = P.projection(o3);
   Vector u;
@@ -109,5 +112,11 @@ IGL_INLINE void igl::copyleft::cgal::half_space_box(
 }
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::half_space_box<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3>&, Eigen::Matrix<int, 12, 3, 0, 12, 3>&);
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::half_space_box<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4> > const&, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3>&, Eigen::Matrix<int, 12, 3, 0, 12, 3>&);
+// generated by autoexplicit.sh
+template void igl::copyleft::cgal::half_space_box<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3>&, Eigen::Matrix<int, 12, 3, 0, 12, 3>&);
 template void igl::copyleft::cgal::half_space_box<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3>&, Eigen::Matrix<int, 12, 3, 0, 12, 3>&);
 #endif

include/igl/copyleft/cgal/intersect_with_half_space.cpp

@@ -19,11 +19,11 @@ IGL_INLINE bool igl::copyleft::cgal::intersect_with_half_space(
   Eigen::PlainObjectBase<DerivedFC > & FC,
   Eigen::PlainObjectBase<DerivedJ > & J)
 {
-  Eigen::Matrix<CGAL::Epeck::FT,8,3> BV;
-  Eigen::Matrix<int,12,3> BF;
-  half_space_box(p,n,V,BV,BF);
-  // Disturbingly, (BV,BF) must be first argument
-  return mesh_boolean(BV,BF,V,F,MESH_BOOLEAN_TYPE_INTERSECT,VC,FC,J);
+  typedef CGAL::Plane_3<CGAL::Epeck> Plane;
+  typedef CGAL::Point_3<CGAL::Epeck> Point;
+  typedef CGAL::Vector_3<CGAL::Epeck> Vector;
+  Plane P(Point(p(0),p(1),p(2)),Vector(n(0),n(1),n(2)));
+  return intersect_with_half_space(V,F,P,VC,FC,J);
 }
 
 template <
@@ -40,15 +40,48 @@ IGL_INLINE bool igl::copyleft::cgal::intersect_with_half_space(
   Eigen::PlainObjectBase<DerivedVC > & VC,
   Eigen::PlainObjectBase<DerivedFC > & FC,
   Eigen::PlainObjectBase<DerivedJ > & J)
+{
+  typedef CGAL::Plane_3<CGAL::Epeck> Plane;
+  Plane P(equ(0),equ(1),equ(2),equ(3));
+  return intersect_with_half_space(V,F,P,VC,FC,J);
+}
+
+template <
+  typename DerivedV,
+  typename DerivedF,
+  typename DerivedVC,
+  typename DerivedFC,
+  typename DerivedJ>
+IGL_INLINE bool igl::copyleft::cgal::intersect_with_half_space(
+  const Eigen::PlainObjectBase<DerivedV > & V,
+  const Eigen::PlainObjectBase<DerivedF > & F,
+  const CGAL::Plane_3<CGAL::Epeck> & P,
+  Eigen::PlainObjectBase<DerivedVC > & VC,
+  Eigen::PlainObjectBase<DerivedFC > & FC,
+  Eigen::PlainObjectBase<DerivedJ > & J)
 {
   Eigen::Matrix<CGAL::Epeck::FT,8,3> BV;
   Eigen::Matrix<int,12,3> BF;
-  half_space_box(equ,V,BV,BF);
+  half_space_box(P,V,BV,BF);
   // Disturbingly, (BV,BF) must be first argument
-  return mesh_boolean(BV,BF,V,F,MESH_BOOLEAN_TYPE_INTERSECT,VC,FC,J);
+  const bool ret = mesh_boolean(BV,BF,V,F,MESH_BOOLEAN_TYPE_INTERSECT,VC,FC,J);
+  // But now J is wrong...
+  std::for_each(
+    J.data(),
+    J.data()+J.size(),
+    [&BF,&F](typename DerivedJ::Scalar & j)
+      {j = (j<BF.rows()?F.rows()+j:j-BF.rows());}
+    );
+  return ret;
 }
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template bool igl::copyleft::cgal::intersect_with_half_space<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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> >&);
+// generated by autoexplicit.sh
+template bool igl::copyleft::cgal::intersect_with_half_space<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 4, 0, -1, 4> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+// generated by autoexplicit.sh
+template bool igl::copyleft::cgal::intersect_with_half_space<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4>, 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<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 1, 4, 1, 1, 4> > 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 bool igl::copyleft::cgal::intersect_with_half_space<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, 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, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > 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/copyleft/cgal/intersect_with_half_space.h

@@ -2,6 +2,8 @@
 #define IGL_COPYLEFT_CGAL_INTERSECT_WITH_HALF_SPACE_H
 #include "../../igl_inline.h"
 #include <Eigen/Core>
+#include <CGAL/Exact_predicates_exact_constructions_kernel.h>
+#include <CGAL/Plane_3.h>
 
 namespace igl
 {
@@ -19,7 +21,8 @@ namespace igl
       // Outputs:
       //   VC  #VC by 3 list of vertex positions of boolean result mesh
       //   FC  #FC by 3 list of triangle indices into VC
-      //   J  #FC list of indices into [F;F.rows()+[1;2]] revealing "birth" facet
+      //   J  #FC list of indices into [F;F.rows()+[1;2]] revealing "birth"
+      //     facet
       template <
         typename DerivedV,
         typename DerivedF,
@@ -36,12 +39,12 @@ namespace igl
         Eigen::PlainObjectBase<DerivedVC > & VC,
         Eigen::PlainObjectBase<DerivedFC > & FC,
         Eigen::PlainObjectBase<DerivedJ > & J);
-
       // Intersect a PWN mesh with a half-space. Plane equation.
       //
       // Inputs:
       //   V  #V by 3 list of mesh vertex positions
-      //   equ  plane equation: a*x+b*y+c*z + d = 0
+      //   equ  plane equation: P(x,y,z) = a*x+b*y+c*z + d = 0, P(x,y,z) < 0 is
+      //     _inside_.
       // Outputs:
       //   VC  #VC by 3 list of vertex positions of boolean result mesh
       //   FC  #FC by 3 list of triangle indices into VC
@@ -60,6 +63,28 @@ namespace igl
         Eigen::PlainObjectBase<DerivedVC > & VC,
         Eigen::PlainObjectBase<DerivedFC > & FC,
         Eigen::PlainObjectBase<DerivedJ > & J);
+      // Intersect a PWN mesh with a half-space. CGAL Plane.
+      //
+      // Inputs:
+      //   V  #V by 3 list of mesh vertex positions
+      //   P  plane 
+      // Outputs:
+      //   VC  #VC by 3 list of vertex positions of boolean result mesh
+      //   FC  #FC by 3 list of triangle indices into VC
+      //   J  #FC list of indices into [F;F.rows()+[1;2]] revealing "birth" facet
+      template <
+        typename DerivedV,
+        typename DerivedF,
+        typename DerivedVC,
+        typename DerivedFC,
+        typename DerivedJ>
+      IGL_INLINE bool intersect_with_half_space(
+        const Eigen::PlainObjectBase<DerivedV > & V,
+        const Eigen::PlainObjectBase<DerivedF > & F,
+        const CGAL::Plane_3<CGAL::Epeck> & P,
+        Eigen::PlainObjectBase<DerivedVC > & VC,
+        Eigen::PlainObjectBase<DerivedFC > & FC,
+        Eigen::PlainObjectBase<DerivedJ > & J);
     }
   }
 }

include/igl/copyleft/cgal/mesh_boolean.cpp

@@ -443,6 +443,8 @@ template bool igl::copyleft::cgal::mesh_boolean<Eigen::Matrix<CGAL::Lazy_exact_n
 template bool igl::copyleft::cgal::mesh_boolean<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> >(std::vector<Eigen::Matrix<double, -1, -1, 0, -1, -1>, std::allocator<Eigen::Matrix<double, -1, -1, 0, -1, -1> > > const&, std::vector<Eigen::Matrix<int, -1, -1, 0, -1, -1>, std::allocator<Eigen::Matrix<int, -1, -1, 0, -1, -1> > > const&, std::function<int (Eigen::Matrix<int, 1, -1, 1, 1, -1>)> const&, std::function<int (int, int)> 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 bool igl::copyleft::cgal::mesh_boolean<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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 bool igl::copyleft::cgal::mesh_boolean<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3>, Eigen::Matrix<int, 12, 3, 0, 12, 3>, 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<CGAL::Lazy_exact_nt<CGAL::Gmpq>, 8, 3, 0, 8, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, 12, 3, 0, 12, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::MeshBooleanType 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 bool igl::copyleft::cgal::mesh_boolean<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<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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, std::basic_string<char, std::char_traits<char>, std::allocator<char> > 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 bool igl::copyleft::cgal::mesh_boolean<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<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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::MeshBooleanType 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> >&);
 #undef IGL_STATIC_LIBRARY
 #include "../../remove_unreferenced.cpp"
 template void igl::remove_unreferenced<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<CGAL::Lazy_exact_nt<CGAL::Gmpq>, -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> >&);

include/igl/copyleft/cgal/propagate_winding_numbers.cpp

@@ -119,6 +119,7 @@ IGL_INLINE bool igl::copyleft::cgal::propagate_winding_numbers(
   bool valid = true;
   if (!piecewise_constant_winding_number(F, uE, uE2E)) 
   {
+    assert(false && "Input mesh is not orientable");
     std::cerr << "Input mesh is not orientable!" << std::endl;
     valid = false;
   }

include/igl/copyleft/cgal/string_to_mesh_boolean_type.cpp

@@ -1,4 +1,4 @@
-#include "string_to_mesh_boolean_type.h"
+#include "string_to_mesh_boolean_type.h"
 #include <algorithm>
 #include <cassert>
 #include <vector>

include/igl/copyleft/swept_volume.cpp

@@ -0,0 +1,49 @@
+#include "swept_volume.h"
+#include "../swept_volume_bounding_box.h"
+#include "../swept_volume_signed_distance.h"
+#include "../voxel_grid.h"
+#include "marching_cubes.h"
+#include <iostream>
+
+IGL_INLINE void igl::copyleft::swept_volume(
+  const Eigen::MatrixXd & V,
+  const Eigen::MatrixXi & F,
+  const std::function<Eigen::Affine3d(const double t)> & transform,
+  const size_t steps,
+  const size_t grid_res,
+  const size_t isolevel_grid,
+  Eigen::MatrixXd & SV,
+  Eigen::MatrixXi & SF)
+{
+  using namespace std;
+  using namespace Eigen;
+  using namespace igl;
+  using namespace igl::copyleft;
+
+  const auto & Vtransform = 
+    [&V,&transform](const size_t vi,const double t)->RowVector3d
+  {
+    Vector3d Vvi = V.row(vi).transpose();
+    return (transform(t)*Vvi).transpose();
+  };
+  AlignedBox3d Mbox;
+  swept_volume_bounding_box(V.rows(),Vtransform,steps,Mbox);
+
+  // Amount of padding: pad*h should be >= isolevel
+  const int pad = isolevel_grid+1;
+  // number of vertices on the largest side
+  const int s = grid_res+2*pad;
+  const double h = Mbox.diagonal().maxCoeff()/(double)(s-2.*pad-1.);
+  const double isolevel = isolevel_grid*h;
+
+  // create grid
+  RowVector3i res;
+  MatrixXd GV;
+  voxel_grid(Mbox,s,pad,GV,res);
+
+  // compute values
+  VectorXd S;
+  swept_volume_signed_distance(V,F,transform,steps,GV,res,h,isolevel,S);
+  S.array()-=isolevel;
+  marching_cubes(S,GV,res(0),res(1),res(2),SV,SF);
+}

include/igl/copyleft/swept_volume.h

@@ -0,0 +1,41 @@
+#ifndef IGL_COPYLEFT_SWEPT_VOLUME_H
+#define IGL_COPYLEFT_SWEPT_VOLUME_H
+#include "../igl_inline.h"
+#include <Eigen/Core>
+#include <Eigen/Geometry>
+namespace igl
+{
+  namespace copyleft
+  {
+    // Compute the surface of the swept volume of a solid object with surface
+    // (V,F) mesh under going rigid motion.
+    // 
+    // Inputs:
+    //   V  #V by 3 list of mesh positions in reference pose
+    //   F  #F by 3 list of mesh indices into V
+    //   transform  function handle so that transform(t) returns the rigid
+    //     transformation at time t∈[0,1]
+    //   steps  number of time steps: steps=3 --> t∈{0,0.5,1}
+    //   grid_res  number of grid cells on the longest side containing the
+    //     motion (isolevel+1 cells will also be added on each side as padding)
+    //   isolevel  distance level to be contoured as swept volume
+    // Outputs:
+    //   SV  #SV by 3 list of mesh positions of the swept surface
+    //   SF  #SF by 3 list of mesh faces into SV
+    IGL_INLINE void swept_volume(
+      const Eigen::MatrixXd & V,
+      const Eigen::MatrixXi & F,
+      const std::function<Eigen::Affine3d(const double t)> & transform,
+      const size_t steps,
+      const size_t grid_res,
+      const size_t isolevel,
+      Eigen::MatrixXd & SV,
+      Eigen::MatrixXi & SF);
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "swept_volume.cpp"
+#endif
+
+#endif

include/igl/doublearea.cpp

@@ -7,6 +7,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "doublearea.h"
 #include "edge_lengths.h"
+#include "parallel_for.h"
 #include "sort.h"
 #include <cassert>
 #include <iostream>
@@ -148,28 +149,21 @@ IGL_INLINE void igl::doublearea(
   assert((Index)s.rows() == m);
   // resize output
   dblA.resize(l.rows(),1);
-  // Minimum number of iterms per openmp thread
-  #ifndef IGL_OMP_MIN_VALUE
-  #  define IGL_OMP_MIN_VALUE 1000
-  #endif
-  #pragma omp parallel for if (m>IGL_OMP_MIN_VALUE)
-  for(Index i = 0;i<m;i++)
-  {
-    //// Heron's formula for area
-    //const typename Derivedl::Scalar arg =
-    //  s(i)*(s(i)-l(i,0))*(s(i)-l(i,1))*(s(i)-l(i,2));
-    //assert(arg>=0);
-    //dblA(i) = 2.0*sqrt(arg);
-    // Kahan's Heron's formula
-    const typename Derivedl::Scalar arg =
-      (l(i,0)+(l(i,1)+l(i,2)))*
-      (l(i,2)-(l(i,0)-l(i,1)))*
-      (l(i,2)+(l(i,0)-l(i,1)))*
-      (l(i,0)+(l(i,1)-l(i,2)));
-    dblA(i) = 2.0*0.25*sqrt(arg);
-    assert( l(i,2) - (l(i,0)-l(i,1)) && "FAILED KAHAN'S ASSERTION");
-    assert(dblA(i) == dblA(i) && "DOUBLEAREA() PRODUCED NaN");
-  }
+  parallel_for(
+    m,
+    [&l,&dblA](const int i)
+    {
+      // Kahan's Heron's formula
+      const typename Derivedl::Scalar arg =
+        (l(i,0)+(l(i,1)+l(i,2)))*
+        (l(i,2)-(l(i,0)-l(i,1)))*
+        (l(i,2)+(l(i,0)-l(i,1)))*
+        (l(i,0)+(l(i,1)-l(i,2)));
+      dblA(i) = 2.0*0.25*sqrt(arg);
+      assert( l(i,2) - (l(i,0)-l(i,1)) && "FAILED KAHAN'S ASSERTION");
+      assert(dblA(i) == dblA(i) && "DOUBLEAREA() PRODUCED NaN");
+    },
+    1000l);
 }
 
 template <typename DerivedV, typename DerivedF, typename DeriveddblA>

include/igl/edge_lengths.cpp

@@ -6,6 +6,7 @@
 // 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 "edge_lengths.h"
+#include "parallel_for.h"
 #include <iostream>
 
 template <typename DerivedV, typename DerivedF, typename DerivedL>
@@ -16,10 +17,6 @@ IGL_INLINE void igl::edge_lengths(
 {
   using namespace std;
   const int m = F.rows();
-  // Minimum number of iterms per openmp thread
-#ifndef IGL_OMP_MIN_VALUE
-#  define IGL_OMP_MIN_VALUE 1000
-#endif
   switch(F.cols())
   {
     case 2:
@@ -35,29 +32,33 @@ IGL_INLINE void igl::edge_lengths(
     {
       L.resize(m,3);
       // loop over faces
-      #pragma omp parallel for if (m>IGL_OMP_MIN_VALUE)
-      for(int i = 0;i<m;i++)
-      {
-        L(i,0) = (V.row(F(i,1))-V.row(F(i,2))).norm();
-        L(i,1) = (V.row(F(i,2))-V.row(F(i,0))).norm();
-        L(i,2) = (V.row(F(i,0))-V.row(F(i,1))).norm();
-      }
+      parallel_for(
+        m,
+        [&V,&F,&L](const int i)
+        {
+          L(i,0) = (V.row(F(i,1))-V.row(F(i,2))).norm();
+          L(i,1) = (V.row(F(i,2))-V.row(F(i,0))).norm();
+          L(i,2) = (V.row(F(i,0))-V.row(F(i,1))).norm();
+        },
+        1000);
       break;
     }
     case 4:
     {
       L.resize(m,6);
       // loop over faces
-      #pragma omp parallel for if (m>IGL_OMP_MIN_VALUE)
-      for(int i = 0;i<m;i++)
-      {
-        L(i,0) = (V.row(F(i,3))-V.row(F(i,0))).norm();
-        L(i,1) = (V.row(F(i,3))-V.row(F(i,1))).norm();
-        L(i,2) = (V.row(F(i,3))-V.row(F(i,2))).norm();
-        L(i,3) = (V.row(F(i,1))-V.row(F(i,2))).norm();
-        L(i,4) = (V.row(F(i,2))-V.row(F(i,0))).norm();
-        L(i,5) = (V.row(F(i,0))-V.row(F(i,1))).norm();
-      }
+      parallel_for(
+        m,
+        [&V,&F,&L](const int i)
+        {
+          L(i,0) = (V.row(F(i,3))-V.row(F(i,0))).norm();
+          L(i,1) = (V.row(F(i,3))-V.row(F(i,1))).norm();
+          L(i,2) = (V.row(F(i,3))-V.row(F(i,2))).norm();
+          L(i,3) = (V.row(F(i,1))-V.row(F(i,2))).norm();
+          L(i,4) = (V.row(F(i,2))-V.row(F(i,0))).norm();
+          L(i,5) = (V.row(F(i,0))-V.row(F(i,1))).norm();
+        },
+        1000);
       break;
     }
     default:

include/igl/embree/ambient_occlusion.cpp

@@ -58,4 +58,5 @@ template void igl::embree::ambient_occlusion<Eigen::Matrix<double, -1, -1, 0, -1
 template void igl::embree::ambient_occlusion<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, -1, 1, 0, -1, 1> >(igl::embree::EmbreeIntersector const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
 template void igl::embree::ambient_occlusion<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<double, -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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
 template void igl::embree::ambient_occlusion<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<double, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -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<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
+template void igl::embree::ambient_occlusion<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<double, -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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/embree/ambient_occlusion.h

@@ -7,7 +7,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #ifndef IGL_EMBREE_AMBIENT_OCCLUSION_H
 #define IGL_EMBREE_AMBIENT_OCCLUSION_H
-#include <igl/igl_inline.h>
+#include "../igl_inline.h"
 #include <Eigen/Core>
 namespace igl
 {

include/igl/embree/reorient_facets_raycast.cpp

@@ -1,9 +1,9 @@
 // This file is part of libigl, a simple c++ geometry processing library.
-// 
+//
 // Copyright (C) 2013 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 "reorient_facets_raycast.h"
 #include "../per_face_normals.h"
@@ -17,8 +17,8 @@
 #include <limits>
 
 template <
-  typename DerivedV, 
-  typename DerivedF, 
+  typename DerivedV,
+  typename DerivedF,
   typename DerivedI,
   typename DerivedC>
 IGL_INLINE void igl::embree::reorient_facets_raycast(
@@ -36,37 +36,37 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
   using namespace std;
   assert(F.cols() == 3);
   assert(V.cols() == 3);
-  
+
   // number of faces
   const int m = F.rows();
-  
+
   MatrixXi FF = F;
   if (facet_wise) {
     C.resize(m);
     for (int i = 0; i < m; ++i) C(i) = i;
-  
+
   } else {
     if (is_verbose) cout << "extracting patches... ";
     bfs_orient(F,FF,C);
   }
   if (is_verbose) cout << (C.maxCoeff() + 1)  << " components. ";
-  
+
   // number of patches
   const int num_cc = C.maxCoeff()+1;
-  
+
   // Init Embree
   EmbreeIntersector ei;
   ei.init(V.template cast<float>(),FF);
-  
+
   // face normal
   MatrixXd N;
   per_face_normals(V,FF,N);
-  
+
   // face area
   Matrix<typename DerivedV::Scalar,Dynamic,1> A;
   doublearea(V,FF,A);
   double area_total = A.sum();
-  
+
   // determine number of rays per component according to its area
   VectorXd area_per_component;
   area_per_component.setZero(num_cc);
@@ -80,7 +80,7 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
     num_rays_per_component(c) = max<int>(static_cast<int>(rays_total * area_per_component(c) / area_total), rays_minimum);
   }
   rays_total = num_rays_per_component.sum();
-  
+
   // generate all the rays
   if (is_verbose) cout << "generating rays... ";
   uniform_real_distribution<float> rdist;
@@ -147,12 +147,12 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
     }
   }
   if (is_verbose) cout << ray_face.size()  << " rays. ";
-  
+
   // per component voting: first=front, second=back
   vector<pair<float, float>> C_vote_distance(num_cc, make_pair(0, 0));      // sum of distance between ray origin and intersection
   vector<pair<int  , int  >> C_vote_infinity(num_cc, make_pair(0, 0));      // number of rays reaching infinity
   vector<pair<int  , int  >> C_vote_parity(num_cc, make_pair(0, 0));        // sum of parity count for each ray
-  
+
   if (is_verbose) cout << "shooting rays... ";
 #pragma omp parallel for
   for (int i = 0; i < (int)ray_face.size(); ++i)
@@ -161,7 +161,7 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
     Vector3f o = ray_ori [i];
     Vector3f d = ray_dir [i];
     int c = C(f);
-    
+
     // shoot ray toward front & back
     vector<Hit> hits_front;
     vector<Hit> hits_back;
@@ -171,13 +171,13 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
     ei.intersectRay(o, -d, hits_back , num_rays_back );
     if (!hits_front.empty() && hits_front[0].id == f) hits_front.erase(hits_front.begin());
     if (!hits_back .empty() && hits_back [0].id == f) hits_back .erase(hits_back .begin());
-    
+
     if (use_parity) {
 #pragma omp atomic
       C_vote_parity[c].first  += hits_front.size() % 2;
 #pragma omp atomic
       C_vote_parity[c].second += hits_back .size() % 2;
-    
+
     } else {
       if (hits_front.empty())
       {
@@ -187,7 +187,7 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
 #pragma omp atomic
         C_vote_distance[c].first += hits_front[0].t;
       }
-    
+
       if (hits_back.empty())
       {
 #pragma omp atomic
@@ -198,14 +198,14 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
       }
     }
   }
-  
+
   I.resize(m);
   for(int f = 0; f < m; ++f)
   {
     int c = C(f);
     if (use_parity) {
       I(f) = C_vote_parity[c].first > C_vote_parity[c].second ? 1 : 0;      // Ideally, parity for the front/back side should be 1/0 (i.e., parity sum for all rays should be smaller on the front side)
-    
+
     } else {
       I(f) = (C_vote_infinity[c].first == C_vote_infinity[c].second && C_vote_distance[c].first <  C_vote_distance[c].second) ||
               C_vote_infinity[c].first <  C_vote_infinity[c].second
@@ -219,8 +219,8 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
 }
 
 template <
-  typename DerivedV, 
-  typename DerivedF, 
+  typename DerivedV,
+  typename DerivedF,
   typename DerivedFF,
   typename DerivedI>
 IGL_INLINE void igl::embree::reorient_facets_raycast(
@@ -255,4 +255,5 @@ IGL_INLINE void igl::embree::reorient_facets_raycast(
 // Explicit template specialization
 template void igl::embree::reorient_facets_raycast<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::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<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 template void igl::embree::reorient_facets_raycast<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<bool, -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&, int, int, bool, bool, bool, Eigen::PlainObjectBase<Eigen::Matrix<bool, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+template void igl::embree::reorient_facets_raycast<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::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&, int, int, bool, bool, bool, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 #endif

include/igl/fit_rotations.cpp

@@ -71,7 +71,7 @@ IGL_INLINE void igl::fit_rotations_planar(
   using namespace std;
   const int dim = S.cols();
   const int nr = S.rows()/dim;
-  assert(dim == 2 && "_planar input should be 2D");
+  //assert(dim == 2 && "_planar input should be 2D");
   assert(nr * dim == S.rows());
 
   // resize output
@@ -105,6 +105,7 @@ IGL_INLINE void igl::fit_rotations_planar(
 #endif  
 
     // Not sure why polar_dec computes transpose...
+    R.block(0,r*dim,dim,dim).setIdentity();
     R.block(0,r*dim,2,2) = ri.transpose();
   }
 }

include/igl/flood_fill.cpp

@@ -0,0 +1,81 @@
+#include "flood_fill.h"
+#include <limits>
+
+IGL_INLINE void igl::flood_fill(
+  const Eigen::RowVector3i & res,
+  Eigen::VectorXd & S)
+{
+  using namespace Eigen;
+  using namespace std;
+  const auto flood = [&res,&S] (
+     const int xi, 
+     const int yi, 
+     const int zi,
+     const int signed_xi, 
+     const int signed_yi, 
+     const int signed_zi,
+     const double s)
+    {
+      // flood fill this value back on this row
+      for(int bxi = xi;signed_xi<--bxi;)
+      {
+        S(bxi+res(0)*(yi + res(1)*zi)) = s;
+      }
+      // flood fill this value back on any previous rows
+      for(int byi = yi;signed_yi<--byi;)
+      {
+        for(int xi = 0;xi<res(0);xi++)
+        {
+          S(xi+res(0)*(byi + res(1)*zi)) = s;
+        }
+      }
+      // flood fill this value back on any previous "sheets"
+      for(int bzi = zi;signed_zi<--bzi;)
+      {
+        for(int yi = 0;yi<res(1);yi++)
+        {
+          for(int xi = 0;xi<res(0);xi++)
+          {
+            S(xi+res(0)*(yi + res(1)*bzi)) = s;
+          }
+        }
+      }
+    };
+  int signed_zi = -1;
+  double s = numeric_limits<double>::quiet_NaN();
+  for(int zi = 0;zi<res(2);zi++)
+  {
+    int signed_yi = -1;
+    if(zi != 0)
+    {
+      s = S(0+res(0)*(0 + res(1)*(zi-1)));
+    }
+    for(int yi = 0;yi<res(1);yi++)
+    {
+      // index of last signed item on this row
+      int signed_xi = -1;
+      if(yi != 0)
+      {
+        s = S(0+res(0)*(yi-1 + res(1)*zi));
+      }
+      for(int xi = 0;xi<res(0);xi++)
+      {
+        int i = xi+res(0)*(yi + res(1)*zi);
+        if(S(i)!=S(i))
+        {
+          if(s == s)
+          {
+            S(i) = s;
+          }
+          continue;
+        }
+        s = S(i);
+        flood(xi,yi,zi,signed_xi,signed_yi,signed_zi,s);
+        // remember this position
+        signed_xi = xi;
+        signed_yi = yi;
+        signed_zi = zi;
+      }
+    }
+  }
+}

include/igl/flood_fill.h

@@ -0,0 +1,23 @@
+#ifndef IGL_FLOOD_FILL_H
+#define IGL_FLOOD_FILL_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+namespace igl
+{
+  // Given a 3D array with sparse non-nan (number?) data fill in the NaNs via
+  // flood fill. This should ensure that, e.g., if data near 0 always has a band
+  // (surface) of numbered and signed data, then components of nans will be
+  // correctly signed.
+  //
+  // Inputs:
+  //   res  3-long list of dimensions of grid
+  //   S  res(0)*res(1)*res(2)  list of scalar values (with (many) nans), see
+  //     output
+  // Outputs:
+  //   S  flood fill data in place
+  IGL_INLINE void flood_fill(const Eigen::RowVector3i& res, Eigen::VectorXd& S);
+}
+#ifndef IGL_STATIC_LIBRARY
+#  include "flood_fill.cpp"
+#endif
+#endif

include/igl/grid.cpp

@@ -0,0 +1,23 @@
+#include "grid.h"
+
+IGL_INLINE void igl::grid(const Eigen::RowVector3i & res, Eigen::MatrixXd & GV)
+{
+  using namespace Eigen;
+  GV.resize(res(0)*res(1)*res(2),3);
+  for(int zi = 0;zi<res(2);zi++)
+  {
+    const auto lerp = 
+      [&](const double di, const int d)->double{return di/(double)(res(d)-1);};
+    const double z = lerp(zi,2);
+    for(int yi = 0;yi<res(1);yi++)
+    {
+      const double y = lerp(yi,1);
+      for(int xi = 0;xi<res(0);xi++)
+      {
+        const double x = lerp(xi,0);
+        GV.row(xi+res(0)*(yi + res(1)*zi)) = RowVector3d(x,y,z);
+      }
+    }
+  }
+}
+

include/igl/grid.h

@@ -0,0 +1,20 @@
+#ifndef IGL_GRID_H
+#define IGL_GRID_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+namespace igl
+{
+  // Construct vertices of a regular grid, suitable for input to
+  // `igl::marching_cubes`
+  //
+  // Inputs:
+  //   res  3-long list of number of vertices along the x y and z dimensions
+  // Outputs:
+  //   GV  res(0)*res(1)*res(2) by 3 list of mesh vertex positions.
+  //   
+  IGL_INLINE void grid(const Eigen::RowVector3i & res, Eigen::MatrixXd & GV);
+}
+#ifndef IGL_STATIC_LIBRARY
+#  include "grid.cpp"
+#endif
+#endif 

include/igl/infinite_cost_stopping_condition.cpp

@@ -55,7 +55,7 @@ IGL_INLINE void igl::infinite_cost_stopping_condition(
       Eigen::RowVectorXd p;
       double cost;
       cost_and_placement(e,V,F,E,EMAP,EF,EI,cost,p);
-      return isinf(cost);
+      return std::isinf(cost);
     };
 }
 

include/igl/internal_angles.cpp

@@ -8,6 +8,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "internal_angles.h"
 #include "edge_lengths.h"
+#include "parallel_for.h"
 #include "get_seconds.h"
 
 template <typename DerivedV, typename DerivedF, typename DerivedK>
@@ -70,28 +71,19 @@ IGL_INLINE void igl::internal_angles(
   assert(L.cols() == 3 && "Edge-lengths should come from triangles");
   const Index m = L.rows();
   K.resize(m,3);
-  //for(int d = 0;d<3;d++)
-  //{
-  //  const auto & s1 = L.col(d).array();
-  //  const auto & s2 = L.col((d+1)%3).array();
-  //  const auto & s3 = L.col((d+2)%3).array();
-  //  K.col(d) = ((s3.square() + s2.square() - s1.square())/(2.*s3*s2)).acos();
-  //}
-  // Minimum number of iterms per openmp thread
-  #ifndef IGL_OMP_MIN_VALUE
-  #  define IGL_OMP_MIN_VALUE 1000
-  #endif
-  #pragma omp parallel for if (m>IGL_OMP_MIN_VALUE)
-  for(Index f = 0;f<m;f++)
-  {
-    for(size_t d = 0;d<3;d++)
+  parallel_for(
+    m,
+    [&L,&K](const Index f)
     {
-      const auto & s1 = L(f,d);
-      const auto & s2 = L(f,(d+1)%3);
-      const auto & s3 = L(f,(d+2)%3);
-      K(f,d) = acos((s3*s3 + s2*s2 - s1*s1)/(2.*s3*s2));
-    }
-  }
+      for(size_t d = 0;d<3;d++)
+      {
+        const auto & s1 = L(f,d);
+        const auto & s2 = L(f,(d+1)%3);
+        const auto & s3 = L(f,(d+2)%3);
+        K(f,d) = acos((s3*s3 + s2*s2 - s1*s1)/(2.*s3*s2));
+      }
+    },
+    1000l);
 }
 
 #ifdef IGL_STATIC_LIBRARY

include/igl/list_to_matrix.cpp

@@ -118,6 +118,16 @@ IGL_INLINE bool igl::list_to_matrix(const std::vector<T > & V,Eigen::PlainObject
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template bool igl::list_to_matrix<int, Eigen::Matrix<int, -1, 2, 0, -1, 2> >(std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 2, 0, -1, 2> >&);
+// generated by autoexplicit.sh
+template bool igl::list_to_matrix<double, Eigen::Matrix<double, -1, 2, 0, -1, 2> >(std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> >&);
+// generated by autoexplicit.sh
+template bool igl::list_to_matrix<double, Eigen::Matrix<double, 4, 1, 0, 4, 1> >(std::vector<double, std::allocator<double> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 4, 1, 0, 4, 1> >&);
+// generated by autoexplicit.sh
+template bool igl::list_to_matrix<double, Eigen::Matrix<double, 3, 1, 0, 3, 1> >(std::vector<double, std::allocator<double> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 3, 1, 0, 3, 1> >&);
+// generated by autoexplicit.sh
+template bool igl::list_to_matrix<double, Eigen::Matrix<double, 2, 1, 0, 2, 1> >(std::vector<double, std::allocator<double> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> >&);
 template bool igl::list_to_matrix<double, Eigen::Matrix<double, -1, 3, 0, -1, 3> >(std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&);
 template bool igl::list_to_matrix<double, Eigen::Matrix<double, -1, 3, 1, -1, 3> >(std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 1, -1, 3> >&);
 template bool igl::list_to_matrix<double, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(std::vector<std::vector<double, std::allocator<double> >, std::allocator<std::vector<double, std::allocator<double> > > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);

include/igl/parallel_for.h

@@ -0,0 +1,180 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+// 
+// Copyright (C) 2016 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 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_PARALLEL_FOR_H
+#define IGL_PARALLEL_FOR_H
+#include "igl_inline.h"
+#include <functional>
+
+namespace igl
+{
+  // PARALLEL_FOR Functional implementation of a basic, open-mp style, parallel
+  // for loop. If the inner block of a for-loop can be rewritten/encapsulated in
+  // a single (anonymous/lambda) function call `func` so that the serial code
+  // looks like:
+  // 
+  //     for(int i = 0;i<loop_size;i++)
+  //     {
+  //       func(i);
+  //     }
+  //
+  // then `parallel_for(loop_size,func,min_parallel)` will use as many threads as
+  // available on the current hardware to parallelize this for loop so long as
+  // loop_size<min_parallel, otherwise it will just use a serial for loop.
+  //
+  // Inputs:
+  //   loop_size  number of iterations. I.e. for(int i = 0;i<loop_size;i++) ...
+  //   func  function handle taking iteration index as only arguement to compute
+  //     inner block of for loop I.e. for(int i ...){ func(i); }
+  //   min_parallel  min size of loop_size such that parallel (non-serial)
+  //     thread pooling should be attempted {0}
+  // Returns true iff thread pool was invoked
+  template<typename Index, typename FunctionType >
+  inline bool parallel_for(
+    const Index loop_size, 
+    const FunctionType & func,
+    const size_t min_parallel=0);
+  // PARALLEL_FOR Functional implementation of an open-mp style, parallel for
+  // loop with accumulation. For example, serial code separated into n chunks
+  // (each to be parallelized with a thread) might look like:
+  //     
+  //     Eigen::VectorXd S;
+  //     const auto & prep_func = [&S](int n){ S = Eigen:VectorXd::Zero(n); };
+  //     const auto & func = [&X,&S](int i, int t){ S(t) += X(i); };
+  //     const auto & accum_func = [&S,&sum](int t){ sum += S(t); };
+  //     prep_func(n);
+  //     for(int i = 0;i<loop_size;i++)
+  //     {
+  //       func(i,i%n);
+  //     }
+  //     double sum = 0;
+  //     for(int t = 0;t<n;t++)
+  //     {
+  //       accum_func(t);
+  //     }
+  // 
+  // Inputs:
+  //   loop_size  number of iterations. I.e. for(int i = 0;i<loop_size;i++) ...
+  //   prep_func function handle taking n >= number of threads as only
+  //     argument 
+  //   func  function handle taking iteration index i and thread id t as only
+  //     arguements to compute inner block of for loop I.e. 
+  //     for(int i ...){ func(i,t); }
+  //   accum_func  function handle taking thread index as only argument, to be
+  //     called after all calls of func, e.g., for serial accumulation across
+  //     all n (potential) threads, see n in description of prep_func.
+  //   min_parallel  min size of loop_size such that parallel (non-serial)
+  //     thread pooling should be attempted {0}
+  // Returns true iff thread pool was invoked
+  template<
+    typename Index, 
+    typename PrepFunctionType, 
+    typename FunctionType, 
+    typename AccumFunctionType 
+    >
+  inline bool parallel_for(
+    const Index loop_size, 
+    const PrepFunctionType & prep_func,
+    const FunctionType & func,
+    const AccumFunctionType & accum_func,
+    const size_t min_parallel=0);
+}
+
+// Implementation
+
+#include <cmath>
+#include <cassert>
+#include <thread>
+#include <vector>
+#include <algorithm>
+
+template<typename Index, typename FunctionType >
+inline bool igl::parallel_for(
+  const Index loop_size, 
+  const FunctionType & func,
+  const size_t min_parallel)
+{
+  using namespace std;
+  // no op preparation/accumulation
+  const auto & no_op = [](const size_t /*n/t*/){};
+  // two-parameter wrapper ignoring thread id
+  const auto & wrapper = [&func](Index i,size_t /*t*/){ func(i); };
+  return parallel_for(loop_size,no_op,wrapper,no_op,min_parallel);
+}
+
+template<
+  typename Index, 
+  typename PreFunctionType,
+  typename FunctionType, 
+  typename AccumFunctionType>
+inline bool igl::parallel_for(
+  const Index loop_size, 
+  const PreFunctionType & prep_func,
+  const FunctionType & func,
+  const AccumFunctionType & accum_func,
+  const size_t min_parallel)
+{
+  assert(loop_size>=0);
+  if(loop_size==0) return false;
+  // Estimate number of threads in the pool
+  // http://ideone.com/Z7zldb
+  const static size_t sthc = std::thread::hardware_concurrency();
+  const size_t nthreads = loop_size<min_parallel?0:(sthc==0?8:sthc);
+  if(nthreads==0)
+  {
+    // serial
+    prep_func(1);
+    for(Index i = 0;i<loop_size;i++) func(i,0);
+    accum_func(0);
+    return false;
+  }else
+  {
+    // Size of a slice for the range functions
+    Index slice = 
+      std::max(
+        (Index)std::round((loop_size+1)/static_cast<double>(nthreads)),(Index)1);
+ 
+    // [Helper] Inner loop
+    const auto & range = [&func](const Index k1, const Index k2, const size_t t)
+    {
+      for(Index k = k1; k < k2; k++) func(k,t);
+    };
+    prep_func(nthreads);
+    // Create pool and launch jobs
+    std::vector<std::thread> pool;
+    pool.reserve(nthreads);
+    // Inner range extents
+    Index i1 = 0;
+    Index i2 = std::min(0 + slice, loop_size);
+    {
+      size_t t = 0;
+      for (; t+1 < nthreads && i1 < loop_size; ++t)
+      {
+        pool.emplace_back(range, i1, i2, t);
+        i1 = i2;
+        i2 = std::min(i2 + slice, loop_size);
+      }
+      if (i1 < loop_size)
+      {
+        pool.emplace_back(range, i1, loop_size, t);
+      }
+    }
+    // Wait for jobs to finish
+    for (std::thread &t : pool) if (t.joinable()) t.join();
+    // Accumulate across threads
+    for(size_t t = 0;t<nthreads;t++)
+    {
+      accum_func(t);
+    }
+    return true;
+  }
+}
+ 
+//#ifndef IGL_STATIC_LIBRARY
+//#include "parallel_for.cpp"
+//#endif
+#endif

include/igl/per_vertex_normals.cpp

@@ -10,6 +10,7 @@
 #include "get_seconds.h"
 #include "per_face_normals.h"
 #include "doublearea.h"
+#include "parallel_for.h"
 #include "internal_angles.h"
 
 template <typename DerivedV, typename DerivedF>
@@ -68,23 +69,34 @@ IGL_INLINE void igl::per_vertex_normals(
   }
 
   // loop over faces
-  const int Frows = F.rows();
-//// Minimum number of iterms per openmp thread
-//#ifndef IGL_OMP_MIN_VALUE
-//#  define IGL_OMP_MIN_VALUE 1000
-//#endif
-//#pragma omp parallel for if (Frows>IGL_OMP_MIN_VALUE)
-  for(int i = 0; i < Frows;i++)
+  for(int i = 0;i<F.rows();i++)
   {
     // throw normal at each corner
     for(int j = 0; j < 3;j++)
     {
-      // Q: Does this need to be critical?
-      // A: Yes. Different (i,j)'s could produce the same F(i,j)
-//#pragma omp critical
-      N.row(F(i,j)) += W(i,j)*FN.row(i);
+      N.row(F(i,j)) += W(i,j) * FN.row(i);
     }
   }
+
+  //// loop over faces
+  //std::mutex critical;
+  //std::vector<DerivedN> NN;
+  //parallel_for(
+  //  F.rows(),
+  //  [&NN,&N](const size_t n){ NN.resize(n,DerivedN::Zero(N.rows(),3));},
+  //  [&F,&W,&FN,&NN,&critical](const int i, const size_t t)
+  //  {
+  //    // throw normal at each corner
+  //    for(int j = 0; j < 3;j++)
+  //    {
+  //      // Q: Does this need to be critical?
+  //      // A: Yes. Different (i,j)'s could produce the same F(i,j)
+  //      NN[t].row(F(i,j)) += W(i,j) * FN.row(i);
+  //    }
+  //  },
+  //  [&N,&NN](const size_t t){ N += NN[t]; },
+  //  1000l);
+
   // take average via normalization
   N.rowwise().normalize();
 }

include/igl/png/readPNG.cpp

@@ -0,0 +1,47 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2016 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can
+// obtain one at http://mozilla.org/MPL/2.0/.
+#include "readPNG.h"
+#include <stb_image.h>
+
+IGL_INLINE bool igl::png::readPNG(
+  const std::string png_file,
+  Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& R,
+  Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& G,
+  Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& B,
+  Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& A
+)
+{
+  int cols,rows,n;
+  unsigned char *data = stbi_load(png_file.c_str(), &cols, &rows, &n, 4);
+  if(data == NULL) {
+    return false;
+  }
+
+  R.resize(cols,rows);
+  G.resize(cols,rows);
+  B.resize(cols,rows);
+  A.resize(cols,rows);
+
+  for (unsigned i=0; i<rows; ++i) {
+    for (unsigned j=0; j<cols; ++j) {
+      R(j,rows-1-i) = data[4*(j + cols * i) + 0];
+      G(j,rows-1-i) = data[4*(j + cols * i) + 1];
+      B(j,rows-1-i) = data[4*(j + cols * i) + 2];
+      A(j,rows-1-i) = data[4*(j + cols * i) + 3];
+    }
+  }
+
+  stbi_image_free(data);
+
+  return true;
+}
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template specialization
+// generated by autoexplicit.sh
+#endif

include/igl/png/readPNG.h

@@ -0,0 +1,39 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2016 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_PNG_READ_PNG_H
+#define IGL_PNG_READ_PNG_H
+#include "../igl_inline.h"
+#include <Eigen/Core>
+#include <string>
+
+namespace igl
+{
+  namespace png
+  {
+    // Read an image from a .png file into 4 memory buffers
+    //
+    // Input:
+    //  png_file  path to .png file
+    // Output:
+    //  R,G,B,A texture channels
+    // Returns true on success, false on failure
+    //
+    IGL_INLINE bool readPNG(const std::string png_file,
+    Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& R,
+    Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& G,
+    Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& B,
+    Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& A
+    );
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "readPNG.cpp"
+#endif
+
+#endif

include/igl/png/writePNG.cpp

@@ -0,0 +1,46 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2016 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can
+// obtain one at http://mozilla.org/MPL/2.0/.
+#include "writePNG.h"
+#include <stb_image_write.h>
+#include <vector>
+
+IGL_INLINE bool igl::png::writePNG(
+  const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& R,
+  const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& G,
+  const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& B,
+  const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& A,
+  const std::string png_file
+)
+{
+  assert((R.rows() == G.rows()) && (G.rows() == B.rows()) && (B.rows() == A.rows()));
+  assert((R.cols() == G.cols()) && (G.cols() == B.cols()) && (B.cols() == A.cols()));
+
+  const int comp = 4;                                  // 4 Channels Red, Green, Blue, Alpha
+  const int stride_in_bytes = R.rows()*comp;           // Lenght of one row in bytes
+  std::vector<unsigned char> data(R.size()*comp,0);     // The image itself;
+
+  for (unsigned i = 0; i<R.rows();++i)
+  {
+    for (unsigned j = 0; j < R.cols(); ++j)
+    {
+        data[(j * R.rows() * comp) + (i * comp) + 0] = R(i,R.cols()-1-j);
+        data[(j * R.rows() * comp) + (i * comp) + 1] = G(i,R.cols()-1-j);
+        data[(j * R.rows() * comp) + (i * comp) + 2] = B(i,R.cols()-1-j);
+        data[(j * R.rows() * comp) + (i * comp) + 3] = A(i,R.cols()-1-j);
+    }
+  }
+
+  stbi_write_png(png_file.c_str(), R.rows(), R.cols(), comp, data.data(), stride_in_bytes);
+
+  return true;
+}
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template specialization
+// generated by autoexplicit.sh
+#endif

include/igl/png/writePNG.h

@@ -0,0 +1,41 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2016 Daniele Panozzo <daniele.panozzo@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can
+// obtain one at http://mozilla.org/MPL/2.0/.
+#ifndef IGL_PNG_WRITE_PNG_H
+#define IGL_PNG_WRITE_PNG_H
+#include "../igl_inline.h"
+#include <Eigen/Core>
+#include <string>
+
+namespace igl
+{
+  namespace png
+  {
+    // Writes an image to a png file
+    //
+    // Input:
+    //  R,G,B,A texture channels
+    // Output:
+    //  png_file  path to .png file
+    // Returns true on success, false on failure
+    //
+    IGL_INLINE bool writePNG
+    (
+      const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& R,
+      const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& G,
+      const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& B,
+      const Eigen::Matrix<unsigned char,Eigen::Dynamic,Eigen::Dynamic>& A,
+      const std::string png_file
+    );
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "writePNG.cpp"
+#endif
+
+#endif

include/igl/project_to_line.cpp

@@ -120,6 +120,8 @@ IGL_INLINE void igl::project_to_line(
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template void igl::project_to_line<Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 1, 0, 1, 1>, Eigen::Matrix<double, 1, 1, 0, 1, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&);
 template void igl::project_to_line<double>(double, double, double, double, double, double, double, double, double, double&, double&);
 template void igl::project_to_line<double>(double, double, double, double, double, double, double, double, double, double&, double&,double&,double&, double&);
 template void igl::project_to_line<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, 1, 0, 1, 1>, Eigen::Matrix<double, 1, 1, 0, 1, 1> >(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, 1, 0, 1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&);

include/igl/project_to_line_segment.cpp

@@ -43,5 +43,7 @@ IGL_INLINE void igl::project_to_line_segment(
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template void igl::project_to_line_segment<Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, 1, 1, 0, 1, 1>, Eigen::Matrix<double, 1, 1, 0, 1, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&);
 template void igl::project_to_line_segment<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, 1, 0, 1, 1>, Eigen::Matrix<double, 1, 1, 0, 1, 1> >(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, 1, 0, 1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 1, 0, 1, 1> >&);
 #endif

include/igl/random_dir.cpp

@@ -1,9 +1,9 @@
 // This file is part of libigl, a simple c++ geometry processing library.
-// 
+//
 // Copyright (C) 2013 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 "random_dir.h"
 #include <igl/PI.h>
@@ -18,14 +18,14 @@ IGL_INLINE Eigen::Vector3d igl::random_dir()
   double r = sqrt(1.0-z*z);
   double x = r * cos(t);
   double y = r * sin(t);
-  return Vector3d(x,y,z); 
+  return Vector3d(x,y,z);
 }
 
 IGL_INLINE Eigen::MatrixXd igl::random_dir_stratified(const int n)
 {
   using namespace Eigen;
   using namespace std;
-  const double m = floor(sqrt(double(n)));
+  const double m = std::floor(sqrt(double(n)));
   MatrixXd N(n,3);
   int row = 0;
   for(int i = 0;i<m;i++)

include/igl/repdiag.cpp

@@ -18,7 +18,7 @@ IGL_INLINE void igl::repdiag(
   using namespace Eigen;
   int m = A.rows();
   int n = A.cols();
-
+#if false
   vector<Triplet<T> > IJV;
   IJV.reserve(A.nonZeros()*d);
   // Loop outer level
@@ -35,31 +35,32 @@ IGL_INLINE void igl::repdiag(
   }
   B.resize(m*d,n*d);
   B.setFromTriplets(IJV.begin(),IJV.end());
-  
-
-  // Q: Why is this **Very** slow?
-
-  //int m = A.rows();
-  //int n = A.cols();
-
-  //B.resize(m*d,n*d);
-  //// Reserve enough space for new non zeros
-  //B.reserve(d*A.nonZeros());
-
-  //// loop over reps
-  //for(int i=0;i<d;i++)
-  //{
-  //  // Loop outer level
-  //  for (int k=0; k<A.outerSize(); ++k)
-  //  {
-  //    // loop inner level
-  //    for (typename Eigen::SparseMatrix<T>::InnerIterator it(A,k); it; ++it)
-  //    {
-  //      B.insert(i*m+it.row(),i*n+it.col()) = it.value();
-  //    }
-  //  }
-  //}
-  //B.makeCompressed();
+#else
+  // This will not work for RowMajor
+  B.resize(m*d,n*d);
+  Eigen::VectorXi per_col = Eigen::VectorXi::Zero(n*d);
+  for (int k=0; k<A.outerSize(); ++k)
+  {
+    for (typename Eigen::SparseMatrix<T>::InnerIterator it(A,k); it; ++it)
+    {
+      for(int r = 0;r<d;r++) per_col(n*r + k)++;
+    }
+  }
+  B.reserve(per_col);
+  for(int r = 0;r<d;r++)
+  {
+    const int mr = m*r;
+    const int nr = n*r;
+    for (int k=0; k<A.outerSize(); ++k)
+    {
+      for (typename Eigen::SparseMatrix<T>::InnerIterator it(A,k); it; ++it)
+      {
+        B.insert(it.row()+mr,k+nr) = it.value();
+      }
+    }
+  }
+  B.makeCompressed();
+#endif
 }
 
 template <typename T>

include/igl/resolve_duplicated_faces.cpp

@@ -86,4 +86,5 @@ IGL_INLINE void igl::resolve_duplicated_faces(
 #ifdef IGL_STATIC_LIBRARY
 template void igl::resolve_duplicated_faces<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> >&);
 template void igl::resolve_duplicated_faces<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<long, -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<long, -1, 1, 0, -1, 1> >&);
+template void igl::resolve_duplicated_faces<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> >&);
 #endif

include/igl/slice.cpp

@@ -279,4 +279,6 @@ template Eigen::Matrix<double, -1, -1, 0, -1, -1> igl::slice<Eigen::Matrix<doubl
 template void igl::slice<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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 template void igl::slice<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::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::slice<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::Matrix<int, -1, 1, 0, -1, 1> const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::slice<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> const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, int, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
+template void igl::slice<Eigen::PlainObjectBase<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> > >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, int, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/slice_mask.cpp

@@ -55,8 +55,8 @@ IGL_INLINE void igl::slice_mask(
     case 1:
     {
       const int ym = R.count();
-      Y.resize(ym,X.cols());
       assert(X.rows() == R.size() && "X.rows() should match R.size()");
+      Y.resize(ym,X.cols());
       {
         int yi = 0;
         for(int i = 0;i<X.rows();i++)

include/igl/sort.cpp

@@ -11,12 +11,11 @@
 #include "reorder.h"
 #include "IndexComparison.h"
 #include "colon.h"
+#include "parallel_for.h"
 
 #include <cassert>
 #include <algorithm>
 #include <iostream>
-#include <thread>
-#include <functional>
 
 template <typename DerivedX, typename DerivedY, typename DerivedIX>
 IGL_INLINE void igl::sort(
@@ -223,98 +222,86 @@ IGL_INLINE void igl::sort3(
     IX.col(2).setConstant(2);// = Eigen::PlainObjectBase<DerivedIX>::Ones (IX.rows(),1);
   }
 
-  const int n = num_outer;
-  const size_t nthreads = n<8000?1:std::thread::hardware_concurrency();
+
+  const auto & inner = [&IX,&Y,&dim,&ascending](const Index & i)
   {
-    std::vector<std::thread> threads(nthreads);
-    for(int t = 0;t<nthreads;t++)
+    YScalar & a = (dim==1 ? Y(0,i) : Y(i,0));
+    YScalar & b = (dim==1 ? Y(1,i) : Y(i,1));
+    YScalar & c = (dim==1 ? Y(2,i) : Y(i,2));
+    Index & ai = (dim==1 ? IX(0,i) : IX(i,0));
+    Index & bi = (dim==1 ? IX(1,i) : IX(i,1));
+    Index & ci = (dim==1 ? IX(2,i) : IX(i,2));
+    if(ascending)
+    {
+      // 123 132 213 231 312 321
+      if(a > b)
+      {
+        std::swap(a,b);
+        std::swap(ai,bi);
+      }
+      // 123 132 123 231 132 231
+      if(b > c)
+      {
+        std::swap(b,c);
+        std::swap(bi,ci);
+        // 123 123 123 213 123 213
+        if(a > b)
+        {
+          std::swap(a,b);
+          std::swap(ai,bi);
+        }
+        // 123 123 123 123 123 123
+      }
+    }else
     {
-      threads[t] = std::thread(std::bind(
-        [&X,&Y,&IX,&dim,&ascending](const int bi, const int ei, const int t)
+      // 123 132 213 231 312 321
+      if(a < b)
       {
-        // loop over columns (or rows)
-        for(int i = bi;i<ei;i++)
+        std::swap(a,b);
+        std::swap(ai,bi);
+      }
+      // 213 312 213 321 312 321
+      if(b < c)
+      {
+        std::swap(b,c);
+        std::swap(bi,ci);
+        // 231 321 231 321 321 321
+        if(a < b)
         {
-          YScalar & a = (dim==1 ? Y(0,i) : Y(i,0));
-          YScalar & b = (dim==1 ? Y(1,i) : Y(i,1));
-          YScalar & c = (dim==1 ? Y(2,i) : Y(i,2));
-          Index & ai = (dim==1 ? IX(0,i) : IX(i,0));
-          Index & bi = (dim==1 ? IX(1,i) : IX(i,1));
-          Index & ci = (dim==1 ? IX(2,i) : IX(i,2));
-          if(ascending)
-          {
-            // 123 132 213 231 312 321
-            if(a > b)
-            {
-              std::swap(a,b);
-              std::swap(ai,bi);
-            }
-            // 123 132 123 231 132 231
-            if(b > c)
-            {
-              std::swap(b,c);
-              std::swap(bi,ci);
-              // 123 123 123 213 123 213
-              if(a > b)
-              {
-                std::swap(a,b);
-                std::swap(ai,bi);
-              }
-              // 123 123 123 123 123 123
-            }
-          }else
-          {
-            // 123 132 213 231 312 321
-            if(a < b)
-            {
-              std::swap(a,b);
-              std::swap(ai,bi);
-            }
-            // 213 312 213 321 312 321
-            if(b < c)
-            {
-              std::swap(b,c);
-              std::swap(bi,ci);
-              // 231 321 231 321 321 321
-              if(a < b)
-              {
-                std::swap(a,b);
-                std::swap(ai,bi);
-              }
-              // 321 321 321 321 321 321
-            }
-          }
+          std::swap(a,b);
+          std::swap(ai,bi);
         }
-      }, t*n/nthreads, (t+1)==nthreads?n:(t+1)*n/nthreads,t));
+        // 321 321 321 321 321 321
+      }
     }
-    std::for_each(threads.begin(),threads.end(),[](std::thread& x){x.join();});
-  }
+  };
+  parallel_for(num_outer,inner,16000);
 }
 
 template <class T>
 IGL_INLINE void igl::sort(
-  const std::vector<T> & unsorted,
-  const bool ascending,
-  std::vector<T> & sorted,
-  std::vector<size_t> & index_map)
+const std::vector<T> & unsorted,
+const bool ascending,
+std::vector<T> & sorted,
+std::vector<size_t> & index_map)
 {
-  // Original unsorted index map
-  index_map.resize(unsorted.size());
-  for(size_t i=0;i<unsorted.size();i++)
-  {
-    index_map[i] = i;
-  }
-  // Sort the index map, using unsorted for comparison
-  std::sort(
-    index_map.begin(),
-    index_map.end(),
-    igl::IndexLessThan<const std::vector<T>& >(unsorted));
+// Original unsorted index map
+index_map.resize(unsorted.size());
+for(size_t i=0;i<unsorted.size();i++)
+{
+  index_map[i] = i;
+}
+// Sort the index map, using unsorted for comparison
+std::sort(
+  index_map.begin(),
+  index_map.end(),
+  igl::IndexLessThan<const std::vector<T>& >(unsorted));
 
-  // if not ascending then reverse
-  if(!ascending)
-  {
-    std::reverse(index_map.begin(),index_map.end());
-  }
+// if not ascending then reverse
+if(!ascending)
+{
+  std::reverse(index_map.begin(),index_map.end());
+}
   // make space for output without clobbering
   sorted.resize(unsorted.size());
   // reorder unsorted into sorted using index map

include/igl/swept_volume_bounding_box.cpp

@@ -0,0 +1,20 @@
+#include "swept_volume_bounding_box.h"
+
+IGL_INLINE void igl::swept_volume_bounding_box(
+  const size_t & n,
+  const std::function<Eigen::RowVector3d(const size_t vi, const double t)> & V,
+  const size_t & steps,
+  Eigen::AlignedBox3d & box)
+{
+  using namespace Eigen;
+  box.setEmpty();
+  const VectorXd t = VectorXd::LinSpaced(steps,0,1);
+  // Find extent over all time steps
+  for(int ti = 0;ti<t.size();ti++)
+  {
+    for(size_t vi = 0;vi<n;vi++)
+    {
+      box.extend(V(vi,t(ti)).transpose());
+    }
+  }
+}

include/igl/swept_volume_bounding_box.h

@@ -0,0 +1,31 @@
+#ifndef IGL_SWEPT_VOLUME_BOUNDING_BOX_H
+#define IGL_SWEPT_VOLUME_BOUNDING_BOX_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+#include <Eigen/Geometry>
+#include <functional>
+namespace igl
+{
+  // Construct an axis-aligned bounding box containing a shape undergoing a
+  // motion sampled at `steps` discrete momements.
+  //
+  // Inputs:
+  //   n  number of mesh vertices
+  //   V  function handle so that V(i,t) returns the 3d position of vertex
+  //     i at time t, for t∈[0,1]
+  //   steps  number of time steps: steps=3 --> t∈{0,0.5,1}
+  // Outputs:
+  //   box  box containing mesh under motion
+  IGL_INLINE void swept_volume_bounding_box(
+    const size_t & n,
+    const std::function<
+      Eigen::RowVector3d(const size_t vi, const double t)> & V,
+    const size_t & steps,
+    Eigen::AlignedBox3d & box);
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "swept_volume_bounding_box.cpp"
+#endif
+
+#endif 

include/igl/swept_volume_signed_distance.cpp

@@ -0,0 +1,113 @@
+#include "swept_volume_signed_distance.h"
+#include "flood_fill.h"
+#include "signed_distance.h"
+#include "AABB.h"
+#include "pseudonormal_test.h"
+#include "per_face_normals.h"
+#include "per_vertex_normals.h"
+#include "per_edge_normals.h"
+#include <Eigen/Geometry>
+#include <cmath>
+
+IGL_INLINE void igl::swept_volume_signed_distance(
+  const Eigen::MatrixXd & V,
+  const Eigen::MatrixXi & F,
+  const std::function<Eigen::Affine3d(const double t)> & transform,
+  const size_t & steps,
+  const Eigen::MatrixXd & GV,
+  const Eigen::RowVector3i & res,
+  const double h,
+  const double isolevel,
+  const Eigen::VectorXd & S0,
+  Eigen::VectorXd & S)
+{
+  using namespace std;
+  using namespace igl;
+  using namespace Eigen;
+  S = S0;
+  const VectorXd t = VectorXd::LinSpaced(steps,0,1);
+  const bool finite_iso = isfinite(isolevel);
+  const double extension = (finite_iso ? isolevel : 0) + sqrt(3.0)*h;
+  Eigen::AlignedBox3d box(
+    V.colwise().minCoeff().array()-extension,
+    V.colwise().maxCoeff().array()+extension);
+  // Precomputation
+  Eigen::MatrixXd FN,VN,EN;
+  Eigen::MatrixXi E;
+  Eigen::VectorXi EMAP;
+  per_face_normals(V,F,FN);
+  per_vertex_normals(V,F,PER_VERTEX_NORMALS_WEIGHTING_TYPE_ANGLE,FN,VN);
+  per_edge_normals(
+    V,F,PER_EDGE_NORMALS_WEIGHTING_TYPE_UNIFORM,FN,EN,E,EMAP);
+  AABB<MatrixXd,3> tree;
+  tree.init(V,F);
+  for(int ti = 0;ti<t.size();ti++)
+  {
+    const Affine3d At = transform(t(ti));
+    for(int g = 0;g<GV.rows();g++)
+    {
+      // Don't bother finding out how deep inside points are.
+      if(finite_iso && S(g)==S(g) && S(g)<isolevel-sqrt(3.0)*h)
+      {
+        continue;
+      }
+      const RowVector3d gv = 
+        (GV.row(g) - At.translation().transpose())*At.linear();
+      // If outside of extended box, then consider it "far away enough"
+      if(finite_iso && !box.contains(gv.transpose()))
+      {
+        continue;
+      }
+      RowVector3d c,n;
+      int i;
+      double sqrd,s;
+      //signed_distance_pseudonormal(tree,V,F,FN,VN,EN,EMAP,gv,s,sqrd,i,c,n);
+      const double min_sqrd = 
+        finite_iso ? 
+        pow(sqrt(3.)*h+isolevel,2) : 
+        numeric_limits<double>::infinity();
+      sqrd = tree.squared_distance(V,F,gv,min_sqrd,i,c);
+      if(sqrd<min_sqrd)
+      {
+        pseudonormal_test(V,F,FN,VN,EN,EMAP,gv,i,c,s,n);
+        if(S(g) == S(g))
+        {
+          S(g) = min(S(g),s*sqrt(sqrd));
+        }else
+        {
+          S(g) = s*sqrt(sqrd);
+        }
+      }
+    }
+  }
+
+  if(finite_iso)
+  {
+    flood_fill(res,S);
+  }else
+  {
+#ifndef NDEBUG
+    // Check for nans
+    for_each(S.data(),S.data()+S.size(),[](const double s){assert(s==s);});
+#endif
+  }
+}
+
+IGL_INLINE void igl::swept_volume_signed_distance(
+  const Eigen::MatrixXd & V,
+  const Eigen::MatrixXi & F,
+  const std::function<Eigen::Affine3d(const double t)> & transform,
+  const size_t & steps,
+  const Eigen::MatrixXd & GV,
+  const Eigen::RowVector3i & res,
+  const double h,
+  const double isolevel,
+  Eigen::VectorXd & S)
+{
+  using namespace std;
+  using namespace igl;
+  using namespace Eigen;
+  S = VectorXd::Constant(GV.rows(),1,numeric_limits<double>::quiet_NaN());
+  return 
+    swept_volume_signed_distance(V,F,transform,steps,GV,res,h,isolevel,S,S);
+}

include/igl/swept_volume_signed_distance.h

@@ -0,0 +1,57 @@
+#ifndef IGL_SWEPT_VOLUME_SIGNED_DISTANCE_H
+#define IGL_SWEPT_VOLUME_SIGNED_DISTANCE_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+#include <Eigen/Geometry>
+#include <functional>
+namespace igl
+{
+  // Compute the signed distance to a sweep surface of a mesh under-going
+  // an arbitrary motion V(t) discretely sampled at `steps`-many moments in
+  // time at a grid.
+  //
+  // Inputs:
+  //   V  #V by 3 list of mesh positions in reference pose
+  //   F  #F by 3 list of triangle indices [0,n)
+  //   transform  function handle so that transform(t) returns the rigid
+  //     transformation at time t∈[0,1]
+  //   steps  number of time steps: steps=3 --> t∈{0,0.5,1}
+  //   GV  #GV by 3 list of evaluation point grid positions
+  //   res  3-long resolution of GV grid
+  //   h  edge-length of grid
+  //   isolevel  isolevel to "focus" on; grid positions far enough away from
+  //     isolevel (based on h) will get approximate values). Set
+  //     isolevel=infinity to get good values everywhere (slow and
+  //     unnecessary if just trying to extract isolevel-level set).
+  //   S0  #GV initial values (will take minimum with these), can be same
+  //     as S)
+  // Outputs:
+  //   S  #GV list of signed distances
+  IGL_INLINE void swept_volume_signed_distance(
+    const Eigen::MatrixXd & V,
+    const Eigen::MatrixXi & F,
+    const std::function<Eigen::Affine3d(const double t)> & transform,
+    const size_t & steps,
+    const Eigen::MatrixXd & GV,
+    const Eigen::RowVector3i & res,
+    const double h,
+    const double isolevel,
+    const Eigen::VectorXd & S0,
+    Eigen::VectorXd & S);
+  IGL_INLINE void swept_volume_signed_distance(
+    const Eigen::MatrixXd & V,
+    const Eigen::MatrixXi & F,
+    const std::function<Eigen::Affine3d(const double t)> & transform,
+    const size_t & steps,
+    const Eigen::MatrixXd & GV,
+    const Eigen::RowVector3i & res,
+    const double h,
+    const double isolevel,
+    Eigen::VectorXd & S);
+  }
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "swept_volume_signed_distance.cpp"
+#endif 
+
+#endif

include/igl/triangle/triangulate.h

@@ -7,7 +7,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #ifndef IGL_TRIANGLE_TRIANGULATE_H
 #define IGL_TRIANGLE_TRIANGULATE_H
-#include <igl/igl_inline.h>
+#include "../igl_inline.h"
 #include <string>
 #include <Eigen/Core>
 

include/igl/triangle_triangle_adjacency.cpp

@@ -9,6 +9,7 @@
 #include "is_edge_manifold.h"
 #include "all_edges.h"
 #include "unique_simplices.h"
+#include "parallel_for.h"
 #include "unique_edge_map.h"
 #include <algorithm>
 #include <iostream>
@@ -172,41 +173,46 @@ template <
   }
 
   // No race conditions because TT*[f][c]'s are in bijection with e's
-  // Minimum number of iterms per openmp thread
+  // Minimum number of items per thread
   //const size_t num_e = E.rows();
-# ifndef IGL_OMP_MIN_VALUE
-#   define IGL_OMP_MIN_VALUE 1000
-# endif
-# pragma omp parallel for if (m>IGL_OMP_MIN_VALUE)
   // Slightly better memory access than loop over E
-  for(Index f = 0;f<(Index)m;f++)
-  {
-    for(Index c = 0;c<3;c++)
+  igl::parallel_for(
+    m,
+    [&](const Index & f)
     {
-      const Index e = f + m*c;
-      //const Index c = e/m;
-      const vector<uE2EType> & N = uE2E[EMAP(e)];
-      for(const auto & ne : N)
+      for(Index c = 0;c<3;c++)
       {
-        const Index nf = ne%m;
-        // don't add self
-        if(nf != f)
+        const Index e = f + m*c;
+        //const Index c = e/m;
+        const vector<uE2EType> & N = uE2E[EMAP(e)];
+        for(const auto & ne : N)
         {
-          TT[f][c].push_back(nf);
-          if(construct_TTi)
+          const Index nf = ne%m;
+          // don't add self
+          if(nf != f)
           {
-            const Index nc = ne/m;
-            TTi[f][c].push_back(nc);
+            TT[f][c].push_back(nf);
+            if(construct_TTi)
+            {
+              const Index nc = ne/m;
+              TTi[f][c].push_back(nc);
+            }
           }
         }
       }
-    }
-  }
+    },
+    1000ul);
+
+
 }
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 // generated by autoexplicit.sh
+template void igl::triangle_triangle_adjacency<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, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+// generated by autoexplicit.sh
+template void igl::triangle_triangle_adjacency<Eigen::Matrix<int, -1, 3, 0, -1, 3>, int>(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, std::vector<std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >, std::allocator<std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > > >&);
+// generated by autoexplicit.sh
 template void igl::triangle_triangle_adjacency<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> >&);
 // generated by autoexplicit.sh
 template void igl::triangle_triangle_adjacency<Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);

include/igl/uniformly_sample_two_manifold.h

@@ -11,8 +11,10 @@
 #include <Eigen/Dense>
 namespace igl
 {
-  // UNIFORMLY_SAMPLE_TWO_MANIFOLD Attempt to sample a mesh uniformly by furthest
-  // point relaxation as described in "Fast Automatic Skinning Transformations"
+  // UNIFORMLY_SAMPLE_TWO_MANIFOLD Attempt to sample a mesh uniformly by
+  // furthest point relaxation as described in "Fast Automatic Skinning
+  // Transformations"
+  //
   // [Jacobson et al. 12] Section 3.3.
   //
   // Inputs:
@@ -37,6 +39,6 @@ namespace igl
     Eigen::VectorXi & S);
 }
 #ifndef IGL_STATIC_LIBRARY
-#  include "uniformly_sample_two_manifold.h"
+#  include "uniformly_sample_two_manifold.cpp"
 #endif
 #endif

include/igl/unique.cpp

@@ -268,4 +268,9 @@ template void igl::unique<double>(std::vector<double, std::allocator<double> > c
 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> >&);
+
+#ifdef WIN32
+template void __cdecl igl::unique_rows<class Eigen::Matrix<int, -1, -1, 0, -1, -1>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> >(class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1> > const &, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> > &);
+#endif
+
 #endif

include/igl/unique_edge_map.cpp

@@ -52,4 +52,9 @@ template void igl::unique_edge_map<Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::M
 template void igl::unique_edge_map<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 2, 0, -1, 2>, Eigen::Matrix<double, -1, 2, 0, -1, 2>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, long>(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&, std::vector<std::vector<long, std::allocator<long> >, std::allocator<std::vector<long, std::allocator<long> > > >&);
 template void igl::unique_edge_map<Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 2, 0, -1, 2>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, int>(Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 2, 0, -1, 2> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&);
 template void igl::unique_edge_map<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>, unsigned long>(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> >&, std::vector<std::vector<unsigned long, std::allocator<unsigned long> >, std::allocator<std::vector<unsigned long, std::allocator<unsigned long> > > >&);
+
+#ifdef WIN32
+template void __cdecl igl::unique_edge_map<class Eigen::Matrix<int, -1, 3, 0, -1, 3>, class Eigen::Matrix<int, -1, 2, 0, -1, 2>, class Eigen::Matrix<int, -1, 2, 0, -1, 2>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, __int64>(class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, 3, 0, -1, 3> > const &, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, 2, 0, -1, 2> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, 2, 0, -1, 2> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> > &, class std::vector<class std::vector<__int64, class std::allocator<__int64> >, class std::allocator<class std::vector<__int64, class std::allocator<__int64> > > > &);
+#endif
+
 #endif 

include/igl/unique_simplices.cpp

@@ -1,6 +1,6 @@
 // This file is part of libigl, a simple c++ geometry processing library.
 // 
-// Copyright (C) 2013 Alec Jacobson <alecjacobson@gmail.com>
+// Copyright (C) 2016 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 
@@ -8,7 +8,7 @@
 #include "unique_simplices.h"
 #include "sort.h"
 #include "unique.h"
-#include "get_seconds.h"
+#include "parallel_for.h"
 
 template <
   typename DerivedF,
@@ -31,16 +31,7 @@ IGL_INLINE void igl::unique_simplices(
   igl::unique_rows(sortF,C,IA,IC);
   FF.resize(IA.size(),F.cols());
   const size_t mff = FF.rows();
-  // Minimum number of iterms per openmp thread
-  #ifndef IGL_OMP_MIN_VALUE
-  #  define IGL_OMP_MIN_VALUE 1000
-  #endif
-  #pragma omp parallel for if (mff>IGL_OMP_MIN_VALUE)
-  // Copy into output
-  for(size_t i = 0;i<mff;i++)
-  {
-    FF.row(i) = F.row(IA(i));
-  }
+  parallel_for(mff,[&F,&IA,&FF](size_t & i){FF.row(i) = F.row(IA(i));},1000ul);
 }
 
 template <
@@ -64,4 +55,7 @@ template void igl::unique_simplices<Eigen::Matrix<int, -1, 3, 0, -1, 3>, Eigen::
 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> >&);
+#ifdef WIN32
+template void __cdecl igl::unique_simplices<class Eigen::Matrix<int, -1, 2, 0, -1, 2>, class Eigen::Matrix<int, -1, 2, 0, -1, 2>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> >(class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, 2, 0, -1, 2> > const &, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, 2, 0, -1, 2> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> > &, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> > &);
+#endif
 #endif

include/igl/unproject_onto_mesh.cpp

@@ -73,5 +73,6 @@ IGL_INLINE bool igl::unproject_onto_mesh(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 template bool igl::unproject_onto_mesh<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<float, 3, 1, 0, 3, 1> >(Eigen::Matrix<float, 2, 1, 0, 2, 1> const&, Eigen::Matrix<float, 4, 4, 0, 4, 4> const&, Eigen::Matrix<float, 4, 4, 0, 4, 4> const&, Eigen::Matrix<float, 4, 1, 0, 4, 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&, Eigen::PlainObjectBase<Eigen::Matrix<float, 3, 1, 0, 3, 1> >&);
+template bool igl::unproject_onto_mesh<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<float, 2, 1, 0, 2, 1> const&, Eigen::Matrix<float, 4, 4, 0, 4, 4> const&, Eigen::Matrix<float, 4, 4, 0, 4, 4> const&, Eigen::Matrix<float, 4, 1, 0, 4, 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&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
 #endif
 

include/igl/upsample.cpp

@@ -105,5 +105,7 @@ IGL_INLINE void igl::upsample(
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template void igl::upsample<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<int, -1, 3, 0, -1, 3> >(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> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
 template void igl::upsample<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>&);
 #endif

include/igl/viewer/Viewer.cpp

@@ -353,7 +353,7 @@ namespace viewer
   :       Toggle face labels)"
 #endif
 );
-    std::cout<<usage;
+    std::cout<<usage<<std::endl;
 #endif
   }
 

include/igl/viewer/ViewerData.cpp

@@ -11,6 +11,7 @@
 #include <iostream>
 
 #include <igl/per_face_normals.h>
+#include <igl/material_colors.h>
 #include <igl/per_vertex_normals.h>
 
 #ifdef ENABLE_SERIALIZATION
@@ -104,9 +105,10 @@ IGL_INLINE void igl::viewer::ViewerData::set_mesh(const Eigen::MatrixXd& _V, con
     F = _F;
 
     compute_normals();
-    uniform_colors(Eigen::Vector3d(51.0/255.0,43.0/255.0,33.3/255.0),
-                   Eigen::Vector3d(255.0/255.0,228.0/255.0,58.0/255.0),
-                   Eigen::Vector3d(255.0/255.0,235.0/255.0,80.0/255.0));
+    uniform_colors(
+      Eigen::Vector3d(GOLD_AMBIENT[0], GOLD_AMBIENT[1], GOLD_AMBIENT[2]),
+      Eigen::Vector3d(GOLD_DIFFUSE[0], GOLD_DIFFUSE[1], GOLD_DIFFUSE[2]),
+      Eigen::Vector3d(GOLD_SPECULAR[0], GOLD_SPECULAR[1], GOLD_SPECULAR[2]));
 
     grid_texture();
   }

include/igl/viewer/ViewerData.h

@@ -108,7 +108,10 @@ public:
   IGL_INLINE void compute_normals();
 
   // Assigns uniform colors to all faces/vertices
-  IGL_INLINE void uniform_colors(Eigen::Vector3d ambient, Eigen::Vector3d diffuse, Eigen::Vector3d specular);
+  IGL_INLINE void uniform_colors(
+    Eigen::Vector3d ambient, 
+    Eigen::Vector3d diffuse, 
+    Eigen::Vector3d specular);
 
   // Generates a default grid texture
   IGL_INLINE void grid_texture();

include/igl/voxel_grid.cpp

@@ -0,0 +1,51 @@
+#include "voxel_grid.h"
+#include "grid.h"
+
+IGL_INLINE void igl::voxel_grid(
+  const Eigen::AlignedBox3d & box, 
+  const int in_s,
+  const int pad_count,
+  Eigen::MatrixXd & GV,
+  Eigen::RowVector3i & side)
+{
+  using namespace Eigen;
+  using namespace std;
+  MatrixXd::Index si = -1;
+  box.diagonal().maxCoeff(&si);
+  //MatrixXd::Index si = 0;
+  //assert(si>=0);
+  const double s_len = box.diagonal()(si);
+  assert(in_s>(pad_count*2+1) && "s should be > 2*pad_count+1");
+  const double s = in_s - 2*pad_count;
+  side(si) = s;
+  for(int i = 0;i<3;i++)
+  {
+    if(i!=si)
+    {
+      side(i) = ceil(s * (box.max()(i)-box.min()(i))/s_len);
+    }
+  }
+  side.array() += 2*pad_count;
+  grid(side,GV);
+  // A *    p/s  + B = min
+  // A * (1-p/s) + B = max
+  // B = min - A * p/s
+  // A * (1-p/s) + min - A * p/s = max
+  // A * (1-p/s) - A * p/s = max-min
+  // A * (1-2p/s) = max-min
+  // A  = (max-min)/(1-2p/s)
+  const Array<double,3,1> ps= 
+    (double)(pad_count)/(side.transpose().cast<double>().array()-1.);
+  const Array<double,3,1> A = box.diagonal().array()/(1.0-2.*ps);
+  //// This would result in an "anamorphic", but perfectly fit grid:
+  //const Array<double,3,1> B = box.min().array() - A.array()*ps;
+  //GV.array().rowwise() *= A.transpose();
+  //GV.array().rowwise() += B.transpose();
+  // Instead scale by largest factor and move to match center
+  Array<double,3,1>::Index ai = -1;
+  double a = A.maxCoeff(&ai);
+  const Array<double,1,3> ratio = 
+    a*(side.cast<double>().array()-1.0)/(double)(side(ai)-1.0);
+  GV.array().rowwise() *= ratio;
+  GV.rowwise() += (box.center().transpose()-GV.colwise().mean()).eval();
+}

include/igl/voxel_grid.h

@@ -0,0 +1,28 @@
+#ifndef IGL_VOXEL_GRID_H
+#define IGL_VOXEL_GRID_H
+#include "igl_inline.h"
+#include <Eigen/Core>
+#include <Eigen/Geometry>
+namespace igl
+{
+  // Construct the cell center positions of a regular voxel grid (lattice) made
+  // of perfectly square voxels.
+  // 
+  // Inputs:
+  //   box  bounding box to enclose by grid
+  //   s  number of cell centers on largest side (including 2*pad_count)
+  //   pad_count  number of cells beyond box
+  // Outputs:
+  //   GV  res(0)*res(1)*res(2) by 3 list of cell center positions
+  //   res  3-long list of dimension of voxel grid
+  IGL_INLINE void voxel_grid(
+    const Eigen::AlignedBox3d & box, 
+    const int s,
+    const int pad_count,
+    Eigen::MatrixXd & GV,
+    Eigen::RowVector3i & side);
+}
+#ifndef IGL_STATIC_LIBRARY
+#  include "voxel_grid.cpp"
+#endif
+#endif

include/igl/writeMESH.cpp

@@ -137,6 +137,8 @@ IGL_INLINE bool igl::writeMESH(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 // generated by autoexplicit.sh
+template bool igl::writeMESH<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&);
+// generated by autoexplicit.sh
 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&);

include/igl/writeSTL.cpp

@@ -110,6 +110,8 @@ IGL_INLINE bool igl::writeSTL(
 
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
+// generated by autoexplicit.sh
+template bool igl::writeSTL<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, bool);
 template bool igl::writeSTL<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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, bool);
 template bool igl::writeSTL<Eigen::Matrix<double, 8, 3, 0, 8, 3>, Eigen::Matrix<int, 12, 3, 0, 12, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 8, 3, 0, 8, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, 12, 3, 0, 12, 3> > const&, bool);
 #endif

include/igl/writeWRL.cpp

@@ -56,6 +56,8 @@ ccw TRUE
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 // generated by autoexplicit.sh
+template bool igl::writeWRL<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&);
+// generated by autoexplicit.sh
 template bool igl::writeWRL<Eigen::Matrix<double, 8, 3, 0, 8, 3>, Eigen::Matrix<int, 12, 3, 0, 12, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 8, 3, 0, 8, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, 12, 3, 0, 12, 3> > const&);
 template bool igl::writeWRL<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> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&);
 #endif

include/igl/write_triangle_mesh.cpp

@@ -63,6 +63,8 @@ IGL_INLINE bool igl::write_triangle_mesh(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template specialization
 // generated by autoexplicit.sh
+template bool igl::write_triangle_mesh<Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(std::basic_string<char, std::char_traits<char>, std::allocator<char> >, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, bool);
+// generated by autoexplicit.sh
 template bool igl::write_triangle_mesh<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<int, -1, -1, 0, -1, -1> > const&, const bool);
 template bool igl::write_triangle_mesh<Eigen::Matrix<double, 8, 3, 0, 8, 3>, 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, 12, 3, 0, 12, 3> > const&, bool);
 #endif

index.html

@@ -214,6 +214,7 @@ few labs/companies/institutions using libigl:</p>
 
 <ul>
 <li><a href="http://www.adobe.com/technology/">Adobe Research</a></li>
+<li><a href="http://www.ea.com">Electronic Arts, Inc</a></li>
 <li><a href="http://meshconsultants.ca/">Mesh</a>, consultants, Canada</li>
 <li><a href="http://graphics.pixar.com/research/">Pixar Research</a></li>
 <li><a href="http://esotericsoftware.com/">Spine by Esoteric Software</a> is an animation tool dedicated to 2D characters.</li>
@@ -263,7 +264,7 @@ page</a>.</p>
 <h2 id="copyright">Copyright</h2>
 
 <p>2016 Alec Jacobson, Daniele Panozzo, Christian Schüller, Olga Diamanti, Qingnan
-Zhou, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
+Zhou, Sebastian Koch, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
 Giorgis, Luigi Rocca, Leonardo Sacht, Kevin Walliman, Olga Sorkine-Hornung, and others.</p>
 
 <p>Please see individual files for appropriate copyright notices.</p>

python/CMakeLists.txt

@@ -47,7 +47,7 @@ endif()
 
 include_directories(${PYTHON_INCLUDE_DIR} include)
 
-## include pybing
+## include pybind
 include_directories(${PROJECT_SOURCE_DIR}/../external/pybind11/include)
 
 ## include libigl
@@ -57,12 +57,12 @@ option(LIBIGL_WITH_NANOGUI          "Use Nanogui menu"   OFF)
 option(LIBIGL_WITH_CGAL             "Use CGAL"           OFF)
 option(LIBIGL_WITH_BOOLEAN          "Use Cork boolean"   OFF)
 option(LIBIGL_WITH_COMISO           "Use CoMiso"         ON)
-option(LIBIGL_WITH_EMBREE           "Use Embree"         OFF)
+option(LIBIGL_WITH_EMBREE           "Use Embree"         ON)
 option(LIBIGL_WITH_LIM              "Use LIM"            ON)
 option(LIBIGL_WITH_MATLAB           "Use Matlab"         OFF)
 option(LIBIGL_WITH_MOSEK            "Use MOSEK"          OFF)
 option(LIBIGL_WITH_BBW              "Use BBW"            ON)
-option(LIBIGL_WITH_OPENGL_AND_PNG   "Use OpenGL"         ON)
+option(LIBIGL_WITH_PNG              "Use PNG"            ON)
 option(LIBIGL_WITH_TETGEN           "Use Tetgen"         ON)
 option(LIBIGL_WITH_TRIANGLE         "Use Triangle"       ON)
 option(LIBIGL_WITH_XML              "Use XML"            ON)
@@ -84,19 +84,44 @@ add_definitions(${LIBIGL_DEFINITIONS})
 ## Optional modules
 if (LIBIGL_WITH_VIEWER)
   add_definitions(-DPY_VIEWER)
-  list(APPEND SHARED_SOURCES "py_igl_viewer.cpp")
+  list(APPEND SHARED_SOURCES "modules/py_igl_viewer.cpp")
 endif ()
 
 if (LIBIGL_WITH_COMISO)
   add_definitions(-DPY_COMISO)
-  list(APPEND SHARED_SOURCES "copyleft/py_igl_comiso.cpp")
+  list(APPEND SHARED_SOURCES "modules/copyleft/py_igl_comiso.cpp")
+endif ()
+
+if (LIBIGL_WITH_TETGEN)
+  add_definitions(-DPY_TETGEN)
+  list(APPEND SHARED_SOURCES "modules/copyleft/py_igl_tetgen.cpp")
+endif ()
+
+if (LIBIGL_WITH_EMBREE)
+  add_definitions(-DPY_EMBREE)
+  list(APPEND SHARED_SOURCES "modules/py_igl_embree.cpp")
+endif ()
+
+if (LIBIGL_WITH_TRIANGLE)
+  add_definitions(-DPY_TRIANGLE)
+  list(APPEND SHARED_SOURCES "modules/py_igl_triangle.cpp")
+endif ()
+
+if (LIBIGL_WITH_CGAL)
+  add_definitions(-DPY_CGAL)
+  list(APPEND SHARED_SOURCES "modules/copyleft/py_igl_cgal.cpp")
+endif ()
+
+if (LIBIGL_WITH_PNG)
+  add_definitions(-DPY_PNG)
+  list(APPEND SHARED_SOURCES "modules/py_igl_png.cpp")
 endif ()
 
 
 ## Prepare the python library
 add_library(pyigl SHARED
   python_shared.cpp
-  py_vector.cpp
+  modules/py_vector.cpp
   py_igl.cpp
   py_doc.cpp
   ${SHARED_SOURCES}
@@ -142,10 +167,10 @@ elseif (UNIX)
   #Enable flag if undefined symbols appear on pyigl module import to get notified about the missing symbols at link time
   option(CHECK_UNDEFINED        "Check for undefined symbols"    OFF)
 
-  # Strip unnecessary sections of the binary on Linux/Mac OS 
+  # Strip unnecessary sections of the binary on Linux/Mac OS
   if(APPLE)
     set_target_properties(pyigl PROPERTIES MACOSX_RPATH ".")
-    
+
     if (NOT CHECK_UNDEFINED)
       set_target_properties(pyigl PROPERTIES LINK_FLAGS "-undefined dynamic_lookup -dead_strip")
     endif()

python/README.md

@@ -5,7 +5,7 @@
 <span style="color:#F62217">
 Everything in this folder is currently being developed and it is likely to be
 changed radically in the next couple of months, breaking compatibility between
-different version. We plan to stabilize the python API by the end of 2015.
+different version. We plan to stabilize the python API by the end of 2016.
 </span>
 
 ## Introduction
@@ -153,7 +153,7 @@ page](https://github.com/libigl/libigl/issues).
 
 ## Copyright
 2015 Alec Jacobson, Daniele Panozzo, Christian Schüller, Olga Diamanti, Qingnan
-Zhou, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
+Zhou, Sebastian Koch, Nico Pietroni, Stefan Brugger, Kenshi Takayama, Wenzel Jakob, Nikolas De
 Giorgis, Luigi Rocca, Leonardo Sacht, Olga Sorkine-Hornung, and others.
 
 Please see individual files for appropriate copyright notices.

python/modules/copyleft/py_igl_cgal.cpp

@@ -0,0 +1,18 @@
+//#include <Eigen/Geometry>
+//#include <Eigen/Dense>
+//#include <Eigen/Sparse>
+
+
+#include "../../python_shared.h"
+
+#include <igl/copyleft/cgal/mesh_boolean.h>
+
+
+void python_export_igl_cgal(py::module &me) {
+
+  py::module m = me.def_submodule(
+    "cgal", "Wrappers for libigl functions that use cgal");
+
+  #include "../../py_igl/copyleft/cgal/py_mesh_boolean.cpp"
+
+}

python/copyleft/py_igl_comiso.cpp → python/modules/copyleft/py_igl_comiso.cpp

@@ -3,7 +3,7 @@
 #include <Eigen/Sparse>
 
 
-#include "../python_shared.h"
+#include "../../python_shared.h"
 
 #include <igl/copyleft/comiso/nrosy.h>
 #include <igl/copyleft/comiso/miq.h>
@@ -13,7 +13,7 @@ void python_export_igl_comiso(py::module &me) {
   py::module m = me.def_submodule(
     "comiso", "Wrappers for libigl functions that use comiso");
 
-  #include "../py_igl/copyleft/comiso/py_nrosy.cpp"
-  #include "../py_igl/copyleft/comiso/py_miq.cpp"
+  #include "../../py_igl/copyleft/comiso/py_nrosy.cpp"
+  #include "../../py_igl/copyleft/comiso/py_miq.cpp"
 
 }

python/modules/copyleft/py_igl_tetgen.cpp

@@ -0,0 +1,18 @@
+//#include <Eigen/Geometry>
+//#include <Eigen/Dense>
+//#include <Eigen/Sparse>
+
+
+#include "../../python_shared.h"
+
+#include <igl/copyleft/tetgen/tetrahedralize.h>
+
+
+void python_export_igl_tetgen(py::module &me) {
+
+  py::module m = me.def_submodule(
+    "tetgen", "Wrappers for libigl functions that use tetgen");
+
+  #include "../../py_igl/copyleft/tetgen/py_tetrahedralize.cpp"
+
+}

python/modules/py_igl_embree.cpp

@@ -0,0 +1,18 @@
+//#include <Eigen/Geometry>
+//#include <Eigen/Dense>
+//#include <Eigen/Sparse>
+
+
+#include "../python_shared.h"
+
+#include <igl/embree/ambient_occlusion.h>
+
+
+void python_export_igl_embree(py::module &me) {
+
+  py::module m = me.def_submodule(
+    "embree", "Wrappers for libigl functions that use embree");
+
+  #include "../py_igl/embree/py_ambient_occlusion.cpp"
+
+}

python/modules/py_igl_png.cpp

@@ -0,0 +1,16 @@
+
+#include "../python_shared.h"
+
+#include <igl/png/readPNG.h>
+#include <igl/png/writePNG.h>
+
+
+void python_export_igl_png(py::module &me) {
+
+  py::module m = me.def_submodule(
+    "png", "Wrappers for libigl functions that use png");
+
+  #include "../py_igl/png/py_readPNG.cpp"
+  #include "../py_igl/png/py_writePNG.cpp"
+
+}

python/modules/py_igl_triangle.cpp

@@ -0,0 +1,18 @@
+//#include <Eigen/Geometry>
+//#include <Eigen/Dense>
+//#include <Eigen/Sparse>
+
+
+#include "../python_shared.h"
+
+#include <igl/triangle/triangulate.h>
+
+
+void python_export_igl_triangle(py::module &me) {
+
+  py::module m = me.def_submodule(
+    "triangle", "Wrappers for libigl functions that use triangle");
+
+  #include "../py_igl/triangle/py_triangulate.cpp"
+
+}

python/py_igl_viewer.cpp → python/modules/py_igl_viewer.cpp

@@ -1,11 +1,12 @@
 #include <Eigen/Dense>
 #include <Eigen/Sparse>
 
-#include "python_shared.h"
+#include "../python_shared.h"
 #define ENABLE_SERIALIZATION
 #include <igl/viewer/Viewer.h>
 #include <igl/viewer/ViewerCore.h>
 #include <igl/viewer/ViewerData.h>
+#include <igl/viewer/OpenGL_state.h>
 #include <igl/serialize.h>
 
 void python_export_igl_viewer(py::module &m)
@@ -118,6 +119,16 @@ py::enum_<igl::viewer::ViewerData::DirtyFlags>(viewerdata_class, "DirtyFlags")
 
     ;
 
+//////////////////////// OPENGL_State
+
+py::class_<igl::viewer::OpenGL_state> opengl_state_class(me, "OpenGL_state");
+
+    opengl_state_class
+    .def(py::init<>())
+    .def("init", &igl::viewer::OpenGL_state::init)
+
+    ;
+
 //////////////////////// CORE
 
 py::class_<igl::viewer::ViewerCore> viewercore_class(me, "ViewerCore");
@@ -311,6 +322,7 @@ py::class_<igl::viewer::ViewerCore> viewercore_class(me, "ViewerCore");
     .def(py::init<>())
     .def_readwrite("data", &igl::viewer::Viewer::data)
     .def_readwrite("core", &igl::viewer::Viewer::core)
+    .def_readwrite("opengl", &igl::viewer::Viewer::opengl)
     .def("launch", &igl::viewer::Viewer::launch, py::arg("resizable") = true, py::arg("fullscreen") = false)
     .def("launch_init", &igl::viewer::Viewer::launch_init, py::arg("resizable") = true, py::arg("fullscreen") = false)
     .def("launch_rendering", &igl::viewer::Viewer::launch_rendering, py::arg("loop") = true)
@@ -356,6 +368,10 @@ py::class_<igl::viewer::ViewerCore> viewercore_class(me, "ViewerCore");
     .def("open_dialog_load_mesh", &igl::viewer::Viewer::open_dialog_load_mesh)
     .def("open_dialog_save_mesh", &igl::viewer::Viewer::open_dialog_save_mesh)
 
+    // Input handling
+    .def_readwrite("current_mouse_x", &igl::viewer::Viewer::current_mouse_x)
+    .def_readwrite("current_mouse_y", &igl::viewer::Viewer::current_mouse_y)
+
     // Callbacks
     .def_readwrite("callback_init", &igl::viewer::Viewer::callback_init)
     .def_readwrite("callback_pre_draw", &igl::viewer::Viewer::callback_pre_draw)

python/py_vector.cpp → python/modules/py_vector.cpp

@@ -3,7 +3,7 @@
 #include <Eigen/Sparse>
 
 
-#include "python_shared.h"
+#include "../python_shared.h"
 
 /// Creates Python bindings for a dynamic Eigen matrix
 template <typename Type>
@@ -126,6 +126,9 @@ py::class_<Type> bind_eigen_2(py::module &m, const char *name,
         .def("rightCols", [](Type &m, const int& k) { return Type(m.rightCols(k)); })
         .def("leftCols", [](Type &m, const int& k) { return Type(m.leftCols(k)); })
 
+        .def("setLeftCols", [](Type &m, const int& k, const Type& v) { return Type(m.leftCols(k) = v); })
+        .def("setRightCols", [](Type &m, const int& k, const Type& v) { return Type(m.rightCols(k) = v); })
+
         .def("topRows", [](Type &m, const int& k) { return Type(m.topRows(k)); })
         .def("bottomRows", [](Type &m, const int& k) { return Type(m.bottomRows(k)); })
 
@@ -170,6 +173,7 @@ py::class_<Type> bind_eigen_2(py::module &m, const char *name,
         .def("cwiseQuotient", [](const Type &m1, const Type &m2) -> Type { return m1.cwiseQuotient(m2); })
 
         /* Row and column-wise operations */
+        .def("rowwiseSet", [](Type &m, const Type &m2) {return Type(m.rowwise() = Eigen::Matrix<Scalar, 1, Eigen::Dynamic>(m2));} )
         .def("rowwiseSum", [](const Type &m) {return Type(m.rowwise().sum());} )
         .def("rowwiseProd", [](const Type &m) {return Type(m.rowwise().prod());} )
         .def("rowwiseMean", [](const Type &m) {return Type(m.rowwise().mean());} )
@@ -178,6 +182,7 @@ py::class_<Type> bind_eigen_2(py::module &m, const char *name,
         .def("rowwiseMinCoeff", [](const Type &m) {return Type(m.rowwise().minCoeff());} )
         .def("rowwiseMaxCoeff", [](const Type &m) {return Type(m.rowwise().maxCoeff());} )
 
+        .def("colwiseSet", [](Type &m, const Type &m2) {return Type(m.colwise() = Eigen::Matrix<Scalar, Eigen::Dynamic, 1>(m2));} )
         .def("colwiseSum", [](const Type &m) {return Type(m.colwise().sum());} )
         .def("colwiseProd", [](const Type &m) {return Type(m.colwise().prod());} )
         .def("colwiseMean", [](const Type &m) {return Type(m.colwise().mean());} )

python/py_doc.h

@@ -1,71 +1,80 @@
-extern const char *__doc_igl_principal_curvature;
-extern const char *__doc_igl_local_basis;
-extern const char *__doc_igl_signed_distance;
-extern const char *__doc_igl_signed_distance_pseudonormal;
-extern const char *__doc_igl_signed_distance_winding_number;
-extern const char *__doc_igl_cotmatrix;
-extern const char *__doc_igl_floor;
-extern const char *__doc_igl_slice;
-extern const char *__doc_igl_per_face_normals;
-extern const char *__doc_igl_per_face_normals_stable;
-extern const char *__doc_igl_readOFF;
-extern const char *__doc_igl_per_vertex_normals;
-extern const char *__doc_igl_sortrows;
+extern const char *__doc_igl_active_set;
+extern const char *__doc_igl_arap_precomputation;
+extern const char *__doc_igl_arap_solve;
+extern const char *__doc_igl_avg_edge_length;
 extern const char *__doc_igl_barycenter;
-extern const char *__doc_igl_jet;
+extern const char *__doc_igl_barycentric_coordinates;
+extern const char *__doc_igl_boundary_facets;
+extern const char *__doc_igl_boundary_loop;
 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_slice_mask;
-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_cross_field;
 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_copyleft_cgal_mesh_boolean;
+extern const char *__doc_igl_copyleft_comiso_miq;
+extern const char *__doc_igl_copyleft_comiso_nrosy;
+extern const char *__doc_igl_copyleft_tetgen_tetrahedralize;
+extern const char *__doc_igl_cotmatrix;
+extern const char *__doc_igl_covariance_scatter_matrix;
+extern const char *__doc_igl_cross_field_missmatch;
 extern const char *__doc_igl_cut_mesh_from_singularities;
-extern const char *__doc_igl_readDMAT;
 extern const char *__doc_igl_doublearea;
 extern const char *__doc_igl_doublearea_single;
 extern const char *__doc_igl_doublearea_quad;
+extern const char *__doc_igl_edge_lengths;
+extern const char *__doc_igl_eigs;
+extern const char *__doc_igl_embree_ambient_occlusion;
+extern const char *__doc_igl_find_cross_field_singularities;
+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_floor;
+extern const char *__doc_igl_gaussian_curvature;
+extern const char *__doc_igl_grad;
+extern const char *__doc_igl_harmonic;
+extern const char *__doc_igl_invert_diag;
+extern const char *__doc_igl_jet;
+extern const char *__doc_igl_local_basis;
+extern const char *__doc_igl_lscm;
+extern const char *__doc_igl_map_vertices_to_circle;
+extern const char *__doc_igl_massmatrix;
 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;
-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_parula;
+extern const char *__doc_igl_per_corner_normals;
+extern const char *__doc_igl_per_edge_normals;
+extern const char *__doc_igl_per_face_normals;
+extern const char *__doc_igl_per_face_normals_stable;
+extern const char *__doc_igl_per_vertex_normals;
+extern const char *__doc_igl_planarize_quad_mesh;
+extern const char *__doc_igl_png_readPNG;
+extern const char *__doc_igl_png_writePNG;
+extern const char *__doc_igl_point_mesh_squared_distance;
 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_principal_curvature;
+extern const char *__doc_igl_quad_planarity;
+extern const char *__doc_igl_readDMAT;
 extern const char *__doc_igl_readMESH;
-extern const char *__doc_igl_copyleft_comiso_miq;
-extern const char *__doc_igl_copyleft_comiso_nrosy;
+extern const char *__doc_igl_readOBJ;
+extern const char *__doc_igl_readOFF;
+extern const char *__doc_igl_read_triangle_mesh;
+extern const char *__doc_igl_rotate_vectors;
+extern const char *__doc_igl_setdiff;
+extern const char *__doc_igl_signed_distance;
+extern const char *__doc_igl_signed_distance_pseudonormal;
+extern const char *__doc_igl_signed_distance_winding_number;
+extern const char *__doc_igl_slice;
+extern const char *__doc_igl_slice_into;
+extern const char *__doc_igl_slice_mask;
+extern const char *__doc_igl_slice_tets;
+extern const char *__doc_igl_sortrows;
+extern const char *__doc_igl_triangle_triangulate;
+extern const char *__doc_igl_unique;
+extern const char *__doc_igl_unique_rows;
+extern const char *__doc_igl_unproject_onto_mesh;
+extern const char *__doc_igl_upsample;
+extern const char *__doc_igl_writeMESH;
+extern const char *__doc_igl_writeOBJ;

python/py_igl.cpp

@@ -2,136 +2,139 @@
 
 #include "python_shared.h"
 
-#include <igl/readOFF.h>
-#include <igl/writeOBJ.h>
-#include <igl/per_face_normals.h>
-#include <igl/per_corner_normals.h>
-#include <igl/per_vertex_normals.h>
-#include <igl/gaussian_curvature.h>
-#include <igl/jet.h>
-#include <igl/read_triangle_mesh.h>
-#include <igl/cotmatrix.h>
-#include <igl/massmatrix.h>
-#include <igl/invert_diag.h>
-#include <igl/principal_curvature.h>
-#include <igl/parula.h>
-#include <igl/readDMAT.h>
-#include <igl/grad.h>
-#include <igl/avg_edge_length.h>
-#include <igl/barycenter.h>
-#include <igl/doublearea.h>
-#include <igl/floor.h>
-#include <igl/slice.h>
-#include <igl/slice_into.h>
-#include <igl/sortrows.h>
-#include <igl/colon.h>
-#include <igl/boundary_facets.h>
-#include <igl/unique.h>
-#include <igl/setdiff.h>
-#include <igl/min_quad_with_fixed.h>
+#include <igl/AABB.h>
+#include <igl/ARAPEnergyType.h>
+#include <igl/MeshBooleanType.h>
 #include <igl/SolverStatus.h>
 #include <igl/active_set.h>
-#include <igl/eigs.h>
-#include <igl/readOBJ.h>
-#include <igl/harmonic.h>
 #include <igl/arap.h>
-#include <igl/ARAPEnergyType.h>
+#include <igl/avg_edge_length.h>
+#include <igl/barycenter.h>
+#include <igl/barycentric_coordinates.h>
+#include <igl/boundary_facets.h>
 #include <igl/boundary_loop.h>
-#include <igl/map_vertices_to_circle.h>
-#include <igl/lscm.h>
-#include <igl/local_basis.h>
-#include <igl/rotate_vectors.h>
-#include <igl/compute_frame_field_bisectors.h>
+#include <igl/cat.h>
+#include <igl/colon.h>
 #include <igl/comb_cross_field.h>
+#include <igl/comb_frame_field.h>
+#include <igl/compute_frame_field_bisectors.h>
+#include <igl/cotmatrix.h>
+#include <igl/covariance_scatter_matrix.h>
 #include <igl/cross_field_missmatch.h>
-#include <igl/find_cross_field_singularities.h>
 #include <igl/cut_mesh_from_singularities.h>
-#include <igl/comb_frame_field.h>
+#include <igl/doublearea.h>
+#include <igl/edge_lengths.h>
+#include <igl/eigs.h>
+#include <igl/find_cross_field_singularities.h>
+#include <igl/fit_rotations.h>
+#include <igl/floor.h>
+#include <igl/gaussian_curvature.h>
+#include <igl/grad.h>
+#include <igl/harmonic.h>
+#include <igl/invert_diag.h>
+#include <igl/jet.h>
+#include <igl/local_basis.h>
+#include <igl/lscm.h>
+#include <igl/map_vertices_to_circle.h>
+#include <igl/massmatrix.h>
+#include <igl/min_quad_with_fixed.h>
 #include <igl/n_polyvector.h>
-
+#include <igl/parula.h>
+#include <igl/per_corner_normals.h>
+#include <igl/per_edge_normals.h>
+#include <igl/per_face_normals.h>
+#include <igl/per_vertex_normals.h>
+#include <igl/planarize_quad_mesh.h>
 #include <igl/point_mesh_squared_distance.h>
-#include <igl/AABB.h>
+#include <igl/polar_svd.h>
+#include <igl/principal_curvature.h>
+#include <igl/quad_planarity.h>
+#include <igl/readDMAT.h>
 #include <igl/readMESH.h>
-#include <igl/writeMESH.h>
+#include <igl/readOBJ.h>
+#include <igl/readOFF.h>
+#include <igl/read_triangle_mesh.h>
+#include <igl/rotate_vectors.h>
+#include <igl/setdiff.h>
+#include <igl/signed_distance.h>
+#include <igl/slice.h>
+#include <igl/slice_into.h>
+#include <igl/slice_mask.h>
 #include <igl/slice_tets.h>
-#include <igl/edge_lengths.h>
+#include <igl/sortrows.h>
+#include <igl/unique.h>
+#include <igl/unproject_onto_mesh.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>
-#include <igl/slice_mask.h>
-#include <igl/signed_distance.h>
-//#include <igl/.h>
+#include <igl/writeMESH.h>
+#include <igl/writeOBJ.h>
 
 
 void python_export_igl(py::module &m)
 {
-#include "py_igl/py_readOFF.cpp"
-#include "py_igl/py_writeOBJ.cpp"
-#include "py_igl/py_per_face_normals.cpp"
-#include "py_igl/py_per_corner_normals.cpp"
-#include "py_igl/py_per_vertex_normals.cpp"
-#include "py_igl/py_gaussian_curvature.cpp"
-#include "py_igl/py_jet.cpp"
-#include "py_igl/py_read_triangle_mesh.cpp"
-#include "py_igl/py_cotmatrix.cpp"
-#include "py_igl/py_massmatrix.cpp"
-#include "py_igl/py_invert_diag.cpp"
-#include "py_igl/py_principal_curvature.cpp"
-#include "py_igl/py_parula.cpp"
-#include "py_igl/py_readDMAT.cpp"
-#include "py_igl/py_grad.cpp"
-#include "py_igl/py_avg_edge_length.cpp"
-#include "py_igl/py_barycenter.cpp"
-#include "py_igl/py_doublearea.cpp"
-#include "py_igl/py_floor.cpp"
-#include "py_igl/py_slice.cpp"
-#include "py_igl/py_slice_into.cpp"
-#include "py_igl/py_sortrows.cpp"
-#include "py_igl/py_colon.cpp"
-#include "py_igl/py_boundary_facets.cpp"
-#include "py_igl/py_unique.cpp"
-#include "py_igl/py_setdiff.cpp"
-#include "py_igl/py_min_quad_with_fixed.cpp"
+#include "py_igl/py_AABB.cpp"
+#include "py_igl/py_ARAPEnergyType.cpp"
+#include "py_igl/py_MeshBooleanType.cpp"
 #include "py_igl/py_SolverStatus.cpp"
 #include "py_igl/py_active_set.cpp"
-#include "py_igl/py_eigs.cpp"
-#include "py_igl/py_readOBJ.cpp"
-#include "py_igl/py_harmonic.cpp"
-#include "py_igl/py_ARAPEnergyType.cpp"
 #include "py_igl/py_arap.cpp"
+#include "py_igl/py_avg_edge_length.cpp"
+#include "py_igl/py_barycenter.cpp"
+#include "py_igl/py_barycentric_coordinates.cpp"
+#include "py_igl/py_boundary_facets.cpp"
 #include "py_igl/py_boundary_loop.cpp"
-#include "py_igl/py_map_vertices_to_circle.cpp"
-#include "py_igl/py_lscm.cpp"
-#include "py_igl/py_local_basis.cpp"
-#include "py_igl/py_rotate_vectors.cpp"
-#include "py_igl/py_compute_frame_field_bisectors.cpp"
+#include "py_igl/py_cat.cpp"
+#include "py_igl/py_colon.cpp"
 #include "py_igl/py_comb_cross_field.cpp"
+#include "py_igl/py_comb_frame_field.cpp"
+#include "py_igl/py_compute_frame_field_bisectors.cpp"
+#include "py_igl/py_cotmatrix.cpp"
+#include "py_igl/py_covariance_scatter_matrix.cpp"
 #include "py_igl/py_cross_field_missmatch.cpp"
-#include "py_igl/py_find_cross_field_singularities.cpp"
 #include "py_igl/py_cut_mesh_from_singularities.cpp"
-#include "py_igl/py_comb_frame_field.cpp"
+#include "py_igl/py_doublearea.cpp"
+#include "py_igl/py_edge_lengths.cpp"
+#include "py_igl/py_eigs.cpp"
+#include "py_igl/py_find_cross_field_singularities.cpp"
+#include "py_igl/py_fit_rotations.cpp"
+#include "py_igl/py_floor.cpp"
+#include "py_igl/py_gaussian_curvature.cpp"
+#include "py_igl/py_grad.cpp"
+#include "py_igl/py_harmonic.cpp"
+#include "py_igl/py_invert_diag.cpp"
+#include "py_igl/py_jet.cpp"
+#include "py_igl/py_local_basis.cpp"
+#include "py_igl/py_lscm.cpp"
+#include "py_igl/py_map_vertices_to_circle.cpp"
+#include "py_igl/py_massmatrix.cpp"
+#include "py_igl/py_min_quad_with_fixed.cpp"
 #include "py_igl/py_n_polyvector.cpp"
-
+#include "py_igl/py_parula.cpp"
+#include "py_igl/py_per_corner_normals.cpp"
+#include "py_igl/py_per_edge_normals.cpp"
+#include "py_igl/py_per_face_normals.cpp"
+#include "py_igl/py_per_vertex_normals.cpp"
+#include "py_igl/py_planarize_quad_mesh.cpp"
 #include "py_igl/py_point_mesh_squared_distance.cpp"
-#include "py_igl/py_AABB.cpp"
+#include "py_igl/py_polar_svd.cpp"
+#include "py_igl/py_principal_curvature.cpp"
+#include "py_igl/py_quad_planarity.cpp"
+#include "py_igl/py_readDMAT.cpp"
 #include "py_igl/py_readMESH.cpp"
-#include "py_igl/py_writeMESH.cpp"
+#include "py_igl/py_readOBJ.cpp"
+#include "py_igl/py_readOFF.cpp"
+#include "py_igl/py_read_triangle_mesh.cpp"
+#include "py_igl/py_rotate_vectors.cpp"
+#include "py_igl/py_setdiff.cpp"
+#include "py_igl/py_signed_distance.cpp"
+#include "py_igl/py_slice.cpp"
+#include "py_igl/py_slice_into.cpp"
+#include "py_igl/py_slice_mask.cpp"
 #include "py_igl/py_slice_tets.cpp"
-#include "py_igl/py_edge_lengths.cpp"
+#include "py_igl/py_sortrows.cpp"
+#include "py_igl/py_unique.cpp"
+#include "py_igl/py_unproject_onto_mesh.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"
-#include "py_igl/py_slice_mask.cpp"
-#include "py_igl/py_signed_distance.cpp"
-
-//#include "py_igl/py_.cpp"
+#include "py_igl/py_writeMESH.cpp"
+#include "py_igl/py_writeOBJ.cpp"
 
 }

python/py_igl/copyleft/cgal/py_mesh_boolean.cpp

@@ -0,0 +1,107 @@
+// COMPLETE BINDINGS ========================
+
+m.def("mesh_boolean", []
+(
+  const Eigen::MatrixXd& VA,
+  const Eigen::MatrixXi& FA,
+  const Eigen::MatrixXd& VB,
+  const Eigen::MatrixXi& FB,
+  igl::MeshBooleanType & type,
+  Eigen::MatrixXd& VC,
+  Eigen::MatrixXi& FC,
+  Eigen::MatrixXi& J
+)
+{
+  return igl::copyleft::cgal::mesh_boolean(VA, FA, VB, FB, type, VC, FC, J);
+}, __doc_igl_copyleft_cgal_mesh_boolean,
+py::arg("VA"), py::arg("FA"), py::arg("VB"), py::arg("FB"), py::arg("type"), py::arg("VC"), py::arg("FC"), py::arg("J"));
+
+
+m.def("mesh_boolean", []
+(
+  const Eigen::MatrixXd& VA,
+  const Eigen::MatrixXi& FA,
+  const Eigen::MatrixXd& VB,
+  const Eigen::MatrixXi& FB,
+  const std::string & type_str,
+  Eigen::MatrixXd& VC,
+  Eigen::MatrixXi& FC,
+  Eigen::MatrixXi& J
+)
+{
+  return igl::copyleft::cgal::mesh_boolean(VA, FA, VB, FB, type_str, VC, FC, J);
+}, __doc_igl_copyleft_cgal_mesh_boolean,
+py::arg("VA"), py::arg("FA"), py::arg("VB"), py::arg("FB"), py::arg("type_str"), py::arg("VC"), py::arg("FC"), py::arg("J"));
+
+m.def("mesh_boolean", []
+(
+  const Eigen::MatrixXd& VA,
+  const Eigen::MatrixXi& FA,
+  const Eigen::MatrixXd& VB,
+  const Eigen::MatrixXi& FB,
+  const igl::MeshBooleanType & type,
+  Eigen::MatrixXd& VC,
+  Eigen::MatrixXi& FC
+)
+{
+  return igl::copyleft::cgal::mesh_boolean(VA, FA, VB, FB, type, VC, FC);
+}, __doc_igl_copyleft_cgal_mesh_boolean,
+py::arg("VA"), py::arg("FA"), py::arg("VB"), py::arg("FB"), py::arg("type"), py::arg("VC"), py::arg("FC"));
+
+
+
+// INCOMPLETE BINDINGS ========================
+
+
+
+
+//m.def("mesh_boolean", []
+//(
+//  const Eigen::MatrixXd& VA,
+//  const Eigen::MatrixXd& FA,
+//  const Eigen::MatrixXd& VB,
+//  const Eigen::MatrixXd& FB,
+//  std::function<int (const Eigen::Matrix<int, 1, Eigen::Dynamic>)> & wind_num_op,
+//  std::function<int (const int, const int)> & keep,
+//  Eigen::MatrixXd& VC,
+//  Eigen::MatrixXd& FC,
+//  Eigen::MatrixXd& J
+//)
+//{
+//  return igl::copyleft::cgal::mesh_boolean(VA, FA, VB, FB, wind_num_op, keep, VC, FC, J);
+//}, __doc_igl_copyleft_cgal_mesh_boolean,
+//py::arg("VA"), py::arg("FA"), py::arg("VB"), py::arg("FB"), py::arg("wind_num_op"), py::arg("keep"), py::arg("VC"), py::arg("FC"), py::arg("J"));
+
+//m.def("mesh_boolean", []
+//(
+//  std::vector<DerivedV> & Vlist,
+//  std::vector<DerivedF> & Flist,
+//  std::function<int (const Eigen::Matrix<int, 1, Eigen::Dynamic>)> & wind_num_op,
+//  std::function<int (const int, const int)> & keep,
+//  Eigen::MatrixXd& VC,
+//  Eigen::MatrixXd& FC,
+//  Eigen::MatrixXd& J
+//)
+//{
+//  return igl::copyleft::cgal::mesh_boolean(Vlist, Flist, wind_num_op, keep, VC, FC, J);
+//}, __doc_igl_copyleft_cgal_mesh_boolean,
+//py::arg("Vlist"), py::arg("Flist"), py::arg("wind_num_op"), py::arg("keep"), py::arg("VC"), py::arg("FC"), py::arg("J"));
+
+//m.def("mesh_boolean", []
+//(
+//  const Eigen::MatrixXd& VV,
+//  const Eigen::MatrixXd& FF,
+//  const Eigen::MatrixXd& sizes,
+//  std::function<int (const Eigen::Matrix<int, 1, Eigen::Dynamic>)> & wind_num_op,
+//  std::function<int (const int, const int)> & keep,
+//  Eigen::MatrixXd& VC,
+//  Eigen::MatrixXd& FC,
+//  Eigen::MatrixXd& J
+//)
+//{
+//  return igl::copyleft::cgal::mesh_boolean(VV, FF, sizes, wind_num_op, keep, VC, FC, J);
+//}, __doc_igl_copyleft_cgal_mesh_boolean,
+//py::arg("VV"), py::arg("FF"), py::arg("sizes"), py::arg("wind_num_op"), py::arg("keep"), py::arg("VC"), py::arg("FC"), py::arg("J"));
+
+
+

python/py_igl/copyleft/tetgen/py_tetrahedralize.cpp

@@ -0,0 +1,16 @@
+
+m.def("tetrahedralize", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const std::string switches,
+  Eigen::MatrixXd& TV,
+  Eigen::MatrixXi& TT,
+  Eigen::MatrixXi& TF
+)
+{
+  return igl::copyleft::tetgen::tetrahedralize(V, F, switches, TV, TT, TF);
+}, __doc_igl_copyleft_tetgen_tetrahedralize,
+py::arg("V"), py::arg("F"), py::arg("switches"), py::arg("TV"), py::arg("TT"), py::arg("TF"));
+
+

python/py_igl/embree/py_ambient_occlusion.cpp

@@ -0,0 +1,16 @@
+
+
+m.def("ambient_occlusion", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  const Eigen::MatrixXd& P,
+  const Eigen::MatrixXd& N,
+  const int num_samples,
+  Eigen::MatrixXd& S
+)
+{
+  return igl::embree::ambient_occlusion(V, F, P, N, num_samples, S);
+}, __doc_igl_embree_ambient_occlusion,
+py::arg("V"), py::arg("F"), py::arg("P"), py::arg("N"), py::arg("num_samples"), py::arg("S"));
+

python/py_igl/png/py_readPNG.cpp

@@ -0,0 +1,14 @@
+
+m.def("readPNG", []
+(
+  const std::string png_file,
+  Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & R,
+  Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & G,
+  Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & B,
+  Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & A
+)
+{
+  return igl::png::readPNG(png_file, R, G, B, A);
+}, __doc_igl_png_readPNG,
+py::arg("png_file"), py::arg("R"), py::arg("G"), py::arg("B"), py::arg("A"));
+

python/py_igl/png/py_writePNG.cpp

@@ -0,0 +1,15 @@
+
+
+m.def("writePNG", []
+(
+  const Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & R,
+  const Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & G,
+  const Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & B,
+  const Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic> & A,
+  const std::string png_file
+)
+{
+  return igl::png::writePNG(R, G, B, A, png_file);
+}, __doc_igl_png_writePNG,
+py::arg("R"), py::arg("G"), py::arg("B"), py::arg("A"), py::arg("png_file"));
+

python/py_igl/py_MeshBooleanType.cpp

@@ -0,0 +1,10 @@
+py::enum_<igl::MeshBooleanType>(m, "MeshBooleanType")
+    .value("MESH_BOOLEAN_TYPE_UNION", igl::MESH_BOOLEAN_TYPE_UNION)
+    .value("MESH_BOOLEAN_TYPE_INTERSECT", igl::MESH_BOOLEAN_TYPE_INTERSECT)
+    .value("MESH_BOOLEAN_TYPE_MINUS", igl::MESH_BOOLEAN_TYPE_MINUS)
+    .value("MESH_BOOLEAN_TYPE_XOR", igl::MESH_BOOLEAN_TYPE_XOR)
+    .value("MESH_BOOLEAN_TYPE_RESOLVE", igl::MESH_BOOLEAN_TYPE_RESOLVE)
+    .value("NUM_MESH_BOOLEAN_TYPES", igl::NUM_MESH_BOOLEAN_TYPES)
+    .export_values();
+
+

python/py_igl/py_planarize_quad_mesh.cpp

@@ -0,0 +1,14 @@
+
+m.def("planarize_quad_mesh", []
+(
+  const Eigen::MatrixXd& Vin,
+  const Eigen::MatrixXi& F,
+  const int maxIter,
+  const double & threshold,
+  Eigen::MatrixXd& Vout
+)
+{
+  return igl::planarize_quad_mesh(Vin, F, maxIter, threshold, Vout);
+}, __doc_igl_planarize_quad_mesh,
+py::arg("Vin"), py::arg("F"), py::arg("maxIter"), py::arg("threshold"), py::arg("Vout"));
+

python/py_igl/py_quad_planarity.cpp

@@ -0,0 +1,15 @@
+
+
+m.def("quad_planarity", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  Eigen::MatrixXd& P
+)
+{
+  Eigen::VectorXd Pv;
+  igl::quad_planarity(V, F, Pv);
+  P = Pv;
+}, __doc_igl_quad_planarity,
+py::arg("V"), py::arg("F"), py::arg("P"));
+

python/py_igl/py_slice.cpp

@@ -41,6 +41,19 @@ m.def("slice", []
 }, __doc_igl_slice,
 py::arg("X"), py::arg("R"), py::arg("C"), py::arg("Y"));
 
+m.def("slice", []
+(
+  const Eigen::MatrixXd& X,
+  const Eigen::MatrixXi& R,
+  const int dim,
+  Eigen::MatrixXd& Y
+)
+{
+  assert_is_VectorX("R",R);
+  return igl::slice(X,R,dim,Y);
+}, __doc_igl_slice,
+py::arg("X"), py::arg("R"), py::arg("dim"), py::arg("Y"));
+
 m.def("slice", []
 (
   const Eigen::MatrixXd& X,

python/py_igl/py_unproject_onto_mesh.cpp

@@ -0,0 +1,57 @@
+// COMPLETE BINDINGS ========================
+
+m.def("unproject_onto_mesh", []
+(
+  const Eigen::MatrixXd & pos,
+  const Eigen::MatrixXd & model,
+  const Eigen::MatrixXd & proj,
+  const Eigen::MatrixXd & viewport,
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& F,
+  Eigen::MatrixXi& fid, // TODO: Can we replace this with integer object reference?
+  Eigen::MatrixXd& bc
+)
+{
+  assert_is_Vector2("pos", pos);
+  Eigen::Vector2f posv;
+  if (pos.size() != 0)
+    posv = Eigen::Vector2f(pos.cast<float>());
+  assert_is_Matrix4("model", model);
+  Eigen::Matrix4f modelm;
+  if (model.size() != 0)
+    modelm = model.cast<float>();
+  assert_is_Matrix4("proj", proj);
+  Eigen::Matrix4f projm;
+  if (proj.size() != 0)
+    projm = proj.cast<float>();
+  assert_is_Vector4("viewport", viewport);
+  Eigen::Vector4f viewportv;
+  if (viewport.size() != 0)
+    viewportv = Eigen::Vector4f(viewport.cast<float>());
+
+  Eigen::VectorXd bcv;
+  int fidi;
+  bool ret = igl::unproject_onto_mesh(posv, modelm, projm, viewportv, V, F, fidi, bcv);
+  fid(0, 0) = fidi;
+  bc = bcv;
+  return ret;
+}, __doc_igl_unproject_onto_mesh,
+py::arg("pos"), py::arg("model"), py::arg("proj"), py::arg("viewport"), py::arg("V"), py::arg("F"), py::arg("fid"), py::arg("bc"));
+
+// INCOMPLETE BINDINGS ========================
+
+//m.def("unproject_onto_mesh", []
+//(
+//  Eigen::Vector2f & pos,
+//  Eigen::Matrix4f & model,
+//  Eigen::Matrix4f & proj,
+//  Eigen::Vector4f & viewport,
+//  std::function<bool (const Eigen::Vector3f &, const Eigen::Vector3f &, igl::Hit &)> & shoot_ray,
+//  int & fid,
+//  Eigen::MatrixXd& bc
+//)
+//{
+//  return igl::unproject_onto_mesh(pos, model, proj, viewport, shoot_ray, fid, bc);
+//}, __doc_igl_unproject_onto_mesh,
+//py::arg("pos"), py::arg("model"), py::arg("proj"), py::arg("viewport"), py::arg("shoot_ray"), py::arg("fid"), py::arg("bc"));
+

python/py_igl/triangle/py_triangulate.cpp

@@ -0,0 +1,16 @@
+
+
+m.def("triangulate", []
+(
+  const Eigen::MatrixXd& V,
+  const Eigen::MatrixXi& E,
+  const Eigen::MatrixXd& H,
+  const std::string flags,
+  Eigen::MatrixXd& V2,
+  Eigen::MatrixXi& F2
+)
+{
+  return igl::triangle::triangulate(V, E, H, flags, V2, F2);
+}, __doc_igl_triangle_triangulate,
+py::arg("V"), py::arg("E"), py::arg("H"), py::arg("flags"), py::arg("V2"), py::arg("F2"));
+

python/python_shared.cpp

@@ -14,6 +14,26 @@ extern void python_export_igl_viewer(py::module &);
 extern void python_export_igl_comiso(py::module &);
 #endif
 
+#ifdef PY_TETGEN
+extern void python_export_igl_tetgen(py::module &);
+#endif
+
+#ifdef PY_EMBREE
+extern void python_export_igl_embree(py::module &);
+#endif
+
+#ifdef PY_TRIANGLE
+extern void python_export_igl_triangle(py::module &);
+#endif
+
+#ifdef PY_CGAL
+extern void python_export_igl_cgal(py::module &);
+#endif
+
+#ifdef PY_PNG
+extern void python_export_igl_png(py::module &);
+#endif
+
 PYBIND11_PLUGIN(pyigl) {
     py::module m("pyigl", R"pyigldoc(
         Python wrappers for libigl
@@ -24,6 +44,80 @@ PYBIND11_PLUGIN(pyigl) {
         .. autosummary::
            :toctree: _generate
 
+           AABB
+           ARAPEnergyType
+           MeshBooleanType
+           SolverStatus
+           active_set
+           arap
+           avg_edge_length
+           barycenter
+           barycentric_coordinates
+           boundary_facets
+           boundary_loop
+           cat
+           colon
+           comb_cross_field
+           comb_frame_field
+           compute_frame_field_bisectors
+           copyleft_cgal_mesh_boolean
+           copyleft_comiso_miq
+           copyleft_comiso_nrosy
+           copyleft_tetgen_tetrahedralize
+           cotmatrix
+           covariance_scatter_matrix
+           cross_field_missmatch
+           cut_mesh_from_singularities
+           doublearea
+           edge_lengths
+           eigs
+           embree_ambient_occlusion
+           find_cross_field_singularities
+           fit_rotations
+           floor
+           gaussian_curvature
+           grad
+           harmonic
+           invert_diag
+           jet
+           local_basis
+           lscm
+           map_vertices_to_circle
+           massmatrix
+           min_quad_with_fixed
+           n_polyvector
+           parula
+           per_corner_normals
+           per_edge_normals
+           per_face_normals
+           per_vertex_normals
+           planarize_quad_mesh
+           png_readPNG
+           png_writePNG
+           point_mesh_squared_distance
+           polar_svd
+           principal_curvature
+           quad_planarity
+           readDMAT
+           readMESH
+           readOBJ
+           readOFF
+           read_triangle_mesh
+           rotate_vectors
+           setdiff
+           signed_distance
+           slice
+           slice_into
+           slice_mask
+           slice_tets
+           sortrows
+           triangle_triangulate
+           unique
+           unproject_onto_mesh
+           upsample
+           writeMESH
+           writeOBJ
+
     )pyigldoc");
 
     python_export_vector(m);
@@ -38,5 +132,25 @@ PYBIND11_PLUGIN(pyigl) {
     python_export_igl_comiso(m);
     #endif
 
+    #ifdef PY_TETGEN
+    python_export_igl_tetgen(m);
+    #endif
+
+    #ifdef PY_EMBREE
+    python_export_igl_embree(m);
+    #endif
+
+    #ifdef PY_TRIANGLE
+    python_export_igl_triangle(m);
+    #endif
+
+    #ifdef PY_CGAL
+    python_export_igl_cgal(m);
+    #endif
+
+    #ifdef PY_PNG
+    python_export_igl_png(m);
+    #endif
+
     return m.ptr();
 }

python/python_shared.h

@@ -31,6 +31,26 @@ void assert_is_RowVectorX(const std::string name, const Eigen::PlainObjectBase<S
     throw std::runtime_error(name + " must be a row vector.");
 }
 
+template<typename Scalar>
+void assert_is_Vector2(const std::string name, const Eigen::PlainObjectBase<Scalar>& v)
+{
+  if (v.size() == 0)
+    return;
+
+  if ((v.cols() != 1) || (v.rows() != 2))
+    throw std::runtime_error(name + " must be a column vector with 2 entries.");
+}
+
+template<typename Scalar>
+void assert_is_RowVector2(const std::string name, const Eigen::PlainObjectBase<Scalar>& v)
+{
+  if (v.size() == 0)
+    return;
+
+  if ((v.cols() != 2) || (v.rows() != 1))
+    throw std::runtime_error(name + " must be a row vector with 2 entries.");
+}
+
 template<typename Scalar>
 void assert_is_Vector3(const std::string name, const Eigen::PlainObjectBase<Scalar>& v)
 {

python/scripts/generate_bindings.py

@@ -63,6 +63,9 @@ def map_parameter_types(name, cpp_type, parsed_types, errors, enum_types):
     if cpp_type.startswith("MatY"):
         result.append("Eigen::SparseMatrix<double>&")
         skip_parsing = True
+    if cpp_type.startswith("Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic>"):
+        result.append("Eigen::Matrix<unsigned char, Eigen::Dynamic, Eigen::Dynamic>")
+        skip_parsing = True
     if cpp_type == "std::vector<std::vector<Scalar> > &":
         result.append("std::vector<std::vector<double> > &")
         skip_parsing = True
@@ -298,4 +301,4 @@ if __name__ == '__main__':
         for k in l:
             fs.write("%s: %i \n" %(k, len(files[k])))
             fs.writelines("\n".join(files[k]))
-            fs.write("\n\n\n")
+            fs.write("\n\n\n")

python/scripts/generate_docstrings.py

@@ -21,6 +21,7 @@ def get_filepaths(directory):
     it yields a 3-tuple (dirpath, dirnames, filenames).
     """
     file_paths = []  # List which will store all of the full filepaths.
+    root_file_paths = []
 
     # Walk the tree.
     for root, directories, files in os.walk(directory):
@@ -29,7 +30,10 @@ def get_filepaths(directory):
             filepath = os.path.join(root, filename)
             file_paths.append(filepath)  # Add it to the list.
 
-    return file_paths  # Self-explanatory.
+            if root.endswith(directory): # Add only the files in the root directory
+                root_file_paths.append(filepath)
+
+    return file_paths, root_file_paths  # file_paths contains all file paths, core_file_paths only the ones in <directory>
 
 
 def get_name_from_path(path, basepath, prefix, postfix):
@@ -53,8 +57,8 @@ if __name__ == '__main__':
     # List all files in the given folder and subfolders
     cpp_base_path = sys.argv[1]
     py_base_path = sys.argv[2]
-    cpp_file_paths = get_filepaths(cpp_base_path)
-    py_file_paths = get_filepaths(py_base_path)
+    cpp_file_paths, cpp_root_file_paths = get_filepaths(cpp_base_path)
+    py_file_paths, py_root_file_paths = get_filepaths(py_base_path)
 
     # Add all the .h filepaths to a dict
     mapping = {}
@@ -65,10 +69,16 @@ if __name__ == '__main__':
 
     # Add all python binding files to a list
     implemented_names = []
+    core_implemented_names = []
     for f in py_file_paths:
         if f.endswith(".cpp"):
             name = get_name_from_path(f, py_base_path, "py_", ".cpp")
             implemented_names.append(name)
+            if f in py_root_file_paths:
+                core_implemented_names.append(name)
+
+    implemented_names.sort()
+    core_implemented_names.sort()
 
     # Create a list of cpp header files for which a python binding file exists
     files_to_parse = []
@@ -124,3 +134,9 @@ if __name__ == '__main__':
     rendered = tpl.render(functions=implemented_names)
     with open("../python_shared.cpp", 'w') as fs:
         fs.write(rendered)
+
+    # Write py_igl_cpp file with all core library files
+    tpl = Template(filename='py_igl.mako')
+    rendered = tpl.render(functions=core_implemented_names)
+    with open("../py_igl.cpp", 'w') as fs:
+        fs.write(rendered)

python/scripts/py_igl.mako

@@ -0,0 +1,16 @@
+#include <Eigen/Dense>
+
+#include "python_shared.h"
+
+% for f in functions:
+#include <igl/${f}.h>
+% endfor
+
+
+void python_export_igl(py::module &m)
+{
+% for f in functions:
+#include "py_igl/py_${f}.cpp"
+% endfor
+
+}

python/scripts/python_shared.mako

@@ -14,6 +14,26 @@ extern void python_export_igl_viewer(py::module &);
 extern void python_export_igl_comiso(py::module &);
 #endif
 
+#ifdef PY_TETGEN
+extern void python_export_igl_tetgen(py::module &);
+#endif
+
+#ifdef PY_EMBREE
+extern void python_export_igl_embree(py::module &);
+#endif
+
+#ifdef PY_TRIANGLE
+extern void python_export_igl_triangle(py::module &);
+#endif
+
+#ifdef PY_CGAL
+extern void python_export_igl_cgal(py::module &);
+#endif
+
+#ifdef PY_PNG
+extern void python_export_igl_png(py::module &);
+#endif
+
 PYBIND11_PLUGIN(pyigl) {
     py::module m("pyigl", R"pyigldoc(
         Python wrappers for libigl
@@ -42,5 +62,25 @@ PYBIND11_PLUGIN(pyigl) {
     python_export_igl_comiso(m);
     #endif
 
+    #ifdef PY_TETGEN
+    python_export_igl_tetgen(m);
+    #endif
+
+    #ifdef PY_EMBREE
+    python_export_igl_embree(m);
+    #endif
+
+    #ifdef PY_TRIANGLE
+    python_export_igl_triangle(m);
+    #endif
+
+    #ifdef PY_CGAL
+    python_export_igl_cgal(m);
+    #endif
+
+    #ifdef PY_PNG
+    python_export_igl_png(m);
+    #endif
+
     return m.ptr();
 }

python/tutorial/001_BasicTypes.py

@@ -4,7 +4,7 @@ from iglhelpers import *
 
 # Create a numpy dense array
 # 2 types are supported by the wrappers: float64 and int64
-dense_matrix = np.array( [ (1,2,3), (4,5,6) , (7,8,9) ], dtype='float64')
+dense_matrix = np.array([(1, 2, 3), (4, 5, 6), (7, 8, 9)], dtype='float64')
 
 # libigl wrappers uses Eigen as a matrix type, you can easily convert between numpy and Eigen using
 # the helper function p2e. This operation duplicates the data.
@@ -18,7 +18,7 @@ dense_matrix_eigen_2 = dense_matrix_eigen * dense_matrix_eigen
 print("Eigen Matrix: \n", dense_matrix_eigen_2, "\n", sep='')
 
 # and access single elements
-print("Eigen Matrix(0,0): ", dense_matrix_eigen_2[0,0], "\n")
+print("Eigen Matrix(0,0): ", dense_matrix_eigen_2[0, 0], "\n")
 
 # To convert it back to a numpy array, use the helper function e2p
 dense_matrix_2 = e2p(dense_matrix_eigen_2)
@@ -39,10 +39,3 @@ print("Sparse matrix Eigen: ", sparse_matrix_eigen, sep='')
 # And converted back with e2p
 sparse_matrix_2 = e2p(sparse_matrix_eigen)
 print("Sparse matrix Numpy: ", sparse_matrix_2.todense(), sep='')
-
-
-
-
-
-
-

python/tutorial/101_FileIO.py

@@ -1,14 +1,16 @@
-from __future__ import print_function
-# Add the igl library to the modules search path
 import sys, os
-sys.path.insert(0, os.getcwd() + "/../")
 
+# Add the igl library to the modules search path
+sys.path.insert(0, os.getcwd() + "/../")
 import pyigl as igl
 
+from shared import TUTORIAL_SHARED_PATH
+
+
 # Load a mesh in OFF format
 V = igl.eigen.MatrixXd()
 F = igl.eigen.MatrixXi()
-igl.readOFF("../../tutorial/shared/cube.off", V, F)
+igl.readOFF(TUTORIAL_SHARED_PATH + "cube.off", V, F)
 
 # Print the vertices and faces matrices
 print("Vertices: \n", V, sep='')

python/tutorial/102_DrawMesh.py

@@ -1,15 +1,21 @@
-# Add the igl library to the modules search path
 import sys, os
-sys.path.insert(0, os.getcwd() + "/../")
 
+# Add the igl library to the modules search path
+sys.path.insert(0, os.getcwd() + "/../")
 import pyigl as igl
 
+from shared import TUTORIAL_SHARED_PATH, check_dependencies
+
+dependencies = ["viewer"]
+check_dependencies(dependencies)
+
+
 # Load a mesh in OFF format
 V = igl.eigen.MatrixXd()
 F = igl.eigen.MatrixXi()
-igl.readOFF("../../tutorial/shared/beetle.off", V, F)
+igl.readOFF(TUTORIAL_SHARED_PATH + "beetle.off", V, F)
 
 # Plot the mesh
-viewer = igl.viewer.Viewer();
-viewer.data.set_mesh(V, F);
-viewer.launch();
+viewer = igl.viewer.Viewer()
+viewer.data.set_mesh(V, F)
+viewer.launch()

python/tutorial/102_DrawMesh_TCP.py

@@ -1,11 +1,15 @@
-## This is a test application for the TCPViewer
-
-# Add the igl library to the modules search path
 import sys, os
+import time
+# Add the igl library to the modules search path
 sys.path.insert(0, os.getcwd() + "/../")
+import pyigl as igl
 
-import os
-import time
+import tcpviewer
+
+from shared import TUTORIAL_SHARED_PATH
+
+
+## This is a test application for the TCPViewer
 
 # Launch the tcp viewer
 os.system("python ../tcpviewer.py&")
@@ -13,13 +17,10 @@ os.system("python ../tcpviewer.py&")
 # Wait for it to set up the socket
 time.sleep(1)
 
-import pyigl as igl
-import tcpviewer
-
 # Read a mesh
 V = igl.eigen.MatrixXd()
 F = igl.eigen.MatrixXi()
-igl.readOFF('../../tutorial/shared/beetle.off', V, F)
+igl.readOFF(TUTORIAL_SHARED_PATH + "beetle.off", V, F)
 
 # Send it to the viewer
 viewer = tcpviewer.TCPViewer()

python/tutorial/103_Events.py

@@ -1,9 +1,15 @@
-# Add the igl library to the modules search path
 import sys, os
-sys.path.insert(0, os.getcwd() + "/../")
 
+# Add the igl library to the modules search path
+sys.path.insert(0, os.getcwd() + "/../")
 import pyigl as igl
 
+from shared import TUTORIAL_SHARED_PATH, check_dependencies
+
+dependencies = ["viewer"]
+check_dependencies(dependencies)
+
+
 V1 = igl.eigen.MatrixXd()
 F1 = igl.eigen.MatrixXi()
 
@@ -15,22 +21,22 @@ def key_pressed(viewer, key, modifier):
 
     if key == ord('1'):
         # # Clear should be called before drawing the mesh
-        viewer.data.clear();
+        viewer.data.clear()
         # # Draw_mesh creates or updates the vertices and faces of the displayed mesh.
         # # If a mesh is already displayed, draw_mesh returns an error if the given V and
         # # F have size different than the current ones
-        viewer.data.set_mesh(V1, F1);
-        viewer.core.align_camera_center(V1,F1);
+        viewer.data.set_mesh(V1, F1)
+        viewer.core.align_camera_center(V1,F1)
     elif key == ord('2'):
-        viewer.data.clear();
-        viewer.data.set_mesh(V2, F2);
-        viewer.core.align_camera_center(V2,F2);
+        viewer.data.clear()
+        viewer.data.set_mesh(V2, F2)
+        viewer.core.align_camera_center(V2,F2)
     return False
 
 
 #  Load two meshes
-igl.readOFF("../../tutorial/shared/bumpy.off", V1, F1);
-igl.readOFF("../../tutorial/shared/fertility.off", V2, F2);
+igl.readOFF(TUTORIAL_SHARED_PATH + "bumpy.off", V1, F1)
+igl.readOFF(TUTORIAL_SHARED_PATH + "fertility.off", V2, F2)
 
 print("1 Switch to bump mesh")
 print("2 Switch to fertility mesh")

python/tutorial/104_Colors.py

@@ -1,25 +1,31 @@
-# Add the igl library to the modules search path
 import sys, os
-sys.path.insert(0, os.getcwd() + "/../")
 
+# Add the igl library to the modules search path
+sys.path.insert(0, os.getcwd() + "/../")
 import pyigl as igl
 
+from shared import TUTORIAL_SHARED_PATH, check_dependencies
+
+dependencies = ["viewer"]
+check_dependencies(dependencies)
+
+
 V = igl.eigen.MatrixXd()
 F = igl.eigen.MatrixXi()
 C = igl.eigen.MatrixXd()
 
 # Load a mesh in OFF format
-igl.readOFF("../../tutorial/shared/screwdriver.off", V, F)
+igl.readOFF(TUTORIAL_SHARED_PATH + "screwdriver.off", V, F)
 
 # Plot the mesh
 viewer = igl.viewer.Viewer()
 viewer.data.set_mesh(V, F)
 
 # Use the z coordinate as a scalar field over the surface
-Z = V.col(2);
+Z = V.col(2)
 
 # Compute per-vertex colors
-igl.jet(Z,True,C)
+igl.jet(Z, True, C)
 
 # Add per-vertex colors
 viewer.data.set_colors(C)