newton.solvers.SolverImplicitMPM#

class newton.solvers.SolverImplicitMPM(model, config, temporary_store=None, verbose=None, enable_timers=False)[source]#

Bases: SolverBase

Implicit MPM solver for granular and elasto-plastic materials.

Implements an implicit Material Point Method (MPM) algorithm roughly following [1], extended with a GPU-friendly rheology solver supporting pressure-dependent yield (Drucker-Prager), viscosity, dilatancy, and isotropic hardening/softening.

This variant is particularly well-suited for very stiff materials and the fully inelastic limit. It is less versatile than traditional explicit MPM but offers unconditional stability with respect to the time step.

Call register_custom_attributes() on your ModelBuilder before building the model to enable the MPM-specific per-particle material parameters and state variables (e.g. mpm:young_modulus, mpm:friction, mpm:particle_elastic_strain).

[1] https://doi.org/10.1145/2897824.2925877

Parameters:

  • model (Model) – The model to simulate.

  • config (Config) – Solver configuration. See SolverImplicitMPM.Config.

  • temporary_store (TemporaryStore | None) – Optional Warp FEM temporary store for reusing scratch allocations across steps.

  • verbose (bool | None) – Enable verbose solver output. Defaults to wp.config.verbose.

  • enable_timers (bool) – Enable per-section wall-clock timings.

class Config(max_iterations=250, tolerance=0.0001, solver='gauss-seidel', warmstart_mode='auto', collider_velocity_mode='forward', voxel_size=0.1, grid_type='sparse', grid_padding=0, max_active_cell_count=-1, transfer_scheme='apic', integration_scheme='pic', critical_fraction=0.0, air_drag=1.0, collider_normal_from_sdf_gradient=False, collider_basis='Q1', strain_basis='P0', velocity_basis='Q1')#

Bases: object

Configuration for SolverImplicitMPM.

Per-particle properties can be configured using custom attributes on the Model. See SolverImplicitMPM.register_custom_attributes() for details.

__init__(max_iterations=250, tolerance=0.0001, solver='gauss-seidel', warmstart_mode='auto', collider_velocity_mode='forward', voxel_size=0.1, grid_type='sparse', grid_padding=0, max_active_cell_count=-1, transfer_scheme='apic', integration_scheme='pic', critical_fraction=0.0, air_drag=1.0, collider_normal_from_sdf_gradient=False, collider_basis='Q1', strain_basis='P0', velocity_basis='Q1')#
air_drag: float = 1.0#

Numerical drag for the background air.

collider_basis: str = 'Q1'#

P0 (piecewise constant), Q1 (trilinear), S2 (quadratic serendipity), pic8 (particle-based with max 8 points per cell)

Type:

Collider basis function string. Examples

collider_normal_from_sdf_gradient: bool = False#

Compute collider normals from sdf gradient rather than closest point

collider_velocity_mode: str = 'forward'#

Collider velocity computation mode, may be one of ‘forward’ or ‘backward’. ‘forward’ uses the current velocity, ‘backward’ uses the previous timestep position. Deprecated aliases ‘instantaneous’ (=’forward’) and ‘finite_difference’ (=’backward’) are also accepted.

critical_fraction: float = 0.0#

Fraction for particles under which the yield surface collapses.

grid_padding: int = 0#

Number of empty cells to add around particles when allocating the grid.

grid_type: str = 'sparse'#

Type of grid to use. May be one of sparse, dense, fixed.

integration_scheme: str = 'pic'#

Integration scheme controlling shape-function support. May be one of pic, gimp.

max_active_cell_count: int = -1#

Maximum number of active cells to use for active subsets of dense grids. -1 means unlimited.

max_iterations: int = 250#

Maximum number of iterations for the rheology solver.

solver: str = 'gauss-seidel'#

Solver to use for the rheology solver. May be one of gauss-seidel, jacobi, cg.

strain_basis: str = 'P0'#

Strain basis functions. May be one of P0, P1d, Q1, Q1d, or pic[n].

tolerance: float = 0.0001#

Tolerance for the rheology solver.

transfer_scheme: str = 'apic'#

Transfer scheme to use for particle-grid transfers. May be one of apic, pic.

velocity_basis: str = 'Q1'#

B2 (quadratic b-spline), Q1 (trilinear)

Type:

Velocity basis function string. Examples

voxel_size: float = 0.1#

Size of the grid voxels.

warmstart_mode: str = 'auto'#

Warmstart mode to use for the rheology solver. May be one of none, auto, particles, grid, smoothed.

classmethod register_custom_attributes(builder)#

Register MPM-specific custom attributes in the ‘mpm’ namespace.

This method registers per-particle material parameters and state variables for the implicit MPM solver.

