NanoVDB Math: refactor onto MatBase/VecBase, fix correctness bugs, add tests

[swahtz]

Consolidates nanovdb/math/Math.h by lifting Mat / Vec arithmetic onto two shared base classes (MatBase<T, R, C> and a new VecBase<T, N>). The five Mat and three Vec derived classes become thin delegating wrappers. The pass also fixes several long-standing correctness bugs and fills inconsistent / partial APIs.

Files changed: Math.h (the refactor), NanoVDB.h (Map API now constexpr — depends on the matMult rework), tools/cuda/PointsToGrid.cuh (three call-site casts from the explicit ctor tightening), and three test files (additive: 11 new TEST_F blocks plus mechanical cast updates), and one python-binding kernel.

Bug fixes

• math::Min / Max for int32_t / uint32_t no longer route through fminf / fmaxf — old form lost precision above 2^24 and could be UB near INT_MAX. Heavily used by HDDA.h.
• Removed dimensionally-invalid operator*(Mat3, Mat2x3) -> Mat2x3. No callers in the tree.
• Mat2::inverse() now returns an explicit zero on singular input (was uninitialised) and uses T(1)/det for T=double.
• Vec*::operator/(T) and Mat*::operator/(T) now per-element-divide. The old (T(1)/s) * (*this) did integer division for integer T, so Vec3i(6,8,10) / 2 returned (0,0,0).
• Vec2's operator+/-/+=/-= against Coord now take Coord2. Old signature silently dropped z.
• math::Round(Vec3<float>) and Round(Vec3<double>) both use floor(x + 0.5) now. The float overload previously used rintf (round-half-to-even), so the two disagreed at half-integers.
• Corrected floor() / ceil() doc swap on Vec2 / Vec3. Implementations were always right.

Refactor

• MatBase<T, R, C>: shared element-wise / scalar / equality / transpose / matmat / matvec implementations, all __hostdev__ constexpr. rows() / cols() / size() promoted public.
• VecBase<T, N> (new): shared element-wise / scalar / equality / reductions / mutating min/max / integer rounding. floorAs<Result> / ceilAs<Result> / roundAs<Result> unify Round semantics across float / double / long double.
• All five Mat and three Vec classes become thin wrappers: ctors + per-type extras (Mat2::inverse, Vec3::cross/outer) + one-line operators delegating to the base.
• Mat3 / Mat3x2 / Mat4 scalar ctors switched from template<typename Source> to T, matching Mat2 / Mat2x3. Mixed-literal construction (Mat3<double>(1, 2.0, 3, ...)) now works uniformly. Array ctors keep template<typename Source> for cross-type interop.

API completeness

• Every Mat type has the same operator surface: unary -, binary + / -, += / -=, *(T) / /(T), *=(T) / /=(T), == / !=. Previously scattered.
• New mat*mat: Mat3 * Mat3x2 standalone, Mat4 * Mat4 member.
• New mat*vec members: Mat2 * Vec2, Mat2x3 * Vec3 -> Vec2, Mat3x2 * Vec2 -> Vec3. Mat4 * Vec4 promoted to member.
• New scalar * MatN standalones for Mat2x3 / Mat3x2 / Mat3 / Mat4.
• Vec4 gains parity with Vec3: min, max, asPointer, floor, ceil, round (returns Vec4<int32_t> — no Coord4).

constexpr sweep

Every function in Math.h and the Map API in NanoVDB.h that C++17 allows is now constexpr. Carve-outs: length / normalize (std::sqrt, C++26), the float / double overloads of Min / Max / Floor / Ceil / Round / Abs (C++23), and CUDA-device-only *Atomic helpers.

