better cm

.github/workflows/continuous.yml

@@ -25,14 +25,14 @@ jobs:
     strategy:
       fail-fast: false
       matrix:
-        os: [ubuntu-24.04, macos-15]
+        os: [ubuntu-24.04, macos-latest]
         config: [Release]
         build-params: [ {static:  ON, tutorials:  ON, tests:  ON }, {static: OFF, tutorials: OFF, tests:  ON }, {static: OFF, tutorials:  ON, tests: OFF }]
     env:
       IGL_NUM_THREADS: 1  # See https://github.com/libigl/libigl/pull/996
     steps:
       - name: Checkout repository
-        uses: actions/checkout@v4
+        uses: actions/checkout@v1
         with:
           fetch-depth: 10
 
@@ -56,7 +56,7 @@ jobs:
         if: runner.os == 'macOS'
         uses: maxim-lobanov/setup-xcode@v1
         with:
-          xcode-version: latest-stable
+          xcode-version: '16.3.0'
 
       - name: Cache Build
         id: cache-build
@@ -118,15 +118,31 @@ jobs:
       ]
     steps:
       - name: Checkout repository
-        uses: actions/checkout@v4
+        uses: actions/checkout@v1
         with:
           fetch-depth: 10
 
       - name: Install Ninja
         uses: seanmiddleditch/gha-setup-ninja@master
 
-      - name: Install sccache
-        uses: mozilla-actions/sccache-action@v0.0.9
+      - name: Set env variable for sccache
+        run: |
+          echo "appdata=$env:LOCALAPPDATA" >> ${env:GITHUB_ENV}
+
+      - name: Cache build
+        id: cache-build
+        uses: actions/cache@v4
+        with:
+          path: ${{ env.appdata }}\Mozilla\sccache
+          key: ${{ runner.os }}-${{ matrix.config }}-${{ matrix.build-params.static }}-cache
+
+      - name: Prepare sccache
+        run: |
+          iwr -useb 'https://raw.githubusercontent.com/scoopinstaller/install/master/install.ps1' -outfile 'install.ps1'
+          .\install.ps1 -RunAsAdmin
+          scoop install sccache --global
+          # Scoop modifies the PATH so we make it available for the next steps of the job
+          echo "${env:PATH}" >> ${env:GITHUB_PATH}
 
       - name: Configure and build
         shell: cmd
@@ -158,3 +174,4 @@ jobs:
 
       - name: Tests
         run: cd build; ctest --show-only; ctest --verbose -j2
+

.gitignore

@@ -51,4 +51,3 @@ LibiglOptions.cmake
 dox/
 latex/
 scripts/
-CLAUDE.md

cmake/recipes/external/eigen.cmake

@@ -8,7 +8,7 @@ include(FetchContent)
 FetchContent_Declare(
     eigen
     GIT_REPOSITORY https://gitlab.com/libeigen/eigen.git
-    GIT_TAG tags/5.0.1
+    GIT_TAG tags/3.4.0
     GIT_SHALLOW TRUE
 )
 FetchContent_GetProperties(eigen)

include/igl/FastWindingNumberForSoups.h