Attributes registered on Model (per-particle):

  • mpm:young_modulus: Young’s modulus in Pa

  • mpm:poisson_ratio: Poisson’s ratio for elasticity

  • mpm:damping: Elastic damping relaxation time in seconds

  • mpm:friction: Friction coefficient

  • mpm:yield_pressure: Yield pressure in Pa

  • mpm:tensile_yield_ratio: Tensile yield ratio

  • mpm:yield_stress: Deviatoric yield stress in Pa

  • mpm:hardening: Hardening factor for plasticity

  • mpm:hardening_rate: Hardening rate for plasticity

  • mpm:softening_rate: Softening rate for plasticity

  • mpm:dilatancy: Dilatancy factor for plasticity

  • mpm:viscosity: Viscosity for plasticity [Pa·s]

Attributes registered on State (per-particle):

  • mpm:particle_qd_grad: Velocity gradient for APIC transfer

  • mpm:particle_elastic_strain: Elastic deformation gradient

  • mpm:particle_Jp: Determinant of plastic deformation gradient

  • mpm:particle_stress: Cauchy stress tensor [Pa]

  • mpm:particle_transform: Overall deformation gradient for rendering

__init__(model, config, temporary_store=None, verbose=None, enable_timers=False)#
collect_collider_impulses(state)#

Collect current collider impulses and their application positions.

Returns a tuple of 3 arrays:

  • Impulse values in world units.

  • Collider positions in world units.

  • Collider id, that can be mapped back to the model’s body ids using the collider_body_index property.

notify_model_changed(flags)#
project_outside(state_in, state_out, dt, gap=None)#

Project particles outside of colliders, and adjust their velocity and velocity gradients

Parameters:

  • state_in (State) – The input state.

  • state_out (State) – The output state. Only particle_q, particle_qd, and particle_qd_grad are written.

  • dt (float) – The time step, for extrapolating the collider end-of-step positions from its current position and velocity.

  • gap (float | None) – Maximum distance for closest-point queries. If None, the default is the voxel size times sqrt(3).

sample_render_grains(state, grains_per_particle)#

Generate per-particle point samples used for high-resolution rendering.

Parameters:

  • state (State) – Current Newton state providing particle positions.

  • grains_per_particle (int) – Number of grains to sample per particle.

Returns:

A wp.array with shape (num_particles, grains_per_particle) of type wp.vec3 containing grain positions.

Return type:

array

setup_collider(collider_meshes=None, collider_body_ids=None, collider_margins=None, collider_friction=None, collider_adhesion=None, collider_projection_threshold=None, model=None, body_com=None, body_mass=None, body_inv_inertia=None, body_q=None)#

Configure collider geometry and material properties.

By default, collisions are set up against all shapes in the model with newton.ShapeFlags.COLLIDE_PARTICLES. Use this method to customize collider sources, materials, or to read colliders from a different model.

Parameters:

  • collider_meshes (list[Mesh] | None) – Warp triangular meshes used as colliders.

  • collider_body_ids (list[int] | None) – For dynamic colliders, per-mesh body ids.

  • collider_margins (list[float] | None) – Per-mesh signed distance offsets (m).

  • collider_friction (list[float] | None) – Per-mesh Coulomb friction coefficients.

  • collider_adhesion (list[float] | None) – Per-mesh adhesion (Pa).

  • collider_projection_threshold (list[float] | None) – Per-mesh projection threshold (m).

  • model (Model | None) – The model to read collider properties from. Default to solver’s model.

  • body_com (array | None) – For dynamic colliders, per-body center of mass.

  • body_mass (array | None) – For dynamic colliders, per-body mass. Pass zeros for kinematic bodies.

  • body_inv_inertia (array | None) – For dynamic colliders, per-body inverse inertia.

  • body_q (array | None) – For dynamic colliders, per-body initial transform.

step(state_in, state_out, control, contacts, dt)#

Advance the simulation by one time step.

Transfers particle data to the grid, solves the implicit rheology system, and transfers the result back to update particle positions, velocities, and stress.

Parameters:

  • state_in (State) – Input state at the start of the step.

  • state_out (State) – Output state written with updated particle data. May be the same object as state_in for in-place stepping.

  • control (Control) – Control input (unused; material parameters come from the model).

  • contacts (Contacts) – Contact information (unused; collisions are handled internally).

  • dt (float) – Time step duration [s].

update_particle_frames(state_prev, state, dt, min_stretch=0.25, max_stretch=2.0)#

Update per-particle deformation frames for rendering and projection.

Integrates the particle deformation gradient using the velocity gradient and clamps its principal stretches to the provided bounds for robustness.

update_render_grains(state_prev, state, grains, dt)#

Advect grain samples with the grid velocity and keep them inside the deformed particle.

Parameters:

  • state_prev (State) – Previous state (t_n).

  • state (State) – Current state (t_{n+1}).

  • grains (array) – 2D array of grain positions per particle to be updated in place. See sample_render_grains.

  • dt (float) – Time step duration.

property collider_body_index: array#

Array mapping collider indices to body indices.

Returns:

Per-collider body index array. Value is -1 for colliders that are not bodies.

property voxel_size: float#

Grid voxel size used by the solver.