merge

.appveyor.yml

@@ -25,7 +25,11 @@ build_script:
   - set PATH=C:\Python%PYTHON%-x64;C:\Python%PYTHON%-x64\Scripts;%PATH%
   - mkdir build
   - cd build
-  - cmake -DCMAKE_BUILD_TYPE=%CONFIG% -DLIBIGL_WITH_CGAL=ON -DLIBIGL_WITH_COMISO=OFF -DEMBREE_ISA_AVX=OFF -DEMBREE_ISA_AVX2=OFF -DEMBREE_ISA_AVX512SKX=OFF -G "Visual Studio 15 2017 Win64" ../
+  - cmake -DCMAKE_BUILD_TYPE=%CONFIG%
+      -DLIBIGL_WITH_CGAL=ON
+      -DLIBIGL_WITH_COMISO=OFF
+      -G "Visual Studio 15 2017 Win64"
+      ../
   - set MSBuildLogger="C:\Program Files\AppVeyor\BuildAgent\Appveyor.MSBuildLogger.dll"
   - set MSBuildOptions=/v:m /m /p:BuildInParallel=true /p:Configuration=%CONFIG% /logger:%MSBuildLogger%
   - msbuild %MSBuildOptions% libigl.sln

.travis.yml

@@ -32,18 +32,26 @@ matrix:
     compiler: gcc-7
     env:
     - MATRIX_EVAL="export CC=gcc-7 CXX=g++-7 CONFIG=Release PYTHON=python3"
+  - os: linux # same config like above but with -DLIBIGL_USE_STATIC_LIBRARY=OFF to test static and header-only builds
+    compiler: gcc-7
+    env:
+    - MATRIX_EVAL="export CC=gcc-7 CXX=g++-7 CONFIG=Release PYTHON=python3 CMAKE_EXTRA='-DLIBIGL_USE_STATIC_LIBRARY=OFF'"
   - os: linux
     compiler: gcc-7
     env:
-    - MATRIX_EVAL="export CC=gcc-7 CXX=g++-7 CONFIG=Release PYTHON=python3 CMAKE_EXTRA='-DLIBIGL_EIGEN_VERSION=3.3.7 -DLIBIGL_EIGEN_MD5=f2a417d083fe8ca4b8ed2bc613d20f07'"
+    - MATRIX_EVAL="export CC=gcc-7 CXX=g++-7 CONFIG=Release PYTHON=python3 CMAKE_EXTRA='-DLIBIGL_EIGEN_VERSION=3.3.7'"
   - os: osx
     compiler: clang
     env:
     - MATRIX_EVAL="export CONFIG=Debug PYTHON=python3 LIBIGL_NUM_THREADS=1"
+  - os: osx # same config like above but with -DLIBIGL_USE_STATIC_LIBRARY=OFF to test static and header-only builds
+    compiler: clang
+    env:
+    - MATRIX_EVAL="export CONFIG=Debug PYTHON=python3 LIBIGL_NUM_THREADS=1 CMAKE_EXTRA='-DLIBIGL_USE_STATIC_LIBRARY=OFF'"
   - os: osx
     compiler: clang
     env:
-    - MATRIX_EVAL="export CONFIG=Debug PYTHON=python3 LIBIGL_NUM_THREADS=1 CMAKE_EXTRA='-DLIBIGL_EIGEN_VERSION=3.3.7 -DLIBIGL_EIGEN_MD5=f2a417d083fe8ca4b8ed2bc613d20f07'""
+    - MATRIX_EVAL="export CONFIG=Debug PYTHON=python3 LIBIGL_NUM_THREADS=1 CMAKE_EXTRA='-DLIBIGL_EIGEN_VERSION=3.3.7'""
 
 install:
 - if [[ "$TRAVIS_OS_NAME" == "osx" ]]; then export PATH="/usr/local/opt/ccache/libexec:$PATH"; fi
@@ -55,7 +63,11 @@ script:
 # Tutorials and tests
 - mkdir build
 - pushd build
-- cmake ${CMAKE_EXTRA} -DCMAKE_BUILD_TYPE=$CONFIG -DLIBIGL_CHECK_UNDEFINED=ON -DLIBIGL_WITH_CGAL=ON -DEMBREE_ISA_AVX=OFF -DEMBREE_ISA_AVX2=OFF -DEMBREE_ISA_AVX512SKX=OFF ../
+- cmake ${CMAKE_EXTRA}
+    -DCMAKE_BUILD_TYPE=$CONFIG
+    -DLIBIGL_CHECK_UNDEFINED=ON
+    -DLIBIGL_WITH_CGAL=ON
+    ../
 - make -j 2
 - ctest --verbose
 - popd

CMakeLists.txt

@@ -2,7 +2,8 @@ cmake_minimum_required(VERSION 3.1)
 project(libigl)
 
 # Detects whether this is a top-level project
-if(CMAKE_SOURCE_DIR STREQUAL CMAKE_CURRENT_SOURCE_DIR)
+get_directory_property(LIBIGL_PARENT_DIR PARENT_DIRECTORY)
+if(NOT LIBIGL_PARENT_DIR)
 	set(LIBIGL_TOPLEVEL_PROJECT ON)
 else()
 	set(LIBIGL_TOPLEVEL_PROJECT OFF)
@@ -29,6 +30,7 @@ option(LIBIGL_WITH_OPENGL_GLFW_IMGUI "Use ImGui"                    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_PREDICATES        "Use exact predicates"         ON)
 option(LIBIGL_WITH_XML               "Use XML"                      ON)
 option(LIBIGL_WITH_PYTHON            "Use Python"                   ${LIBIGL_BUILD_PYTHON})
 ### End

cmake/LibiglDownloadExternal.cmake

@@ -2,12 +2,21 @@
 include(DownloadProject)
 
 # With CMake 3.8 and above, we can hide warnings about git being in a
-# detached head by passing an extra GIT_CONFIG option
+# detached head by passing an extra GIT_CONFIG option.
 set(LIBIGL_EXTRA_OPTIONS TLS_VERIFY OFF)
 if(NOT (${CMAKE_VERSION} VERSION_LESS "3.8.0"))
 	list(APPEND LIBIGL_EXTRA_OPTIONS GIT_CONFIG advice.detachedHead=false)
 endif()
 
+# On CMake 3.6.3 and above, there is an option to use shallow clones of git repositories.
+# The shallow clone option only works with real tags, not SHA1, so we use a separate option.
+set(LIBIGL_BRANCH_OPTIONS)
+if(NOT (${CMAKE_VERSION} VERSION_LESS "3.6.3"))
+	# Disabled for now until we can make sure that it has no adverse effects
+	# (Downside is that the eigen mirror is huge again)
+	# list(APPEND LIBIGL_BRANCH_OPTIONS GIT_SHALLOW 1)
+endif()
+
 option(LIBIGL_SKIP_DOWNLOAD "Skip downloading external libraries" OFF)
 
 # Shortcut functions
@@ -55,24 +64,21 @@ function(igl_download_cork)
 endfunction()
 
 ## Eigen
-# 3.2.10 8ad10ac703a78143a4062c9bda9d8fd3
-# 3.3.7  f2a417d083fe8ca4b8ed2bc613d20f07
-set(LIBIGL_EIGEN_VERSION 3.2.10                           CACHE STRING "Default version of Eigen used by libigl.")
-set(LIBIGL_EIGEN_MD5     8ad10ac703a78143a4062c9bda9d8fd3 CACHE STRING "md5sum of the tar.gz archive corresponding to this version.")
+set(LIBIGL_EIGEN_VERSION 3.2.10 CACHE STRING "Default version of Eigen used by libigl.")
 function(igl_download_eigen)
 	igl_download_project(eigen
-		URL           http://bitbucket.org/eigen/eigen/get/${LIBIGL_EIGEN_VERSION}.tar.gz
-		URL_MD5       ${LIBIGL_EIGEN_MD5}
+		GIT_REPOSITORY https://github.com/eigenteam/eigen-git-mirror.git
+		GIT_TAG        ${LIBIGL_EIGEN_VERSION}
+		${LIBIGL_BRANCH_OPTIONS}
 	)
 endfunction()
 
 ## Embree
 function(igl_download_embree)
 	igl_download_project(embree
-		URL            https://github.com/embree/embree/archive/v3.2.3.tar.gz
-		URL_MD5        1868cda1c97d83d7a0b67b0b64b18cef
-		# GIT_REPOSITORY https://github.com/embree/embree.git
-		# GIT_TAG        cb61322db3bb7082caed21913ad14869b436fe78
+		GIT_REPOSITORY https://github.com/embree/embree.git
+		GIT_TAG        v3.5.2
+		${LIBIGL_BRANCH_OPTIONS}
 	)
 endfunction()
 
@@ -88,7 +94,8 @@ endfunction()
 function(igl_download_glfw)
 	igl_download_project(glfw
 		GIT_REPOSITORY https://github.com/glfw/glfw.git
-		GIT_TAG        53c8c72c676ca97c10aedfe3d0eb4271c5b23dba
+		GIT_TAG        3.3
+		${LIBIGL_BRANCH_OPTIONS}
 	)
 endfunction()
 
