Using Implicit MPM

Using Implicit MPM#

Newton’s implicit Material Point Method (MPM) solver models particle materials such as granular media. MPM support and rigid-MPM coupling are experimental. Start with the compact scripts/demos/mpm/newton_mpm_granular.py example; snowball_smash.py adds coupling and teapot_fill.py adds cavity sampling.

Train and Regenerate Franka Pour#

The IsaacContrib-Franka-Pour task restores episodes from a reset artifact containing the connected 14-phase reset distribution. The canonical 20,000-row artifact downloads from the standard Isaac Lab asset root on first use, so training needs no artifact setup:

uv run isaaclab train --rl_library rsl_rl --task IsaacContrib-Franka-Pour \
  --num_envs 2048 --device cuda:0

The checked-in generator remains the executable reference for reproducing or customizing the distribution. It takes about two minutes on an L40S-class GPU and writes a local artifact that can be selected explicitly:

uv run python scripts/tools/generate_franka_pour_reset_dataset.py --device cuda:0
uv run isaaclab train --rl_library rsl_rl --task IsaacContrib-Franka-Pour \
  --num_envs 2048 --device cuda:0 \
  env.reset_dataset_path=datasets/franka_pour/reset_dataset.pt

The task validates the payload’s stored content digest automatically. Setting ISAACSIM_ASSET_ROOT redirects the canonical artifact to a compatible local or self-hosted asset tree. Digest pinning remains available for custom reproducible experiments.

The source uses a non-colliding analytic fill volume whose height is controlled by env.source_fill_level in (0, 1]. The default 0.70 produces a 735-particle jittered lattice up to roughly 70% of the cup height. env.pour_target_frac independently controls the fraction of that live payload that must reach the receiver.

Play a checkpoint in Kit with the canonical task configuration; no external callback or particle override is required:

uv run isaaclab play --rl_library rsl_rl --task IsaacContrib-Franka-Pour \
  --checkpoint /path/to/model.pt --num_envs 1 --device cuda:0 --visualizer kit

Minimal Setup#

Use the same voxel size for the solver grid and particle generator. Add the generated object to an InteractiveSceneCfg like any other declarative asset.

import isaaclab.sim as sim_utils
from isaaclab_newton.assets import MPMObjectCfg
from isaaclab_newton.physics import MPMSolverCfg, NewtonCfg
from isaaclab_newton.sim.spawners.mpm import MPMGridCfg

voxel_size = 0.02

sim_cfg = sim_utils.SimulationCfg(
    dt=1.0 / 100.0,
    physics=NewtonCfg(
        solver_cfg=MPMSolverCfg(
            voxel_size=voxel_size,
            max_iterations=100,
            tolerance=1.0e-4,
        ),
        num_substeps=2,
    ),
)

media = MPMObjectCfg(
    prim_path="{ENV_REGEX_NS}/Media",
    spawn=MPMGridCfg(
        lower=(-0.1, -0.1, 0.0),
        upper=(0.1, 0.1, 0.2),
        voxel_size=voxel_size,
        particles_per_cell=2.0,
        particle_placement="cell_center",
    ),
)

Tune the particle material separately through MPMParticleMaterialCfg. The implicit solve already resolves collider contact. project_outside_colliders adds a hard post-step correction for particles that remain inside colliders; use it only when that geometric correction is intentional, and not as a substitute for valid initial states, collision geometry, or a stable timestep. Coupled MPM entries do not support this manager-level projection pass.

Render a Particle Surface#

MPM simulation state remains a set of particles. Surface reconstruction is an optional visualization pass: it does not change particle motion, collisions, or material behavior. Run the teapot example to compare the available modes:

# Reconstructed surface (default)
uv run python scripts/demos/mpm/teapot_fill.py --device cuda:0 \
  --visualizer newton_gl --fluid_render_mode surface
# Surface and source particles together
uv run python scripts/demos/mpm/teapot_fill.py --device cuda:0 \
  --visualizer newton_gl --fluid_render_mode both
# Path-traced translucent surface
uv run --extra ovrtx python scripts/demos/mpm/teapot_fill.py --device cuda:0 \
  --visualizer newton_rtx --fluid_render_mode surface

Surface rendering is available in the Newton GL and Newton RTX visualizers. The Kit visualizer continues to render the MPM particles directly.

To reconstruct a surface in another Newton MPM script, create one reusable newton.geometry.ParticleSurface after sim.reset(). On each render update, extract from the current Newton particle positions, radii, flags, and world indices, then pass the returned vertex, triangle-index, and normal arrays to NewtonGLVisualizer.log_mesh() or NewtonRTXVisualizer.log_mesh() with dynamic=True before calling sim.render(). The visualizer stages the latest mesh by name and submits it inside Newton’s required viewer-frame lifecycle.

Tune reconstruction independently from the simulation:

  • voxel_size controls surface detail and memory use. It can be smaller than the MPM solver voxel size.

  • kernel_radius controls how far each particle contributes to the surface. Start near three times the particle spacing.

  • max_grid_cells provides fixed-capacity storage for CUDA graph capture. Increase it if the reconstructed domain outgrows the reserved grid.

  • Anisotropic kernels preserve sheets and stretched fluid features better, but cost more than isotropic kernels.

The demo handles CUDA graph capture, empty surfaces, inactive particles, and dynamic topology in one reusable helper:

