Add SIMD batch support to the distance, barrier, and friction functions - #251
Merged
Conversation
Specialize `Eigen::NumTraits` for `xsimd::batch` and add `ipc::SimdBatch<T>` in the new `ipc/utils/simd.hpp`, so `Eigen::Vector3<ipc::SimdBatch<double>>` is a valid argument to the distance functions. A caller holding a structure-of-arrays layout can then evaluate one independent problem per SIMD lane per call. The scalar templating this builds on is what makes it possible; the values, gradients, and Hessians of point-point, point-line, point-edge, line-line, edge-edge, and point-triangle are instantiated for `SimdBatch<float>` and `SimdBatch<double>`. The generated `autogen` kernels needed instantiating for the batch types too, since they are where the derivatives actually land. Verified lane-by-lane against the scalar path: the values agree to one ulp, and every gradient and Hessian entry agrees to 1e-9 relative to the magnitude of the result. The derivatives need the looser, norm-scaled bound because they are ill-conditioned for near-parallel edges -- a random configuration routinely produces entries of order 1e6 whose small differences are cancellations of large terms -- while the batch and scalar paths contract multiply-adds differently. AUTO and the `*_distance_type` predicates are unavailable for batch scalars and throw. This is semantic, not an omission: the distance type is a per-lane property, but the predicates return a single enum, so two lanes cannot report different closest features. Resolve the distance types scalar-side and group problems by type before batching, which is the natural SoA layout anyway. [Breaking] `xsimd` and `SIMD_CXX_FLAGS` become PUBLIC rather than PRIVATE. `xsimd::default_arch` is selected from each translation unit's own compiler flags: with the flags PRIVATE, this library's TUs resolve it to `i8mm+neon64` while a consumer's resolve it to `arm64+neon`, a genuinely different type, so the instantiations here would not match what a caller names. The cost is that consumers are now compiled with the detected SIMD flags (typically `-march=native`), which makes their binaries non-portable; disable IPC_TOOLKIT_WITH_SIMD if that is unwanted. Moving the edge-edge and point-triangle kernel definitions into a `.tpp` so callers instantiate in their own TU would avoid this, at the cost of a larger refactor. A caveat on the payoff: on an Apple M-series (NEON, 2 lanes per double) a batched point-line sweep measured 1.06-1.13x over the scalar loop, and 1.13-1.18x with float at 4 lanes -- well short of the lane count, because the compiler already auto-vectorizes a loop over independent problems and these kernels are largely memory-bound. Wider ISAs may do better; that is unmeasured. Only `SimdBatch<double>` is covered by the tests; the `float` instantiations compile and link but are not numerically verified. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
`sqr` was defined four times across the library: as `Math<T>::sqr` in math.hpp and as three separate anonymous-namespace copies, two of which were `double`-only and so unusable from templated code. - Add free `ipc::sqr` and `ipc::cubic` to a new `math/scalar_math.hpp`. - `Math<T>::sqr`/`cubic` forward to them, so the public API is unchanged and there is a single implementation. - Drop the duplicates in normal_collisions.cpp and voxel_size_heuristic.cpp. Split math.hpp along its dependencies so consumers that only want the scalar helpers do not pay for the rest: - `math/scalar_math.hpp` -- no includes. - `math/heaviside_type.hpp` -- `HeavisideType`, `OrientationTypes`. - `math/math.hpp` -- the parts needing eigen_ext, and math.tpp which pulls in the autodiff scalars. It still includes both new headers, so existing `#include <ipc/math/math.hpp>` keeps working. normal_collisions.cpp and voxel_size_heuristic.cpp no longer pull in TinyAD as a result (verified with -H). Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Instantiate the remaining scalar-templated geometry for `SimdBatch<float>` and `SimdBatch<double>`: the barrier functions and classes, the tangent bases and their Jacobians, the closest-point Jacobians/Hessians, the unnormalized-normal Jacobians/Hessians, the signed-distance Hessians, and the triangle-area gradient. Three things blocked this and are fixed here: - Qualified `std::sqrt` cannot find `xsimd::sqrt`, which is only reachable by ADL. `ipc::sqrt` does the lookup in one place via a block-scope `using std::sqrt`, so it also covers the autodiff scalars. - Eigen's `normalized()` guards a zero-length vector with an `if`, which a batch cannot answer with one bool. `ipc::normalized` applies that rule per-lane and defers to Eigen for scalars. - Functions that pick a case from the values themselves (a barrier clamping at dhat, the 3D point-point tangent basis choosing a reference axis) have no single answer for a batch. `select_lazy` expresses these as one first-match-wins cascade: a scalar evaluates only the case it lands in, while a batch evaluates every case and blends per-lane. Each barrier is now written once rather than as separate scalar and batch bodies. The `*_distance_type` predicates remain scalar-only by design -- they return an enum, which has no per-lane blend. Also fixes TwoStageBarrier at penetration. It had no `d <= 0` guard, so its log returned NaN for d < 0 and its first derivative flipped sign, reporting an attractive force. It now matches the other log barriers: +inf for the value, 0 for the derivatives. Tests cover the batch/scalar agreement lane-wise, sweeping lane offsets so every branch is exercised even when a batch holds fewer lanes than there are cases, plus the barrier penetration and stage-boundary conventions. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
Time the per-collision barrier value, gradient, and Hessian with double, float, SimdBatch<double>, and SimdBatch<float> on the assembly benchmark scenes, so the four are compared on identical work. Collisions are grouped by kind and distance type and packed one per lane; a single templated loop evaluates the chain d -> f(d) -> derivatives for every scalar type, and the library path is timed alongside as a reference. A compute-only Hessian column sums the entries instead of storing them, separating the arithmetic from the 144 stores per collision, which are memory-bound once threads share the bus. Every variant is checked against the double path, and the double path against the library on the unmollified collisions. Two caveats the report states explicitly: the edge-edge mollifier has no batch implementation and is excluded, and raw float overflows the generated line-line and point-plane Hessians at scene scale (edge lengths ~1e-4), so "rescaled" float variants re-center each stencil and divide by dhat before conversion, folding the scale back into kappa. On puffer-ball (512k collisions, AVX2, single-threaded) SimdBatch<double> is 3.1x faster than the scalar loop on the gradient and 2.2x on the Hessian without stores; SimdBatch<float> 7.1x and 4.5x. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
The normalized normals (`point_line_normal`, `triangle_normal`, `line_line_normal`) called Eigen's `normalized()`, whose zero-length guard is an `if` on a value that a batch cannot answer with one bool. They now use `ipc::normalized`, which applies that rule per-lane, as the tangent bases already do. That was the only thing blocking the values and gradients of the `point_line`, `line_line`, and `point_plane` signed distances, which route through those normals; previously only their Hessians were instantiated for batch scalars. All of these are header-only templates, so no new instantiations were needed. `edge_length_gradient` asserted `(e1 - e0).norm() != 0`, which likewise has no single answer for a batch. That assert is now scalar-only, matching how `tangent_basis.hpp` handles the same pattern. The relative-velocity functions turned out to already be batch-compatible; they are noted in the release notes and covered by tests separately. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
…ions Six SIMD test files each carried a private copy of the same lane-packing and lane-comparison helpers. They now share `tests/simd_utils.hpp`, which owns `Batch`, `L`, `Lanes`, `Points<dim>`, `pack`, `random_points`, `lane_cases`, `approx`, `check_lanes`, `check_scalar_lanes`, and `check_swept_lanes`. The four pre-existing files shrank by 37-81 lines each. Three changes beyond a straight extraction: - `Lanes` is a `std::array<double, Batch::size>` aligned to the architecture's requirement, so packing uses `load_aligned`/`store_aligned` rather than staging through an unaligned buffer. A `static_assert` checks that the type actually meets the alignment its load claims. - Comparisons go through `Catch::Approx` instead of a hand-rolled `close`. `epsilon(0)` is not optional there: `Approx` ors its margin test with an epsilon test whose default is float-grade (~1e-5 relative), which would silently pass a wrong double. `Approx` also compares without subtracting, so it already gives the sign-respecting infinity matching the barriers need. - The barrier's region sweep and the mollifier's threshold sweep were the same loop, and are now one `check_swept_lanes`. New coverage for the batch functions that had none: the normals and the signed distances, the triangle area and edge length, and the relative velocities, which were already batch-compatible but untested. The relative velocities are held to the value tolerance rather than the derivative one, since they are linear combinations with no cancellation of large terms. `test_simd_utils.cpp` pins down the shared helper's own contract, six files now depending on it: the epsilon leak above, infinity matching by sign, and NaN never passing. It caught `tol * scale` going infinite when the expected value is `+inf` -- the case the barrier tests rely on -- which would have made those lanes accept anything. The signed-distance values needed a looser bound than the unsigned ones, and structurally so: a signed distance is `normal.dot(p - x0)` with the normal a normalized cross product, so the cross cancels and the sqrt and division that follow amplify the two paths' differing multiply-add contraction. The reason is recorded at the call site rather than hidden in a per-file constant.
The packed buffers were plain `std::vector<R>`, guaranteeing only `alignof(R)`, so the batch loads and stores had to be unaligned. They now use `xsimd::default_allocator<R>`, which falls back to `std::allocator` on an architecture with no alignment requirement. Every batch access is a whole number of lanes into one of these buffers, and one lane's width times the lane count is the architecture's alignment, so an aligned base makes all of them aligned. This matters here rather than in the tests: an unaligned load can cost throughput on some instruction sets, which is the very thing the benchmark measures. The traits assert `xsimd::is_aligned<A>` on every load and store instead of leaving that reasoning as a comment. An aligned load off a misaligned address happens to work on NEON and faults on AVX, so it cannot be validated by running on one machine; the assert checks it on every access in a debug build, including after any future change to the packing layout. Also drops the note that a batch cannot evaluate the edge-edge mollifier, which stopped being true when the mollifier gained batch instantiations. The benchmark still does not apply it, so that the four variants stay on the same arithmetic the original numbers were taken with. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
- Replace `A.ldlt().solve()` with a branchless closed form, since pivoting branches on matrix values and cannot vectorize - Add `scalar_of_t` and `all_of` to `ipc/utils/simd.hpp`, fixing the residual tolerance for batches and keeping its assert alive per-lane Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
These were the last hardcoded-`double` potentials, which is what ruled out a batch-capable friction or adhesion path. - Move the friction mollifier, smooth-mu, and adhesion functions into header-only templates, deleting `smooth_friction_mollifier.cpp` and `adhesion.cpp` - Rewrite their piecewise branches as `select_lazy` cascades, so one batch may carry lanes on either side of a threshold. Several inactive branches divide by `y`, which a batch evaluates even where `y == 0`; the per-lane blend discards the infinity rather than propagating a NaN - Split `dihedral_angle` and its derivatives into `ipc::detail` kernels instantiated for `float`, `double`, and both batch types, behind front ends that deduce `T` - Add `ipc::abs` and `ipc::atan2`. `Math<T>::abs` picks the sign with a ternary, which asks a batch for one `bool` its lanes may disagree on - Keep the anisotropic friction helpers `double`-only: their branches test the material rather than a per-problem speed, so there is nothing per-lane Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
- Return `decltype(auto)` from `ipc::abs` so it stops hijacking Eigen's array `abs` under `using namespace ipc` - Zero the 2x2 closest-point solve on a singular Gram matrix instead of dividing by zero, and scale its residual assert to the system's magnitude instead of an absolute bound - Add `ipc::literal<T>` and use it wherever a fractional constant would narrow for a `SimdBatch<float>` - Collapse the generic `std::` math forwarders into one variadic macro over a single name list - Delegate `Math<T>::abs` to `ipc::abs` instead of a batch-incompatible ternary - Fix the point-edge distance test's inverted transition-point guard and add the degenerate-edge guard its finite-difference comparison needs - Factor the dihedral-angle Jacobian assembly out of the gradient and Hessian, consolidate the SIMD test helpers, and correct several stale or misleading comments - Decouple config.hpp.in from <Eigen/Core> by mirroring the StorageOptions values it needs Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
- Fix single bracket array initializer in tests
Suppress -Warray-bounds only on GCC 16+, the one compiler that trips it. In a translation unit instantiating tbb::enumerable_thread_specific for more than one type, GCC 16 speculatively devirtualizes the type-erased construct callback in combine() to an override for the wrong type and reports the placement-new into the smaller ets_element as out of bounds. GCC <= 15 and Clang are clean across all 182 translation units, so they keep the warning instead of losing it project-wide. Drop -Wmaybe-uninitialized from the non-GNU branch: it is GCC-only, so check_cxx_compiler_flag filters it out on Clang while the -Wno- form disables it on GCC. It was enabled nowhere. Remove <type_traits> and the EIGEN_USING_STD(abs) workaround from eigen_ext.hpp; neither the header nor its .tpp uses either one, and the latter injected a using-declaration into the global namespace from a public header.
zfergus
marked this pull request as draft
September 6, 2026 01:53
- Add MeshFEMSparse, which is on by default but was missing - Drop nlohmann/json and Tracy, which are off by default - Dash the edges of optional dependencies and note it in the legend - Correct several node comments that named the wrong target
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #251 +/- ##
==========================================
+ Coverage 96.41% 96.70% +0.29%
==========================================
Files 179 191 +12
Lines 16808 17291 +483
Branches 962 928 -34
==========================================
+ Hits 16205 16721 +516
+ Misses 603 570 -33
Flags with carried forward coverage won't be shown. Click here to find out more. ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Contributor
There was a problem hiding this comment.
🟡 Changes recommended
A debug-only assertion in solve_spd_2x2 uses || on a potentially per-lane SIMD mask, which is unsafe for batch types and should be replaced with a per-lane OR to ensure correct compilation/behavior.
Once you've addressed the issues Copilot identified, you can request another Copilot review.
Pull request overview
Extends IPC Toolkit’s geometry-layer scalar templating to support xsimd::batch (ipc::SimdBatch<T>) as a valid “scalar” throughout distance, barrier, tangent-basis, closest-point, friction, adhesion, and dihedral-angle codepaths, enabling structure-of-arrays callers to evaluate one independent problem per SIMD lane.
Changes:
- Add SIMD scalar infrastructure (
ipc/utils/simd.hpp) and ADL-friendly scalar math wrappers (ipc/math/scalar_math.hpp), then refactor value-dependent branches toselect_lazy. - Template/instantiate the remaining geometry-layer hardcoded-
doublepaths (friction/adhesion/dihedral angle) and add SIMD-focused unit tests + a SIMD barrier benchmark. - Update build/docs to reflect that SIMD flags and
xsimdare now applied/linkedPUBLICwhenIPC_TOOLKIT_WITH_SIMD=ON.
File summaries
| File | Description |
|---|---|
| tests/src/tests/utils/test_simd_utils.cpp | Adds unit tests pinning SIMD test helper semantics (Approx behavior, lane packing). |
| tests/src/tests/utils/CMakeLists.txt | Registers test_simd_utils.cpp. |
| tests/src/tests/tangent/test_simd_tangent_basis.cpp | Adds lane-wise SIMD vs scalar tangent-basis tests (including per-lane reference-axis blend). |
| tests/src/tests/tangent/test_simd_relative_velocity.cpp | Adds SIMD coverage for relative velocity values/Jacobians/dx_dbeta in 2D/3D. |
| tests/src/tests/tangent/test_simd_closest_point.cpp | Adds SIMD closest-point tests including near-parallel edge conditioning. |
| tests/src/tests/tangent/test_closest_point.cpp | Tightens closest-point scalar assertions (zero epsilon, tiny margin). |
| tests/src/tests/tangent/CMakeLists.txt | Registers tangent SIMD test files. |
| tests/src/tests/simd_utils.hpp | New shared SIMD test utilities (lane packing, comparisons, tolerances). |
| tests/src/tests/potential/CMakeLists.txt | Registers SIMD barrier potential benchmark source. |
| tests/src/tests/geometry/test_simd_geometry.cpp | Adds SIMD tests for normals, normalization, area, and dihedral angle (value/grad/Hess). |
| tests/src/tests/geometry/CMakeLists.txt | Registers test_simd_geometry.cpp. |
| tests/src/tests/friction/test_simd_friction.cpp | Adds SIMD tests for friction mollifier and smooth-μ families across thresholds. |
| tests/src/tests/friction/CMakeLists.txt | Registers test_simd_friction.cpp. |
| tests/src/tests/distance/test_simd_signed_distance.cpp | Adds SIMD tests for signed distances (value/grad/Hess) with adjusted tolerance. |
| tests/src/tests/distance/test_simd_edge_edge_mollifier.cpp | Adds SIMD tests for edge-edge mollifier + per-lane threshold blending. |
| tests/src/tests/distance/test_simd_distance.cpp | Adds SIMD tests for core distances plus explicit rejection of AUTO. |
| tests/src/tests/distance/test_point_edge.cpp | Makes finite-difference derivative checks robust by skipping degenerate edges and fixing branch condition. |
| tests/src/tests/distance/CMakeLists.txt | Registers distance SIMD test files. |
| tests/src/tests/CMakeLists.txt | Adds simd_utils.hpp to test target sources. |
| tests/src/tests/broad_phase/test_lbvh.cpp | Removes unused iostream include. |
| tests/src/tests/barrier/test_simd_barrier.cpp | Adds SIMD barrier tests over all piecewise regions. |
| tests/src/tests/barrier/test_barrier.cpp | Adds regression tests for penetration convention + stage boundaries; cleans up abs usage. |
| tests/src/tests/barrier/CMakeLists.txt | Registers SIMD barrier tests. |
| tests/src/tests/adhesion/test_simd_adhesion.cpp | Adds SIMD tests for normal/tangential adhesion and smooth-μ adhesion variants. |
| tests/src/tests/adhesion/CMakeLists.txt | Registers test_simd_adhesion.cpp. |
| src/ipc/utils/simd.hpp | New SIMD glue: SimdBatch, NumTraits specialization, select_lazy, literal, all_of, normalized, etc. |
| src/ipc/utils/eigen_ext.hpp | Removes conditional EIGEN_USING_STD(abs) workaround. |
| src/ipc/utils/CMakeLists.txt | Registers new simd.hpp. |
| src/ipc/tangent/tangent_basis.hpp | Switches to ipc::normalized and adds per-lane selection for SIMD reference-axis choice. |
| src/ipc/tangent/closest_point.hpp | Replaces LDLT solve with branchless 2×2 solver for SIMD compatibility. |
| src/ipc/tangent/closest_point.cpp | Instantiates closest-point autogen for SIMD batch types. |
| src/ipc/smooth_contact/common.hpp | Qualifies abs as std::abs in parameter validation. |
| src/ipc/math/scalar_math.hpp | New ADL-friendly wrappers for abs/atan2/fma/log/sqrt + shared sqr/cubic. |
| src/ipc/math/math.tpp | Uses Math<T>::abs consistently. |
| src/ipc/math/math.hpp | Splits heaviside types + scalar helpers out; routes Math<T> helpers through ipc:: wrappers. |
| src/ipc/math/heaviside_type.hpp | New header extracting HeavisideType and OrientationTypes. |
| src/ipc/math/CMakeLists.txt | Registers new math headers. |
| src/ipc/geometry/normal.hpp | Switches to ipc::normalized for batch-safe normalization. |
| src/ipc/geometry/normal.cpp | Uses ipc::sqrt and instantiates for SIMD batch types. |
| src/ipc/geometry/area.hpp | Avoids scalar-only degeneracy assert when templated on SIMD batch. |
| src/ipc/geometry/area.cpp | Uses ipc::sqrt and instantiates for SIMD batch types. |
| src/ipc/geometry/angle.hpp | Converts dihedral-angle API to templated front ends calling detail kernels. |
| src/ipc/friction/smooth_mu.hpp | Makes smooth-μ family templated and SIMD-safe via select_lazy. |
| src/ipc/friction/smooth_mu.cpp | Removes scalar smooth-μ implementations (now header-only templates); keeps anisotropic helpers double-only. |
| src/ipc/friction/smooth_friction_mollifier.hpp | Makes friction mollifier functions templated and SIMD-safe via select_lazy. |
| src/ipc/friction/smooth_friction_mollifier.cpp | Deletes old scalar implementations (now header-only templates). |
| src/ipc/friction/CMakeLists.txt | Stops compiling deleted friction mollifier .cpp. |
| src/ipc/distance/signed/point_plane.cpp | Uses fixed-size reshaped(fix<…>) and instantiates Hessian for SIMD batch types. |
| src/ipc/distance/signed/point_line.cpp | Uses fixed-size reshaped(fix<…>), guards scalar-only asserts, instantiates for SIMD batch types. |
| src/ipc/distance/signed/line_line.cpp | Uses fixed-size reshaped(fix<…>) and instantiates for SIMD batch types. |
| src/ipc/distance/point_triangle.cpp | Instantiates point-triangle distance family for SIMD batch types (explicit type only). |
| src/ipc/distance/point_plane.cpp | Instantiates point-plane autogen for SIMD batch types. |
| src/ipc/distance/point_line.cpp | Instantiates point-line autogen for SIMD batch types. |
| src/ipc/distance/point_edge.cpp | Instantiates point-edge Hessian kernels for SIMD batch types. |
| src/ipc/distance/line_line.cpp | Instantiates line-line autogen for SIMD batch types. |
| src/ipc/distance/edge_edge.cpp | Instantiates edge-edge distance family for SIMD batch types. |
| src/ipc/distance/edge_edge_mollifier.hpp | Refactors mollifier branches to select_lazy, uses ipc::fma, adds scalar fast-path guards. |
| src/ipc/distance/edge_edge_mollifier.cpp | Adds scalar fast-path guards + instantiates mollifier templates/autogen for SIMD batch types. |
| src/ipc/config.hpp.in | Removes Eigen include by mirroring storage option constants and substituting layout token directly. |
| src/ipc/collisions/normal/normal_collisions.cpp | Reuses shared ipc::sqr helper by including scalar_math.hpp. |
| src/ipc/broad_phase/voxel_size_heuristic.cpp | Reuses shared ipc::sqr helper by including scalar_math.hpp. |
| src/ipc/barrier/barrier.hpp | Updates TwoStageBarrier documentation to include d <= 0 convention. |
| src/ipc/barrier/barrier.cpp | Refactors barriers to select_lazy, adds SIMD instantiations, and uses ipc::log/sqr/infinity. |
| src/ipc/adhesion/CMakeLists.txt | Stops compiling deleted adhesion .cpp. |
| src/ipc/adhesion/adhesion.cpp | Deletes old scalar adhesion implementation (now header-only templates). |
| python/src/geometry/angle.cpp | Binds Python angle APIs to detail::…<double> template instantiations. |
| python/src/friction/smooth_mu.cpp | Binds Python smooth-μ APIs to templated <double> instantiations. |
| python/src/friction/smooth_friction_mollifier.cpp | Binds Python friction mollifier APIs to templated <double> instantiations. |
| python/src/adhesion/adhesion.cpp | Binds Python adhesion APIs to templated <double> instantiations. |
| docs/source/about/release_notes.rst | Documents new templating/SIMD support and the PUBLIC SIMD flag/link behavior. |
| docs/source/about/dependencies.rst | Adds warning about PUBLIC SIMD flags and binary portability implications. |
| docs/source/_static/graphviz/dependencies.dot | Updates dependency graph legend and indicates optional deps via dashed edges and xsimd scope. |
| CMakeLists.txt | Makes SIMD compile options and xsimd linkage PUBLIC when SIMD is enabled. |
| cmake/recipes/meshfem_sparse.cmake | Switches generated export header creation to file(CONFIGURE CONTENT ...). |
| cmake/ipc_toolkit/ipc_toolkit_warnings.cmake | Adds targeted warning suppressions for specific GCC diagnostics. |
| CLAUDE.md | Adds benchmark/test-running guidance (notably excluding [!benchmark] when filtering [simd]). |
| .gitignore | Ignores __cmake_systeminformation. |
Review details
- Files reviewed: 80/82 changed files
- Comments generated: 1
- Review effort level: Lite
💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.
The nearly-parallel sweep stopped at 1e-3, so the singular branch of solve_spd_2x2 was never taken per-lane. Slide the parallel pair along its shared direction so the zeroed lane leaves a residual of order ‖b‖; neutering the det <= 0 guard now trips the assert.
Random segments are unconstrained in direction and offset, so the closest points routinely land far off both ends and the coordinates amplify the inputs. On an 8-lane runner one draw drifted 1.2e-14 relative, just over VALUE_TOL, while the same test passes at 2 lanes. The sibling near-parallel test already uses this bound for the same quantity.
The --coverage compile options were PRIVATE, so only the library's own translation units were instrumented. Much of the library is header-only templates instantiated in the caller's TU, so a header exercised solely by tests recorded no coverage and read as untested: adhesion.hpp at 34.6% despite scalar and batch tests covering every function in it. The link options were already PUBLIC, and the workflow strips '*tests/*' from the report, which only has an effect if tests are instrumented.
A single "*" does not cross a "/", so "tests/*" only matched files sitting directly in tests/ and everything under tests/src/tests/ was reported. Confirmed with codecov's validator: the old pattern compiles to (?s:tests/[^/]*)\Z, the new one to (?s:tests/.*)\Z.
The closest-point Hessians had no test at all: the existing cases finite-difference each Jacobian against the value, then stop. Check each Hessian against a finite difference of the corresponding Jacobian, in the same cases. The front ends also have a runtime branch on size() for argument types that do not carry their dimension. Every test passed a fixed-size vector, so only the if constexpr arms ran. Passing the same points as a VectorMax3d takes the fallback; both reach the same kernel, so the results must be identical. Takes tangent_basis.hpp from 73.3% to 100% and closest_point.hpp from 65.4% to 100%.
Coverage flags on ipc_toolkit go back to PRIVATE; the test target sets the same flags on itself. Instrumentation still reaches the tests, which it must because a template is only compiled where it is instantiated, but consumers of ipc::toolkit no longer inherit --coverage.
The naive a00*a11 - a01*a10 cancels catastrophically as the edges approach parallel, which is where the closest-point solve is used. Recovering the dropped rounding error with fma makes the determinant good to a couple of ulp however badly the products cancel. Measured against the exact determinant of the same double inputs, the relative error at 1e-4 rad drops from 2.5e-9 to 4.5e-17, and the full solve improves by the same factor: the determinant, not the numerators, was the dominant error for this geometry. Well-conditioned angles are unchanged. Degrades to the naive expression, rather than misbehaving, where xsimd has no hardware fma; NEON and FMA3 both provide a fused one.
The probe needs -mavx2 -mfma and <immintrin.h>, so it can only ever answer for x86. On AArch64 it failed while the chip does have hardware FMA, which reads like a missing capability. Only the log line changes; FMA_FOUND and FMA_FLAGS are untouched, and an empty FMA_FLAGS stays correct on AArch64 because fmadd and vfmaq_* are mandatory there with no flag to opt into.
The old note said accuracy was unchanged, which was measured before the determinant used Kahan's fused form. Measured against the exact solution of the same double inputs, the closed form now beats LDLT by up to 5.6e7x near parallel, because LDLT's second pivot suffers the same cancellation the naive determinant did.
abs, atan2, fma, log and sqrt carry the most collision-prone names in C++. At ipc scope they join the overload set of anyone writing 'using namespace ipc;', where ipc::abs was already found to hijack Eigen's array abs. They are mechanism, not API: each exists only so a template can reach the xsimd and autodiff overloads by ADL. Following Eigen::numext, they now live in their own namespace. Not ipc::detail, since the kernels live there and would pick these up by ordinary lookup ahead of ADL, reproducing the same hijack one level down. sqr, cubic and MOLLIFIER_THRESHOLD_EPS stay at ipc scope; they collide with nothing and Math<T>::sqr forwards to them as public API.
zfergus
added a commit
that referenced
this pull request
Sep 7, 2026
Main's SIMD batch support (#251) instantiates the distance, geometry, tangent, and barrier templates for xsimd batches. Sixteen of those files are shared device sources here, and a batch has no device-callable operations, so nvcc would have compiled host-only calls into device code for every one of them. Rather than split each file, we let a shared device source be compiled twice under CUDA and split the instantiations between the passes: - config.hpp gains IPC_TOOLKIT_INSTANTIATE_DEVICE_SCALARS and IPC_TOOLKIT_INSTANTIATE_HOST_SCALARS, both 1 when CUDA is off. - ipc_toolkit_target_shared_device_sources() now adds the .cpp alongside its generated .cu wrapper, defining one switch per pass, so float and double are emitted by nvcc and the autodiff and batch scalars by the host compiler. Every symbol is still emitted exactly once. This replaces the edge_edge/point_triangle _impl.hpp split from the previous merge, so both files go back to main's layout. Two of main's changes subsume workarounds this branch was carrying: - solve_spd_2x2() replaces Eigen's pivoting LDLT with Cramer's rule for the same reason we needed a device path, so the __CUDA_ARCH__ branches in closest_point.hpp are gone. That also drops a stale IPC_TOOLKIT_HOST_DEVICE on a constexpr variable, which is not valid CUDA. - scalar_math.hpp routes the std math functions through the global namespace under nvcc, which is what our unqualified fma and log calls were doing by hand. Its forwarders, sqr, cubic, and the scalar helpers in simd.hpp are annotated, since device code reaches them through the barrier functions. Kept from this branch where main's version is not device-safe: the Maps standing in for .reshaped() in the signed-distance contractions, and the explicit column scatter in place of Eigen::all index slicing in angle.cpp. SIMD compile flags stay PUBLIC, as main made them, but remain restricted to C++ sources so nvcc does not read them as input files.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Description
Makes an
xsimdbatch a valid scalar type throughout the geometry layer, so a caller holding a structure-of-arrays layout evaluatesSimdBatch<T>::sizeindependent collisions per call instead of one. Builds on the scalar templating of #249 and covers the distance, barrier, tangent-basis, closest-point, normal, mollifier, friction, adhesion, and dihedral-angle functions.The mechanism is two new headers plus a rewrite of every value-dependent branch:
ipc/utils/simd.hpp(new) — specializesEigen::NumTraitsforxsimd::batchsoEigen::Vector3<ipc::SimdBatch<double>>is a valid argument, and supplies the scalar/batch bridges the kernels need:scalar_of_t,all_of,select,select_lazy,literal<T>,infinity<T>, andipc::normalized.ipc/math/scalar_math.hpp(new) —ipc::abs,atan2,fma,log, andsqrtresolve through a block-scopeusing std::…, so ADL reaches thexsimd::overloads a qualifiedstd::sqrtcannot find. It also holds the freeipc::sqr/ipc::cubicthat replace four separate definitions, and splits the scalar helpers out ofmath.hppso consumers stop pulling in Eigen and TinyAD for them.select_lazycascades replace theif/elsein the barriers, the friction and edge-edge mollifiers, smooth-μ, adhesion, and the tangent basis' reference-axis pick. A scalar still evaluates only the case it lands in; a batch evaluates every case and blends per lane, so one batch may carry lanes on either side of a threshold.T— the last hardcodeddoubles in the geometry layer, which is what ruled out a batch friction or adhesion path.smooth_friction_mollifier.cppandadhesion.cppare gone; both headers are header-only templates now, followingedge_edge_mollifier.hpp.Impact:
SimdBatch<double>and 1.80–4.32× withSimdBatch<float>on NEON, rising to 1.01–3.09× and 2.09–7.36× on AVX2; see Performance.doubleinstantiations keep their original early-return form.API changes
ipc::SimdBatch<T>aliasesxsimd::batch<T, xsimd::default_arch>. PassingEigen::Vector3<ipc::SimdBatch<double>>to a distance function evaluates one independent problem per lane.SimdBatch<float>andSimdBatch<double>: the values/gradients/Hessians of the point-point, point-line, point-edge, line-line, edge-edge, and point-triangle distances; the barrier functions and classes (barrier,ClampedLogBarrier,ClampedLogSqBarrier,CubicBarrier,TwoStageBarrier); the tangent bases and their Jacobians; the closest-point Jacobians/Hessians; the unnormalized-normal Jacobians/Hessians; the signed-distance Hessians; the triangle-area gradient; the edge-edge mollifier family; and the dihedral angle with its gradient and Hessian.edge_edge_closest_pointandpoint_triangle_closest_pointsolve their 2×2 symmetric positive-definite system with Cramer's rule instead ofA.ldlt().solve(). Pivoting requires branching on matrix values, which breaks vectorization since different lanes would need different control flow. Accuracy is limited by the Gram matrix's conditioning rather than the algorithm, so the branchless form wins. On a singular matrix it returns zero branchlessly instead of dividing by it, matching Eigen's LDLT pseudo-inverse behavior.point_line_normal,triangle_normal,line_line_normal) useipc::normalizedinstead of Eigen'snormalized(), whose zero-length guard is anifa batch cannot answer with onebool. That was the only thing blocking the values and gradients of thepoint_line,line_line, andpoint_planesigned distances, which route through those normals.AUTOand the*_distance_typepredicates are not available for batch scalars and throwstd::invalid_argument. This is semantic, not an omission: the distance type is a per-lane property, but the predicates return a single enum, so two lanes cannot report different closest features. Resolve the distance types scalar-side and group problems by type before batching, which is the natural structure-of-arrays layout anyway.🐛 Fixed in passing
TwoStageBarrierhad nod <= 0guard, so its log returned NaN ford < 0and its first derivative flipped sign, reporting an attractive force at penetration. It now matches the other log barriers:+inffor the value,0for the derivatives.💥 Breaking changes
xsimdandSIMD_CXX_FLAGSare linked/appliedPUBLICrather thanPRIVATE.xsimd::default_archis selected from each translation unit's own compiler flags: with the flagsPRIVATE, this library's TUs resolve it toi8mm+neon64while a consumer's resolve toarm64+neon— a genuinely different type, so the instantiations here would not match what a caller names and the link fails. The cost is that consumers are now compiled with the detected SIMD flags (typically-march=native), making their binaries non-portable; disableIPC_TOOLKIT_WITH_SIMDif that is unwanted.smooth_mu(float_y, 0.5, 0.3, 0.001)must becomesmooth_mu(float_y, 0.5f, 0.3f, 0.001f).anisotropic_mu_eff_f,anisotropic_x_from_tau_aniso,anisotropic_mu_eff_from_tau_aniso) staydouble-only. Their branches test the material rather than a per-problem speed, so there is nothing per-lane to blend yet.eigen_ext.hppno longer leaks a globalusing std::absunderEIGEN_DONT_VECTORIZE. Callers relying on unqualifiedabsat namespace scope must qualify it; the affected tests in this repo are updated.Performance
Per-collision barrier potential — value, gradient, and Hessian — evaluated with five scalar types on the same collision set across the assembly benchmark scenes, so every variant does identical work. The distance type is resolved scalar-side and passed in, which a batch requires; the edge-edge mollifier is excluded because no batch caller applies it yet. "Hessian (compute only)" sums the entries instead of storing them, separating the arithmetic from the 144 stores per collision.
Apple M3 Pro — NEON, 128-bit (2 lanes per
double, 4 perfloat)SimdBatch<double>runs 0.97–2.07× andSimdBatch<float>1.80–4.32× over the scalardoubleloop, across all six scenes and all four quantities.Intel Core Ultra 9 285K — AVX2+FMA (4 lanes per
double, 8 perfloat)SimdBatch<double>reaches 1.01–3.09× andSimdBatch<float>2.09–7.36×, against 2.07× and 4.32× on NEON.What the Hessian columns show. Speedup falls off exactly where the stores dominate. On the two largest scenes on AVX2,
SimdBatch<double>'s stored Hessian collapses to parity — 1.01× on rod-twist and 1.02× on puffer-ball — while the compute-only Hessian on the same data holds 1.71× and 1.82×. Each collision writes 144 doubles, and that traffic does not shrink with lane width; 24 threads saturate the bus before the arithmetic is the limit. NEON, with half the lanes and fewer threads, never reaches that point: 1.37× stored against 1.51× compute-only on puffer-ball.floatneeds rescaling to be usable. Rawfloatoverflows the generated line-line and point-plane Hessians at scene scale (edge lengths ~1e-4), and does so identically for the scalar and batch float paths — a precision property, not a SIMD one. The plottedfloattimings are from the rescaled inputs, which re-center each stencil on its centroid and divide by d̂ before conversion, folding the scale back into κ; the instructions are the same either way, only the packing differs.Accuracy.
SimdBatch<double>reproduces the scalar path: on puffer-ball's 512k collisions the value agrees to 6.4e-16 relative and the gradient to 1.8e-16 median. The 9.0e-11 gradient max is the near-parallel edge-edge configurations where the derivative is ill-conditioned and the two paths contract multiply-adds differently.New dependencies
None.
xsimdwas already fetched via CPM by #207; this PR only widens its link scope toPUBLIC.Type of change
How Has This Been Tested?
test_simd_{distance,signed_distance,edge_edge_mollifier,barrier,friction,adhesion,geometry,closest_point,tangent_basis,relative_velocity,utils}.cpp, sharing the helpers intests/src/tests/simd_utils.hpp. Each packs independent problems one per lane and checks every lane against the scalar path, sweeping lane offsets so every branch of aselect_lazycascade is exercised even when a batch holds fewer lanes than there are cases. Run withipc_toolkit_tests "[simd] ~[!benchmark]".ipc_toolkit_tests "~[!benchmark]"on a Debug build: 320 cases / 3,510,369 assertions pass, 5 skipped (the spatial-hash and STQ cases that skip themselves in Debug). No CUDA on this host, soIPC_TOOLKIT_WITH_CUDA=OFF.ipc_toolkit_tests "[simd_barrier_potential]", which checks every scalar variant against thedoublepath, and thedoublepath against the library, before reporting any timing. Run on both machines plotted above:-O3, Apple clang 21.0.0, arm64 NEON 128-bit.-O3, GCC 11.5,-march=native(AVX2+FMA).Test Configuration:
-march=native; Debug for the unit tests, Release for the benchmarkChecklist