@@ -97,6 +104,7 @@ function(igl_download_imgui)
 	igl_download_project(imgui
 		GIT_REPOSITORY https://github.com/ocornut/imgui.git
 		GIT_TAG        v1.69
+		${LIBIGL_BRANCH_OPTIONS}
 	)
 	igl_download_project(libigl-imgui
 		GIT_REPOSITORY https://github.com/libigl/libigl-imgui.git
@@ -139,8 +147,8 @@ endfunction()
 ## Triangle
 function(igl_download_triangle)
 	igl_download_project(triangle
-		GIT_REPOSITORY https://github.com/jdumas/triangle.git
-		GIT_TAG        2cd0672ff1f67f9f6bb8e556e84901293e637b76
+		GIT_REPOSITORY https://github.com/libigl/triangle.git
+		GIT_TAG        d284c4a843efac043c310f5fa640b17cf7d96170
 	)
 endfunction()
 
@@ -152,6 +160,14 @@ function(igl_download_catch2)
 	)
 endfunction()
 
+## Predicates
+function(igl_download_predicates)
+	igl_download_project(predicates
+		GIT_REPOSITORY https://github.com/libigl/libigl-predicates.git
+		GIT_TAG        4c57c1d3f31646b010d1d58bfbe201e75c2b2ad8
+	)
+endfunction()
+
 ################################################################################
 
 ## Test data

cmake/libigl.cmake

@@ -32,6 +32,7 @@ option(LIBIGL_WITH_OPENGL_GLFW_IMGUI "Use ImGui"                    OFF)
 option(LIBIGL_WITH_PNG               "Use PNG"                      OFF)
 option(LIBIGL_WITH_TETGEN            "Use Tetgen"                   OFF)
 option(LIBIGL_WITH_TRIANGLE          "Use Triangle"                 OFF)
+option(LIBIGL_WITH_PREDICATES        "Use exact predicates"         OFF)
 option(LIBIGL_WITH_XML               "Use XML"                      OFF)
 option(LIBIGL_WITH_PYTHON            "Use Python"                   OFF)
 option(LIBIGL_WITHOUT_COPYLEFT       "Disable Copyleft libraries"   OFF)
@@ -270,11 +271,11 @@ endif()
 if(LIBIGL_WITH_EMBREE)
   set(EMBREE_DIR "${LIBIGL_EXTERNAL}/embree")
 
-  set(EMBREE_TESTING_INTENSITY 0 CACHE INT "" FORCE)
+  set(EMBREE_TESTING_INTENSITY 0 CACHE STRING "" FORCE)
   set(EMBREE_ISPC_SUPPORT OFF CACHE BOOL " " FORCE)
   set(EMBREE_TASKING_SYSTEM "INTERNAL" CACHE BOOL " " FORCE)
   set(EMBREE_TUTORIALS OFF CACHE BOOL " " FORCE)
-  set(EMBREE_MAX_ISA NONE CACHE STRINGS " " FORCE)
+  set(EMBREE_MAX_ISA "SSE2" CACHE STRING " " FORCE)
   set(EMBREE_STATIC_LIB ON CACHE BOOL " " FORCE)
   if(MSVC)
     set(EMBREE_STATIC_RUNTIME ON CACHE BOOL " " FORCE)
@@ -296,7 +297,7 @@ endif()
 ################################################################################
 ### Compile the matlab part ###
 if(LIBIGL_WITH_MATLAB)
-  find_package(Matlab REQUIRED COMPONENTS MEX_COMPILER MX_LIBRARY ENG_LIBRARY)
+  find_package(Matlab REQUIRED COMPONENTS MEX_COMPILER MX_LIBRARY ENG_LIBRARY MAT_LIBRARY)
   compile_igl_module("matlab")
   target_link_libraries(igl_matlab ${IGL_SCOPE} ${Matlab_LIBRARIES})
   target_include_directories(igl_matlab ${IGL_SCOPE} ${Matlab_INCLUDE_DIRS})
@@ -411,6 +412,20 @@ if(LIBIGL_WITH_TRIANGLE)
   target_include_directories(igl_triangle ${IGL_SCOPE} ${TRIANGLE_DIR})
 endif()
 
+################################################################################
+### Compile the predicates part ###
+if(LIBIGL_WITH_PREDICATES)
+  set(PREDICATES_DIR "${LIBIGL_EXTERNAL}/predicates")
+  if(NOT TARGET predicates)
+    igl_download_predicates()
+    add_subdirectory("${PREDICATES_DIR}" "predicates")
+  endif()
+  compile_igl_module("predicates")
+  target_link_libraries(igl_predicates ${IGL_SCOPE} predicates)
+  target_include_directories(igl_predicates ${IGL_SCOPE} ${PREDICATES_DIR})
+  target_compile_definitions(igl_predicates ${IGL_SCOPE} -DLIBIGL_WITH_PREDICATES)
+endif()
+
 ################################################################################
 ### Compile the xml part ###
 if(LIBIGL_WITH_XML)

include/igl/AABB.cpp

@@ -1073,6 +1073,9 @@ template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 2>::init<Eigen
 template double igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_distance<Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<double, 1, 3, 1, 1, 3> const&, double, int&, Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> >&) const;
 template bool igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::intersect_ray<Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<double, 1, 3, 1, 1, 3> const&, Eigen::Matrix<double, 1, 3, 1, 1, 3> const&, igl::Hit&) const;
 template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_distance<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, 2, 3, 0, 2, 3>, Eigen::Matrix<double, 2, 1, 0, 2, 1>, Eigen::Matrix<int, 2, 1, 0, 2, 1>, Eigen::Matrix<double, 2, 3, 0, 2, 3> >(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<double, 2, 3, 0, 2, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 1, 0, 2, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, 2, 1, 0, 2, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, 2, 3, 0, 2, 3> >&) const;
