Nonlinear CCD (#72)

* Implement nonlinear CCD of linear geometry
* Optional filib

Author: zfergus

CMakeLists.txt

@@ -63,12 +63,13 @@ project(IPCToolkit
 
 option(IPC_TOOLKIT_BUILD_TESTS                "Build unit-tests"  ${IPC_TOOLKIT_TOPLEVEL_PROJECT})
 option(IPC_TOOLKIT_BUILD_PYTHON               "Build Python bindings"                         OFF)
-option(IPC_TOOLKIT_WITH_CORRECT_CCD           "Use the TightInclusion CCD"                     ON)
+option(IPC_TOOLKIT_WITH_CORRECT_CCD           "Use Tight Inclusion CCD"                        ON)
 option(IPC_TOOLKIT_WITH_SIMD                  "Enable SIMD"                                   OFF)
 option(IPC_TOOLKIT_WITH_CUDA                  "Enable CUDA CCD"                               OFF)
 option(IPC_TOOLKIT_WITH_RATIONAL_INTERSECTION "Use rational edge-triangle intersection check" OFF)
 option(IPC_TOOLKIT_WITH_ROBIN_MAP             "Use Tessil's robin-map rather than std maps"    ON)
 option(IPC_TOOLKIT_WITH_ABSEIL                "Use Abseil's hash functions"                    ON)
+option(IPC_TOOLKIT_WITH_FILIB                 "Use filib for interval arithmetic"              ON)
 option(IPC_TOOLKIT_WITH_CODE_COVERAGE         "Enable coverage reporting"                     OFF)
 mark_as_advanced(IPC_TOOLKIT_WITH_CODE_COVERAGE)
 
@@ -181,6 +182,12 @@ if(IPC_TOOLKIT_WITH_ABSEIL)
   target_link_libraries(ipc_toolkit PUBLIC absl::hash)
 endif()
 
+# Intervals
+if(IPC_TOOLKIT_WITH_FILIB)
+  include(filib)
+  target_link_libraries(ipc_toolkit PUBLIC filib::filib)
+endif()
+
 # Extra warnings (link last for highest priority)
 include(ipc_toolkit_warnings)
 target_link_libraries(ipc_toolkit PRIVATE ipc::toolkit::warnings)

IPCToolkitOptions.cmake.sample

@@ -30,8 +30,14 @@
 
 # option(IPC_TOOLKIT_BUILD_TESTS                "Build unit-tests"                               ON)
 # option(IPC_TOOLKIT_BUILD_PYTHON               "Build Python bindings"                         OFF)
-# option(IPC_TOOLKIT_WITH_CORRECT_CCD           "Use the TightInclusion CCD"                     ON)
+# option(IPC_TOOLKIT_WITH_CORRECT_CCD           "Use Tight Inclusion CCD"                        ON)
 # option(IPC_TOOLKIT_WITH_SIMD                  "Enable SIMD"                                   OFF)
 # option(IPC_TOOLKIT_WITH_CUDA                  "Enable CUDA CCD"                               OFF)
 # option(IPC_TOOLKIT_WITH_RATIONAL_INTERSECTION "Use rational edge-triangle intersection check" OFF)
-# option(IPC_TOOLKIT_WITH_ROBIN_MAP             "Use Tessil's robin-map rather than std maps"    ON)
+# option(IPC_TOOLKIT_WITH_ROBIN_MAP             "Use Tessil's robin-map rather than std maps"    ON)
+# option(IPC_TOOLKIT_WITH_ABSEIL                "Use Abseil's hash functions"                    ON)
+# option(IPC_TOOLKIT_WITH_FILIB                 "Use filib for interval arithmetic"              ON)
+# option(IPC_TOOLKIT_WITH_CODE_COVERAGE         "Enable coverage reporting"                     OFF)
+# option(IPC_TOOLKIT_TESTS_CCD_BENCHMARK        "Enable CCD benchmark test"                      ON)
+# set(IPC_TOOLKIT_TESTS_CCD_BENCHMARK_DIR     "" CACHE PATH "Path to the CCD benchmark directory")
+# set(IPC_TOOLKIT_TESTS_NEW_CCD_BENCHMARK_DIR "" CACHE PATH "Path to the new CCD benchmark directory")

cmake/ipc_toolkit/ipc_toolkit_warnings.cmake

@@ -13,7 +13,7 @@ set(IPC_TOOLKIT_WARNING_FLAGS
   -pedantic
 
   # -Wconversion
-  #-Wunsafe-loop-optimizations # broken with C++11 loops
+  # -Wunsafe-loop-optimizations # broken with C++11 loops
   -Wunused
 
   -Wno-long-long
@@ -39,14 +39,13 @@ set(IPC_TOOLKIT_WARNING_FLAGS
   -Werror=int-to-pointer-cast
   -Werror=pointer-to-int-cast
 
-  -Wno-unused-variable
+  -Wunused-variable
   -Wunused-but-set-variable
   -Wno-unused-parameter
 
   #-Weffc++
-  -Wno-old-style-cast
-  # -Wno-sign-conversion
-  #-Wsign-conversion
+  -Wold-style-cast
+  # -Wsign-conversion
 
   # -Wshadow
 
@@ -69,7 +68,7 @@ set(IPC_TOOLKIT_WARNING_FLAGS
 
   -Wcast-align
   -Wdisabled-optimization
-  #-Winline # produces warning on default implicit destructor
+  # -Winline # produces warning on default implicit destructor
   -Winvalid-pch
   # -Wmissing-include-dirs
   -Wpacked

cmake/recipes/filib.cmake

@@ -0,0 +1,10 @@
+# filib (https://github.com/zfergus/filib.git)
+# License: LGPL
+if(TARGET filib::filib)
+  return()
+endif()
+
+message(STATUS "Third-party: creating target 'filib::filib'")
+
+include(CPM)
+CPMAddPackage("gh:zfergus/filib#63acb1f78dcdf7667ca2d214f7a3fa3505fdf339")

src/ipc/broad_phase/voxel_size_heuristic.cpp

@@ -6,6 +6,11 @@
 
 namespace ipc {
 
+namespace {
+    // Avoid unused variable warnings
+    inline void check_success(bool success) { assert(success); }
+} // namespace
+
 double suggest_good_voxel_size(
     const Eigen::MatrixXd& vertices,
     const Eigen::MatrixXi& edges,
@@ -129,16 +134,14 @@ double median_edge_length(
     }
 
     double median = -1;
-    const bool success = igl::median(lengths, median);
-    assert(success);
+    check_success(igl::median(lengths, median));
     return median;
 }
 
 double median_displacement_length(const Eigen::MatrixXd& displacements)
 {
     double median = -1;
-    const bool success = igl::median(displacements.rowwise().norm(), median);
-    assert(success);
+    check_success(igl::median(displacements.rowwise().norm(), median));
     return median;
 }
 

src/ipc/ccd/CMakeLists.txt

@@ -7,6 +7,8 @@ set(SOURCES
   ccd.hpp
   inexact_point_edge.cpp
   inexact_point_edge.hpp
+  nonlinear_ccd.cpp
+  nonlinear_ccd.hpp
   point_static_plane.cpp
   point_static_plane.hpp
 )

src/ipc/ccd/ccd.cpp

@@ -9,7 +9,6 @@
 
 #ifdef IPC_TOOLKIT_WITH_CORRECT_CCD
 #include <tight_inclusion/ccd.hpp>
-#include <tight_inclusion/interval_root_finder.hpp>
 #else
 #include <ipc/ccd/inexact_point_edge.hpp>
 #include <CTCD.h>
@@ -20,39 +19,14 @@
 
 namespace ipc {
 
-namespace {
-    inline Eigen::Vector3d to_3D(const VectorMax3d& v)
-    {
-        assert(v.size() == 2 || v.size() == 3);
-        return v.size() == 2 ? Eigen::Vector3d(v.x(), v.y(), 0) : v.head<3>();
-    }
-
-    inline bool check_initial_distance(
-        const double initial_distance, const double min_distance, double& toi)
-    {
-        if (initial_distance > min_distance) {
-            return false;
-        }
-
-        logger().warn(
-            "Initial distance {} ≤ d_min={}, returning toi=0!",
-            initial_distance, min_distance);
-
-        toi = 0; // Initially touching
-
-        return true;
-    }
-} // namespace
-
-/// @brief Scale the distance tolerance to be at most this fraction of the initial distance.
+/// Scale the distance tolerance to be at most this fraction of the initial
+/// distance.
 static constexpr double INITIAL_DISTANCE_TOLERANCE_SCALE = 0.5;
 
-/// @brief Special value for max_iterations to run tight inclusion without a maximum number of iterations.
+/// Special value for max_iterations to run tight inclusion without a maximum
+/// number of iterations.
 static constexpr long TIGHT_INCLUSION_UNLIMITED_ITERATIONS = -1;
 
-/// @brief Tolerance for small time of impact which triggers rerunning CCD without a minimum separation.
-static constexpr double SMALL_TOI = 1e-6;
-
 bool ccd_strategy(
     const std::function<bool(
         long /*max_iterations*/,
@@ -87,15 +61,15 @@ bool ccd_strategy(
     // #ifdef IPC_TOOLKIT_WITH_CORRECT_CCD
     //     // Tight inclusion will have higher accuracy and better performance
     //     // if we shrink the minimum distance. The value 1e-10 is arbitrary.
-    //     while (is_impacting && toi < SMALL_TOI && min_distance > 1e-10) {
+    //     while (is_impacting && toi < CCD_SMALL_TOI && min_distance > 1e-10) {
     //         min_distance /= 10;
     //         is_impacting =
     //             ccd(max_iterations, min_distance, /*no_zero_toi=*/false,
     //             toi);
     //     }
     // #endif
 
-    if (is_impacting && toi < SMALL_TOI) {
+    if (is_impacting && toi < CCD_SMALL_TOI) {
         is_impacting = ccd(
             /*max_iterations=*/TIGHT_INCLUSION_UNLIMITED_ITERATIONS,
             /*min_distance=*/min_distance, /*no_zero_toi=*/true, toi);
@@ -219,7 +193,7 @@ bool point_edge_ccd_3D(
             max_iterations,               // maximum number of iterations
             output_tolerance,             // delta_actual
             no_zero_toi);
-        if (adjusted_tolerance < output_tolerance && toi < SMALL_TOI) {
+        if (adjusted_tolerance < output_tolerance && toi < CCD_SMALL_TOI) {
             logger().trace(
                 "ticcd::edgeEdgeCCD exceeded iteration limit (min_dist={:g} "
                 "max_iterations={:d} input_tol={:g} output_tol={:g} toi={:g})",
@@ -252,10 +226,9 @@ bool point_edge_ccd(
     const long max_iterations,
     const double conservative_rescaling)
 {
-    int dim = p_t0.size();
-    assert(p_t1.size() == dim);
-    assert(e0_t0.size() == dim && e1_t0.size() == dim);
-    assert(e0_t1.size() == dim && e1_t1.size() == dim);
+    assert(p_t1.size() == p_t0.size());
+    assert(e0_t0.size() == p_t0.size() && e1_t0.size() == p_t0.size());
+    assert(e0_t1.size() == p_t0.size() && e1_t1.size() == p_t0.size());
     return point_edge_ccd_3D(
         to_3D(p_t0), to_3D(e0_t0), to_3D(e1_t0), to_3D(p_t1), to_3D(e0_t1),
         to_3D(e1_t1), toi, min_distance, tmax, tolerance, max_iterations,
@@ -305,7 +278,7 @@ bool edge_edge_ccd(
             max_iterations,               // maximum number of iterations
             output_tolerance,             // delta_actual
             no_zero_toi);
-        if (adjusted_tolerance < output_tolerance && toi < SMALL_TOI) {
+        if (adjusted_tolerance < output_tolerance && toi < CCD_SMALL_TOI) {
             logger().trace(
                 "ticcd::edgeEdgeCCD exceeded iteration limit (min_dist={:g} "
                 "max_iterations={:d} input_tol={:g} output_tol={:g} toi={:g})",
@@ -368,7 +341,7 @@ bool point_triangle_ccd(
             max_iterations,               // maximum number of iterations
             output_tolerance,             // delta_actual
             no_zero_toi);
-        if (adjusted_tolerance < output_tolerance && toi < SMALL_TOI) {
+        if (adjusted_tolerance < output_tolerance && toi < CCD_SMALL_TOI) {
             logger().trace(
                 "ticcd::vertexFaceCCD exceeded iteration limit (min_dist={:g} "
                 "max_iterations={:d} input_tol={:g} output_tol={:g} toi={:g})",
@@ -388,4 +361,22 @@ bool point_triangle_ccd(
         conservative_rescaling, toi);
 }
 
+// -----------------------------------------------------------------------------
+
+bool check_initial_distance(
+    const double initial_distance, const double min_distance, double& toi)
+{
+    if (initial_distance > min_distance) {
+        return false;
+    }
+
+    logger().warn(
+        "Initial distance {} ≤ d_min={}, returning toi=0!", initial_distance,
+        min_distance);
+
+    toi = 0; // Initially touching
+
+    return true;
+}
+
 } // namespace ipc