@@ -281,7 +281,7 @@ static constexpr inline fpreal32 SYSmax(fpreal32 a, fpreal32 b)	{ return h_max(a
 static constexpr inline fpreal64 SYSmin(fpreal64 a, fpreal64 b)	{ return h_min(a,b); }
 static constexpr inline fpreal64 SYSmax(fpreal64 a, fpreal64 b)	{ return h_max(a,b); }
 
-// Some systems have size_t as a separate type from uint.  Some don't.
+// Some systems have size_t as a seperate type from uint.  Some don't.
 #if (defined(LINUX) && defined(IA64)) || defined(MBSD)
 static constexpr inline size_t SYSmin(size_t a, size_t b)		{ return h_min(a,b); }
 static constexpr inline size_t SYSmax(size_t a, size_t b)		{ return h_max(a,b); }
@@ -742,408 +742,408 @@ vm_allbits(const v4si &a)
 
 #endif
 #endif
-#pragma once
-#ifndef __SSE__
-#ifndef __VM_SIMDFunc__
-#define __VM_SIMDFunc__
-
-
-
-#include <cmath>
+#pragma once
+#ifndef __SSE__
+#ifndef __VM_SIMDFunc__
+#define __VM_SIMDFunc__
+
 
+
+#include <cmath>
+
 namespace igl {
   /// @private
-  namespace FastWindingNumber {
-
-struct v4si {
-	int32 v[4];
-};
-
-struct v4sf {
-	float v[4];
-};
-
-static SYS_FORCE_INLINE v4sf V4SF(const v4si &v) {
-	static_assert(sizeof(v4si) == sizeof(v4sf) && alignof(v4si) == alignof(v4sf), "v4si and v4sf must be compatible");
-	return *(const v4sf*)&v;
-}
-
-static SYS_FORCE_INLINE v4si V4SI(const v4sf &v) {
-	static_assert(sizeof(v4si) == sizeof(v4sf) && alignof(v4si) == alignof(v4sf), "v4si and v4sf must be compatible");
-	return *(const v4si*)&v;
-}
-
-static SYS_FORCE_INLINE int32 conditionMask(bool c) {
-	return c ? int32(0xFFFFFFFF) : 0;
-}
-
-static SYS_FORCE_INLINE v4sf
-VM_SPLATS(float f) {
-	return v4sf{{f, f, f, f}};
-}
-
-static SYS_FORCE_INLINE v4si
-VM_SPLATS(uint32 i) {
-	return v4si{{int32(i), int32(i), int32(i), int32(i)}};
-}
-
-static SYS_FORCE_INLINE v4si
-VM_SPLATS(int32 i) {
-	return v4si{{i, i, i, i}};
-}
-
-static SYS_FORCE_INLINE v4sf
-VM_SPLATS(float a, float b, float c, float d) {
-	return v4sf{{a, b, c, d}};
-}
-
-static SYS_FORCE_INLINE v4si
-VM_SPLATS(uint32 a, uint32 b, uint32 c, uint32 d) {
-	return v4si{{int32(a), int32(b), int32(c), int32(d)}};
-}
-
-static SYS_FORCE_INLINE v4si
-VM_SPLATS(int32 a, int32 b, int32 c, int32 d) {
-	return v4si{{a, b, c, d}};
-}
-
-static SYS_FORCE_INLINE v4si
-VM_LOAD(const int32 v[4]) {
-	return v4si{{v[0], v[1], v[2], v[3]}};
-}
-
-static SYS_FORCE_INLINE v4sf
-VM_LOAD(const float v[4]) {
-	return v4sf{{v[0], v[1], v[2], v[3]}};
-}
-
-
-static inline v4si VM_ICMPEQ(v4si a, v4si b) {
-	return v4si{{
-		conditionMask(a.v[0] == b.v[0]),
-		conditionMask(a.v[1] == b.v[1]),
-		conditionMask(a.v[2] == b.v[2]),
-		conditionMask(a.v[3] == b.v[3])
-	}};
-}
-
-static inline v4si VM_ICMPGT(v4si a, v4si b) {
-	return v4si{{
-		conditionMask(a.v[0] > b.v[0]),
-		conditionMask(a.v[1] > b.v[1]),
-		conditionMask(a.v[2] > b.v[2]),
-		conditionMask(a.v[3] > b.v[3])
-	}};
-}
-
-static inline v4si VM_ICMPLT(v4si a, v4si b) {
-	return v4si{{
-		conditionMask(a.v[0] < b.v[0]),
-		conditionMask(a.v[1] < b.v[1]),
-		conditionMask(a.v[2] < b.v[2]),
-		conditionMask(a.v[3] < b.v[3])
-	}};
-}
-
-static inline v4si VM_IADD(v4si a, v4si b) {
-	return v4si{{
-		(a.v[0] + b.v[0]),
-		(a.v[1] + b.v[1]),
-		(a.v[2] + b.v[2]),
-		(a.v[3] + b.v[3])
-	}};
-}
-
-static inline v4si VM_ISUB(v4si a, v4si b) {
-	return v4si{{
-		(a.v[0] - b.v[0]),
-		(a.v[1] - b.v[1]),
-		(a.v[2] - b.v[2]),
-		(a.v[3] - b.v[3])
-	}};
-}
-
-static inline v4si VM_OR(v4si a, v4si b) {
-	return v4si{{
-		(a.v[0] | b.v[0]),
-		(a.v[1] | b.v[1]),
-		(a.v[2] | b.v[2]),
-		(a.v[3] | b.v[3])
-	}};
-}
-
-static inline v4si VM_AND(v4si a, v4si b) {
-	return v4si{{
-		(a.v[0] & b.v[0]),
-		(a.v[1] & b.v[1]),
-		(a.v[2] & b.v[2]),
-		(a.v[3] & b.v[3])
-	}};
-}
-
-static inline v4si VM_ANDNOT(v4si a, v4si b) {
-	return v4si{{
-		((~a.v[0]) & b.v[0]),
-		((~a.v[1]) & b.v[1]),
-		((~a.v[2]) & b.v[2]),
-		((~a.v[3]) & b.v[3])
-	}};
-}
-
-static inline v4si VM_XOR(v4si a, v4si b) {
-	return v4si{{
-		(a.v[0] ^ b.v[0]),
-		(a.v[1] ^ b.v[1]),
-		(a.v[2] ^ b.v[2]),
-		(a.v[3] ^ b.v[3])
-	}};
-}
-
-static SYS_FORCE_INLINE int
-VM_EXTRACT(const v4si v, int index) {
-	return v.v[index];
-}
-
-static SYS_FORCE_INLINE float
-VM_EXTRACT(const v4sf v, int index) {
-	return v.v[index];
-}
-
-static SYS_FORCE_INLINE v4si
-VM_INSERT(v4si v, int32 value, int index) {
-	v.v[index] = value;
-	return v;
-}
-
-static SYS_FORCE_INLINE v4sf
-VM_INSERT(v4sf v, float value, int index) {
-	v.v[index] = value;
-	return v;
-}
-
-static inline v4si VM_CMPEQ(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] == b.v[0]),
-		conditionMask(a.v[1] == b.v[1]),
-		conditionMask(a.v[2] == b.v[2]),
-		conditionMask(a.v[3] == b.v[3])
-	}};
-}
-
-static inline v4si VM_CMPNE(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] != b.v[0]),
-		conditionMask(a.v[1] != b.v[1]),
-		conditionMask(a.v[2] != b.v[2]),
-		conditionMask(a.v[3] != b.v[3])
-	}};
-}
-
-static inline v4si VM_CMPGT(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] > b.v[0]),
-		conditionMask(a.v[1] > b.v[1]),
-		conditionMask(a.v[2] > b.v[2]),
-		conditionMask(a.v[3] > b.v[3])
-	}};
-}
-
-static inline v4si VM_CMPLT(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] < b.v[0]),
-		conditionMask(a.v[1] < b.v[1]),
-		conditionMask(a.v[2] < b.v[2]),
-		conditionMask(a.v[3] < b.v[3])
-	}};
-}
-
-static inline v4si VM_CMPGE(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] >= b.v[0]),
-		conditionMask(a.v[1] >= b.v[1]),
-		conditionMask(a.v[2] >= b.v[2]),
-		conditionMask(a.v[3] >= b.v[3])
-	}};
-}
-
-static inline v4si VM_CMPLE(v4sf a, v4sf b) {
-	return v4si{{
-		conditionMask(a.v[0] <= b.v[0]),
-		conditionMask(a.v[1] <= b.v[1]),
-		conditionMask(a.v[2] <= b.v[2]),
-		conditionMask(a.v[3] <= b.v[3])
-	}};
-}
-
-static inline v4sf VM_ADD(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] + b.v[0]),
-		(a.v[1] + b.v[1]),
-		(a.v[2] + b.v[2]),
-		(a.v[3] + b.v[3])
-	}};
-}
-
-static inline v4sf VM_SUB(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] - b.v[0]),
-		(a.v[1] - b.v[1]),
-		(a.v[2] - b.v[2]),
-		(a.v[3] - b.v[3])
-	}};
-}
-
-static inline v4sf VM_NEG(v4sf a) {
-	return v4sf{{
-		(-a.v[0]),
-		(-a.v[1]),
-		(-a.v[2]),
-		(-a.v[3])
-	}};
-}
-
-static inline v4sf VM_MUL(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] * b.v[0]),
-		(a.v[1] * b.v[1]),
-		(a.v[2] * b.v[2]),
-		(a.v[3] * b.v[3])
-	}};
-}
-
-static inline v4sf VM_DIV(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] / b.v[0]),
-		(a.v[1] / b.v[1]),
-		(a.v[2] / b.v[2]),
-		(a.v[3] / b.v[3])
-	}};
-}
-
-static inline v4sf VM_MADD(v4sf a, v4sf b, v4sf c) {
-	return v4sf{{
-		(a.v[0] * b.v[0]) + c.v[0],
-		(a.v[1] * b.v[1]) + c.v[1],
-		(a.v[2] * b.v[2]) + c.v[2],
-		(a.v[3] * b.v[3]) + c.v[3]
-	}};
-}
-
-static inline v4sf VM_ABS(v4sf a) {
-	return v4sf{{
-		(a.v[0] < 0) ? -a.v[0] : a.v[0],
-		(a.v[1] < 0) ? -a.v[1] : a.v[1],
-		(a.v[2] < 0) ? -a.v[2] : a.v[2],
-		(a.v[3] < 0) ? -a.v[3] : a.v[3]
-	}};
-}
-
-static inline v4sf VM_MAX(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] < b.v[0]) ? b.v[0] : a.v[0],
-		(a.v[1] < b.v[1]) ? b.v[1] : a.v[1],
-		(a.v[2] < b.v[2]) ? b.v[2] : a.v[2],
-		(a.v[3] < b.v[3]) ? b.v[3] : a.v[3]
-	}};
-}
-
-static inline v4sf VM_MIN(v4sf a, v4sf b) {
-	return v4sf{{
-		(a.v[0] > b.v[0]) ? b.v[0] : a.v[0],
-		(a.v[1] > b.v[1]) ? b.v[1] : a.v[1],
-		(a.v[2] > b.v[2]) ? b.v[2] : a.v[2],
-		(a.v[3] > b.v[3]) ? b.v[3] : a.v[3]
-	}};
-}
-
-static inline v4sf VM_INVERT(v4sf a) {
-	return v4sf{{
-		(1.0f/a.v[0]),
-		(1.0f/a.v[1]),
-		(1.0f/a.v[2]),
-		(1.0f/a.v[3])
-	}};
-}
-
-static inline v4sf VM_SQRT(v4sf a) {
-	return v4sf{{
-		std::sqrt(a.v[0]),
-		std::sqrt(a.v[1]),
-		std::sqrt(a.v[2]),
-		std::sqrt(a.v[3])
-	}};
-}
-
-static inline v4si VM_INT(v4sf a) {
-	return v4si{{
-		int32(a.v[0]),
-		int32(a.v[1]),
-		int32(a.v[2]),
-		int32(a.v[3])
-	}};
-}
-
-static inline v4sf VM_IFLOAT(v4si a) {
-	return v4sf{{
-		float(a.v[0]),
-		float(a.v[1]),
-		float(a.v[2]),
-		float(a.v[3])
-	}};
-}
-
-static SYS_FORCE_INLINE void VM_P_FLOOR() {}
-
-static SYS_FORCE_INLINE int32 singleIntFloor(float f) {
-	// Casting to int32 usually truncates toward zero, instead of rounding down,
-	// so subtract one if the result is above f.
-	int32 i = int32(f);
-	i -= (float(i) > f);
-	return i;
-}
-static inline v4si VM_FLOOR(v4sf a) {
-	return v4si{{
-		singleIntFloor(a.v[0]),
-		singleIntFloor(a.v[1]),
-		singleIntFloor(a.v[2]),
-		singleIntFloor(a.v[3])
-	}};
-}
-
-static SYS_FORCE_INLINE void VM_E_FLOOR() {}
-
-static SYS_FORCE_INLINE bool vm_allbits(v4si a) {
-	return (
-		(a.v[0] == -1) && 
-		(a.v[1] == -1) && 
-		(a.v[2] == -1) && 
-		(a.v[3] == -1)
-	);
-}
-
-int SYS_FORCE_INLINE _mm_movemask_ps(const v4si& v) {
-	return (
-		int(v.v[0] < 0) |
-		(int(v.v[1] < 0)<<1) |
-		(int(v.v[2] < 0)<<2) |
-		(int(v.v[3] < 0)<<3)
-	);
-}
-
-int SYS_FORCE_INLINE _mm_movemask_ps(const v4sf& v) {
-	// Use std::signbit just in case it needs to distinguish between +0 and -0
-	// or between positive and negative NaN values (e.g. these could really
-	// be integers instead of floats).
-	return (
-		int(std::signbit(v.v[0])) |
-		(int(std::signbit(v.v[1]))<<1) |
-		(int(std::signbit(v.v[2]))<<2) |
-		(int(std::signbit(v.v[3]))<<3)
-	);
-}
-}}
-#endif
-#endif
+  namespace FastWindingNumber {
+
+struct v4si {
+	int32 v[4];
+};
+
+struct v4sf {
+	float v[4];
+};
+
+static SYS_FORCE_INLINE v4sf V4SF(const v4si &v) {
+	static_assert(sizeof(v4si) == sizeof(v4sf) && alignof(v4si) == alignof(v4sf), "v4si and v4sf must be compatible");
+	return *(const v4sf*)&v;
+}
+
+static SYS_FORCE_INLINE v4si V4SI(const v4sf &v) {
+	static_assert(sizeof(v4si) == sizeof(v4sf) && alignof(v4si) == alignof(v4sf), "v4si and v4sf must be compatible");
+	return *(const v4si*)&v;
+}
+
+static SYS_FORCE_INLINE int32 conditionMask(bool c) {
+	return c ? int32(0xFFFFFFFF) : 0;
+}
+
+static SYS_FORCE_INLINE v4sf
+VM_SPLATS(float f) {
+	return v4sf{{f, f, f, f}};
+}
+
+static SYS_FORCE_INLINE v4si
+VM_SPLATS(uint32 i) {
+	return v4si{{int32(i), int32(i), int32(i), int32(i)}};
+}
+
+static SYS_FORCE_INLINE v4si
+VM_SPLATS(int32 i) {
+	return v4si{{i, i, i, i}};
+}
+
+static SYS_FORCE_INLINE v4sf
+VM_SPLATS(float a, float b, float c, float d) {
+	return v4sf{{a, b, c, d}};
+}
+
+static SYS_FORCE_INLINE v4si
+VM_SPLATS(uint32 a, uint32 b, uint32 c, uint32 d) {
+	return v4si{{int32(a), int32(b), int32(c), int32(d)}};
+}
+
+static SYS_FORCE_INLINE v4si
+VM_SPLATS(int32 a, int32 b, int32 c, int32 d) {
+	return v4si{{a, b, c, d}};
+}
+
+static SYS_FORCE_INLINE v4si
+VM_LOAD(const int32 v[4]) {
+	return v4si{{v[0], v[1], v[2], v[3]}};
+}
+
+static SYS_FORCE_INLINE v4sf
+VM_LOAD(const float v[4]) {
+	return v4sf{{v[0], v[1], v[2], v[3]}};
+}
+
+
+static inline v4si VM_ICMPEQ(v4si a, v4si b) {
+	return v4si{{
+		conditionMask(a.v[0] == b.v[0]),
+		conditionMask(a.v[1] == b.v[1]),
+		conditionMask(a.v[2] == b.v[2]),
+		conditionMask(a.v[3] == b.v[3])
+	}};
+}
+
+static inline v4si VM_ICMPGT(v4si a, v4si b) {
+	return v4si{{
+		conditionMask(a.v[0] > b.v[0]),
+		conditionMask(a.v[1] > b.v[1]),
+		conditionMask(a.v[2] > b.v[2]),
+		conditionMask(a.v[3] > b.v[3])
+	}};
+}
+
+static inline v4si VM_ICMPLT(v4si a, v4si b) {
+	return v4si{{
+		conditionMask(a.v[0] < b.v[0]),
+		conditionMask(a.v[1] < b.v[1]),
+		conditionMask(a.v[2] < b.v[2]),
+		conditionMask(a.v[3] < b.v[3])
+	}};
+}
+
+static inline v4si VM_IADD(v4si a, v4si b) {
+	return v4si{{
+		(a.v[0] + b.v[0]),
+		(a.v[1] + b.v[1]),
+		(a.v[2] + b.v[2]),
+		(a.v[3] + b.v[3])
+	}};
+}
+
+static inline v4si VM_ISUB(v4si a, v4si b) {
+	return v4si{{
+		(a.v[0] - b.v[0]),
+		(a.v[1] - b.v[1]),
+		(a.v[2] - b.v[2]),
+		(a.v[3] - b.v[3])
+	}};
+}
+
+static inline v4si VM_OR(v4si a, v4si b) {
+	return v4si{{
+		(a.v[0] | b.v[0]),
+		(a.v[1] | b.v[1]),
+		(a.v[2] | b.v[2]),
+		(a.v[3] | b.v[3])
+	}};
+}
+
+static inline v4si VM_AND(v4si a, v4si b) {
+	return v4si{{
+		(a.v[0] & b.v[0]),
+		(a.v[1] & b.v[1]),
+		(a.v[2] & b.v[2]),
+		(a.v[3] & b.v[3])
+	}};
+}
+
+static inline v4si VM_ANDNOT(v4si a, v4si b) {
+	return v4si{{
+		((~a.v[0]) & b.v[0]),
+		((~a.v[1]) & b.v[1]),
+		((~a.v[2]) & b.v[2]),
+		((~a.v[3]) & b.v[3])
+	}};
+}
+
+static inline v4si VM_XOR(v4si a, v4si b) {
+	return v4si{{
+		(a.v[0] ^ b.v[0]),
+		(a.v[1] ^ b.v[1]),
+		(a.v[2] ^ b.v[2]),
+		(a.v[3] ^ b.v[3])
+	}};
+}
+
+static SYS_FORCE_INLINE int
+VM_EXTRACT(const v4si v, int index) {
+	return v.v[index];
+}
+
+static SYS_FORCE_INLINE float
+VM_EXTRACT(const v4sf v, int index) {
+	return v.v[index];
+}
+
+static SYS_FORCE_INLINE v4si
+VM_INSERT(v4si v, int32 value, int index) {
+	v.v[index] = value;
+	return v;
+}
+
+static SYS_FORCE_INLINE v4sf
+VM_INSERT(v4sf v, float value, int index) {
+	v.v[index] = value;
+	return v;
+}
+
+static inline v4si VM_CMPEQ(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] == b.v[0]),
+		conditionMask(a.v[1] == b.v[1]),
+		conditionMask(a.v[2] == b.v[2]),
+		conditionMask(a.v[3] == b.v[3])
+	}};
+}
+
+static inline v4si VM_CMPNE(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] != b.v[0]),
+		conditionMask(a.v[1] != b.v[1]),
+		conditionMask(a.v[2] != b.v[2]),
+		conditionMask(a.v[3] != b.v[3])
+	}};
+}
+
+static inline v4si VM_CMPGT(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] > b.v[0]),
+		conditionMask(a.v[1] > b.v[1]),
+		conditionMask(a.v[2] > b.v[2]),
+		conditionMask(a.v[3] > b.v[3])
+	}};
+}
+
+static inline v4si VM_CMPLT(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] < b.v[0]),
+		conditionMask(a.v[1] < b.v[1]),
+		conditionMask(a.v[2] < b.v[2]),
+		conditionMask(a.v[3] < b.v[3])
+	}};
+}
+
+static inline v4si VM_CMPGE(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] >= b.v[0]),
+		conditionMask(a.v[1] >= b.v[1]),
+		conditionMask(a.v[2] >= b.v[2]),
+		conditionMask(a.v[3] >= b.v[3])
+	}};
+}
+
+static inline v4si VM_CMPLE(v4sf a, v4sf b) {
+	return v4si{{
+		conditionMask(a.v[0] <= b.v[0]),
+		conditionMask(a.v[1] <= b.v[1]),
+		conditionMask(a.v[2] <= b.v[2]),
+		conditionMask(a.v[3] <= b.v[3])
+	}};
+}
+
+static inline v4sf VM_ADD(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] + b.v[0]),
+		(a.v[1] + b.v[1]),
+		(a.v[2] + b.v[2]),
+		(a.v[3] + b.v[3])
+	}};
+}
+
+static inline v4sf VM_SUB(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] - b.v[0]),
+		(a.v[1] - b.v[1]),
+		(a.v[2] - b.v[2]),
+		(a.v[3] - b.v[3])
+	}};
+}
+
+static inline v4sf VM_NEG(v4sf a) {
+	return v4sf{{
+		(-a.v[0]),
+		(-a.v[1]),
+		(-a.v[2]),
+		(-a.v[3])
+	}};
+}
+
+static inline v4sf VM_MUL(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] * b.v[0]),
+		(a.v[1] * b.v[1]),
+		(a.v[2] * b.v[2]),
+		(a.v[3] * b.v[3])
+	}};
+}
+
+static inline v4sf VM_DIV(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] / b.v[0]),
+		(a.v[1] / b.v[1]),
+		(a.v[2] / b.v[2]),
+		(a.v[3] / b.v[3])
+	}};
+}
+
+static inline v4sf VM_MADD(v4sf a, v4sf b, v4sf c) {
+	return v4sf{{
+		(a.v[0] * b.v[0]) + c.v[0],
+		(a.v[1] * b.v[1]) + c.v[1],
+		(a.v[2] * b.v[2]) + c.v[2],
+		(a.v[3] * b.v[3]) + c.v[3]
+	}};
+}
+
+static inline v4sf VM_ABS(v4sf a) {
+	return v4sf{{
+		(a.v[0] < 0) ? -a.v[0] : a.v[0],
+		(a.v[1] < 0) ? -a.v[1] : a.v[1],
+		(a.v[2] < 0) ? -a.v[2] : a.v[2],
+		(a.v[3] < 0) ? -a.v[3] : a.v[3]
+	}};
+}
+
+static inline v4sf VM_MAX(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] < b.v[0]) ? b.v[0] : a.v[0],
+		(a.v[1] < b.v[1]) ? b.v[1] : a.v[1],
+		(a.v[2] < b.v[2]) ? b.v[2] : a.v[2],
+		(a.v[3] < b.v[3]) ? b.v[3] : a.v[3]
+	}};
+}
+
+static inline v4sf VM_MIN(v4sf a, v4sf b) {
+	return v4sf{{
+		(a.v[0] > b.v[0]) ? b.v[0] : a.v[0],
+		(a.v[1] > b.v[1]) ? b.v[1] : a.v[1],
+		(a.v[2] > b.v[2]) ? b.v[2] : a.v[2],
+		(a.v[3] > b.v[3]) ? b.v[3] : a.v[3]
+	}};
+}
+
+static inline v4sf VM_INVERT(v4sf a) {
+	return v4sf{{
+		(1.0f/a.v[0]),
+		(1.0f/a.v[1]),
+		(1.0f/a.v[2]),
+		(1.0f/a.v[3])
+	}};
+}
+
+static inline v4sf VM_SQRT(v4sf a) {
+	return v4sf{{
+		std::sqrt(a.v[0]),
+		std::sqrt(a.v[1]),
+		std::sqrt(a.v[2]),
+		std::sqrt(a.v[3])
+	}};
+}
+
+static inline v4si VM_INT(v4sf a) {
+	return v4si{{
+		int32(a.v[0]),
+		int32(a.v[1]),
+		int32(a.v[2]),
+		int32(a.v[3])
+	}};
+}
+
+static inline v4sf VM_IFLOAT(v4si a) {
+	return v4sf{{
+		float(a.v[0]),
+		float(a.v[1]),
+		float(a.v[2]),
+		float(a.v[3])
+	}};
+}
+
+static SYS_FORCE_INLINE void VM_P_FLOOR() {}
+
+static SYS_FORCE_INLINE int32 singleIntFloor(float f) {
+	// Casting to int32 usually truncates toward zero, instead of rounding down,
+	// so subtract one if the result is above f.
+	int32 i = int32(f);
+	i -= (float(i) > f);
+	return i;
+}
+static inline v4si VM_FLOOR(v4sf a) {
+	return v4si{{
+		singleIntFloor(a.v[0]),
+		singleIntFloor(a.v[1]),
+		singleIntFloor(a.v[2]),
+		singleIntFloor(a.v[3])
+	}};
+}
+
+static SYS_FORCE_INLINE void VM_E_FLOOR() {}
+
+static SYS_FORCE_INLINE bool vm_allbits(v4si a) {
+	return (
+		(a.v[0] == -1) && 
+		(a.v[1] == -1) && 
+		(a.v[2] == -1) && 
+		(a.v[3] == -1)
+	);
+}
+
+int SYS_FORCE_INLINE _mm_movemask_ps(const v4si& v) {
+	return (
+		int(v.v[0] < 0) |
+		(int(v.v[1] < 0)<<1) |
+		(int(v.v[2] < 0)<<2) |
+		(int(v.v[3] < 0)<<3)
+	);
+}
+
+int SYS_FORCE_INLINE _mm_movemask_ps(const v4sf& v) {
+	// Use std::signbit just in case it needs to distinguish between +0 and -0
+	// or between positive and negative NaN values (e.g. these could really
+	// be integers instead of floats).
+	return (
+		int(std::signbit(v.v[0])) |
+		(int(std::signbit(v.v[1]))<<1) |
+		(int(std::signbit(v.v[2]))<<2) |
+		(int(std::signbit(v.v[3]))<<3)
+	);
+}
+}}
+#endif
+#endif
 /*
  * Copyright (c) 2018 Side Effects Software Inc.
  *
@@ -3758,7 +3758,7 @@ namespace UT_Thread { inline int getNumProcessors() {
 /////
 ///// Requirements for the Range functor are:
 /////   - the requirements of the tbb Range Concept
-/////   - UT_estimatorNumItems<Range> must return the estimated number of work items
+/////   - UT_estimatorNumItems<Range> must return the the estimated number of work items
 /////     for the range. When Range::size() is not the correct estimate, then a 
 /////     (partial) specialization of UT_estimatorNumItemsimatorRange must be provided
 /////     for the type Range.

include/igl/HalfEdgeIterator.cpp

@@ -8,8 +8,6 @@
 
 #include "HalfEdgeIterator.h"
 
-#include <cassert>
-
 template <typename DerivedF, typename DerivedFF, typename DerivedFFi>
 IGL_INLINE igl::HalfEdgeIterator<DerivedF,DerivedFF,DerivedFFi>::HalfEdgeIterator(
     const Eigen::MatrixBase<DerivedF>& _F,

include/igl/MshLoader.cpp

@@ -114,31 +114,10 @@ IGL_INLINE igl::MshLoader::MshLoader(const std::string &filename) {
     fin.close();
 }
 
-IGL_INLINE int igl::MshLoader::node_dense_index(int node_tag) const {
-    const auto it = m_node_tag_to_dense.find(node_tag);
-    if (it == m_node_tag_to_dense.end()) {
-        std::stringstream err_msg;
-        err_msg << "Unknown node tag: " << node_tag;
-        throw std::runtime_error(err_msg.str());
-    }
-	return it->second;
-}
-
-IGL_INLINE int igl::MshLoader::element_dense_index(int elem_tag) const {
-	const auto it = m_element_tag_to_dense.find(elem_tag);
-    if (it == m_element_tag_to_dense.end()) {
-		std::stringstream err_msg;
-        err_msg << "Unknown element tag: " << elem_tag;
-		throw std::runtime_error(err_msg.str());
-    }
-	return it->second;
-}
-
 IGL_INLINE void igl::MshLoader::parse_nodes(std::ifstream& fin) {
     size_t num_nodes;
     fin >> num_nodes;
     m_nodes.resize(num_nodes*3);
-    m_node_tag_to_dense.clear();
 
     if (m_binary) {
 		size_t stride = (4+3*m_data_size);
@@ -148,37 +127,23 @@ IGL_INLINE void igl::MshLoader::parse_nodes(std::ifstream& fin) {
         fin.read(data, num_bytes);
 
         for (size_t i=0; i<num_nodes; i++) {
-            int node_tag;
-			memcpy(&node_tag, data+i*stride, sizeof(int));
-
-            if (node_tag <= 0) {
-                throw std::runtime_error("Invalid node tag");
-            }
-            if (m_node_tag_to_dense.find(node_tag) != m_node_tag_to_dense.end()) {
-                throw std::runtime_error("Duplicate node tag");
-            }
-            m_node_tag_to_dense[node_tag] = static_cast<int>(i);
-            // directly move into vector storage
-            // this works only when m_data_size==sizeof(Float)==sizeof(double)
-			memcpy(&m_nodes[i*3], data+i*stride + 4, m_data_size*3);
+            int node_idx;
+			memcpy(&node_idx, data+i*stride, sizeof(int));
+			node_idx-=1;
+			// directly move into vector storage
+			// this works only when m_data_size==sizeof(Float)==sizeof(double)
+			memcpy(&m_nodes[node_idx*3], data+i*stride + 4, m_data_size*3);
         }
         delete [] data;
     } else {
-        int node_tag;
+        int node_idx;
         for (size_t i=0; i<num_nodes; i++) {
-            fin >> node_tag;
-
-            if (node_tag <= 0) {
-				throw std::runtime_error("Invalid node tag");
-            }
-            if (m_node_tag_to_dense.find(node_tag) != m_node_tag_to_dense.end()) {
-                throw std::runtime_error("Duplicate node tag");
-            }
-            m_node_tag_to_dense[node_tag] = static_cast<int>(i);
+            fin >> node_idx;
+            node_idx -= 1;
             // here it's 3D node explicitly
-            fin >> m_nodes[i*3]
-                >> m_nodes[i*3+1]
-                >> m_nodes[i*3+2];
+            fin >> m_nodes[node_idx*3]
+                >> m_nodes[node_idx*3+1]
+                >> m_nodes[node_idx*3+2];
         }
     }
 }
@@ -187,7 +152,6 @@ IGL_INLINE void igl::MshLoader::parse_elements(std::ifstream& fin) {
     m_elements_tags.resize(2); //hardcoded to have 2 tags
     size_t num_elements;
     fin >> num_elements;
-    m_element_tag_to_dense.clear();
 
     size_t nodes_per_element;
 
@@ -204,24 +168,15 @@ IGL_INLINE void igl::MshLoader::parse_elements(std::ifstream& fin) {
 
             // store node info
             for (size_t i=0; i<num_elems; i++) {
-                int elem_tag;
+                int elem_idx;
 
                 // all elements in the segment share the same elem_type and number of nodes per element
                 m_elements_types.push_back(elem_type);
                 m_elements_lengths.push_back(nodes_per_element);
 
-                fin.read((char*)&elem_tag, sizeof(int));
-
-                if (elem_tag <= 0) {
-                    throw std::runtime_error("Invalid element tag");
-                }
-                if (m_element_tag_to_dense.find(elem_tag) != m_element_tag_to_dense.end()) {
-                    throw std::runtime_error("Duplicate element tag");
-                }
-                m_element_tag_to_dense[elem_tag] = static_cast<int>(m_elements_ids.size());
-
-                elem_tag -= 1;
-                m_elements_ids.push_back(elem_tag);
+                fin.read((char*)&elem_idx, sizeof(int));
+                elem_idx -= 1;
+                m_elements_ids.push_back(elem_idx);
 
                 // read first two tags
                 for (size_t j=0; j<num_tags; j++) {
@@ -236,10 +191,10 @@ IGL_INLINE void igl::MshLoader::parse_elements(std::ifstream& fin) {
                 m_elements_nodes_idx.push_back(m_elements.size());
                 // Element values.
                 for (size_t j=0; j<nodes_per_element; j++) {
-                    int node_tag;
-                    fin.read((char*)&node_tag, sizeof(int));
+                    int idx;
+                    fin.read((char*)&idx, sizeof(int));
                     
-                    m_elements.push_back(node_dense_index(node_tag));
+                    m_elements.push_back(idx-1);
                 }
             }
             elem_read += num_elems;
@@ -247,16 +202,8 @@ IGL_INLINE void igl::MshLoader::parse_elements(std::ifstream& fin) {
     } else {
         for (size_t i=0; i<num_elements; i++) {
             // Parse per element header
-            int elem_tag, elem_type, num_tags;
-            fin >> elem_tag >> elem_type >> num_tags;
-
-            if (elem_tag <= 0) {
-                throw std::runtime_error("Invalid element tag");
-            }
-            if (m_element_tag_to_dense.find(elem_tag) != m_element_tag_to_dense.end()) {
-                throw std::runtime_error("Duplicate element tag");
-            }
-            m_element_tag_to_dense[elem_tag] = static_cast<int>(m_elements_ids.size());
+            int elem_num, elem_type, num_tags;
+            fin >> elem_num >> elem_type >> num_tags;
 
             // read tags.
             for (size_t j=0; j<num_tags; j++) {
@@ -271,14 +218,14 @@ IGL_INLINE void igl::MshLoader::parse_elements(std::ifstream& fin) {
             m_elements_types.push_back(elem_type);
             m_elements_lengths.push_back(nodes_per_element);
 
-            elem_tag -= 1;
-            m_elements_ids.push_back(elem_tag);
+            elem_num -= 1;
+            m_elements_ids.push_back(elem_num);
             m_elements_nodes_idx.push_back(m_elements.size());
             // Parse node idx.
             for (size_t j=0; j<nodes_per_element; j++) {
-                int node_tag;
-                fin >> node_tag;
-                m_elements.push_back(node_dense_index(node_tag)); // msh index starts from 1.
+                int idx;
+                fin >> idx;
+                m_elements.push_back(idx-1); // msh index starts from 1.
             }
         }
     }
@@ -327,7 +274,7 @@ IGL_INLINE void igl::MshLoader::parse_node_field( std::ifstream& fin ) {
     int num_components    = int_tags[1];
     int num_entries       = int_tags[2];
 
-    std::vector<Float> field((m_nodes.size()/3)*num_components);
+    std::vector<Float> field( num_entries*num_components );
 
     if (m_binary) {
         size_t num_bytes = (num_components * m_data_size + 4) * num_entries;
@@ -335,20 +282,23 @@ IGL_INLINE void igl::MshLoader::parse_node_field( std::ifstream& fin ) {
         igl::_msh_eat_white_space(fin);
         fin.read(data, num_bytes);
         for (size_t i=0; i<num_entries; i++) {
-			int node_tag;
-			memcpy(&node_tag,&data[i*(4+num_components*m_data_size)],4);
-
-            const int node_idx = node_dense_index(node_tag);
+			int node_idx;
+			memcpy(&node_idx,&data[i*(4+num_components*m_data_size)],4);
+			
+            if(node_idx<1) throw std::runtime_error("Negative or zero index");
+            node_idx -= 1;
+			
+            if(node_idx>=num_entries) throw std::runtime_error("Index too big");
             size_t base_idx = i*(4+num_components*m_data_size) + 4;
             // TODO: make this work when m_data_size != sizeof(double) ?
 			memcpy(&field[node_idx*num_components], &data[base_idx], num_components*m_data_size);
         }
         delete [] data;
     } else {
-        int node_tag;
+        int node_idx;
         for (size_t i=0; i<num_entries; i++) {
-            fin >> node_tag;
-            const int node_idx = node_dense_index(node_tag);
+            fin >> node_idx;
+            node_idx -= 1;
             for (size_t j=0; j<num_components; j++) {
                 fin >> field[node_idx*num_components+j];
             }
@@ -396,7 +346,7 @@ IGL_INLINE void igl::MshLoader::parse_element_field(std::ifstream& fin) {
     std::string fieldname = str_tags[0];
     int num_components = int_tags[1];
     int num_entries = int_tags[2];
-    std::vector<Float> field(m_elements_ids.size()*num_components);
+    std::vector<Float> field(num_entries*num_components);
 
     if (m_binary) {
         size_t num_bytes = (num_components * m_data_size + 4) * num_entries;
@@ -404,20 +354,20 @@ IGL_INLINE void igl::MshLoader::parse_element_field(std::ifstream& fin) {
         igl::_msh_eat_white_space(fin);
         fin.read(data, num_bytes);
         for (int i=0; i<num_entries; i++) {
-			int elem_tag;
+			int elem_idx;
 			// works with sizeof(int)==4
-			memcpy(&elem_tag, &data[i*(4+num_components*m_data_size)],4);
-            const int elem_idx = element_dense_index(elem_tag);
+			memcpy(&elem_idx, &data[i*(4+num_components*m_data_size)],4);
+            elem_idx -= 1;
 			
 			// directly copy data into vector storage space
 			memcpy(&field[elem_idx*num_components], &data[i*(4+num_components*m_data_size) + 4], m_data_size*num_components);
         }
         delete [] data;
     } else {
-        int elem_tag;
+        int elem_idx;
         for (size_t i=0; i<num_entries; i++) {
-            fin >> elem_tag;
-            const int elem_idx = element_dense_index(elem_tag);
+            fin >> elem_idx;
+            elem_idx -= 1;
             for (size_t j=0; j<num_components; j++) {
                 fin >> field[elem_idx*num_components+j];
             }

include/igl/MshLoader.h

@@ -15,7 +15,6 @@
 #include <string>
 #include <vector>
 #include <algorithm>
-#include <unordered_map>
 
 namespace igl {
 
@@ -156,15 +155,11 @@ class MshLoader {
         void parse_element_field(std::ifstream& fin);
         void parse_unknown_field(std::ifstream& fin,
                 const std::string& fieldname);
-        int node_dense_index(int node_tag) const;
-        int element_dense_index(int elem_tag) const;
 
     private:
         bool   m_binary;
         size_t m_data_size;
-        std::unordered_map<int, int> m_node_tag_to_dense;
-        std::unordered_map<int, int> m_element_tag_to_dense;
-
+        
         FloatVector m_nodes;    // len x 3 vector 
 
         IndexVector m_elements; // linear array for nodes corresponding to each element 

include/igl/SphereMeshWedge.cpp

@@ -8,7 +8,6 @@
 #include "SphereMeshWedge.h"
 #include "round_cone_signed_distance.h"
 #include "sign.h"
-#include <cassert>
 #include <Eigen/QR>
 #include <Eigen/Geometry>
 

include/igl/WindingNumberAntipodalScene.h

@@ -1,259 +0,0 @@
-// This file is part of libigl, a simple c++ geometry processing library.
-//
-// Copyright (C) 2026 Philip Trettner <trettner@shapedcode.com>, Cedric Martens <cedric.martens@umontreal.ca>
-//
-// 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_WINDINGNUMBERANTIPODALSCENE_H
-#define IGL_WINDINGNUMBERANTIPODALSCENE_H
-
-#include "PI.h"
-#include "parallel_for.h"
-
-#include <Eigen/Core>
-
-#include <cassert>
-#include <cmath>
-#include <cstdint>
-#include <unordered_map>
-#include <vector>
-
-namespace igl
-{
-  /// Precomputed scene for the Antipodal Method generalized winding number
-  ///
-  /// The scene stores only the weighted open-boundary edges of the input
-  /// triangle mesh and a fixed antipodal reference direction `x0`. Closed
-  /// (manifold) meshes have an empty boundary and the fractional term is
-  /// exactly zero.
-  ///
-  /// Querying the winding number additionally requires an `Intersector` that
-  /// returns the signed ray-mesh crossing count along `(p, x0)` over the
-  /// original mesh; the scene itself is intersector-agnostic. See
-  /// `igl::embree::EmbreeIntersector::signedIntersectionsRay` for an
-  /// optimized concrete implementation.
-  ///
-  /// ### Intersector concept
-  /// A type `I` satisfies the concept when it exposes:
-  ///   - `using OriginType    = ...;`     (3D row vector type)
-  ///   - `using DirectionType = ...;`     (3D row vector type)
-  ///   - `int signedIntersectionsRay(
-  ///         OriginType origin, DirectionType direction,
-  ///         /* defaulted tnear, tfar, mask */) const;`
-  ///
-  /// `winding_number` casts query point/direction to the intersector's types
-  /// at the call site, so a `double`-precision scene against a float-only
-  /// intersector works without an adaptor.
-  template <typename Scalar>
-  class WindingNumberAntipodalScene
-  {
-  public:
-    using Point     = Eigen::Matrix<Scalar, 1, 3>;
-    using Direction = Eigen::Matrix<Scalar, 1, 3>;
-
-  private:
-    struct WeightedSeg
-    {
-      Point a;
-      Point b;
-      Scalar w;
-    };
-
-  public:
-    /// Build a scene from a 3D triangle mesh.
-    ///
-    /// @param[in] V   #V by 3 list of vertex positions
-    /// @param[in] F   #F by 3 list of triangle indices
-    /// @param[in] x0  unit reference direction (defaults to a fixed non-axis-aligned vector)
-    template <typename DerivedV, typename DerivedF>
-    WindingNumberAntipodalScene(
-      const Eigen::MatrixBase<DerivedV> & V,
-      const Eigen::MatrixBase<DerivedF> & F,
-      const Direction & x0 = default_x0())
-      : m_x0(x0), m_face_count(static_cast<size_t>(F.rows()))
-    {
-      assert(V.cols() == 3 && "WindingNumberAntipodalScene: only 3D vertex positions are supported");
-      assert(F.cols() == 3 && "WindingNumberAntipodalScene: only triangle meshes are supported");
-      build_boundary_segments(V, F, m_boundary);
-    }
-
-    /// Single-point query: full generalized winding number at `p`
-    /// (fractional + signed integer crossings).
-    template <typename Intersector, typename Derivedp>
-    Scalar winding_number(
-      const Intersector & intersector,
-      const Eigen::MatrixBase<Derivedp> & p) const
-    {
-      const Point pp(static_cast<Scalar>(p(0)),
-                     static_cast<Scalar>(p(1)),
-                     static_cast<Scalar>(p(2)));
-      const Direction x1 = -m_x0;
-      Scalar area = Scalar(0);
-      for (const auto & ws : m_boundary)
-      {
-        const Point v0 = ws.a - pp;
-        const Point v1 = ws.b - pp;
-        area += ws.w * half_solid_angle_unorm(x1, v0, v1);
-      }
-      const Scalar frac = area / (Scalar(2) * Scalar(igl::PI));
-
-      using IO = typename Intersector::OriginType;
-      using ID = typename Intersector::DirectionType;
-      const IO io(static_cast<typename IO::Scalar>(pp(0)),
-                  static_cast<typename IO::Scalar>(pp(1)),
-                  static_cast<typename IO::Scalar>(pp(2)));
-      const ID id(static_cast<typename ID::Scalar>(m_x0(0)),
-                  static_cast<typename ID::Scalar>(m_x0(1)),
-                  static_cast<typename ID::Scalar>(m_x0(2)));
-      const int c = intersector.signedIntersectionsRay(io, id);
-
-      return frac + Scalar(c);
-    }
-
-    /// Batch query, parallelized via `igl::parallel_for`.
-    ///
-    /// @param[in]  intersector  Concept-compatible intersector built over the
-    ///                          same mesh used to construct the scene. Must
-    ///                          be safe to query concurrently.
-    /// @param[in]  O            #O by 3 list of query points
-    /// @param[out] W            #O by 1 list of winding numbers
-    template <typename Intersector, typename DerivedO, typename DerivedW>
-    void winding_number(
-      const Intersector & intersector,
-      const Eigen::MatrixBase<DerivedO> & O,
-      Eigen::PlainObjectBase<DerivedW> & W) const
-    {
-      W.resize(O.rows(), 1);
-
-      // Adaptive parallel-for threshold.
-      //
-      // The libigl thread pool has a roughly fixed ~1 ms TOTAL overhead per
-      // parallel_for invocation (not per iteration). So we only spawn the
-      // pool when the WHOLE batch is expected to take ≥ 1 ms.
-      //
-      // Per-query work heuristic:
-      //   t_q ≈ 50 * B + 100 * sqrt(F)  ns
-      // with B = boundary segment count, F = triangle count. Pool wins once
-      //   t_q * O > 10^6 ns  ⇒  O > 10^6 / t_q
-      // which is exactly parallel_for's `min_parallel` semantics.
-      //
-      // (This is a rough heuristic and should be revisited once parallel_for becomes lower-overhead)
-      const double t_q_ns =
-          50.0  * static_cast<double>(m_boundary.size()) +
-          100.0 * std::sqrt(static_cast<double>(m_face_count));
-      const size_t min_parallel = static_cast<size_t>(
-          std::ceil(1.0e6 / std::max(t_q_ns, 1.0)));
-
-      igl::parallel_for(O.rows(), [&](const int o)
-      {
-        W(o) = winding_number(intersector, O.row(o));
-      }, min_parallel);
-    }
-
-    /// Reference direction `x0` used to evaluate this scene.
-    const Direction & x0() const { return m_x0; }
-    /// Number of weighted boundary edge segments.
-    size_t num_boundary_segments() const { return m_boundary.size(); }
-    /// Number of triangles in the original mesh.
-    size_t num_faces() const { return m_face_count; }
-
-    /// Default reference direction: normalize(1, sqrt(2), sqrt(3)). A fixed
-    /// non-axis-aligned unit vector; any unit vector works per the paper;
-    /// this choice avoids accidental alignment with axis-aligned geometry.
-    static Direction default_x0()
-    {
-      Direction d(Scalar(1),
-                  Scalar(std::sqrt(2.0)),
-                  Scalar(std::sqrt(3.0)));
-      return d / d.norm();
-    }
-
-    /// Half the signed solid angle subtended by the spherical triangle
-    /// (x1, v0, v1) at the origin, via the unnormalized
-    /// Van Oosterom-Strackee formula. `x1` is expected to be a unit vector
-    /// (the antipodal "south pole" `-x0`); `v0`, `v1` need not be normalized.
-    /// Callers accumulate per-edge contributions and divide by 2π.
-    static Scalar half_solid_angle_unorm(
-      const Direction & x1, const Point & v0, const Point & v1)
-    {
-      const Scalar l0 = v0.norm();
-      const Scalar l1 = v1.norm();
-      const Scalar num   = x1.dot(v0.cross(v1));
-      const Scalar denom = l0 * l1
-                         + l1 * x1.dot(v0)
-                         + l0 * x1.dot(v1)
-                         + v0.dot(v1);
-      return std::atan2(num, denom);
-    }
-
-    /// Extract the open boundary as oriented, weighted edge segments.
-    ///
-    /// Each undirected edge {min, max} accumulates +1 for every triangle
-    /// that traverses it as min→max and -1 for every traversal max→min.
-    /// Interior edges of an oriented manifold cancel to zero and are
-    /// dropped. Surviving edges are emitted with positive weight; the
-    /// segment direction is flipped when the net count is negative so the
-    /// weight is always > 0. Non-manifold (≥3 incident triangles per edge)
-    /// is handled the same way. The surviving net count becomes the weight.
-    template <typename DerivedV, typename DerivedF>
-    static void build_boundary_segments(
-      const Eigen::MatrixBase<DerivedV> & V,
-      const Eigen::MatrixBase<DerivedF> & F,
-      std::vector<WeightedSeg> & out)
-    {
-      out.clear();
-      std::unordered_map<std::uint64_t, int> counts;
-      counts.reserve(static_cast<size_t>(F.rows()) * 3);
-
-      const auto pack = [](int lo, int hi) -> std::uint64_t
-      {
-        return (static_cast<std::uint64_t>(static_cast<std::uint32_t>(lo)) << 32)
-             |  static_cast<std::uint64_t>(static_cast<std::uint32_t>(hi));
-      };
-
-      for (Eigen::Index t = 0; t < F.rows(); ++t)
-      {
-        const int tri[3] = {
-          static_cast<int>(F(t, 0)),
-          static_cast<int>(F(t, 1)),
-          static_cast<int>(F(t, 2))
-        };
-        for (int e = 0; e < 3; ++e)
-        {
-          const int u = tri[e];
-          const int v = tri[(e + 1) % 3];
-          const int lo = u < v ? u : v;
-          const int hi = u < v ? v : u;
-          const std::uint64_t k = pack(lo, hi);
-          counts[k] += (u < v) ? +1 : -1;
-        }
-      }
-
-      out.reserve(counts.size());
-      for (const auto & kv : counts)
-      {
-        const int sum = kv.second;
-        if (sum == 0) continue;
-        const std::uint32_t lo = static_cast<std::uint32_t>(kv.first >> 32);
-        const std::uint32_t hi = static_cast<std::uint32_t>(kv.first & 0xFFFFFFFFu);
-        const Point a(static_cast<Scalar>(V(lo, 0)),
-                      static_cast<Scalar>(V(lo, 1)),
-                      static_cast<Scalar>(V(lo, 2)));
-        const Point b(static_cast<Scalar>(V(hi, 0)),
-                      static_cast<Scalar>(V(hi, 1)),
-                      static_cast<Scalar>(V(hi, 2)));
-        WeightedSeg ws;
-        if (sum > 0) { ws.a = a; ws.b = b; ws.w = static_cast<Scalar>(sum); }
-        else         { ws.a = b; ws.b = a; ws.w = static_cast<Scalar>(-sum); }
-        out.push_back(ws);
-      }
-    }
-
-    std::vector<WeightedSeg> m_boundary;
-    Direction m_x0;
-    size_t m_face_count = 0;
-  };
-}
-
-#endif

include/igl/WindingNumberTree.h

@@ -11,7 +11,6 @@
 #include <map>
 #include <Eigen/Dense>
 #include "WindingNumberMethod.h"
-#include <cassert>
 #include <memory>
 
 namespace igl

include/igl/active_set.cpp

@@ -14,7 +14,6 @@
 #include "placeholders.h"
 #include "PlainMatrix.h"
 
-#include <cassert>
 #include <iostream>
 #include <limits>
 #include <algorithm>

include/igl/adjacency_matrix.cpp

@@ -9,7 +9,6 @@
 
 #include "verbose.h"
 
-#include <cassert>
 #include <vector>
 
 template <typename DerivedF, typename T>

include/igl/all_pairs_distances.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "all_pairs_distances.h"
 #include <Eigen/Dense>
-#include <cassert>
 
 template <typename Mat>
 IGL_INLINE void igl::all_pairs_distances(

include/igl/angular_distance.cpp

@@ -8,9 +8,6 @@
 #include "angular_distance.h"
 #include "EPS.h"
 #include "PI.h"
-
-#include <cassert>
-
 IGL_INLINE double igl::angular_distance(
   const Eigen::Quaterniond & A,
   const Eigen::Quaterniond & B)

include/igl/arap_linear_block.cpp

@@ -9,7 +9,6 @@
 #include "verbose.h"
 #include "cotmatrix_entries.h"
 #include <Eigen/Dense>
-#include <cassert>
 
 template <typename MatV, typename MatF, typename MatK>
 IGL_INLINE void igl::arap_linear_block(

include/igl/arap_rhs.cpp

@@ -10,7 +10,6 @@
 #include "verbose.h"
 #include "repdiag.h"
 #include "cat.h"
-#include <cassert>
 #include <iostream>
 
 template<typename DerivedV, typename DerivedF, typename DerivedK>

include/igl/average_from_edges_onto_vertices.cpp

@@ -7,8 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "average_from_edges_onto_vertices.h"
 
-#include <cassert>
-
 template<typename DerivedF,typename DerivedE,typename DerivedoE,
 typename DeriveduE,typename DeriveduV>
 IGL_INLINE void

include/igl/barycentric_coordinates.cpp

@@ -8,8 +8,6 @@
 #include "barycentric_coordinates.h"
 #include "volume.h"
 
-#include <cassert>
-
 template <
   typename DerivedP,
   typename DerivedA,

include/igl/biharmonic_coordinates.cpp

@@ -15,7 +15,6 @@
 #include "normal_derivative.h"
 #include "on_boundary.h"
 #include <Eigen/Sparse>
-#include <cassert>
 
 template <
   typename DerivedV,

include/igl/bijective_composite_harmonic_mapping.cpp

@@ -11,7 +11,6 @@
 #include "harmonic.h"
 #include "placeholders.h"
 //#include "matlab/MatlabWorkspace.h"
-#include <cassert>
 #include <iostream>
 
 template <

include/igl/blue_noise.cpp

@@ -17,7 +17,6 @@
 #include <vector>
 #include <random>
 #include <cstdint>
-#include <cassert>
 
 namespace igl
 {

include/igl/boundary_conditions.cpp

@@ -11,7 +11,6 @@
 #include "EPS.h"
 #include "project_to_line.h"
 
-#include <cassert>
 #include <vector>
 #include <map>
 #include <iostream>

include/igl/boundary_facets.cpp

@@ -15,7 +15,6 @@
 
 #include <Eigen/Core>
 
-#include <cassert>
 #include <map>
 #include <iostream>
 

include/igl/bounding_box.cpp

@@ -6,7 +6,6 @@
 // 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 "bounding_box.h"
-#include <cassert>
 #include <iostream>
 
 template <typename DerivedV, typename DerivedBV, typename DerivedBF>

include/igl/cat.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "cat.h"
 
-#include <cassert>
 #include <cstdio>
 
 // Bug in unsupported/Eigen/SparseExtra needs iostream first

include/igl/centroid.cpp

@@ -8,8 +8,6 @@
 #include "centroid.h"
 #include <Eigen/Geometry>
 
-#include <cassert>
-
 template <
   typename DerivedV, 
   typename DerivedF, 

include/igl/circulation.h

@@ -47,7 +47,7 @@ namespace igl
   ///     F(f,:) opposite the vth corner, where EI(e,0)=v. Similarly EF(e,1) "
   ///     e=(j->i)
   /// @param[in] EI  #E by 2 list of edge flap corners (see above).
-  /// @param[out] #vN list of faces touched by circulation (in cyclically order).
+  /// @param[out] #vN list of of faces touched by circulation (in cyclically order).
   ///   
   /// \see edge_flaps
   template <typename DerivedEMAP, typename DerivedEF, typename DerivedEI, typename DerivedvN>

include/igl/circumradius.cpp

@@ -10,8 +10,6 @@
 #include "doublearea.h"
 #include "placeholders.h"
 #include <Eigen/QR>
-#include <cassert>
-
 template <
   typename DerivedV, 
   typename DerivedF,

include/igl/collapse_edge.cpp

@@ -8,7 +8,6 @@
 #include "collapse_edge.h"
 #include "circulation.h"
 #include "edge_collapse_is_valid.h"
-#include <cassert>
 #include <vector>
 
 template <

include/igl/collapse_edge.h

@@ -36,10 +36,10 @@ namespace igl
   ///     e=(j->i)
   /// @param[in,out] EI  #E by 2 list of edge flap corners (see above).
   /// [mesh inputs]
-  /// @param[out] e1  index into E of edge collapsed on left
-  /// @param[out] e2  index into E of edge collapsed on right
-  /// @param[out] f1  index into F of face collapsed on left
-  /// @param[out] f2  index into F of face collapsed on right
+  /// @param[out] e1  index into E of edge collpased on left
+  /// @param[out] e2  index into E of edge collpased on right
+  /// @param[out] f1  index into F of face collpased on left
+  /// @param[out] f2  index into F of face collpased on right
   /// @return true if edge was collapsed
   ///
   ///

include/igl/collapse_least_cost_edge.cpp

@@ -9,8 +9,6 @@
 #include "collapse_edge.h"
 #include "circulation.h"
 
-#include <cassert>
-
 IGL_INLINE bool igl::collapse_least_cost_edge(
   const decimate_cost_and_placement_callback & cost_and_placement,
   const decimate_pre_collapse_callback       & pre_collapse,

include/igl/collapse_least_cost_edge.h

@@ -47,10 +47,10 @@ namespace igl
   /// @param[in] C  #E by dim list of stored placements
   /// @param[out] e  index into E of attempted collapsed edge. Set to -1 if Q is empty or
   ///               contains only infinite cost edges.
-  /// @param[out] e1  index into E of edge collapsed on left.
-  /// @param[out] e2  index into E of edge collapsed on right.
-  /// @param[out] f1  index into F of face collapsed on left.
-  /// @param[out] f2  index into F of face collapsed on right.
+  /// @param[out] e1  index into E of edge collpased on left.
+  /// @param[out] e2  index into E of edge collpased on right.
+  /// @param[out] f1  index into F of face collpased on left.
+  /// @param[out] f2  index into F of face collpased on right.
   ///
   /// \bug This function is not templated nicely and refactoring it and its
   /// dependencies to do so is non-trivial, see

include/igl/collapse_small_triangles.cpp

@@ -13,8 +13,6 @@
 #include "colon.h"
 #include "faces_first.h"
 
-#include <cassert>
-
 #include <limits>
 
 #include <iostream>

include/igl/comb_cross_field.cpp

@@ -18,8 +18,6 @@
 
 #include "triangle_triangle_adjacency.h"
 
-#include <cassert>
-
 namespace igl {
   template <typename DerivedV, typename DerivedF>
   class Comb

include/igl/comb_frame_field.cpp

@@ -16,8 +16,6 @@
 #include "PI.h"
 #include "PlainMatrix.h"
 
-#include <cassert>
-
 template <typename DerivedV, typename DerivedF, typename DerivedP>
 IGL_INLINE void igl::comb_frame_field(
   const Eigen::MatrixBase<DerivedV> &V,

include/igl/connected_components.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 
 #include "connected_components.h"
-#include <cassert>
 #include <queue>
 
 template < typename Atype, typename DerivedC, typename DerivedK>

include/igl/copyleft/cgal/BinaryWindingNumberOperations.h

@@ -13,7 +13,6 @@
 #include "../../igl_inline.h"
 #include "../../MeshBooleanType.h"
 #include <Eigen/Core>
-#include <cassert>
 
 // TODO: This is not written according to libigl style. These should be
 // function handles.

include/igl/copyleft/cgal/extract_cells_single_component.cpp

@@ -15,7 +15,6 @@
 #include <set>
 #include <map>
 #include <queue>
-#include <cassert>
 
 //#define EXTRACT_CELLS_SINGLE_COMPONENT_TIMING
 

include/igl/copyleft/cgal/is_self_intersecting.cpp

@@ -7,8 +7,6 @@
 #include "remesh_self_intersections.h"
 #include "../../collapse_edge.h"
 
-#include <cassert>
-
 template <
   typename DerivedV,
   typename DerivedF>

include/igl/copyleft/cgal/mesh_boolean_type_to_funcs.cpp

@@ -1,8 +1,6 @@
 #include "mesh_boolean_type_to_funcs.h"
 #include "BinaryWindingNumberOperations.h"
 
-#include <cassert>
-
 IGL_INLINE void igl::copyleft::cgal::mesh_boolean_type_to_funcs(
   const MeshBooleanType & type,
   std::function<int(const Eigen::Matrix<int,1,Eigen::Dynamic>) >& wind_num_op,

include/igl/copyleft/cgal/oriented_bounding_box.cpp

@@ -1,5 +1,4 @@
 #include "oriented_bounding_box.h"
-#include <boost/range/has_range_iterator.hpp>
 #include <CGAL/Aff_transformation_3.h>
 //#include <CGAL/Simple_cartesian.h>
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>

include/igl/copyleft/cgal/outer_edge.cpp

@@ -8,7 +8,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "outer_edge.h"
 #include "outer_vertex.h"
-#include <cassert>
 #include <iostream>
 #include <vector>
 

include/igl/copyleft/cgal/outer_vertex.cpp

@@ -8,7 +8,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "outer_vertex.h"
 #include "outer_edge.h"
-#include <cassert>
 #include <iostream>
 #include <vector>
 

include/igl/copyleft/cgal/peel_outer_hull_layers.cpp

@@ -10,7 +10,6 @@
 #include "../../LinSpaced.h"
 #include <vector>
 #include <iostream>
-#include <cassert>
 //#define IGL_PEEL_OUTER_HULL_LAYERS_DEBUG
 #ifdef IGL_PEEL_OUTER_HULL_LAYERS_DEBUG
 #include "../../writePLY.h"

include/igl/copyleft/cgal/trim_with_solid.h

@@ -29,11 +29,11 @@ namespace igl
         /// output and strictly fewer FLOPS than CHECK_EACH_FACE. There will be
         /// many tiny patches along the intersection of A and B. 
         CHECK_EACH_PATCH = 2,
-        /// Merge A and B into the same mesh and resolve all self-intersections.
+        /// Merge A and B into the same mesh and resolve all self-itnersections.
         /// Then "undo" remeshing on faces of A not involved in intersections
-        /// with B (i.e., self-intersections in A). Then separate into patches
+        /// with B (i.e., self-intersections in A). Then seperate into patches
         /// based on connected faces --- where connected means sharing a
-        /// _manifold edge_ --- then label each patch as inside or outside.
+        /// _manifold edge_ --- then label each aptch as inside or outside.
         /// Results in fewer patches than CHECK_EACH_PATCH but finding, meshing and
         /// mesh-undoing self-intersections in A can be costly. This
         /// could result in different output from CHECK_EACH_PATCH because of

include/igl/copyleft/cgal/wire_mesh.cpp

@@ -7,7 +7,6 @@
 #include "coplanar.h"
 #include "mesh_boolean.h"
 #include <Eigen/Geometry>
-#include <cassert>
 #include <vector>
 
 template <

include/igl/copyleft/marching_cubes.cpp

@@ -26,7 +26,6 @@
 #include "marching_cubes.h"
 #include "marching_cubes_tables.h"
 
-#include <cassert>
 #include <unordered_map>
 
 

include/igl/copyleft/tetgen/cdt.cpp

@@ -10,8 +10,6 @@
 #include "../../PlainMatrix.h"
 #include "tetrahedralize.h"
 
-#include <cassert>
-
 template <
   typename DerivedV, 
   typename DerivedF, 

include/igl/copyleft/tetgen/mesh_with_skeleton.cpp

@@ -11,7 +11,6 @@
 #include "../../sample_edges.h"
 #include "../../cat.h"
 
-#include <cassert>
 #include <iostream>
 // Default settings pq2Y tell tetgen to mesh interior of triangle mesh and
 // to produce a graded tet mesh

include/igl/copyleft/tetgen/tetgenio_to_tetmesh.cpp

@@ -11,7 +11,6 @@
 #include "../../list_to_matrix.h"
 
 // STL includes
-#include <cassert>
 #include <iostream>
 
 template <

include/igl/cotmatrix.cpp

@@ -15,8 +15,6 @@
 // Bug in unsupported/Eigen/SparseExtra needs iostream first
 #include <iostream>
 
-#include <cassert>
-
 template <typename DerivedV, typename DerivedF, typename Scalar>
 IGL_INLINE void igl::cotmatrix(
   const Eigen::MatrixBase<DerivedV> & V, 

include/igl/cotmatrix_entries.cpp

@@ -17,7 +17,6 @@
 
 #include "verbose.h"
 
-#include <cassert>
 
 template <typename DerivedV, typename DerivedF, typename DerivedC>
 IGL_INLINE void igl::cotmatrix_entries(

include/igl/covariance_scatter_matrix.cpp

@@ -13,7 +13,6 @@
 #include "verbose.h"
 #include "cat.h"
 #include "PI.h"
-#include <cassert>
 
 template <
   typename DerivedV, 

include/igl/cr_vector_curvature_correction.cpp

@@ -18,7 +18,6 @@
 
 #include "PI.h"
 
-#include <cassert>
 
 template <typename DerivedV, typename DerivedF, typename DerivedE,
 typename DerivedOE, typename ScalarK>

include/igl/cr_vector_laplacian.cpp

@@ -14,7 +14,6 @@
 #include "doublearea.h"
 #include "squared_edge_lengths.h"
 
-#include <cassert>
 
 template <typename DerivedV, typename DerivedF, typename DerivedE,
 typename DerivedOE, typename ScalarL>

include/igl/cr_vector_mass.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "cr_vector_mass.h"
 
-#include <cassert>
 #include <vector>
 
 #include "orient_halfedges.h"

include/igl/cross.cpp

@@ -6,7 +6,6 @@
 // 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 "cross.h"
-#include <cassert>
 
 // http://www.antisphere.com/Wiki/tools:anttweakbar
 IGL_INLINE void igl::cross(

include/igl/crouzeix_raviart_cotmatrix.cpp

@@ -10,7 +10,6 @@
 #include "oriented_facets.h"
 #include "is_edge_manifold.h"
 #include "cotmatrix_entries.h"
-#include <cassert>
 
 template <typename DerivedV, typename DerivedF, typename LT, typename DerivedE, typename DerivedEMAP>
 void igl::crouzeix_raviart_cotmatrix(

include/igl/cubic.cpp

@@ -1,6 +1,4 @@
 #include "cubic.h"
-#include <cassert>
-
   template
     <
     typename  DerivedC,

include/igl/cubic_monomial_bases.cpp

@@ -1,5 +1,4 @@
 #include "cubic_monomial_bases.h"
-#include <cassert>
 
 template <
   typename DerivedC,

include/igl/cycodebase/box_cubic.cpp

@@ -8,7 +8,6 @@
 #include "box_cubic.h"
 #include "../cubic.h"
 #include "../parallel_for.h"
-#include "../placeholders.h"
 #include <cyPolynomial.h>
 
 template < 
@@ -24,8 +23,8 @@ IGL_INLINE void igl::cycodebase::box_cubic(
   typedef Eigen::Matrix<Scalar,DerivedC::RowsAtCompileTime,1> VectorSC;
 
   // Using the control points is a simple (but not tight) bound
-  B1 = C({0,3},igl::placeholders::all).colwise().minCoeff();
-  B2 = C({0,3},igl::placeholders::all).colwise().maxCoeff();
+  B1 = C({0,3},Eigen::all).colwise().minCoeff();
+  B2 = C({0,3},Eigen::all).colwise().maxCoeff();
 
   // Better to find each of the t values where dC/dt = 0 and evaluate C there
   for(int d = 0;d<C.cols();d++)
@@ -68,7 +67,7 @@ IGL_INLINE void igl::cycodebase::box_cubic(
   {
     RowVectorP B1_c, B2_c;
     // Eval copies, but is it really good to make a template for the non copy?
-    box_cubic(P(C.row(c),igl::placeholders::all).eval(),B1_c,B2_c);
+    box_cubic(P(C.row(c),Eigen::all).eval(),B1_c,B2_c);
     B1.row(c) = B1_c;
     B2.row(c) = B2_c;
   },1000);

include/igl/cycodebase/point_cubic_squared_distance.cpp

@@ -74,6 +74,7 @@ void igl::cycodebase::point_cubic_squared_distance(
   S.resize(Q.rows());
   K.resize(Q.rows(),dim);
   const int n = Q.rows();
+  for(int i = 0;i<n;i++)
   igl::parallel_for(n,[&](const int i)
   {
     RowVectorSD k_i;

include/igl/cycodebase/point_spline_squared_distance.cpp

@@ -4,7 +4,6 @@
 #include "../eytzinger_aabb_sdf.h"
 #include "../parallel_for.h"
 #include "../cubic_monomial_bases.h"
-#include "../placeholders.h"
 
 template <
   typename DerivedQ,
@@ -70,7 +69,7 @@ void igl::cycodebase::point_spline_squared_distance(
     MatrixS4D M_unused;
     for(int j = 0; j < C.rows(); j++)
     {
-      C_vec[j] = P(C.row(j), igl::placeholders::all);
+      C_vec[j] = P(C.row(j), Eigen::all);
       cubic_monomial_bases( C_vec[j], M_unused, D_vec[j], B_vec[j]);
     }
   }
@@ -111,8 +110,6 @@ void igl::cycodebase::point_spline_squared_distance(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
-template void igl::cycodebase::point_spline_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::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<double, -1, -1, 0, -1, -1>>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<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>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&);
-// generated by autoexplicit.sh
 template void igl::cycodebase::point_spline_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, 2, 1, -1, 2>, 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<double, -1, -1, 0, -1, -1>>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1>> const&, 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<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>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&);
 // generated by autoexplicit.sh
 template void igl::cycodebase::point_spline_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, 1, -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<double, -1, -1, 0, -1, -1>>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -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<double, -1, 1, 0, -1, 1>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&);

include/igl/cycodebase/spline_eytzinger_aabb.cpp

@@ -24,8 +24,6 @@ void igl::cycodebase::spline_eytzinger_aabb(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
-template void igl::cycodebase::spline_eytzinger_aabb<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::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<double, -1, -1, 0, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);
-// generated by autoexplicit.sh
 template void igl::cycodebase::spline_eytzinger_aabb<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, 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, 2, 1, -1, 2>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 1, -1, 2>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);
 // generated by autoexplicit.sh
 template void igl::cycodebase::spline_eytzinger_aabb<Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, -1, 1, -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, 1, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);

include/igl/decimate_callback_types.h

@@ -64,10 +64,10 @@ namespace igl
   /// @param[in] C  #E by dim list of stored placements
   /// @param[in] e  index into E of attempted collapsed edge. Set to -1 if Q is empty or
   ///               contains only infinite cost edges.
-  /// @param[in] e1  index into E of edge collapsed on left.
-  /// @param[in] e2  index into E of edge collapsed on right.
-  /// @param[in] f1  index into F of face collapsed on left.
-  /// @param[in] f2  index into F of face collapsed on right.
+  /// @param[in] e1  index into E of edge collpased on left.
+  /// @param[in] e2  index into E of edge collpased on right.
+  /// @param[in] f1  index into F of face collpased on left.
+  /// @param[in] f2  index into F of face collpased on right.
   /// @return whether to stop
   using decimate_stopping_condition_callback = 
     std::function<bool(
@@ -139,10 +139,10 @@ namespace igl
   /// @param[in] C  #E by dim list of stored placements
   /// @param[in] e  index into E of attempted collapsed edge. Set to -1 if Q is empty or
   ///               contains only infinite cost edges.
-  /// @param[in] e1  index into E of edge collapsed on left.
-  /// @param[in] e2  index into E of edge collapsed on right.
-  /// @param[in] f1  index into F of face collapsed on left.
-  /// @param[in] f2  index into F of face collapsed on right.
+  /// @param[in] e1  index into E of edge collpased on left.
+  /// @param[in] e2  index into E of edge collpased on right.
+  /// @param[in] f1  index into F of face collpased on left.
+  /// @param[in] f2  index into F of face collpased on right.
   /// @param[in] collapsed whether collapse actual took place
   using decimate_post_collapse_callback = 
     std::function<void(

include/igl/deform_skeleton.cpp

@@ -6,8 +6,6 @@
 // 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 "deform_skeleton.h"
-#include <cassert>
-
 void igl::deform_skeleton(
   const Eigen::MatrixXd & C,
   const Eigen::MatrixXi & BE,

include/igl/delaunay_triangulation.cpp

@@ -12,7 +12,6 @@
 #include "unique_edge_map.h"
 #include "is_delaunay.h"
 
-#include <cassert>
 #include <vector>
 #include <sstream>
 

include/igl/direct_delta_mush.cpp

@@ -9,8 +9,6 @@
 #include "cotmatrix.h"
 #include "PlainMatrix.h"
 
-#include <cassert>
-
 template <
   typename DerivedV,
   typename DerivedOmega,

include/igl/edge_collapse_is_valid.cpp

@@ -11,7 +11,6 @@
 #include "intersect.h"
 #include "unique.h"
 #include "list_to_matrix.h"
-#include <cassert>
 #include <vector>
 
 template <

include/igl/edge_exists_near.cpp

@@ -6,7 +6,6 @@
 // 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_exists_near.h"
-#include <cassert>
 
 template <
   typename DeriveduE,

include/igl/edge_midpoints.cpp

@@ -7,8 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "edge_midpoints.h"
 
-#include <cassert>
-
 template<typename DerivedV,typename DerivedF,typename DerivedE,
 typename DerivedoE, typename Derivedmps>
 IGL_INLINE void

include/igl/edge_topology.cpp

@@ -8,7 +8,6 @@
 #include "edge_topology.h"
 #include "is_edge_manifold.h"
 #include <algorithm>
-#include <cassert>
 
 template<typename DerivedV, typename DerivedF, typename DerivedE>
 IGL_INLINE void igl::edge_topology(

include/igl/edge_vectors.cpp

@@ -13,7 +13,6 @@
 
 #include "PI.h"
 
-#include <cassert>
 
 template<typename DerivedV,typename DerivedF,typename DerivedE,
 typename DerivedoE, typename Derivedvec>

include/igl/edges.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 #include "edges.h"
 #include "adjacency_matrix.h"
-#include <cassert>
 #include <iostream>
 
 template <typename DerivedF, typename DerivedE>

include/igl/edges_to_path.cpp

@@ -5,7 +5,6 @@
 #include "unique.h"
 #include "adjacency_list.h"
 #include "PlainMatrix.h"
-#include <cassert>
 
 template <
   typename DerivedE,

include/igl/eigs.cpp

@@ -12,7 +12,6 @@
 #include "massmatrix.h"
 #include "placeholders.h"
 #include <iostream>
-#include <cassert>
 
 template <
   typename Atype,

include/igl/embree/EmbreeDevice.h

@@ -12,9 +12,7 @@
 #ifndef IGL_EMBREE_EMBREE_DEVICE_H
 #define IGL_EMBREE_EMBREE_DEVICE_H
 #include <embree4/rtcore.h>
-#include <cassert>
 #include <iostream>
-#include <cassert>
 
 namespace igl
 {

include/igl/embree/EmbreeIntersector.cpp

@@ -386,58 +386,3 @@ igl::embree::EmbreeIntersector
   ray.hit.instID[0] = RTC_INVALID_GEOMETRY_ID;
   ray.hit.primID = RTC_INVALID_GEOMETRY_ID;
 }
-
-IGL_INLINE int
-igl::embree::EmbreeIntersector
-::signedIntersectionsRay(
-  const OriginType    & origin,
-  const DirectionType & direction,
-  float tnear,
-  float tfar,
-  int   mask) const
-{
-  struct query_context
-  {
-    RTCRayQueryContext base; // MUST be first for the reinterpret_cast in the filter
-    int sum;
-  };
-
-  query_context q;
-  rtcInitRayQueryContext(&q.base);
-  q.sum = 0;
-
-  RTCRay ray{};
-  ray.org_x = origin[0];
-  ray.org_y = origin[1];
-  ray.org_z = origin[2];
-  ray.dir_x = direction[0];
-  ray.dir_y = direction[1];
-  ray.dir_z = direction[2];
-  ray.tnear = tnear;
-  ray.tfar  = tfar;
-  ray.mask  = static_cast<unsigned int>(mask);
-  ray.flags = 0;
-
-  RTCOccludedArguments rargs;
-  rtcInitOccludedArguments(&rargs);
-  rargs.flags = (RTCRayQueryFlags)(
-      RTC_RAY_QUERY_FLAG_COHERENT | RTC_RAY_QUERY_FLAG_INVOKE_ARGUMENT_FILTER);
-  rargs.feature_mask = RTC_FEATURE_FLAG_ALL;
-  rargs.context = &q.base;
-  rargs.filter = +[](RTCFilterFunctionNArguments const* fargs)
-  {
-    assert(fargs->N == 1 && fargs->valid[0]);
-    auto const& fray = reinterpret_cast<RTCRay&>(*fargs->ray);
-    auto const& fhit = reinterpret_cast<RTCHit&>(*fargs->hit);
-    auto const d = fray.dir_x * fhit.Ng_x
-                 + fray.dir_y * fhit.Ng_y
-                 + fray.dir_z * fhit.Ng_z;
-    reinterpret_cast<query_context*>(fargs->context)->sum += d > 0.0f ? +1 : -1;
-    // Reject the hit so traversal continues and we visit every crossing.
-    fargs->valid[0] = 0;
-  };
-  rargs.occluded = nullptr;
-
-  rtcOccluded1(scene, &ray, &rargs);
-  return q.sum;
-}

include/igl/embree/EmbreeIntersector.h

@@ -37,8 +37,6 @@ namespace igl
     public:
       typedef Eigen::Matrix<float,Eigen::Dynamic,3> PointMatrixType;
       typedef Eigen::Matrix<int,Eigen::Dynamic,3> FaceMatrixType;
-      typedef Eigen::RowVector3f OriginType;
-      typedef Eigen::RowVector3f DirectionType;
     public:
       EmbreeIntersector();
     private:
@@ -152,23 +150,6 @@ namespace igl
         Hit<float> &hit,
         int mask = 0xFFFFFFFF) const;
 
-      /// Signed ray-mesh crossing count along `(origin, direction)`. Visits
-      /// every hit (occluded ray + argument filter rejecting each) and
-      /// accumulates `sign(direction · Ng)` per hit.
-      ///
-      /// @param[in] origin     ray origin
-      /// @param[in] direction  ray direction (need not be normalized)
-      /// @param[in] tnear      start of ray segment
-      /// @param[in] tfar       end of ray segment
-      /// @param[in] mask       a 32 bit mask to identify active geometries
-      /// @return signed crossing count
-      int signedIntersectionsRay(
-        const OriginType   & origin,
-        const DirectionType& direction,
-        float tnear = 0.0f,
-        float tfar  = std::numeric_limits<float>::infinity(),
-        int   mask  = 0xFFFFFFFF) const;
-
     private:
 
       struct Vertex   {float x,y,z,a;};

include/igl/embree/winding_number_antipodal.cpp

@@ -1,29 +0,0 @@
-// This file is part of libigl, a simple c++ geometry processing library.
-//
-// Copyright (C) 2026 Philip Trettner <trettner@shapedcode.com>, Cedric Martens <cedric.martens@umontreal.ca>
-//
-// 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/.
-//
-// Embree-target satellite that emits explicit template instantiations of
-// `igl::winding_number_antipodal` against `igl::embree::EmbreeIntersector`.
-// igl_core does not link embree, so these instantiations cannot live in
-// include/igl/winding_number_antipodal.cpp. In header-only mode this file
-// is empty after preprocessing.
-#include "../winding_number_antipodal.h"
-
-#ifdef IGL_STATIC_LIBRARY
-// Pull in the function body so we can instantiate it here. In static-lib
-// mode the header does not transitively include the body.
-#include "../winding_number_antipodal.cpp"
-#include "EmbreeIntersector.h"
-
-// Mirror of the (V, F, O, W) instantiations in include/igl/winding_number.cpp,
-// adapted to the antipodal API (3D triangle meshes only, with EmbreeIntersector
-// as the intersector).
-template void igl::winding_number_antipodal<igl::embree::EmbreeIntersector, 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::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::embree::EmbreeIntersector const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
-template void igl::winding_number_antipodal<igl::embree::EmbreeIntersector, Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 3, 0, -1, 3>, Eigen::Matrix<double, -1, 1, 0, -1, 1> >(Eigen::MatrixBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::embree::EmbreeIntersector const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 1, 0, -1, 1> >&);
-template void igl::winding_number_antipodal<igl::embree::EmbreeIntersector, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, 1, 0, -1, 1> >(Eigen::MatrixBase<Eigen::Matrix<float, -1, -1, 0, -1, -1> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::embree::EmbreeIntersector const&, Eigen::MatrixBase<Eigen::Matrix<float, -1, -1, 0, -1, -1> > const&, Eigen::PlainObjectBase<Eigen::Matrix<float, -1, 1, 0, -1, 1> >&);
-template void igl::winding_number_antipodal<igl::embree::EmbreeIntersector, Eigen::Matrix<float, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<float, -1, 3, 0, -1, 3>, Eigen::Matrix<float, -1, 1, 0, -1, 1> >(Eigen::MatrixBase<Eigen::Matrix<float, -1, 3, 0, -1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, -1, 0, -1, -1> > const&, igl::embree::EmbreeIntersector const&, Eigen::MatrixBase<Eigen::Matrix<float, -1, 3, 0, -1, 3> > const&, Eigen::PlainObjectBase<Eigen::Matrix<float, -1, 1, 0, -1, 1> >&);
-#endif

include/igl/exploded_view.cpp

@@ -8,7 +8,6 @@
 #include "exploded_view.h"
 #include "barycenter.h"
 #include "volume.h"
-#include <cassert>
 
 template <
   typename DerivedV,

include/igl/eytzinger_aabb.cpp

@@ -1,10 +1,8 @@
 #include "eytzinger_aabb.h"
 #include "PlainMatrix.h"
 #include "median.h"
-#include "placeholders.h"
 #include <vector>
 #include <algorithm>
-#include <cassert>
 
 template <
   typename DerivedPB,
@@ -43,8 +41,8 @@ IGL_INLINE void igl::eytzinger_aabb(
   {
     if(I.size() == 0) { return; }
 
-    B1.row(i) = PB1(I, igl::placeholders::all).colwise().minCoeff();
-    B2.row(i) = PB2(I, igl::placeholders::all).colwise().maxCoeff();
+    B1.row(i) = PB1(I, Eigen::all).colwise().minCoeff();
+    B2.row(i) = PB2(I, Eigen::all).colwise().maxCoeff();
 
     if(I.size() == 1) { leaf(i) = I[0]; return; }
 
@@ -87,8 +85,6 @@ IGL_INLINE void igl::eytzinger_aabb(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
-template void igl::eytzinger_aabb<Eigen::Matrix<double, -1, -1, 1, -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, 1, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);
-// generated by autoexplicit.sh
 template void igl::eytzinger_aabb<Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 1, -1, 2>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, 2, 1, -1, 2>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);
 // generated by autoexplicit.sh
 template void igl::eytzinger_aabb<Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<double, -1, -1, 1, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>> const&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<double, -1, -1, 1, -1, -1>>&, Eigen::PlainObjectBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>>&);

include/igl/eytzinger_aabb_sdf.cpp

@@ -103,8 +103,6 @@ IGL_INLINE void igl::eytzinger_aabb_sdf(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
-template void igl::eytzinger_aabb_sdf<true, Eigen::Matrix<double, 1, -1, 1, 1, -1>, std::function<double (int)>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, 1, -1, 1, 1, -1>> const&, std::function<double (int)> 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>> const&, Eigen::Matrix<double, 1, -1, 1, 1, -1>::Scalar&);
-// generated by autoexplicit.sh
 template void igl::eytzinger_aabb_sdf<true, Eigen::Matrix<double, 1, -1, 1, 1, -1>, std::function<double (int)>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, 1, -1, 1, 1, -1>> const&, std::function<double (int)> const&, 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<int, -1, 1, 0, -1, 1>> const&, Eigen::Matrix<double, 1, -1, 1, 1, -1>::Scalar&);
 // generated by autoexplicit.sh
 template void igl::eytzinger_aabb_sdf<false, Eigen::Matrix<double, 1, 2, 1, 1, 2>, std::function<double (int)>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, Eigen::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, 1, 2, 1, 1, 2>> const&, std::function<double (int)> const&, 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<int, -1, 1, 0, -1, 1>> const&, Eigen::Matrix<double, 1, 2, 1, 1, 2>::Scalar&);

include/igl/eytzinger_aabb_winding_number.cpp

@@ -1,6 +1,5 @@
 #include "eytzinger_aabb_winding_number.h"
 #include "PI.h"
-#include <cassert>
 
 // signed_angle helper. Not sure why igl::signed_angle is doing all that
 // normalization (divisions and square roots) when atan2 doesn't need it.
@@ -89,7 +88,7 @@ IGL_INLINE void igl::eytzinger_aabb_winding_number(
         (p(1) < B1(r,1)) || (p(1) > B2(r,1)) )
     {
       // consider the pairs in the sequence at I(C(r)) to I(C(r)+1)
-      assert(((C(r+1) - C(r)) % 2) == 0);
+      assert(C(r+1) - C(r) % 2 == 0);
       for(int i = C(r); i < C(r+1); i+=2)
       {
         wn += ::signed_angle(V,p,I(i),I(i+1)) / Scalar(2.0 * igl::PI);
@@ -107,8 +106,6 @@ IGL_INLINE void igl::eytzinger_aabb_winding_number(
 #ifdef IGL_STATIC_LIBRARY
 // Explicit template instantiation
 // generated by autoexplicit.sh
-template void igl::eytzinger_aabb_winding_number<Eigen::Matrix<double, 1, -1, 1, 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::Matrix<int, -1, 1, 0, -1, 1>>(Eigen::MatrixBase<Eigen::Matrix<double, 1, -1, 1, 1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, std::function<Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar (int)> 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>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, 1, 0, -1, 1>> const&, Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar&);
-// generated by autoexplicit.sh
 template void igl::eytzinger_aabb_winding_number<Eigen::Matrix<double, 1, -1, 1, 1, -1>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, 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, 1, 1, -1>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, -1, 0, -1, -1>> const&, std::function<Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar (int)> const&, 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<int, -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&, Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar&);
 template void igl::eytzinger_aabb_winding_number<Eigen::Matrix<double, 1, 2, 1, 1, 2>, Eigen::Matrix<double, -1, -1, 0, -1, -1>, Eigen::Matrix<int, -1, -1, 0, -1, -1>, Eigen::Matrix<double, -1, 2, 1, -1, 2>, 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, 2, 1, 1, 2>> const&, 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, 2, 1, -1, 2>> const&, Eigen::MatrixBase<Eigen::Matrix<double, -1, 2, 1, -1, 2>> const&, Eigen::MatrixBase<Eigen::Matrix<int, -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&, Eigen::Matrix<double, -1, -1, 0, -1, -1>::Scalar&);
 #endif

include/igl/eytzinger_aabb_winding_number_tree.cpp

@@ -2,7 +2,6 @@
 #include "matlab_format.h"
 
 #include <Eigen/Sparse>
-#include <cassert>
 #include <iostream>
 #include <vector>
 #include <queue>
@@ -100,7 +99,6 @@ IGL_INLINE void igl::eytzinger_aabb_winding_number_tree(
   C(0) = 0;
   for(int r = 0;r<leaf.size();r++)
   {
-    assert(vI[r].size() % 2 == 0);
     C(r+1) = C(r) + vI[r].size();
   }
   I.resize(C(leaf.size()));

include/igl/face_areas.cpp

@@ -8,7 +8,6 @@
 #include "face_areas.h"
 #include "edge_lengths.h"
 #include "doublearea.h"
-#include <cassert>
 
 template <typename DerivedV, typename DerivedT, typename DerivedA>
 IGL_INLINE void igl::face_areas(

include/igl/faces_first.cpp

@@ -10,8 +10,6 @@
 #include <vector>
 #include <Eigen/Dense>
 
-#include <cassert>
-
 template <typename MatV, typename MatF, typename VecI>
 IGL_INLINE void igl::faces_first(
   const MatV & V, 

include/igl/facet_adjacency_matrix.cpp

@@ -8,8 +8,6 @@
 #include "facet_adjacency_matrix.h"
 #include "unique_edge_map.h"
 
-#include <cassert>
-
 template <typename DerivedF, typename Atype>
 IGL_INLINE void igl::facet_adjacency_matrix(
   const Eigen::MatrixBase<DerivedF> & F, Eigen::SparseMatrix<Atype> & A)

include/igl/facet_components.cpp

@@ -10,7 +10,6 @@
 #include "triangle_triangle_adjacency.h"
 #include "facet_adjacency_matrix.h"
 #include "connected_components.h"
-#include <cassert>
 #include <vector>
 #include <queue>
 

include/igl/fast_winding_number.cpp

@@ -473,4 +473,5 @@ template void igl::fast_winding_number<Eigen::Matrix<float, -1, -1, 0, -1, -1>,
 template Eigen::Matrix<float, 1, 3, 1, 1, 3>::Scalar igl::fast_winding_number<Eigen::Matrix<float, 1, 3, 1, 1, 3> >(igl::FastWindingNumberBVH const&, float, Eigen::MatrixBase<Eigen::Matrix<float, 1, 3, 1, 1, 3> > const&);
 template void igl::fast_winding_number<Eigen::Matrix<float, -1, 3, 0, -1, 3>, Eigen::Matrix<int, -1, 3, 0, -1, 3> >(Eigen::MatrixBase<Eigen::Matrix<float, -1, 3, 0, -1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, 3, 0, -1, 3> > const&, int, igl::FastWindingNumberBVH&);
 template void igl::fast_winding_number<Eigen::Matrix<float, -1, 3, 1, -1, 3>, Eigen::Matrix<int, -1, 3, 1, -1, 3> >(Eigen::MatrixBase<Eigen::Matrix<float, -1, 3, 1, -1, 3> > const&, Eigen::MatrixBase<Eigen::Matrix<int, -1, 3, 1, -1, 3> > const&, int, igl::FastWindingNumberBVH&);
+template Eigen::CwiseUnaryOp<Eigen::internal::scalar_cast_op<double, float>, Eigen::Matrix<double, 1, 3, 1, 1, 3> const>::Scalar igl::fast_winding_number<Eigen::CwiseUnaryOp<Eigen::internal::scalar_cast_op<double, float>, Eigen::Matrix<double, 1, 3, 1, 1, 3> const> >(igl::FastWindingNumberBVH const&, float, Eigen::MatrixBase<Eigen::CwiseUnaryOp<Eigen::internal::scalar_cast_op<double, float>, Eigen::Matrix<double, 1, 3, 1, 1, 3> const> > const&);
 #endif

include/igl/find_cross_field_singularities.cpp

@@ -14,7 +14,6 @@
 #include "vertex_triangle_adjacency.h"
 #include "is_border_vertex.h"
 
-#include <cassert>
 
 template <typename DerivedV, typename DerivedF, typename DerivedM, typename DerivedO>
 IGL_INLINE void igl::find_cross_field_singularities(const Eigen::MatrixBase<DerivedV> &V,

include/igl/find_zero.cpp

@@ -2,8 +2,6 @@
 #include "for_each.h"
 #include "any.h"
 
-#include <cassert>
-
 template <typename AType, typename DerivedI>
 IGL_INLINE void igl::find_zero(
   const Eigen::SparseMatrix<AType> & A,

include/igl/fit_plane.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 
 #include "fit_plane.h"
-#include <cassert>
 #include <iostream>
 
 IGL_INLINE void igl::fit_plane(

include/igl/fit_rotations.cpp

@@ -12,7 +12,6 @@
 #include "polar_dec.h"
 #include "polar_svd.h"
 #include "C_STR.h"
-#include <cassert>
 #include <iostream>
 
 template <typename DerivedS, typename DerivedD>

include/igl/flip_avoiding_line_search.cpp

@@ -10,7 +10,6 @@
 #include "PI.h"
 
 #include <Eigen/Dense>
-#include <cassert>
 #include <vector>
 
 namespace igl

include/igl/flip_edge.cpp

@@ -7,7 +7,6 @@
 // obtain one at http://mozilla.org/MPL/2.0/.
 
 #include "flip_edge.h"
-#include <cassert>
 
 template <
   typename DerivedF,

include/igl/flipped_triangles.cpp

@@ -9,8 +9,6 @@
 
 #include "list_to_matrix.h"
 #include <vector>
-#include <cassert>
-
 template <typename DerivedV, typename DerivedF, typename DerivedX>
 IGL_INLINE void igl::flipped_triangles(
   const Eigen::MatrixBase<DerivedV> & V,

include/igl/forward_kinematics.cpp

@@ -8,7 +8,6 @@
 #include "forward_kinematics.h"
 #include <functional>
 #include <iostream>
-#include <cassert>
 
 IGL_INLINE void igl::forward_kinematics(
   const Eigen::MatrixXd & C,

include/igl/frame_field_deformer.cpp

@@ -9,7 +9,6 @@
 
 #include <Eigen/Dense>
 #include <Eigen/Sparse>
-#include <cassert>
 #include <vector>
 
 #include "cotmatrix_entries.h"

include/igl/gaussian_curvature.cpp

@@ -8,7 +8,6 @@
 #include "gaussian_curvature.h"
 #include "internal_angles.h"
 #include "PI.h"
-#include <cassert>
 #include <iostream>
 template <typename DerivedV, typename DerivedF, typename DerivedK>
 IGL_INLINE void igl::gaussian_curvature(