+template void igl::AABB<Eigen::Matrix<double, -1, -1, 0, -1, -1>, 3>::squared_distance<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(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<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>, 2>::squared_distance<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(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<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<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(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<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;
 #ifdef WIN32
 template void igl::AABB<class Eigen::Matrix<double,-1,-1,0,-1,-1>,2>::squared_distance<class Eigen::Matrix<int,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,1,0,-1,1>,class Eigen::Matrix<__int64,-1,1,0,-1,1>,class Eigen::Matrix<double,-1,3,0,-1,3> >(class Eigen::MatrixBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::MatrixBase<class Eigen::Matrix<int,-1,-1,0,-1,-1> > const &,class Eigen::MatrixBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<__int64,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,3,0,-1,3> > &)const;
 template void igl::AABB<class Eigen::Matrix<double,-1,-1,0,-1,-1>,3>::squared_distance<class Eigen::Matrix<int,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,1,0,-1,1>,class Eigen::Matrix<__int64,-1,1,0,-1,1>,class Eigen::Matrix<double,-1,3,0,-1,3> >(class Eigen::MatrixBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::MatrixBase<class Eigen::Matrix<int,-1,-1,0,-1,-1> > const &,class Eigen::MatrixBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<__int64,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,3,0,-1,3> > &)const;

include/igl/adjacency_list.cpp

@@ -132,8 +132,18 @@ IGL_INLINE void igl::adjacency_list(
   const std::vector<std::vector<Index> > & F,
   std::vector<std::vector<Index> >& A)
 {
-  A.clear(); 
-  A.resize(F.maxCoeff()+1);
+  A.clear();
+  
+  // Find maxCoeff
+  Index maxCoeff = 0;
+  for(const auto &vec : F)
+  {
+    for(int coeff : vec)
+    {
+      maxCoeff = std::max(coeff, maxCoeff);
+    }
+  }
+  A.resize(maxCoeff + 1);
   
   // Loop over faces
   for(int i = 0;i<F.size();i++)
@@ -142,8 +152,8 @@ IGL_INLINE void igl::adjacency_list(
     for(int j = 0;j<F[i].size();j++)
     {
       // Get indices of edge: s --> d
-      int s = F(i,j);
-      int d = F(i,(j+1)%F[i].size());
+      int s = F[i][j];
+      int d = F[i][(j+1)%F[i].size()];
       A.at(s).push_back(d);
       A.at(d).push_back(s);
     }
@@ -166,4 +176,5 @@ template void igl::adjacency_list<Eigen::Matrix<int, -1, 2, 0, -1, 2>, int>(Eige
 template void igl::adjacency_list<Eigen::Matrix<int, -1, -1, 0, -1, -1>, int>(Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&, bool);
 template void igl::adjacency_list<Eigen::Matrix<int, -1, 3, 0, -1, 3>, int>(Eigen::MatrixBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&, bool);
 template void igl::adjacency_list<class Eigen::Matrix<int, -1, -1, 0, -1, -1>, unsigned int>(class Eigen::MatrixBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1> > const &, class std::vector<class std::vector<unsigned int, class std::allocator<unsigned int> >, class std::allocator<class std::vector<unsigned int, class std::allocator<unsigned int> > > > &, bool);
+template void igl::adjacency_list<int>(std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > > const&, std::vector<std::vector<int, std::allocator<int> >, std::allocator<std::vector<int, std::allocator<int> > > >&);
 #endif

include/igl/adjacency_list.h

@@ -31,7 +31,7 @@ namespace igl
   // See also: edges, cotmatrix, diag
   template <typename Index, typename IndexVector>
   IGL_INLINE void adjacency_list(
-    const Eigen::MatrixBase<Index> & F, 
+    const Eigen::MatrixBase<Index>  & F,
     std::vector<std::vector<IndexVector> >& A,
     bool sorted = false);
 

include/igl/cat.cpp

@@ -251,6 +251,65 @@ IGL_INLINE void igl::cat(const std::vector<std::vector< Mat > > & A, Mat & C)
   }
 }
 
+template <typename T, typename DerivedC>
+IGL_INLINE void igl::cat(const int dim, const std::vector<T> & A, Eigen::PlainObjectBase<DerivedC> & C)
+{
+  assert(dim == 1 || dim == 2);
+  using namespace Eigen;
+
+  const int num_mat = A.size();
+  if(num_mat == 0)
+  {
+    C.resize(0,0);
+    return;
+  }
+
+  if(dim == 1)
+  {
+    const int A_cols = A[0].cols();
+
+    int tot_rows = 0;
+    for(const auto & m : A)
+    {
+      tot_rows += m.rows();
+    }
+
+    C.resize(tot_rows, A_cols);
+
+    int cur_row = 0;
+    for(int i = 0; i < num_mat; i++)
+    {
+      assert(A_cols == A[i].cols());
+      C.block(cur_row,0,A[i].rows(),A_cols) = A[i];
+      cur_row += A[i].rows();
+    }
+  }
+  else if(dim == 2)
+  {
+    const int A_rows = A[0].rows();
+
+    int tot_cols = 0;
+    for(const auto & m : A)
+    {
+      tot_cols += m.cols();
+    }
+
+    C.resize(A_rows,tot_cols);
+
+    int cur_col = 0;
+    for(int i = 0; i < num_mat; i++)
+    {
+      assert(A_rows == A[i].rows());
+      C.block(0,cur_col,A_rows,A[i].cols()) = A[i];
+      cur_col += A[i].cols();
+    }
+  }
+  else
+  {
+    fprintf(stderr,"cat.h: Error: Unsupported dimension %d\n",dim);
+  }
+}
+
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
@@ -265,4 +324,13 @@ template Eigen::Matrix<double, -1, 1, 0, -1, 1> igl::cat<Eigen::Matrix<double, -
 template void igl::cat<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(int, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
 template void igl::cat<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1>&);
 template void igl::cat<Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1> >(int, Eigen::MatrixBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::Matrix<int, -1, 1, 0, -1, 1>&);
+template void igl::cat<Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, -1, -1, 0, -1, -1>, std::allocator<Eigen::Matrix<int, -1, -1, 0, -1, -1> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 1, 4, 1, 1, 4>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 1, 4, 1, 1, 4>, std::allocator<Eigen::Matrix<int, 1, 4, 1, 1, 4> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 1, 15, 1, 1, 15>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 1, 15, 1, 1, 15>, std::allocator<Eigen::Matrix<int, 1, 15, 1, 1, 15> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 1, 2, 1, 1, 2>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 1, 2, 1, 1, 2>, std::allocator<Eigen::Matrix<int, 1, 2, 1, 1, 2> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 1, 27, 1, 1, 27>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 1, 27, 1, 1, 27>, std::allocator<Eigen::Matrix<int, 1, 27, 1, 1, 27> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 1, 3, 1, 1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 1, 3, 1, 1, 3>, std::allocator<Eigen::Matrix<int, 1, 3, 1, 1, 3> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<int, 3, 1, 0, 3, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<int, 3, 1, 0, 3, 1>, std::allocator<Eigen::Matrix<int, 3, 1, 0, 3, 1> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<double, 1, 3, 1, 1, 3>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<double, 1, 3, 1, 1, 3>, std::allocator<Eigen::Matrix<double, 1, 3, 1, 1, 3> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
+template void igl::cat<Eigen::Matrix<double, 3, 1, 0, 3, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(int, std::vector<Eigen::Matrix<double, 3, 1, 0, 3, 1>, std::allocator<Eigen::Matrix<double, 3, 1, 0, 3, 1> > > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/cat.h

@@ -62,6 +62,16 @@ namespace igl
   //   C
   template <class Mat>
   IGL_INLINE void cat(const std::vector<std::vector< Mat > > & A, Mat & C);
+
+  // Concatenate a std::vector of matrices along the specified dimension
+  //
+  // Inputs:
+  //   dim  dimension along which to concatenate, 1 or 2
+  //   A  std::vector of eigen matrices. Must have identical # cols if dim == 1 or rows if dim == 2
+  // Outputs:
+  //   C  output matrix
+  template <typename T, typename DerivedC>
+  IGL_INLINE void cat(const int dim, const std::vector<T> & A, Eigen::PlainObjectBase<DerivedC> & C);
 }
 
 #ifndef IGL_STATIC_LIBRARY

include/igl/colon.cpp

@@ -59,8 +59,10 @@ template void igl::colon<int, double, double, double>(int, double, double, Eigen
 template void igl::colon<double, double, double>(double, double, Eigen::Matrix<double, -1, 1, 0, -1, 1>&);
 template void igl::colon<double, double, double, double>(double, double, double, Eigen::Matrix<double, -1, 1, 0, -1, 1>&);
 template void igl::colon<int, int, long>(int, int, Eigen::Matrix<long, -1, 1, 0, -1, 1>&);
+template void igl::colon<int, long, double>(int, long, Eigen::Matrix<double, -1, 1, 0, -1, 1>&);
 #ifdef WIN32
 template void igl::colon<int, long long,long>(int, long long, class Eigen::Matrix<long,-1,1,0,-1,1> &);
 template void igl::colon<int, __int64, __int64>(int, __int64, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> &);
+template void igl::colon<int,__int64,double>(int,__int64,class Eigen::Matrix<double,-1,1,0,-1,1> &);
 #endif
 #endif

include/igl/copyleft/marching_cubes.cpp

@@ -544,4 +544,6 @@ template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, 1, 0, -1,
 template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, 1, 0, -1, 1>,Eigen::Matrix<double, -1, -1, 0, -1, -1>,Eigen::Matrix<double, -1, -1, 0, -1, -1>,Eigen::Matrix<int, -1, -1, 0, -1, -1>,Eigen::Matrix<int, -1, 3, 0, -1, 3> >(const Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&,const Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&,const Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&,Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&,Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
 template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, 1, 0, -1, 1>,Eigen::Matrix<double, -1, -1, 0, -1, -1>,Eigen::Matrix<double, -1, -1, 0, -1, -1>,Eigen::Matrix<int, -1, 8, 0, -1, 8>,Eigen::Matrix<int, -1, 3, 0, -1, 3> >(const Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&,const Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&,const Eigen::MatrixBase<Eigen::Matrix<int, -1, 8, 0, -1, 8> >&,Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&,Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
 template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, 1, 0, -1, 1>,Eigen::Matrix<double, -1, 3, 0, -1, 3>,Eigen::Matrix<double, -1, 3, 0, -1, 3>,Eigen::Matrix<int, -1, 8, 0, -1, 8>,Eigen::Matrix<int, -1, 3, 0, -1, 3> >(const Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&,const Eigen::MatrixBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&,const Eigen::MatrixBase<Eigen::Matrix<int, -1, 8, 0, -1, 8> >&,Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&,Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> >&);
+template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
+template void igl::copyleft::marching_cubes<Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> > const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, unsigned int, unsigned int, unsigned int, double, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/copyleft/opengl2/render_to_tga.cpp

@@ -72,6 +72,7 @@ IGL_INLINE bool igl::opengl::render_to_tga(
     if(NULL==imgFile)
     {
       printf("IOError: %s could not be opened...\n",tga_file.c_str());
+      free(genericImage);
       return false;
     }
   }

include/igl/copyleft/opengl2/tga.cpp

@@ -373,6 +373,7 @@ igl::opengl::gliReadTGA(FILE *fp, char *name, int /*hflip*/, int vflip)
       sprintf(error, "TGA: error reading colormap (ftell == %ld)\n",
         ftell (fp));
       if (_verbose) printf("%s\n", error);
+      free(cmap);
       return NULL;
     }
 

include/igl/matrix_to_list.cpp

@@ -76,4 +76,6 @@ template void igl::matrix_to_list<Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::D
 template void igl::matrix_to_list<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, std::vector<Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar, std::allocator<Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar> >&);
 template void igl::matrix_to_list<Eigen::Matrix<int, 1, -1, 1, 1, -1> >(Eigen::DenseBase<Eigen::Matrix<int, 1, -1, 1, 1, -1> > const&, std::vector<Eigen::Matrix<int, 1, -1, 1, 1, -1>::Scalar, std::allocator<Eigen::Matrix<int, 1, -1, 1, 1, -1>::Scalar> >&);
 template std::vector<Eigen::Matrix<int, -1, -1, 0, -1, -1>::Scalar, std::allocator<Eigen::Matrix<int, -1, -1, 0, -1, -1>::Scalar> > igl::matrix_to_list<Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&);
+template void igl::matrix_to_list<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false> >(Eigen::DenseBase<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false> > const&, std::vector<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false>::Scalar, std::allocator<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false>::Scalar> >&);
+template std::vector<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false>::Scalar, std::allocator<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false>::Scalar> > igl::matrix_to_list<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false> >(Eigen::DenseBase<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1>, 1, -1, false> > const&);
 #endif

include/igl/octree.cpp

@@ -106,7 +106,7 @@ namespace igl {
         for(ChildrenType i = 0; i < 8; i++){
           children.emplace_back(neg_ones);
           point_indices.emplace_back(std::vector<IndexType>());
-          centers.emplace_back(translate_center(curr_center,h,i));
+          centers.emplace_back(translate_center(curr_center,h/2,i));
           widths.emplace_back(h);
         }
 

include/igl/opengl/ViewerCore.cpp

@@ -7,6 +7,7 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 
 #include "ViewerCore.h"
+#include "ViewerData.h"
 #include "gl.h"
 #include "../quat_to_mat.h"
 #include "../snap_to_fixed_up.h"

include/igl/opengl/ViewerCore.h

@@ -9,7 +9,6 @@
 #define IGL_OPENGL_VIEWERCORE_H
 
 #include <igl/opengl/MeshGL.h>
-#include <igl/opengl/ViewerData.h>
 
 #include <igl/igl_inline.h>
 #include <Eigen/Geometry>
@@ -20,9 +19,11 @@ namespace igl
 namespace opengl
 {
 
+// Forward declaration
+class ViewerData;
+
 // Basic class of the 3D mesh viewer
 // TODO: write documentation
-
 class ViewerCore
 {
 public:

include/igl/opengl/glfw/Viewer.cpp

@@ -124,7 +124,7 @@ namespace glfw
 {
 
   IGL_INLINE int Viewer::launch(bool resizable /*= true*/, bool fullscreen /*= false*/,
-    const std::string &name, int windowWidth /*= 1280*/, int windowHeight /*= 800*/)
+    const std::string &name, int windowWidth /*= 0*/, int windowHeight /*= 0*/)
   {
     // TODO return values are being ignored...
     launch_init(resizable,fullscreen,name,windowWidth,windowHeight);
@@ -153,9 +153,21 @@ namespace glfw
       GLFWmonitor *monitor = glfwGetPrimaryMonitor();
       const GLFWvidmode *mode = glfwGetVideoMode(monitor);
       window = glfwCreateWindow(mode->width,mode->height,name.c_str(),monitor,nullptr);
+      windowWidth = mode->width;
+      windowHeight = mode->height;
     }
     else
     {
+      // Set default windows width
+      if (windowWidth <= 0 & core_list.size() == 1 && core().viewport[2] > 0)
+        windowWidth = core().viewport[2];
+      else if (windowWidth <= 0)
+        windowWidth = 1280;
+      // Set default windows height
+      if (windowHeight <= 0 & core_list.size() == 1 && core().viewport[3] > 0)
+        windowHeight = core().viewport[3];
+      else if (windowHeight <= 0)
+        windowHeight = 800;
       window = glfwCreateWindow(windowWidth,windowHeight,name.c_str(),nullptr,nullptr);
     }
     if (!window)
@@ -988,6 +1000,20 @@ namespace glfw
     return data_list[index];
   }
 
+  IGL_INLINE const ViewerData& Viewer::data(int mesh_id /*= -1*/) const
+  {
+    assert(!data_list.empty() && "data_list should never be empty");
+    int index;
+    if (mesh_id == -1)
+      index = selected_data_index;
+    else
+      index = mesh_index(mesh_id);
+
+    assert((index >= 0 && index < data_list.size()) &&
+      "selected_data_index or mesh_id should be in bounds");
+    return data_list[index];
+  }
+
   IGL_INLINE int Viewer::append_mesh(bool visible /*= true*/)
   {
     assert(data_list.size() >= 1);
@@ -1043,6 +1069,18 @@ namespace glfw
     return core_list[core_index];
   }
 
+  IGL_INLINE const ViewerCore& Viewer::core(unsigned core_id /*= 0*/) const
+  {
+    assert(!core_list.empty() && "core_list should never be empty");
+    int core_index;
+    if (core_id == 0)
+      core_index = selected_core_index;
+    else
+      core_index = this->core_index(core_id);
+    assert((core_index >= 0 && core_index < core_list.size()) && "selected_core_index should be in bounds");
+    return core_list[core_index];
+  }
+
   IGL_INLINE bool Viewer::erase_core(const size_t index)
   {
     assert((index >= 0 && index < core_list.size()) && "index should be in bounds");

include/igl/opengl/glfw/Viewer.h

@@ -45,8 +45,8 @@ namespace glfw
     // UI Enumerations
     enum class MouseButton {Left, Middle, Right};
     enum class MouseMode { None, Rotation, Zoom, Pan, Translation} mouse_mode;
-    IGL_INLINE int launch(bool resizable = true, bool fullscreen = false, const std::string &name = "libigl viewer", int width = 1280, int height = 800);
-    IGL_INLINE int launch_init(bool resizable = true, bool fullscreen = false, const std::string &name = "libigl viewer", int width = 1280, int height = 800);
+    IGL_INLINE int launch(bool resizable = true, bool fullscreen = false, const std::string &name = "libigl viewer", int width = 0, int height = 0);
+    IGL_INLINE int launch_init(bool resizable = true, bool fullscreen = false, const std::string &name = "libigl viewer", int width = 0, int height = 0);
     IGL_INLINE bool launch_rendering(bool loop = true);
     IGL_INLINE void launch_shut();
     IGL_INLINE void init();
@@ -89,6 +89,7 @@ namespace glfw
     // Inputs:
     //   mesh_id  unique identifier associated to the desired mesh (current mesh if -1)
     IGL_INLINE ViewerData& data(int mesh_id = -1);
+    IGL_INLINE const ViewerData& data(int mesh_id = -1) const;
 
     // Append a new "slot" for a mesh (i.e., create empty entries at the end of
     // the data_list and opengl_state_list.
@@ -133,6 +134,7 @@ namespace glfw
     // Inputs:
     //   core_id  unique identifier corresponding to the desired viewport (current viewport if 0)
     IGL_INLINE ViewerCore& core(unsigned core_id = 0);
+    IGL_INLINE const ViewerCore& core(unsigned core_id = 0) const;
 
     // Append a new "slot" for a viewport (i.e., copy properties of the current viewport, only
     // changing the viewport size/position)

include/igl/ply.h

@@ -539,9 +539,12 @@ inline PlyFile *ply_open_for_writing(
 
   fp = fopen (name, "w");
   if (fp == NULL) {
+    free(name);
     return (NULL);
   }
 
+  free(name);
+
   /* create the actual PlyFile structure */
 
   plyfile = ply_write (fp, nelems, elem_names, file_type);
@@ -2253,6 +2256,7 @@ inline char **get_words(FILE *fp, int *nwords, char **orig_line)
   if (result == NULL) {
     *nwords = 0;
     *orig_line = NULL;
+    free(words);
     return (NULL);
   }
 
@@ -2316,7 +2320,15 @@ inline char **get_words(FILE *fp, int *nwords, char **orig_line)
     /* save pointer to beginning of word */
     if (num_words >= max_words) {
       max_words += 10;
-      words = (char **) realloc (words, sizeof (char *) * max_words);
+      char **temp = (char **) realloc (words, sizeof (char *) * max_words);
+
+      if(temp){
+          words = temp;
+      }
+      else{
+          free(words);
+          return NULL;
+      }
     }
     words[num_words++] = ptr;
 
@@ -2330,7 +2342,7 @@ inline char **get_words(FILE *fp, int *nwords, char **orig_line)
 
   /* return the list of words */
   *nwords = num_words;
-  *orig_line = str_copy;
+  *orig_line = str_copy; // ToDo: This looks like UB, returns pointer to local variable on stack.
   return (words);
 }
 

include/igl/point_mesh_squared_distance.cpp

@@ -51,6 +51,7 @@ IGL_INLINE void igl::point_mesh_squared_distance(
 template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 1, 0, -1, 1>, Eigen::Matrix<long, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, 3, 0, -1, 3> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<long, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> >&);
+template void igl::point_mesh_squared_distance<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 #ifdef WIN32
 template void igl::point_mesh_squared_distance<class Eigen::Matrix<double,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,-1,0,-1,-1>,class Eigen::Matrix<double,-1,1,0,-1,1>,class Eigen::Matrix<__int64,-1,1,0,-1,1>,class Eigen::Matrix<double,-1,3,0,-1,3> >(class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,-1,0,-1,-1> > const &,class Eigen::Matrix<int,-1,-1,0,-1,-1> const &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<__int64,-1,1,0,-1,1> > &,class Eigen::PlainObjectBase<class Eigen::Matrix<double,-1,3,0,-1,3> > &);
 #endif

include/igl/polar_svd3x3.cpp

@@ -75,25 +75,6 @@ IGL_INLINE void igl::polar_svd3x3_avx(const Eigen::Matrix<T, 3*8, 3>& A, Eigen::
     R.block(3*k, 0, 3, 3) = U.block(3*k, 0, 3, 3) * Vt.block(3*k, 0, 3, 3).transpose();
   }
 
-  // test:
-  for (int k=0; k<8; k++)
-  {
-    Eigen::Matrix3f Apart = A.block(3*k, 0, 3, 3);
-    Eigen::Matrix3f Rpart;
-    polar_svd3x3(Apart, Rpart);
-
-    Eigen::Matrix3f Rpart_SSE = R.block(3*k, 0, 3, 3);
-    Eigen::Matrix3f diff = Rpart - Rpart_SSE;
-    float diffNorm = diff.norm();
-    if (std::abs(diffNorm) > 0.001)
-    {
-      printf("Huh: diffNorm = %15f (k = %i)\n", diffNorm, k);
-    }
-
-    // Unused
-    //int hu = 1;
-  }
-  // eof test
 }
 #endif
 

include/igl/predicates/predicates.cpp

@@ -0,0 +1,170 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2019 Qingnan Zhou <qnzhou@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can
+// obtain one at http://mozilla.org/MPL/2.0/.
+#include <igl/predicates/predicates.h>
+#include <predicates.h>
+#include <type_traits>
+
+namespace igl {
+namespace predicates {
+
+using REAL = IGL_PREDICATES_REAL;
+
+#ifdef LIBIGL_PREDICATES_USE_FLOAT
+#define IGL_ASSERT_SCALAR(Vector)                        \
+  static_assert(                                         \
+    std::is_same<typename Vector::Scalar, float>::value, \
+    "Shewchuk's exact predicates only support float")
+#else
+#define IGL_ASSERT_SCALAR(Vector)                           \
+  static_assert(                                            \
+    std::is_same<typename Vector::Scalar, double>::value || \
+    std::is_same<typename Vector::Scalar, float>::value,    \
+    "Shewchuk's exact predicates only support float and double")
+#endif
+
+IGL_INLINE void exactinit() {
+  static bool initialized = false;
+  if (! initialized) {
+    ::exactinit();
+    initialized = true;
+  }
+}
+
+template<typename Vector2D>
+IGL_INLINE Orientation orient2d(
+    const Eigen::MatrixBase<Vector2D>& pa,
+    const Eigen::MatrixBase<Vector2D>& pb,
+    const Eigen::MatrixBase<Vector2D>& pc)
+{
+  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2D, 2);
+  IGL_ASSERT_SCALAR(Vector2D);
+
+  using Point = Eigen::Matrix<REAL, 2, 1>;
+  Point a{pa[0], pa[1]};
+  Point b{pb[0], pb[1]};
+  Point c{pc[0], pc[1]};
+
+  const auto r = ::orient2d(a.data(), b.data(), c.data());
+
+  if (r > 0) return Orientation::POSITIVE;
+  else if (r < 0) return Orientation::NEGATIVE;
+  else return Orientation::COLLINEAR;
+}
+
+template<typename Vector3D>
+IGL_INLINE Orientation orient3d(
+    const Eigen::MatrixBase<Vector3D>& pa,
+    const Eigen::MatrixBase<Vector3D>& pb,
+    const Eigen::MatrixBase<Vector3D>& pc,
+    const Eigen::MatrixBase<Vector3D>& pd)
+{
+  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector3D, 3);
+  IGL_ASSERT_SCALAR(Vector3D);
+
+  using Point = Eigen::Matrix<REAL, 3, 1>;
+  Point a{pa[0], pa[1], pa[2]};
+  Point b{pb[0], pb[1], pb[2]};
+  Point c{pc[0], pc[1], pc[2]};
+  Point d{pd[0], pd[1], pd[2]};
+
+  const auto r = ::orient3d(a.data(), b.data(), c.data(), d.data());
+
+  if (r > 0) return Orientation::POSITIVE;
+  else if (r < 0) return Orientation::NEGATIVE;
+  else return Orientation::COPLANAR;
+}
+
+template<typename Vector2D>
+IGL_INLINE Orientation incircle(
+    const Eigen::MatrixBase<Vector2D>& pa,
+    const Eigen::MatrixBase<Vector2D>& pb,
+    const Eigen::MatrixBase<Vector2D>& pc,
+    const Eigen::MatrixBase<Vector2D>& pd)
+{
+  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector2D, 2);
+  IGL_ASSERT_SCALAR(Vector2D);
+
+  using Point = Eigen::Matrix<REAL, 2, 1>;
+  Point a{pa[0], pa[1]};
+  Point b{pb[0], pb[1]};
+  Point c{pc[0], pc[1]};
+  Point d{pd[0], pd[1]};
+
+  const auto r = ::incircle(a.data(), b.data(), c.data(), d.data());
+
+  if (r > 0) return Orientation::INSIDE;
+  else if (r < 0) return Orientation::OUTSIDE;
+  else return Orientation::COCIRCULAR;
+}
+
+template<typename Vector3D>
+IGL_INLINE Orientation insphere(
+    const Eigen::MatrixBase<Vector3D>& pa,
+    const Eigen::MatrixBase<Vector3D>& pb,
+    const Eigen::MatrixBase<Vector3D>& pc,
+    const Eigen::MatrixBase<Vector3D>& pd,
+    const Eigen::MatrixBase<Vector3D>& pe)
+{
+  EIGEN_STATIC_ASSERT_VECTOR_SPECIFIC_SIZE(Vector3D, 3);
+  IGL_ASSERT_SCALAR(Vector3D);
+
+  using Point = Eigen::Matrix<REAL, 3, 1>;
+  Point a{pa[0], pa[1], pa[2]};
+  Point b{pb[0], pb[1], pb[2]};
+  Point c{pc[0], pc[1], pc[2]};
+  Point d{pd[0], pd[1], pd[2]};
+  Point e{pe[0], pe[1], pe[2]};
+
+  const auto r = ::insphere(a.data(), b.data(), c.data(), d.data(), e.data());
+
+  if (r > 0) return Orientation::INSIDE;
+  else if (r < 0) return Orientation::OUTSIDE;
+  else return Orientation::COSPHERICAL;
+}
+
+}
+}
+
+#undef IGL_ASSERT_SCALAR
+
+#ifdef IGL_STATIC_LIBRARY
+// Explicit template instantiation
+
+#define IGL_ORIENT2D(Vector) template igl::predicates::Orientation igl::predicates::orient2d<Vector>(const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&)
+#define IGL_INCIRCLE(Vector) template igl::predicates::Orientation igl::predicates::incircle<Vector>(const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&)
+#define IGL_ORIENT3D(Vector) template igl::predicates::Orientation igl::predicates::orient3d<Vector>(const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&)
+#define IGL_INSPHERE(Vector) template igl::predicates::Orientation igl::predicates::insphere<Vector>(const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&, const Eigen::MatrixBase<Vector>&)
+
+#define IGL_MATRIX(T, R, C) Eigen::Matrix<T, R, C>
+IGL_ORIENT2D(IGL_MATRIX(float, 1, 2));
+IGL_INCIRCLE(IGL_MATRIX(float, 1, 2));
+IGL_ORIENT2D(IGL_MATRIX(float, 2, 1));
+IGL_INCIRCLE(IGL_MATRIX(float, 2, 1));
+IGL_ORIENT3D(IGL_MATRIX(float, 1, 3));
+IGL_INSPHERE(IGL_MATRIX(float, 1, 3));
+IGL_ORIENT3D(IGL_MATRIX(float, 3, 1));
+IGL_INSPHERE(IGL_MATRIX(float, 3, 1));
+
+#ifndef LIBIGL_PREDICATES_USE_FLOAT
+IGL_ORIENT2D(IGL_MATRIX(double, 1, 2));
+IGL_INCIRCLE(IGL_MATRIX(double, 1, 2));
+IGL_ORIENT2D(IGL_MATRIX(double, 2, 1));
+IGL_INCIRCLE(IGL_MATRIX(double, 2, 1));
+IGL_ORIENT3D(IGL_MATRIX(double, 1, 3));
+IGL_INSPHERE(IGL_MATRIX(double, 1, 3));
+IGL_ORIENT3D(IGL_MATRIX(double, 3, 1));
+IGL_INSPHERE(IGL_MATRIX(double, 3, 1));
+#endif
+#undef IGL_MATRIX
+
+#undef IGL_ORIENT2D
+#undef IGL_ORIENT3D
+#undef IGL_INCIRCLE
+#undef IGL_INSPHERE
+
+#endif

include/igl/predicates/predicates.h

@@ -0,0 +1,99 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2019 Qingnan Zhou <qnzhou@gmail.com>
+//
+// This Source Code Form is subject to the terms of the Mozilla Public License 
+// v. 2.0. If a copy of the MPL was not distributed with this file, You can 
+// obtain one at http://mozilla.org/MPL/2.0/.
+#pragma once
+#ifndef IGL_PREDICATES_PREDICATES_H
+#define IGL_PREDICATES_PREDICATES_H
+
+#include <igl/igl_inline.h>
+#include <Eigen/Core>
+
+namespace igl {
+  namespace predicates {
+    enum class Orientation {
+      POSITIVE=1, INSIDE=1,
+      NEGATIVE=-1, OUTSIDE=-1,
+      COLLINEAR=0, COPLANAR=0, COCIRCULAR=0, COSPHERICAL=0, DEGENERATE=0
+    };
+
+    // Initialize internal variable used by predciates.  Must be called before
+    // using exact predicates.
+    IGL_INLINE void exactinit();
+
+    // Compute the orientation of the triangle formed by pa, pb, pc.
+    //
+    // Input:
+    //   pa, pb, pc  2D points.
+    //
+    // Output:
+    //   Return POSITIVE if pa, pb, pc are counterclockwise oriented.
+    //          NEGATIVE if they are clockwise oriented.
+    //          COLLINEAR if they are collinear.
+    template<typename Vector2D>
+    IGL_INLINE Orientation orient2d(
+        const Eigen::MatrixBase<Vector2D>& pa,
+        const Eigen::MatrixBase<Vector2D>& pb,
+        const Eigen::MatrixBase<Vector2D>& pc);
+
+    // Compute the orientation of the tetrahedron formed by pa, pb, pc, pd.
+    //
+    // Input:
+    //   pa, pb, pc, pd  3D points.
+    //
+    // Output:
+    //   Return POSITIVE if pd is "below" the oriented plane formed by pa, pb and pc.
+    //          NEGATIVE if pd is "above" the plane.
+    //          COPLANAR if pd is on the plane.
+    template<typename Vector3D>
+    IGL_INLINE Orientation orient3d(
+        const Eigen::MatrixBase<Vector3D>& pa,
+        const Eigen::MatrixBase<Vector3D>& pb,
+        const Eigen::MatrixBase<Vector3D>& pc,
+        const Eigen::MatrixBase<Vector3D>& pd);
+
+    // Decide whether a point is inside/outside/on a circle.
+    //
+    // Input:
+    //   pa, pb, pc  2D points that defines an oriented circle.
+    //   pd          2D query point.
+    //
+    // Output:
+    //   Return INSIDE if pd is inside of the circle defined by pa, pb and pc.
+    //          OUSIDE if pd is outside of the circle.
+    //          COCIRCULAR pd is exactly on the circle.
+    template<typename Vector2D>
+    IGL_INLINE Orientation incircle(
+        const Eigen::MatrixBase<Vector2D>& pa,
+        const Eigen::MatrixBase<Vector2D>& pb,
+        const Eigen::MatrixBase<Vector2D>& pc,
+        const Eigen::MatrixBase<Vector2D>& pd);
+
+    // Decide whether a point is inside/outside/on a sphere.
+    //
+    // Input:
+    //   pa, pb, pc, pd  3D points that defines an oriented sphere.
+    //   pe              3D query point.
+    //
+    // Output:
+    //   Return INSIDE if pe is inside of the sphere defined by pa, pb, pc and pd.
+    //          OUSIDE if pe is outside of the sphere.
+    //          COSPHERICAL pd is exactly on the sphere.
+    template<typename Vector3D>
+    IGL_INLINE Orientation insphere(
+        const Eigen::MatrixBase<Vector3D>& pa,
+        const Eigen::MatrixBase<Vector3D>& pb,
+        const Eigen::MatrixBase<Vector3D>& pc,
+        const Eigen::MatrixBase<Vector3D>& pd,
+        const Eigen::MatrixBase<Vector3D>& pe);
+  }
+}
+
+#ifndef IGL_STATIC_LIBRARY
+#  include "predicates.cpp"
+#endif
+
+#endif

include/igl/principal_curvature.cpp

@@ -448,7 +448,7 @@ IGL_INLINE void CurvatureCalculator::getKRing(const int start, const double r, s
   int bufsize=vertices.rows();
   vv.reserve(bufsize);
   std::list<std::pair<int,int> > queue;
-  bool* visited = (bool*)calloc(bufsize,sizeof(bool));
+  std::vector<bool> visited(bufsize, false);
   queue.push_back(std::pair<int,int>(start,0));
   visited[start]=true;
   while (!queue.empty())
@@ -470,8 +470,6 @@ IGL_INLINE void CurvatureCalculator::getKRing(const int start, const double r, s
       }
     }
   }
-  free(visited);
-  return;
 }
 
 
@@ -479,16 +477,16 @@ IGL_INLINE void CurvatureCalculator::getSphere(const int start, const double r,
 {
   int bufsize=vertices.rows();
   vv.reserve(bufsize);
-  std::list<int>* queue= new std::list<int>();
-  bool* visited = (bool*)calloc(bufsize,sizeof(bool));
-  queue->push_back(start);
+  std::list<int> queue;
+  std::vector<bool> visited(bufsize, false);
+  queue.push_back(start);
   visited[start]=true;
   Eigen::Vector3d me=vertices.row(start);
-  std::priority_queue<std::pair<int, double>, std::vector<std::pair<int, double> >, comparer >* extra_candidates= new  std::priority_queue<std::pair<int, double>, std::vector<std::pair<int, double> >, comparer >();
-  while (!queue->empty())
+  std::priority_queue<std::pair<int, double>, std::vector<std::pair<int, double> >, comparer > extra_candidates;
+  while (!queue.empty())
   {
-    int toVisit=queue->front();
-    queue->pop_front();
+    int toVisit=queue.front();
+    queue.pop_front();
     vv.push_back(toVisit);
     for (unsigned int i=0; i<vertex_to_vertices[toVisit].size(); ++i)
     {
@@ -498,17 +496,17 @@ IGL_INLINE void CurvatureCalculator::getSphere(const int start, const double r,
         Eigen::Vector3d neigh=vertices.row(neighbor);
         double distance=(me-neigh).norm();
         if (distance<r)
-          queue->push_back(neighbor);
+          queue.push_back(neighbor);
         else if ((int)vv.size()<min)
-          extra_candidates->push(std::pair<int,double>(neighbor,distance));
+          extra_candidates.push(std::pair<int,double>(neighbor,distance));
         visited[neighbor]=true;
       }
     }
   }
-  while (!extra_candidates->empty() && (int)vv.size()<min)
+  while (!extra_candidates.empty() && (int)vv.size()<min)
   {
-    std::pair<int, double> cand=extra_candidates->top();
-    extra_candidates->pop();
+    std::pair<int, double> cand=extra_candidates.top();
+    extra_candidates.pop();
     vv.push_back(cand.first);
     for (unsigned int i=0; i<vertex_to_vertices[cand.first].size(); ++i)
     {
@@ -517,14 +515,11 @@ IGL_INLINE void CurvatureCalculator::getSphere(const int start, const double r,
       {
         Eigen::Vector3d neigh=vertices.row(neighbor);
         double distance=(me-neigh).norm();
-        extra_candidates->push(std::pair<int,double>(neighbor,distance));
+        extra_candidates.push(std::pair<int,double>(neighbor,distance));
         visited[neighbor]=true;
       }
     }
   }
-  free(extra_candidates);
-  free(queue);
-  free(visited);
 }
 
 IGL_INLINE Eigen::Vector3d CurvatureCalculator::project(const Eigen::Vector3d& v, const Eigen::Vector3d& vp, const Eigen::Vector3d& ppn)

include/igl/readDMAT.cpp

@@ -195,8 +195,8 @@ IGL_INLINE bool igl::readDMAT(
     assert(W.size() == 0);
     // Resize for output
     W.resize(num_rows,typename std::vector<Scalar>(num_cols));
-    double * Wraw = new double[num_rows*num_cols];
-    fread(Wraw, sizeof(double), num_cols*num_rows, fp);
+    std::unique_ptr<double[]> Wraw(new double[num_rows*num_cols]);
+    fread(Wraw.get(), sizeof(double), num_cols*num_rows, fp);
     // Loop over columns slowly
     for(int j = 0;j < num_cols;j++)
     {

include/igl/readMSH.cpp

@@ -328,7 +328,13 @@ IGL_INLINE bool igl::readMSH(
     for (size_t i=0; i<num_int_tags; i++)
         fin >> int_tags[i];
 
-    if (num_string_tags <= 0 || num_int_tags <= 2) { assert(false && "invalid format"); return; }
+    if (num_string_tags <= 0 || num_int_tags <= 2) {
+        delete[] str_tags;
+        delete[] real_tags;
+        delete[] int_tags;
+        assert(false && "invalid format");
+        return;
+    }
     std::string fieldname = str_tags[0];
     int num_components = int_tags[1];
     int num_entries = int_tags[2];
@@ -397,7 +403,13 @@ IGL_INLINE bool igl::readMSH(
     for (size_t i=0; i<num_int_tags; i++)
         fin >> int_tags[i];
 
-    if (num_string_tags <= 0 || num_int_tags <= 2) { assert(false && "invalid format"); return; }
+    if (num_string_tags <= 0 || num_int_tags <= 2) {
+        delete[] str_tags;
+        delete[] real_tags;
+        delete[] int_tags;
+        assert(false && "invalid format");
+        return;
+    }
     std::string fieldname = str_tags[0];
     int num_components = int_tags[1];
     int num_entries = int_tags[2];

include/igl/remesh_along_isoline.h

@@ -28,7 +28,7 @@ namespace igl
   //  U  #U by dim list of mesh vertex positions #U>=#V
   //  G  #G by 3 list of mesh triangle indices into U, #G>=#F
   //  SU  #U list of scalar field values over new mesh
-  //  J  #G list of indices into G revealing birth triangles
+  //  J  #G list of indices into F revealing birth triangles
   //  BC  #U by #V sparse matrix of barycentric coordinates so that U = BC*V
   //  L  #G list of bools whether scalar field in triangle below or above val
   template <

include/igl/remove_unreferenced.cpp

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

include/igl/slice.cpp

@@ -363,6 +363,9 @@ template void igl::slice<unsigned int, unsigned int>(Eigen::SparseMatrix<unsigne
 template void igl::slice<Eigen::SparseMatrix<double, 0, int>, Eigen::Array<int, -1, 1, 0, -1, 1>, Eigen::SparseMatrix<double, 0, int> >(Eigen::SparseMatrix<double, 0, int> const&, Eigen::DenseBase<Eigen::Array<int, -1, 1, 0, -1, 1> > const&, int, Eigen::SparseMatrix<double, 0, int>&);
 template void igl::slice<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1> const, -1, 1, true>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::Matrix<double, -1, -1, 0, -1, -1> const&, Eigen::DenseBase<Eigen::Block<Eigen::Matrix<int, -1, -1, 0, -1, -1> const, -1, 1, true> > const&, int, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
 template void igl::slice<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Block<Eigen::Matrix<int, -1, 1, 0, -1, 1> const, -1, 1, true>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::Matrix<double, -1, -1, 0, -1, -1> const&, Eigen::DenseBase<Eigen::Block<Eigen::Matrix<int, -1, 1, 0, -1, 1> const, -1, 1, true> > const&, int, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
+template void igl::slice<Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > >(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::DenseBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, int, 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::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::Matrix<int, -1, -1, 0, -1, -1> const&, Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, int, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
+template void igl::slice<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1>, -1, 1, true>, Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::Matrix<double, -1, -1, 0, -1, -1> const&, Eigen::DenseBase<Eigen::Block<Eigen::Matrix<double, -1, -1, 0, -1, -1>, -1, 1, true> > const&, int, Eigen::Matrix<double, -1, -1, 0, -1, -1>&);
 #ifdef WIN32
 template void igl::slice<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>>>(class Eigen::Matrix<__int64, -1, 1, 0, -1, 1> const &, class Eigen::DenseBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> const &, int, class Eigen::PlainObjectBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> &);
 template void igl::slice<class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>>, class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>>>(class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>> const &, class Eigen::DenseBase<class Eigen::Matrix<__int64, -1, 1, 0, -1, 1>> const &, int, class Eigen::PlainObjectBase<class Eigen::Matrix<int, -1, -1, 0, -1, -1>> &);

include/igl/slice_mask.cpp

@@ -165,4 +165,5 @@ template void igl::slice_mask<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::
 template void igl::slice_mask<Eigen::Matrix<double, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::Array<bool, -1, 1, 0, -1, 1> const&, int, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> >&);
 template void igl::slice_mask<Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, Eigen::Array<bool, -1, 1, 0, -1, 1> const&, Eigen::Array<bool, -1, 1, 0, -1, 1> const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 template void igl::slice_mask<Eigen::Matrix<int, -1, 1, 0, -1, 1> >(Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::Array<bool, -1, 1, 0, -1, 1> const&, int, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+template void igl::slice_mask<Eigen::Matrix<int, -1, 1, 0, -1, 1>, Eigen::Matrix<int, -1, -1, 0, -1, -1> >(Eigen::DenseBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> > const&, Eigen::Array<bool, -1, 1, 0, -1, 1> const&, int, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> >&);
 #endif

include/igl/snap_points.cpp

@@ -85,5 +85,6 @@ IGL_INLINE void igl::snap_points(
 // Explicit template instantiation
 template void igl::snap_points<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<double, -1, 1, 0, -1, 1> >(Eigen::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
 template void igl::snap_points<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::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
+template void igl::snap_points<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::PlainObjectBase<Eigen::Matrix<double, 1, 3, 1, 1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1> >&);
 #endif
 

include/igl/writeOBJ.cpp

@@ -118,6 +118,39 @@ IGL_INLINE bool igl::writeOBJ(
   return true;
 }
 
+template <typename DerivedV, typename T>
+IGL_INLINE bool igl::writeOBJ(
+  const std::string &str,
+  const Eigen::MatrixBase<DerivedV>& V,
+  const std::vector<std::vector<T> >& F)
+{
+  using namespace std;
+  using namespace Eigen;
+  assert(V.cols() == 3 && "V should have 3 columns");
+  ofstream s(str);
+  if(!s.is_open())
+  {
+    fprintf(stderr,"IOError: writeOBJ() could not open %s\n",str.c_str());
+    return false;
+  }
+  s<<V.format(IOFormat(FullPrecision,DontAlignCols," ","\n","v ","","","\n"));
+  
+  for(const auto& face : F)
+  {
+    int face_size = face.size();
+    assert(face_size != 0);
+
+    s << (face_size == 2 ? "l" : "f");
+
+    for(const auto& vi : face)
+    {
+      s<<" "<<vi; 
+    }
+    s<<"\n";
+  }
+  return true;
+}
+
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh

include/igl/writeOBJ.h

@@ -13,6 +13,7 @@
 
 #include <Eigen/Core>
 #include <string>
+#include <vector>
 
 namespace igl 
 {
@@ -50,6 +51,18 @@ namespace igl
     const Eigen::MatrixBase<DerivedV>& V,
     const Eigen::MatrixBase<DerivedF>& F);
 
+  // Write a mesh of mixed tris and quads to an ascii obj file
+  // Inputs:
+  //   str  path to outputfile
+  //   V  #V by 3 mesh vertex positions
+  //   F  #F std::vector of std::vector<Index> of size 3 or 4 mesh indices into V
+  // Returns true on success, false on error
+  template <typename DerivedV, typename T>
+  IGL_INLINE bool writeOBJ(
+    const std::string &str,
+    const Eigen::MatrixBase<DerivedV>& V,
+    const std::vector<std::vector<T> >& F);
+
 }
 
 #ifndef IGL_STATIC_LIBRARY

include/igl/xml/serialize_xml.cpp

@@ -124,7 +124,7 @@ namespace igl
       {
         std::cerr << "File not found!" << std::endl;
         doc->PrintError();
-        doc = NULL;
+        delete doc;
       }
       else
       {

python/modules/py_igl_opengl_glfw.cpp

@@ -353,14 +353,6 @@ py::class_<igl::opengl::ViewerCore> viewercore_class(me, "ViewerCore");
 // UI Enumerations
     py::class_<igl::opengl::glfw::Viewer> viewer_class(me, "Viewer");
 
-//    #ifdef IGL_VIEWER_WITH_NANOGUI
-//    py::object fh = (py::object) py::module::import("nanogui").attr("FormHelper");
-//    py::class_<nanogui::FormHelper> form_helper_class(me, "FormHelper", fh);
-
-//    py::object screen = (py::object) py::module::import("nanogui").attr("Screen");
-//    py::class_<nanogui::Screen> screen_class(me, "Screen", screen);
-//    #endif
-
     py::enum_<igl::opengl::glfw::Viewer::MouseButton>(viewer_class, "MouseButton")
         .value("Left", igl::opengl::glfw::Viewer::MouseButton::Left)
         .value("Middle", igl::opengl::glfw::Viewer::MouseButton::Middle)
@@ -378,16 +370,12 @@ py::class_<igl::opengl::ViewerCore> viewercore_class(me, "ViewerCore");
     //   viewer.data() = data;
     // })
 
-    .def("data", &igl::opengl::glfw::Viewer::data,pybind11::return_value_policy::reference)
+    .def("data", (igl::opengl::ViewerData & (igl::opengl::glfw::Viewer::*)(int)) &igl::opengl::glfw::Viewer::data,pybind11::return_value_policy::reference)
+    // .def("data", (const igl::opengl::ViewerData & (igl::opengl::glfw::Viewer::*)(int) const) &igl::opengl::glfw::Viewer::data,pybind11::return_value_policy::reference)
 
     //.def_readwrite("core", &igl::opengl::glfw::Viewer::core)
     //.def_readwrite("opengl", &igl::opengl::glfw::Viewer::opengl)
 
-    #ifdef IGL_VIEWER_WITH_NANOGUI
-    .def_readwrite("ngui", &igl::opengl::glfw::Viewer::ngui)
-    .def_readwrite("screen", &igl::opengl::glfw::Viewer::screen)
-    #endif
-
     .def("launch", &igl::opengl::glfw::Viewer::launch, py::arg("resizable") = true,
       py::arg("fullscreen") = false, py::arg("name") = "libigl viewer",
       py::arg("windowWidth") = 1280, py::arg("windowHeight") = 800)

python/py_doc.cpp

@@ -615,6 +615,26 @@ const char *__doc_igl_exact_geodesic = R"igl_Qu8mg5v7(
     // Note:  
     //      Specifying a face as target/source means its center.  
     //)igl_Qu8mg5v7";
+const char *__doc_igl_heat_geodesics_precompute = R"igl_Qu8mg5v7(
+    // Precompute factorized solvers for computing a fast approximation of
+    // geodesic distances on a mesh (V,F). [Crane et al. 2013]
+    //
+    // Inputs:
+    //   V  #V by dim list of mesh vertex positions
+    //   F  #F by 3 list of mesh face indices into V
+    // Outputs:
+    //   data  precomputation data (see heat_geodesics_solve)
+    //)igl_Qu8mg5v7";
+const char *__doc_igl_heat_geodesics_solve = R"igl_Qu8mg5v7(
+    // Compute fast approximate geodesic distances using precomputed data from a
+    // set of selected source vertices (gamma)
+    //
+    // Inputs:
+    //   data  precomputation data (see heat_geodesics_precompute)
+    //   gamma  #gamma list of indices into V of source vertices
+    // Outputs:
+    //   D  #V list of distances to gamma
+    //)igl_Qu8mg5v7";
 const char *__doc_igl_find_cross_field_singularities = R"igl_Qu8mg5v7(// Inputs:
   //   V                #V by 3 eigen Matrix of mesh vertex 3D positions
   //   F                #F by 3 eigen Matrix of face (quad) indices

python/py_doc.h

@@ -49,6 +49,8 @@ extern const char *__doc_igl_embree_ambient_occlusion;
 extern const char *__doc_igl_embree_line_mesh_intersection;
 extern const char *__doc_igl_embree_reorient_facets_raycast;
 extern const char *__doc_igl_exact_geodesic;
+extern const char *__doc_igl_heat_geodesics_precompute;
+extern const char *__doc_igl_heat_geodesics_solve;
 extern const char *__doc_igl_find_cross_field_singularities;
 extern const char *__doc_igl_fit_rotations;
 extern const char *__doc_igl_fit_rotations_planar;

python/py_igl.cpp

@@ -46,6 +46,7 @@
 #include <igl/edge_topology.h>
 #include <igl/eigs.h>
 #include <igl/exact_geodesic.h>
+#include <igl/heat_geodesics.h>
 #include <igl/find_cross_field_singularities.h>
 #include <igl/fit_rotations.h>
 #include <igl/floor.h>
@@ -145,6 +146,7 @@ void python_export_igl(py::module &m)
 #include "py_igl/py_edge_topology.cpp"
 #include "py_igl/py_eigs.cpp"
 #include "py_igl/py_exact_geodesic.cpp"
+#include "py_igl/py_heat_geodesics.cpp"
 #include "py_igl/py_find_cross_field_singularities.cpp"
 #include "py_igl/py_fit_rotations.cpp"
 #include "py_igl/py_floor.cpp"

python/py_igl/py_heat_geodesics.cpp

@@ -0,0 +1,45 @@
+// This file is part of libigl, a simple c++ geometry processing library.
+//
+// Copyright (C) 2018 Amrollah Seifoddini <a.seif67@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/.
+
+// Wrap the data class, no properties are exposed since it is not necessary
+py::class_<igl::HeatGeodesicsData<double> > heat_geodesics_data(m, "heat_geodesics_data");
+
+heat_geodesics_data.def(py::init<>());
+
+m.def("heat_geodesics_precompute", []
+(
+    const Eigen::MatrixXd &V,
+    const Eigen::MatrixXi &F,
+    igl::HeatGeodesicsData<double> &data
+)
+{
+  return igl::heat_geodesics_precompute(V, F, data);
+}, __doc_igl_heat_geodesics_precompute,
+py::arg("V"), py::arg("F"), py::arg("data"));
+
+
+m.def("heat_geodesics_solve", []
+(
+    const igl::HeatGeodesicsData<double> &data,
+    const Eigen::MatrixXi &gamma,
+    Eigen::MatrixXd &D
+)
+{
+assert_is_VectorX("D", D);
+assert_is_VectorX("gamma", gamma);
+Eigen::VectorXd vD;
+igl::heat_geodesics_solve(data, Eigen::VectorXi(gamma), vD);
+D.resize(vD.size(), 1);
+for (int i = 0; i < vD.size(); i++)
+{
+    D(i, 0) = vD[i];
+}
+return true;
+}, __doc_igl_heat_geodesics_solve,
+py::arg("data"), py::arg("gamma"), py::arg("D"));
+

tests/CMakeLists.txt

@@ -76,6 +76,14 @@ if(LIBIGL_WITH_EMBREE)
   target_link_libraries(libigl_tests PUBLIC igl::embree)
 endif()
 
+if(LIBIGL_WITH_PREDICATES)
+  file(GLOB TEST_SRC_FILES ./include/igl/predicates/*.cpp)
+  file(GLOB TEST_INC_FILES ./include/igl/predicates/*.h ./include/igl/predicates/*.inl)
+  target_sources(libigl_tests PRIVATE ${TEST_SRC_FILES} ${TEST_INC_FILES})
+
+  target_link_libraries(libigl_tests PUBLIC igl::predicates igl::triangle)
+endif()
+
 file(GLOB TEST_SRC_FILES ./include/igl/*.cpp)
 file(GLOB TEST_INC_FILES ./include/igl/*.h ./include/igl/*.inl)
 target_sources(libigl_tests PRIVATE ${TEST_SRC_FILES} ${TEST_INC_FILES})

tests/include/igl/cat.cpp

@@ -0,0 +1,43 @@
+#include <test_common.h>
+#include <igl/cat.h>
+
+#include <vector>
+
+TEST_CASE("cat: rows", "[igl]")
+{
+  std::vector<Eigen::RowVector3i> rows = {
+     Eigen::RowVector3i(1, 2, 3),
+     Eigen::RowVector3i(4, 5, 6),
+     Eigen::RowVector3i(7, 8, 9)
+  };
+
+  Eigen::MatrixXi actual;
+  igl::cat(1,rows,actual);
+
+  Eigen::Matrix3i expected;
+  expected << 1, 2, 3,
+              4, 5, 6,
+              7, 8, 9;  
+
+  test_common::assert_eq(actual,expected);
+}
+
+TEST_CASE("cat: cols", "[igl]")
+{
+  std::vector<Eigen::Vector3i> cols = {
+     Eigen::Vector3i(1, 4, 7),
+     Eigen::Vector3i(2, 5, 8),
+     Eigen::Vector3i(3, 6, 9)
+  };
+
+  Eigen::MatrixXi actual;
+  igl::cat(2,cols,actual);
+
+  Eigen::Matrix3i expected;
+  expected << 1, 2, 3,
+              4, 5, 6,
+              7, 8, 9;  
+
+  test_common::assert_eq(actual,expected);
+}
+

tests/include/igl/predicates/predicates.cpp

@@ -0,0 +1,77 @@
+#include <test_common.h>
+#include <igl/predicates/predicates.h>
+#include <limits>
+#include <igl/triangle/triangulate.h>
+
+TEST_CASE("predicates", "[igl][predicates]") {
+    using namespace igl::predicates;
+    using Scalar = double;
+    igl::predicates::exactinit();
+
+    SECTION("2D") {
+        using Point = Eigen::Matrix<Scalar, 2, 1>;
+        Point a(2,1),b(2,1),c(2,1),d(2,1),e(2,1),f(2,1);
+        a << 0.0, 0.0;
+        b << 1.0, 0.0;
+        c << 1.0, 1.0;
+        d << 2.0, 0.0;
+        e << 0.0, 1.0;
+        f << 0.5, 0.5;
+
+        REQUIRE(orient2d(a, b, c) == Orientation::POSITIVE);
+        REQUIRE(orient2d(a, c, b) == Orientation::NEGATIVE);
+        REQUIRE(orient2d(a, b, b) == Orientation::COLLINEAR);
+        REQUIRE(orient2d(a, a, a) == Orientation::COLLINEAR);
+        REQUIRE(orient2d(a, b, d) == Orientation::COLLINEAR);
+        REQUIRE(orient2d(a, f, c) == Orientation::COLLINEAR);
+
+        REQUIRE(incircle(a,b,c,e) == Orientation::COCIRCULAR);
+        REQUIRE(incircle(a,b,c,a) == Orientation::COCIRCULAR);
+        REQUIRE(incircle(a,b,c,d) == Orientation::OUTSIDE);
+        REQUIRE(incircle(a,b,c,f) == Orientation::INSIDE);
+    }
+
+    SECTION("3D") {
+        using Point = Eigen::Matrix<Scalar, 3, 1>;
+        Point a(3,1),b(3,1),c(3,1),d(3,1),e(3,1),f(3,1);
+        a << 0.0, 0.0, 0.0;
+        b << 1.0, 0.0, 0.0;
+        c << 0.0, 1.0, 0.0;
+        d << 0.0, 0.0, 1.0;
+        e << 1.0, 1.0, 1.0;
+        f << std::numeric_limits<Scalar>::epsilon(), 0.0, 0.0;
+
+        REQUIRE(orient3d(a, b, c, d) == Orientation::NEGATIVE);
+        REQUIRE(orient3d(a, b, d, c) == Orientation::POSITIVE);
+        REQUIRE(orient3d(a, b, d, d) == Orientation::COPLANAR);
+        REQUIRE(orient3d(a, a, a, a) == Orientation::COPLANAR);
+        REQUIRE(orient3d(a, b, f, c) == Orientation::COPLANAR);
+
+        REQUIRE(insphere(a, b, c, d, e) == Orientation::COSPHERICAL);
+        REQUIRE(insphere(a, b, d, e, c) == Orientation::COSPHERICAL);
+        REQUIRE(insphere(b, c, e, d, ((a+b)*0.5).eval()) == Orientation::INSIDE);
+        REQUIRE(insphere(b, c, e, d, (-f).eval()) == Orientation::OUTSIDE);
+        REQUIRE(insphere(f, b, d, c, e) == Orientation::INSIDE);
+    }
+
+    SECTION("Predicate and triangle") {
+        Eigen::Matrix<double, -1, -1> vertices(4, 2);
+        Eigen::Matrix<double, -1, -1> holes;
+        Eigen::Matrix<int, -1, -1> edges;
+        vertices << 0.0, 0.0,
+                    1.0, 0.0,
+                    0.0, 1.0,
+                    1.0, 1.0;
+
+        Eigen::Matrix<double, -1, -1> out_vertices;
+        Eigen::Matrix<int, -1, -1> out_faces;
+
+        // Run constrained Delaunay.
+        igl::triangle::triangulate(vertices, edges, holes, "QcYY",
+                out_vertices, out_faces);
+        REQUIRE(out_vertices.rows() == 4);
+        REQUIRE(out_vertices.cols() == 2);
+        REQUIRE(out_faces.rows() == 2);
+        REQUIRE(out_faces.cols() == 3);
+    }
+}