Rewriting matMult / matMultT from fma(a, b, c) to a * b + c was required (fma isn't constexpr until C++23, P1383R2). GPU codegen unchanged — nvcc -fmad=true contracts back to hardware FMA. Host codegen unchanged under default -ffp-contract=on. Only compile-time constexpr evaluation observes the double-rounded result. Worst-case ≤1 ulp drift, well below NanoVDB's geometric precision.

Other quality passes

• Bounds checks on every operator[] in Math.h via NANOVDB_ASSERT (twelve sites). No-op under NDEBUG; const-index accesses fold at compile time even in debug.
• Doxygen sweep: every previously-undocumented public method gains an @brief; several // @brief typos (missing leading slash — silently invisible to Doxygen) fixed.
• Hidden friends: scalar*Vec / scalar/Vec / scalar*Mat converted from namespace-scope free templates to in-class friend definitions — ADL-found only, no longer pollute the enclosing namespace.

Tests

Eleven additive TEST_F blocks in TestNanoVDB.cc:

Test	Covers
IntegerMinMax	2^24+1 precision; INT_MAX / UINT_MAX overflow safety.
Round	Half-integer agreement between float and double (especially -1.5).
Vec2 / Vec3Ops / Vec4Ops	Full operator coverage; new Vec4 methods; integer-T division.
Mat2 / Mat3 / Mat4 / Mat2x3_Mat3x2	Constructors, every operator, transpose, matmat, matvec; Mat2::inverse non-singular roundtrip + singular-returns-zero; mixed-literal construction.
MatMul	Every dimensionally-valid matmat plus every member matvec across the five matrix types.
MatVecIntrospection	static_asserts for rows/cols/size/SIZE/ValueType; data()/asPointer() ordering.

Risk and compatibility

• Behaviour preserved for every previously-correct code path.
• Behaviour changes are exactly the bug fixes above — each moves an incorrect / undefined output to a defined, correct one.
• matMult rounding: ≤1 ulp drift only in pure abstract-machine evaluation. GPU and host codegen identical to before.
• API additions are pure. The only removed overload is the broken Mat3 * Mat2x3 (no callers).
• One signature tightening: Vec2's Coord overloads now require Coord2. Code passing a 3D Coord to Vec2 will fail to compile — but it was silently dropping z.
• Cross-template ctor tightened to explicit: Vec*(const Vec*T<T2>&) now matches the same-class ctor's convention. Audit found 18 affected sites across the nanovdb + openvdb trees, all updated here (3 in PointsToGrid.cuh, 3 in TestNanoVDB.cu, 12 in TestOpenVDB.cc, 1 in PySampleFromVoxels.cu). All mechanical and value-preserving. Downstream consumers that relied on the implicit conversion will need the same one-line cast — typically nanovdb::Vec3f(yourVec) or an explicit template argument to a returning helper.

[Copilot]

Pull request overview

This PR refactors NanoVDB’s math/Math.h by consolidating shared vector/matrix arithmetic into new VecBase<T, N> and enhanced MatBase<T, R, C> base classes, while also fixing several long-standing correctness issues and adding unit tests to lock in the corrected behavior.

Changes:

• Refactor vector and matrix implementations onto VecBase/MatBase helpers to remove duplicated per-element loops and normalize APIs across types.
• Fix correctness issues (e.g., integer Min/Max, consistent Round semantics, singular Mat2::inverse, scalar division for integer Vec*/Mat*, invalid mat-mul overload removal).
• Add 11 new unit tests covering the refactor, new APIs, and each bug fix with numeric checks and compile-time introspection.

Reviewed changes

Copilot reviewed 2 out of 2 changed files in this pull request and generated 1 comment.

File	Description
nanovdb/nanovdb/math/Math.h	Introduces VecBase, expands MatBase, rewrites Vec*/Mat* as thin wrappers, fixes math correctness bugs, and rationalizes overload sets.
nanovdb/nanovdb/unittest/TestNanoVDB.cc	Adds new TEST_F blocks that exercise the refactored math types and pin down the corrected semantics.

---

💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.

[harrism]

This is a really great cleanup and improvement. I have comments in 4 areas.

1. constexpr -- might as well constexpr all the things
2. Alignas -- on the classes where it makes sense, let's guarantee alignment so the compiler can optimize / coalesce loads (this allows each row of a 4x4 matrix to be loaded with a single 128-bit global load, for example. Without it, the compiler will generate 4 32-bit loads).
3. Other modernization / best practices topics: [[nodiscard]], noexcept, consistent explicit ctors, hidden friends.
4. Documentation -- didn't comment on this much inline, but make sure all methods have doxygen docs. Looks a bit inconsistent?

[harrism]

I also asked an agent to review with modern C++ best practices in mind and then had it prune its list to not include what I already mentioned in my review.

This first one would be good to implement.

 Body-assignment in derived constructors. The new derived ctors do this->mVec[0] = x;
  this->mVec[1] = y; instead of the pre-PR : mVec{x, y} {}. Invisible for fundamental T, a
  regression for any future Vec<UserType>. Fix is a protected variadic constructor on VecBase:

  template<typename... Args>
  __hostdev__ explicit VecBase(Args... args) : mVec{T(args)...} {
      static_assert(sizeof...(Args) == N, "wrong number of args");
  }

  Derived ctors delegate: Vec2(T x, T y) : Base(x, y) {}.

The rest are all speculative / nice to have or depend on C++20 or later.

Defaulted operator== (C++20 only). bool operator==(const Vec3&) const noexcept = default;
  replaces the manual equals() chain; operator!= auto-synthesizes. Only worth raising if NanoVDB's
   C++ baseline allows.

  Tuple-like interface (tuple_size / tuple_element / get) for auto [x, y, z] = v; and std::apply.
  C++17.

  begin()/end() so Vec/Mat work with range-based for and <algorithm>. Six lines per base.

  [[deprecated]] on the lowercase size alias if the goal is to migrate callers to SIZE. Nudges
  without breaking.

  static_cast<T>(x) over functional-style T(x). Identical for arithmetic types but greppable and
  doesn't silently become a constructor call for non-arithmetic T. Minor.

  CTAD deduction guides so Vec3 v(1.0f, 2.0f, 3.0f); deduces Vec3<float>. Three lines per Vec /
  Mat class. C++17.

  C++20 requires clauses on cross-Vec templated functions (dot<V>, mergeMin<V>, the templated
  converting constructor). Replaces static_assert(V::SIZE == N) in the body with constraints at
  the declaration. Better diagnostics, SFINAE-friendly.

  std::array<T, N> instead of T[N] in the base. Free .data(), iterators, structured bindings,
  .size(), .fill(), ==. Layout-compatible with T[N] for fundamental T. Possibly tangles with CUDA
  portability (cuda::std::array exists). Out of scope for this PR; worth knowing.

[swahtz]

Implemented in a1072df. Added a protected variadic ctor to both VecBase and MatBase and converted every derived Vec/Mat ctor (Vec2/3/4, Mat2/2x3/3x2/3/4) to delegate via initializer list — no more body-assignment, no more default-construct-then-overwrite of the base subobject. The static_assert(sizeof...(Args) == N) (and == ROWS*COLS for MatBase) catches arity mismatches at the call site.

Tracking the other suggestions for follow-up but not in this PR — happy to file separate issues for the C++17-feasible ones (tuple-like interface, begin/end, CTAD guides) if useful.

[Copilot]

Pull request overview

Copilot reviewed 6 out of 7 changed files in this pull request and generated no new comments.

[harrism]

Fantastic. Thanks for addressing all my suggestions.

[kmuseth]

In addition to my (old) comments above I'd like to discuss the following (alignment and performance) issues identified by Astra:

Review of OpenVDB PR #2207

Recommendation: request changes. The PR introduces a serialized-format break, a measured CPU normalization slowdown, and several arithmetic regressions. Its source compatibility impact is wider than the constructor tightening highlighted in the description.

Reviewed PR #2207, head f4beed62efca49f4f0a431fde70fdabd0954c8da, against its merge base 7191b429e6bb0547d78d23041515db1b6771867d. Review date: September 15, 2026. The comparison uses the PR's actual branch point, avoiding unrelated changes subsequently added to master.

Findings

1. [P1] Preserve Vec4 alignment to avoid silently changing the file format

Location: Math.h:2007.

The new alignas(alignof(T) * 4) leaves standalone sizeof(Vec4) unchanged but changes padding inside structures containing it. Comparing compiled layouts gives:

Property	Base	PR
alignof(Vec4f)	4	16
alignof(Vec4d)	8	32
Vec4f root background offset	28	32
Vec4f root tile value offset	20	32
Vec4d root tile value offset	24	32
Vec4d leaf minimum offset	80	96
Vec4d leaf voxel-array offset	160	192
Vec4d leaf size	16544	16576

Existing Vec4 grids therefore acquire different interpretations. The major format version is unchanged. PNanoVDB's offset table also remains unchanged, so C++ and PNanoVDB in the same checkout already disagree: a byte buffer populated with background (1,2,3,4) at the PR's C++ offsets reads as (0,1,2,3) using PNanoVDB's offsets. The Vec4d leaf change also affects the stride between adjacent leaves, beyond just statistics.

Requested change: retain the existing alignment for serialized types; use separate aligned types or explicit alignment at allocation/use sites where useful. If a format change is intentional, it needs versioning and coordinated reader/writer changes. Add field-offset and leaf-stride tests for Vec4f/Vec4d, including PNanoVDB. Unchanged standalone sizes are insufficient.

The Vec2/Mat2/Mat4 alignment changes also alter downstream enclosing-structure layouts, although this review's concrete NanoVDB file-format failure is Vec4.

2. [P2] Keep the reciprocal path for floating-point normalization/division

Locations: Math.h:1010, Math.h:1052.

The integer division fix replaces the existing floating-point reciprocal-and-multiply implementation too. For Vec3 normalization, generated vectorized ARM64 code changes from one division plus three multiplies to three divisions. normalize() reaches this path through scalar /=. It is used by the repository's ray tracing and collision examples.

Local measurements, in nanoseconds per vector:

Operation	Base	PR	Change
Vec3f normalization	0.852	1.064	+25%
Vec3d normalization	1.779	2.259	+27%
Vec3d scalar division	1.091	1.216	+11%

These are microbenchmarks on an ARM64 Mac with Clang 22.1.8, C++17, -O3 -DNDEBUG -ffp-contract=on, without fast math. Each process used 1,024 input vectors, varying inputs/divisors, 10,000 repetitions and seven timed samples. The table combines the two process medians per revision from alternating base/head runs, after compilation completed. Kernels were compiled in a separate translation unit without LTO. Normalization was slower in both comparisons.

Performance is operation- and platform-dependent: Vec4f scalar division became faster in this test, and matrix multiplication did not show a comparably large slowdown. These measurements establish a CPU normalization regression, not a blanket slowdown or a GPU result. No CUDA hardware benchmark was run.

Requested change: specialize the integer correction while preserving the existing floating-point path, then benchmark representative CUDA normalization and transform workloads. Also preserve the original reduction evaluation where practical; the new reductions introduce an initial addition to zero.

3. [P2] Preserve runtime FMA semantics and the order of accumulation

Location: Math.h:2144, with analogous changes in the other matMult/matMultT overloads.

The rewrite changes both fusion guarantees and the association of additions: the old nested FMA expression accumulates from the right, while the replacement expression accumulates from the left. Automatic FMA contraction does not restore that order.

Concrete runtime test: use matrix rows (1,1,1), (0,1,0), (0,0,1) and input (1,1e8f,-1e8f). The base returns first component 1, the PR returns 0. This reproduces with -O3 -ffp-contract=on; it is not limited to constant evaluation. The matrix is invertible and all inputs are finite.

The description's claim that drift is at most one ulp and host/device code generation remains identical is therefore unsupported. These helpers feed Map transformations used for world/index coordinates and sampling.

Requested change: retain the original nesting and explicit FMA operations at runtime. Make constant evaluation a separate supported path, or defer constexpr on these operations. Add cancellation-sensitive tests rather than only small, well-conditioned inputs.

4. [P2] Avoid introducing a rounding step before rounding the coordinate

Location: Math.h:317.

Adding 0.5f can round the float before floorf sees it. Two examples:

Input	Base float overload	PR float overload	PR double overload for the same value
8388609.f (exact integer)	8388609	8388610	8388609
0.4999999701976776123046875f (next float below 0.5)	0	1	0

Thus the change affects inputs beyond the intended half-integer tie cases and still disagrees across float/double. The order-zero nearest-neighbor sampler calls this function, so it can select an adjacent voxel for these inputs.

Requested change: implement the desired tie rule without this intermediate rounding error. Test exact large integers and floats immediately adjacent to half-integers. The related member rounding helpers should use the same corrected implementation where appropriate.

5. [P2] Snapshot a scalar divisor that aliases the vector

Location: Math.h:1052–1053.

For Vec3f v(2,4,6); v /= v[0];, the base returns (1,2,3) and the PR returns (1,4,6). The scalar is a reference to the first component; the loop overwrites it with 1 before processing the remaining components. The old reciprocal was computed before any mutation.

Requested change: copy the divisor before the loop or accept it by value. Cover Vec2, Vec3 and Vec4. The newly added matrix /= helper has the same aliasing hazard.

6. [P2] Accumulate mixed-type dot products in the promoted type

Location: Math.h:1067–1069.

The old expression performed arithmetic in the promoted operand type and converted once on return. The new loop stores each intermediate sum in T. Vec3i(1,1,1).dot(Vec3d(.4,.4,.4)) returns 1 in the base and 0 in the PR, because each partial sum is truncated. A float receiver with double operands likewise loses intermediate precision that it previously retained.

Requested change: use the promoted product/sum type for the accumulator and convert at return, seeding from the first product where appropriate.

Source API compatibility inventory

All six small compile probes in the attached bundle compile against the base and fail against the PR. These are observable compatibility changes, regardless of whether the new API policy is preferred.

Change	Previously valid usage affected	Suggested handling
Converting constructors become explicit	Passing an external vector or Vec3d to an argument of type Vec3f	Preserve implicit interop, or clearly release-note the change and provide migration examples. The PR acknowledges this change.
All eight Vec/Mat leaf classes become final	Deriving a small utility type from Vec3f, or extending a matrix type	Remove final for compatibility. These classes have no virtual dispatch to devirtualize. This break is absent from the PR's compatibility summary.
Scalar operators become hidden friends	Qualified nanovdb::math::operator*(2.f, v) calls and explicit function-template use	Preserve a namespace-visible compatibility surface if supported. Ordinary 2.f * v expressions continue to work.
Mat4–Vec4 free multiplication becomes a member	Qualified nanovdb::math::operator*(m, v) calls	Assess or preserve the old callable interface; ordinary m * v still works.
Vec2–Coord arithmetic now takes Coord2	Code intentionally using the first two components of a 3D coordinate	Prefer a deprecation/migration policy if source compatibility is required. The PR acknowledges this change.

The dimensionally invalid Mat3–Mat2x3 overload is also removed. That removal is documented and dimensionally justified. It should still be distinguished from API additions.

Verification and scope

• Compiled both revisions' headers and executed independent layout/arithmetic reproducers.
• Compared C++ offsets with PNanoVDB's own constants using a byte-buffer fixture; this was not a full historical .nvdb file round trip.
• Compiled six source-compatibility probes against both revisions.
• Examined optimized ARM64 assembly and benchmarked eight existing operations.
• Built the PR's TestNanoVDB.cc with GoogleTest v1.17.0. 21 selected tests passed, covering the new math tests, vectors, matrices, Map, coordinates, vector sampling, leaves/roots and the existing CNanoVDB/PNanoVDB checks. These checks do not catch the demonstrated edge cases.
• Inspected the full changed-file list, surrounding math code, call sites, PR discussion and current checks. No full OpenVDB build, CUDA execution, end-to-end rendering benchmark, or downstream application build was performed.
• No review or comments were posted to GitHub.

Reproducing the checks

Unzip pr2207-reproducers.zip. Obtain two OpenVDB source trees at the commits listed above. The included script takes repository roots, each containing the nanovdb/ directory:

python3 reproduce.py --base /path/to/base-checkout --head /path/to/pr-checkout --out results --benchmark

It requires Python 3 and Clang with C++17 support. The standalone probes need no OpenVDB libraries, GoogleTest, CUDA or TBB. Omit --benchmark for just the layout, arithmetic and compilation checks. The bundle also contains the observed outputs, selected GoogleTest results and generated assembly.

Suggested review decision

Resolve the Vec4 format break first. Keep the integer correctness fixes separate from floating-point algorithm changes, preserve existing public API unless a breaking release is intended, and add the failing regression cases above. Then require representative CUDA performance measurements before accepting claims of unchanged GPU performance.

[kmuseth]

Overall this looks good, but I have some concerns regarding alignment ant performance (see review above by Astra)

The MatBase/VecBase consolidation reads cleanly, and the bug fixes are genuine, pre-existing issues worth having caught: Min/Max for integer types routing through fminf/fmaxf, Vec/Mat integer operator/ returning zeros, Mat2::inverse() returning uninitialized data on singular input, and Vec2 arithmetic silently dropping z when given a 3D Coord.

I've left five inline comments. Two I'd like your take on before this merges — the matMult accumulation-order claim, and the voxelToWorld change in the CUDA tests, which I think alters which code path is under test rather than just the conversion. The other three are minor (a stale doc comment, a pre-existing dead assertion, and a dead store).

Verification caveat: locally I built and ran nanovdb_test_nanovdb (163/163 pass) and confirmed sizeof/alignof parity for every Vec/Mat type against master, so the ABI-visible layout is unchanged. However, I have no nvcc on this machine, so I could not compile or run any of the CUDA translation units. My comments on TestNanoVDB.cu are from reading the code only — in particular the claim that voxelToWorld flips those three assertions from the double applyMap path to the float applyMapF path is derived from the if constexpr dispatch in PointsToGrid.cuh, not from observed behaviour. Would appreciate a second pair of eyes there from someone who can build the CUDA targets.

[kmuseth]

The claim that "device-side codegen is unchanged in practice" (and the PR description's "Host codegen unchanged under default -ffp-contract=on" / "only compile-time constexpr evaluation observes the double-rounded result") isn't quite accurate, because this rewrite changed more than just FMA usage — it also reversed the accumulation order:

// before: right-to-left, FMA at each step
fmaf(x, m[0], fmaf(y, m[1], z * m[2]))
// after: left-to-right sum
x * m[0] + y * m[1] + z * m[2]

Contraction can re-fuse the multiplies, but it can't restore the original association, so results differ at runtime, not only under constexpr. Over 500k random inputs, -O2 with contraction enabled:

215222/500000 differ (43.0%), worst absolute err 3.1e-05

To be clear, I don't think this blocks the change — the conclusion holds. Typical drift really is ~1 ulp, and for the case that actually matters here, a uniform-scale Map, results are bit-identical (0 differing out of 500k). It's just the justification that's overstated. Suggest rewording to something like: "accumulation order changes from right-to-left to left-to-right, so results may differ by ~1 ulp on host and device; exact for uniform-scale maps."

[harrism]

We can have our cake and eat it in this case.

Compilers honor parentheses as an association barrier (they can't reassociate across them without -ffast-math/-fassociative-math), so the robust, constexpr-eligible solution is to write each component as

x * mat[0] + (y * mat[1] + z * mat[2])

This restores the right-to-left association, and with contraction enabled the compiler fuses it back into exactly the original fmaf(x, m0, fmaf(y, m1, z*m2)). Measured with clang -O2, 500k random inputs vs the base fmaf form:

┌────────────────────────┬─────────────────────┬────────────────────────────────┐
│     -ffp-contract      │ current PR (a+b+c)  │         parenthesized          │
├────────────────────────┼─────────────────────┼────────────────────────────────┤
│ on (clang/gcc default) │ 215042 differ (43%) │                       0 differ │
├────────────────────────┼─────────────────────┼────────────────────────────────┤
│ fast                   │       215042 differ │                       0 differ │
├────────────────────────┼─────────────────────┼────────────────────────────────┤
│ off                    │       223487 differ │ 178257 differ (~1 ulp, no FMA) │
└────────────────────────┴─────────────────────┴────────────────────────────────┘

Notes (from Fable 5.1):

• Astra's cancellation case mat=(1,1,1)·(1,1e8,-1e8): base=1, current PR=0, parenthesized=1 in all three modes — the catastrophic case is an association bug, so parentheses fix it regardless of FMA.
• Bit-identity to the old code relies on contraction being on. That's the default for nvcc (-fmad=true), GCC (-ffp-contract=fast), and clang (on), but MSVC /fp:precise and anyone building with -ffp-contract=off get a double-rounded result — ~1 ulp drift, same as the PR description already claims. That's an honest version of the "codegen unchanged" statement.
• Same change is needed in all six matMult/matMultT overloads (float/double, and the Mat3-typed variants).

[kmuseth]

This changes which code path is tested, not just the conversion. voxelToWorld is template and dispatches internally:

if constexpr(util::is_same<Vec3T,Vec3d>::value) return map.applyMap(...);   // double
else                                            return map.applyMapF(...);  // float

Vec3T is nanovdb::Vec3f in this test, so adding the explicit argument flips all three sites from the double applyMap path to the float applyMapF path — while the tolerances stay at 1e-9 / 1e-6 / 1e-2.

The 1e-9 here is the concern: the sphere has radius 100 and float carries ~1e-7 relative precision, so this assertion is only satisfiable because createPointSphere(8, 100.0, Vec3d(0.0), 0.5) yields an exactly-representable uniform 0.5 scale. Coverage of the double path is silently dropped, and the test would fail for any non-power-of-two voxel size.

Since the goal was just to satisfy the now-explicit ctor, converting the result preserves the original behaviour:

test = (begin[i] - Vec3T(nanovdb::voxelToWorld(voxel, ijk, grid->map()))).length() < 1e-9;

Same applies at lines 2205 and 2334.

[kmuseth]

The doc comment says "Implicit to preserve foreign-type interop" directly above a ctor declared explicit. Since making this ctor explicit is the PR's one source-breaking change, this is the comment most likely to be read by someone hitting the resulting compile error — worth getting right.

Same contradiction at lines 1181 (Vec2) and 2025 (Vec4).

[kmuseth]

Drive-by, not caused by this PR: test is computed here but never asserted — there's no EXPECT_TRUE(test) after the loop, so this check is dead. Every sibling site pairs them (1581/1583, 1699/1704, 1818/1823, 1942/1947, 2069/2080, 2332/2336); only this Voxel16 one is missing, and it's already like that in master. Since you're touching this line anyway, might be worth fixing.

[kmuseth]

Minor: Derived out{}; value-initializes before the loop overwrites every element, so the zero-fill is a dead store on some of the hottest paths in the library. Compilers usually elide it at -O2, but it's avoidable. Applies to the whole family of helpers here and in MatBase (968, 975, 982, 989, 996, 1003, 1011, 1379, 1387, 1395).

[harrism]

@kmuseth regarding your Astra point 1 (keep vec4 default alignment). This change was based on my original review comment. This is classic CUDA optimization. The way the code is written, with vec4 4-byte-aligned, every load and store of a vec4 will result in 4x 32-bit load/store instructions. If the compiler can see it is 16-byte aligned, it will instead generate a single 128-bit load/store, which is what we want. We can't get close to peak GPU bandwidth with serialized vec4 loads. (Note for a Mat4x4, 4 128-bit loads become 16 32-bit loads, and for all the other types we get some fraction of this vectorization (e.g. vec4d, vec2f, etc.))

My mistake in the May review was not thinking about format layout and data versioning. We can't break that, as you point out.

The good news is there is a cake-and-eat-it solution for this too.

1. At a major release we should consider adding alignas() to the base types as in b2cf24c and updating the PNanoVDB constants table.
2. For now (Claude helped with the following plan):
   A. Leave Vec4/Vec2/Mat4 without alignas (revert that part of b2cf24c), with the same "this is on-disk layout" comment the Vec3 family already has.
   B. Add an aligned load helper and use it in LeafData::getValue, InternalData::getValue/tile reads, and the value-array iterators — the paths that dominate kernel time:

template<typename ValueT>
__hostdev__ inline ValueT loadAligned(const ValueT* p) {
    if constexpr (sizeof(ValueT) == 16 || sizeof(ValueT) == 32)
        return *static_cast<const ValueT*>(__builtin_assume_aligned(p, sizeof(ValueT)));
    else
        return *p;
}

__builtin_assume_aligned is supported by nvcc, clang, and gcc; for MSVC fall back to __assume((uintptr_t)p % 16 == 0) or plain load. With that, nvcc emits LD.E.128 for Vec4f and 2×128 for Vec4d. Alternatively an opt-in Vec4fAligned (alignas(16), same storage) that the accessor reinterprets through — same effect, slightly more type noise (and reinterpret_cast smell).
3. Pin the assumption with static_asserts in LeafData/InternalData: offsetof(..., mValues) % 16 == 0, sizeof(DataType) % 32 == 0 (already there), so nobody breaks the invariant silently.
4. For user-owned arrays (e.g. PointsToGrid inputs, blind data), document that Vec4 is naturally aligned and offer Vec4Aligned as an opt-in if they want vectorized loads — that's my Vec3Aligned idea from the May review generalized.