FluidSurfaceRenderer implementation
class FluidSurfaceRenderer:
    """Extract and display a dynamic water surface in Newton visualizers."""

    def __init__(self, sim) -> None:
        import warp as wp
        from isaaclab_newton.physics import NewtonManager
        from isaaclab_visualizers.newton import NewtonGLVisualizer, NewtonRTXVisualizer
        from newton.geometry import ParticleSurface

        self._wp = wp
        self._visualizers = tuple(
            visualizer
            for visualizer in sim.visualizers
            if isinstance(visualizer, (NewtonGLVisualizer, NewtonRTXVisualizer))
        )
        if not self._visualizers:
            raise RuntimeError("Particle surface rendering requires a Newton GL or RTX visualizer.")

        self._model = NewtonManager.get_model()
        self._state = NewtonManager.get_state_0()
        self._surface = ParticleSurface(
            voxel_size=SURFACE_VOXEL_SIZE,
            max_grid_cells=SURFACE_MAX_GRID_CELLS,
            world_count=max(self._model.world_count, 1),
            kernel_radius=SURFACE_KERNEL_RADIUS,
            threshold=0.4,
            smooth_lambda=0.0,
            anisotropic=True,
            kernel_scale=0.5,
            anisotropy_ratio=16.0,
            anisotropy_scale=1.0,
            anisotropy_min_neighbors=4,
            anisotropy_binning=True,
            anisotropy_strength=0.95,
            field_smooth_iterations=0,
            mesh_smooth_iterations=1,
            device=self._model.device,
        )
        self._empty_points = wp.empty(0, dtype=wp.vec3, device=self._model.device)
        self._empty_indices = wp.empty(0, dtype=wp.int32, device=self._model.device)
        self._empty_normals = wp.empty(0, dtype=wp.vec3, device=self._model.device)
        self._surface_mesh = None
        self._surface_graph = None
        self._capture_surface_extraction()

    def _extract_surface(self):
        """Extract the water surface from the current Newton particle state."""
        return self._surface.extract(
            self._state.particle_q,
            self._model.particle_radius,
            particle_flags=self._model.particle_flags,
            particle_world=self._model.particle_world if self._surface.world_count > 1 else None,
        )

    def _capture_surface_extraction(self) -> None:
        """Capture reconstruction separately from the MPM physics graph."""
        if not self._model.device.is_cuda or args_cli.disable_cuda_graph:
            return
        self._surface_mesh = self._extract_surface()
        with self._wp.ScopedCapture(device=self._model.device) as capture:
            self._surface_mesh = self._extract_surface()
        self._surface_graph = capture.graph

    def update(self) -> int:
        """Reconstruct and publish the current water surface, returning its triangle count."""
        if self._surface_graph is None:
            self._surface_mesh = self._extract_surface()
        else:
            self._wp.capture_launch(self._surface_graph)

        vertices, indices, normals = self._surface_mesh.to_arrays()
        if vertices is None:
            vertices = self._empty_points
            indices = self._empty_indices
            normals = self._empty_normals
            hidden = True
            triangle_count = 0
        else:
            hidden = False
            triangle_count = indices.shape[0] // 3

        for visualizer in self._visualizers:
            visualizer.log_mesh(
                SURFACE_PATH,
                vertices,
                indices,
                normals=normals,
                hidden=hidden,
                backface_culling=False,
                color=WATER_COLOR,
                roughness=0.1,
                metallic=0.0,
                dynamic=True,
                opacity=WATER_OPACITY,
            )
        return triangle_count

Tune Resolution, Time, Then Convergence#

Tune one group at a time in this order:

  1. Voxel and particle resolution. MPMSolverCfg.voxel_size controls the background grid. Smaller voxels resolve thinner geometry but increase active cells and memory. MPMGridCfg.particles_per_cell controls particle density; doubling it along each axis creates about eight times as many particles in 3D. Start coarse, then refine until the measured behavior stops changing.

  2. Timestep and substeps. Each Newton substep uses SimulationCfg.dt / NewtonCfg.num_substeps. Reduce dt or increase num_substeps first when contacts tunnel, jitter, or become unstable. Substeps do not change the policy period, which also includes environment decimation.

  3. Iterations and tolerance. MPMSolverCfg.max_iterations caps the rheology solve; tolerance permits an earlier exit after convergence. Increase the cap only when the solver reaches it, and lower the tolerance only when tighter convergence improves a physical metric. These settings do not repair an unstable timestep, invalid reset, or incorrect collider.

Tune Rigid-MPM Coupling#

For CouplerProxyCfg, first stabilize each solver alone. Then tune the additional controls:

  • CouplerEntryCfg.substeps divides one coupled step for that entry. Increase the MPM entry’s value when only the particle solve needs a smaller timestep.

  • CouplerProxyCfg.iterations repeats the proxy exchange and relaxation; it does not replace smaller physical timesteps.

  • CouplerProxyMappingCfg.mass_scale scales the source body’s effective mass and inertia only in the destination proxy view. It does not change the body’s authored mass in the rigid solver.

Start mass_scale at 1 for a freely moving collider. Increase it when the rigid solver strongly constrains the collider during MPM contact. For example, a cup resting on a table has much greater effective resistance in the supported direction than its free-body mass suggests. Sweep finite values geometrically, such as 1, 10, and 100, and keep the smallest value that prevents unrealistic proxy motion. Newton requires a finite positive value: do not use infinity. An excessively large scalar also suppresses legitimate motion in unsupported directions and can make the interaction effectively one-way.

Validate both the supported and free-moving cases after changing coupling. If the uncoupled systems are unstable, fix their timestep, contacts, and reset states before adjusting mass_scale or coupling iterations.