warp.geometry.sparse_marching_cubes#

warp.geometry.sparse_marching_cubes(
field,
nx,
ny,
nz,
*,
lower=None,
upper=None,
threshold=0.0,
lipschitz_bound=1.0,
device=None,
return_stats=False,
)[source]#

Extract an isosurface from an implicit function using a Lipschitz octree.

Rather than sampling a dense grid, this routine builds a sparse octree that provably brackets the level set of a 1-Lipschitz implicit function (such as a signed distance function) and runs marching cubes only on the near-surface cells, at a cost that scales with the surface area rather than the volume.

The grid is specified exactly as for warp.geometry.IsoSurfaceMarchingCubes.extract(): nx, ny, nz grid nodes over the box [lower, upper], which may be anisotropic. Calling this function and extract() with the same nx, ny, nz and bounds produces the same surface.

Both the pruning pass (at cell centers) and the extraction pass (at cell corners) query field in batches, never one point at a time. If field is implemented entirely with Warp (a kernel launch, with no host round trip), the whole pipeline stays on the GPU; see warp/examples/geometry/example_sparse_marching_cubes.py for a mesh-query implicit function that does this. A callable that wraps a host library (NumPy, PyTorch, …) is equally valid – it just pays a device/host sync on every call, since the values it returns must still land back on the query points’ device.

This function supports backward-mode automatic differentiation: if field’s output array has requires_grad=True (the caller’s responsibility to allocate) and the call is wrapped in a warp.Tape(), gradient flows from the output verts back through the corner values field produced. Cell selection (lipschitz_octree()) is not differentiated – it is a discrete, threshold-based search – and does not print tape warnings either.

Parameters:

  • field (Callable[[array], array]) – The implicit function, as a batched callable with the contract evaluate(points: warp.array[warp.vec3]) -> warp.array[warp.float32], returning the signed distance (or any 1-Lipschitz field whose threshold level set is the surface) at each query point. A bare single-point @warp.func is not accepted directly – wrap it in a kernel that evaluates it over the batch (see above), which keeps evaluation on the GPU instead of hiding a per-point kernel launch behind what looks like a single call.

  • nx (int) – Number of grid nodes in the x-direction.

  • ny (int) – Number of grid nodes in the y-direction.

  • nz (int) – Number of grid nodes in the z-direction.

  • lower (vec3f | tuple[float, float, float] | None) – The 3D coordinate that the grid’s corner at index (0, 0, 0) maps to. Defaults to (0.0, 0.0, 0.0) if None. Anything outside [lower, upper] is never visited, so parts of the level set that leave it are simply missing from the output, leaving the mesh open where it exits.

  • upper (vec3f | tuple[float, float, float] | None) – The 3D coordinate that the grid’s corner at index (nx - 1, ny - 1, nz - 1) maps to. Defaults to align with the grid’s maximal indices if None.

  • threshold (float) – The isovalue defining the surface.

  • lipschitz_bound (float) – An upper bound on the Lipschitz constant of field. Use 1.0 for a true signed distance function. Larger values widen the retained band, trading speed for a stronger guarantee when the field varies faster than unit rate.

  • device (Device | str | None) – The Warp device to run on. Defaults to the current device.

  • return_stats (bool) – If True, also return a dictionary of diagnostics (leaf-cell count, unique-corner count, implicit-function evaluation count, culled-cell count) useful for benchmarking against a dense grid.

Returns:

A tuple (vertices, indices) where vertices is a warp.array[warp.vec3] and indices is a flat warp.array[warp.int32] with three consecutive entries per triangle. If return_stats is True, returns (vertices, indices, stats).

Raises:

  • ValueError – If nx, ny, or nz is less than 2, or lipschitz_bound is negative.

  • TypeError – If field is not callable (including a bare single-point @warp.func, which must be wrapped in a batched kernel first).

Return type:

tuple[array, array] | tuple[array, array, dict]