voxelize

GPU voxelization for triangle meshes and point clouds, loadable through kernels. Meshes rasterize into bit-packed occupancy grids along a conservative surface route and a parity-filled solid route; point clouds bin into occupancy or pooled per-voxel features. The reference baselines are analytic voxel counts, trimesh containment, cuda_voxelizer, and a sequential fp32 pooling reference matched bitwise. All grid outputs are bitwise deterministic: the same inputs produce the same bytes on every card, in any triangle or point order.

Turning geometry into voxels is the entry point for 3D learning, simulation, and spatial indexing, and doing it fast normally costs reproducibility: parallel rasterizers race on atomics, so two runs of the same mesh disagree in the low bits, and CPU tools that are exact are thousands of times slower. This kernel is both: billions of triangles and points per second, and a grid that is a pure function of the geometry, certified by pinned SHA-256 digests across four GPU architectures.

Two million points bin into a trefoil-knot grid at 5 billion points per second; an icosphere's surface voxelizes at 48, 96, and 192 cubed under an orbiting camera (327,680 triangles into 512 cubed in 5.9 ms end to end); a clip plane sweeps the parity-filled solid, exposing the full cross-section disc; and shuffling every triangle reproduces the identical bytes, torch.equal, with digests bit-identical across sm86, sm89, sm90, and sm120.

Usage

from kernels import get_kernel

vox = get_kernel("phanerozoic/voxelize", version=1, trust_remote_code=True)

grid = vox.voxelize_mesh(verts, faces, resolution=512, mode="solid")
grid.occupancy                            # bool [DZ, DY, DX]
grid.origin, grid.voxel_size, grid.dims   # world transform

occ = vox.voxelize_points(points, voxel_size=0.1)
pv = vox.pool_points(points, feats, voxel_size=0.1, reduce="mean")
pv.coords, pv.features, pv.counts         # int32 [M,3], fp32 [M,C], int32 [M]

version selects the release branch; trust_remote_code is required by kernels for publishers without the trusted-publisher mark. mode is 'surface', 'solid', or 'conservative' (their union); pass voxel_size (and optionally origin, dims) for caller-controlled grids; packed=True returns the raw int32 bit grid.

API

Method

The surface pass tests each triangle against every cell of its grid-space bounding box with the Akenine-Möller separating-axis test and sets occupied bits with atomicOr; triangles are bucketed by bounding-box cell count and processed by a thread, a warp, or a block, so large polygons cannot serialize a launch of small ones. The solid pass computes crossing parity per grid column: each triangle XORs one bit at the first voxel center above its intersection with every column inside its 2D projection, and a prefix-XOR converts crossing flags to occupancy. The projected inside test runs in snapped int64 fixed point with a top-left fill rule, so every boundary column is claimed by exactly one of the two triangles sharing an edge and watertight parity cannot leak. Feature pooling stable-sorts points by voxel key and reduces each voxel sequentially in that fixed order. Occupancy is written only through integer OR/XOR atomics, which commute exactly, and every floating-point predicate uses explicit round-per-op intrinsics, so grids are invariant to triangle order, point order, launch configuration, and GPU architecture.

Measured

Event-timed medians of the voxelization kernels (mesh on device, grid pre-allocated, matching how the compared tools report their own time). Surface voxelization against cuda_voxelizer v0.6, itself a deterministic integer-atomic method, so this is a pure speed comparison:

At 1.31M triangles the kernel sustains 7.7 billion triangles per second (11.5 for solid); trimesh (CPU) is three to four orders of magnitude slower. The same kernel on a GeForce RTX 3070 Ti laptop voxelizes the 1.31M-triangle mesh at 512³ in 0.74 ms, bit-identical to the workstation result.

Point pooling into voxels (mean, 32 channels), against the unique + index_add_ recipe:

index_add_ uses float atomics and returns different bytes on every run; making it reproducible runs about 2x slower than voxelize, which matches a sequential fp32 reference to the bit. Occupancy binning runs at 12 to 13 billion points per second.

Correctness

An axis-aligned box placed off voxel boundaries voxelizes to its exact analytic counts. An icosphere's parity fill matches its divergence-theorem volume to 2 parts in 10⁴ at 128³ and agrees with trimesh's ray-cast contains on every sampled voxel; a torus lands within a part in 10³. Pooled features match a sequential fp32 reference bitwise in mean, max, and sum modes. Shuffling faces or points changes no output bit. The suite pins SHA-256 digests of five canonical grids and two pooled outputs; the digests, minted on sm89, reproduce byte-for-byte on L4, H200, RTX PRO 6000, and a GeForce RTX 3070 Ti: four architecture generations.

Requirements and limits

• NVIDIA GPU with compute capability 8.0+; float32 geometry.
• Grids to 8192 voxels per axis, one geometry per call, single device.
• Solid and conservative modes assume a watertight mesh (crossing parity of an open mesh is undefined); surface mode takes any triangle soup, and degenerate triangles voxelize conservatively as their bounding cells.
• Mean and sum pooling are deterministic for a fixed input point order; max pooling is order-invariant outright.

References

Akenine-Möller, "Fast 3D Triangle-Box Overlap Testing" (2001); Schwarz and Seidel, "Fast Parallel Surface and Solid Voxelization on GPUs" (2010); cuda_voxelizer as the GPU reference tool.

License

Apache-2.0.