VBD#
SolverVBD is a unified implicit solver based on
Vertex Block Descent (VBD): cloth, soft bodies, rigid bodies (joints,
contacts, rods), and their interactions are all solved in one VBD loop, with
joint and contact constraints enforced through an augmented-Lagrangian (ALM)
extension of the method. The compliant ALM constraint path is enabled
by default (rigid_compliant_alm=True). Pass rigid_compliant_alm=False
explicitly to retain the deprecated legacy path during its migration period.
Experimental
SolverVBD is experimental. Its public API,
behavior, feature support, performance, and implementation may change
without prior notice.
This page is the VBD backend guide: when to choose the solver, how to set it up correctly, and where to find runnable examples. For symptom-driven diagnosis and parameter tuning, see VBD Tuning and the Simulation Tuning landing page.
When to Choose VBD#
Start from the Supported Features matrix and the
Joint Feature Support tables in the solver overview. The underlying methods are described in the references
cited by SolverVBD:
Anka He Chen, Ziheng Liu, Yin Yang, and Cem Yuksel. 2024. Vertex Block Descent. ACM Trans. Graph. 43, 4. https://doi.org/10.1145/3658179
Chris Giles, Elie Diaz, and Cem Yuksel. 2025. Augmented Vertex Block Descent. ACM Trans. Graph. 44, 4. https://doi.org/10.1145/3731195 (the deprecated legacy rigid constraint path)
Anka H. Chen, Jerry Hsu, Youssef Ayman, and Miles Macklin. 2026. Divide and Truncate: A Penetration and Inversion Free Framework for Coupled Multi-physics System. In ACM SIGGRAPH 2026 Conference Papers. https://doi.org/10.1145/3799902.3811143 (penetration-free truncation, DAT)
Anka He Chen, Jerry Hsu, Ziheng Liu, Miles Macklin, Yin Yang, and Cem Yuksel. 2025. Offset Geometric Contact. ACM Trans. Graph. 44, 4. https://doi.org/10.1145/3731205 (feature regions of full-surface mesh contact, OGC)
API Explanations#
Construction and Stepping#
Colour the model before you finalize it. color()
fills particle_color_groups and
body_color_groups; the solver raises at construction when
particles are present without particle colours, or when rigid bodies are
integrated by VBD without body colours. Pass include_bending=True when the
model has bending edges. Then build the collision pipeline, its contacts buffer
and the solver, in that order:
builder.color(include_bending=True) # the model has bending edges
model = builder.finalize()
collision_pipeline = newton.CollisionPipeline(model)
contacts = collision_pipeline.contacts()
solver = newton.solvers.SolverVBD(model, iterations=10)
state_in, state_out = model.state(), model.state()
control = model.control()
for _ in range(num_steps):
collision_pipeline.collide(state_in, contacts)
solver.step(state_in, state_out, control, contacts, dt)
state_in, state_out = state_out, state_in
For rigid joints the solver captures the structural rest pose from
body_q at construction, so the body transforms must match the joint
angles when the solver is built; ROD rest bend and twist come from
joint_rest_q.
With a solver-owned pipeline, SolverVBD(model, collision_pipeline=...), the
solver allocates its own contacts buffer, seeds the pipeline’s self-contact
configuration from the particle_self_contact_* arguments, runs rigid
detection itself on the schedule of collision_frequency_type, and step
takes contacts=None. Rigid-soft DAT requires this form, and so does a
RIGID slot on ITERATIONS, which also needs contact matching in the
pipeline (below). The self-contact slot works in both
forms; Jitter: where to look first shows the arguments that put it on the
ITERATIONS schedule.
Contact matching between steps and between iterations. Each detection
pass of the pipeline can label every body-body contact with the index of the
same contact in the previous pass. contact_matching="latest" writes the
labels and keeps the fresh geometry; "sticky" also replays the saved
contact points and normal (experimental); the default "disabled" writes no
labels (Contact Matching). The solver uses the labels in two places.
Between steps, rigid_contact_history=True snapshots the solved state of
every contact at the end of the step (multiplier, penalty stiffness, normal)
and restores it onto the matched contacts at the next step’s contact refresh:
the normal multiplier along the new normal and, with "latest" matching,
the tangential multiplier projected onto the friction disk. Unmatched contacts
start cold, and so does every contact when two detection passes run between
two steps. Between iterations, a RIGID slot on the ITERATIONS schedule
re-detects inside the iteration loop and replaces the contact set in the
middle of the solve; the solver snapshots the in-flight state before each pass
and restores it onto the new set in the same way, so the multipliers built up
in the earlier iterations survive the refresh. The first use is opt-in; the
second is not. The constructor therefore raises ValueError when the
RIGID slot is on ITERATIONS and the pipeline has "disabled"
matching, and also when rigid_contact_history=True meets an owned pipeline
without matching.
For CUDA graph capture keep the same order, pipeline, contacts, solver, all
before capture. Body-body contact state is pre-sized from the pipeline’s
published rigid contact maximum and body-particle contact state from the
particle-shape pair count; both grow from Contacts on the first step.
During capture that growth needs Warp’s stream-ordered memory pool, otherwise
the solver raises and asks you to pre-size. rigid_contact_history=True
needs its warm-start buffers allocated before capture, on every device:
construct the CollisionPipeline before the solver, which
publishes rigid_contact_max, or run one uncaptured step first.
Start at dt = 1/600 s with 10 iterations and adjust following
Rules of Thumb. The VBD Tuning parameter table lists the
constructor arguments that matter most for tuning, with their
effect and cost; the signature of SolverVBD is the
complete list.
Model Inputs#
The solver reads three systems from the Model:
Cloth: particles with triangles and bending edges, from
add_cloth_mesh(),add_cloth_grid(), oradd_triangle()andadd_edge(). Material:tri_ke,tri_kaandtri_kdper triangle,edge_keandedge_kdper edge, andparticle_radiusper particle.validate_mesh=Trueruns the element-quality checker.Soft bodies: particles with tetrahedra, from
add_soft_mesh(),add_soft_grid()oradd_tetrahedron(). Material:k_mu,k_lambdaandk_dampper tetrahedron. The surface triangles take part in self-contact and rigid-soft contact.Rigid bodies:
add_body()with shapes, and joints of type BALL, FIXED, FREE, REVOLUTE, PRISMATIC, D6 and ROD (add_joint_rod()for cables). DISTANCE joints are not supported. Shape materialske,kd,muandmargincome fromShapeConfig. Drives and limits readjoint_target_ke,joint_target_kd,joint_limit_lower,joint_limit_upper,joint_limit_keandjoint_limit_kd.
Springs (add_spring()) and particle-body
attachments (add_attachment_particle_body()) connect
the systems. The contact material between particles and shapes mixes
soft_contact_ke, soft_contact_kd and soft_contact_mu on the Model
with the shape values (The energy terms). The Model Parameters table of
VBD Tuning gives the effect and typical values of each of these.
Rigid constraint rows use the compliant ALM formulation by default
(rigid_compliant_alm=True): the authored stiffness is the physical
stiffness of the row and the solver picks the metric \(\rho\) itself
(Rigid bodies and joints). The legacy AVBD path (rigid_compliant_alm=False)
is deprecated since 1.6, together with its controls: the per-joint
vbd:joint_is_hard custom attribute (model.vbd.joint_is_hard, declared
by register_custom_attributes()) and
set_joint_constraint_mode(). Under compliant
ALM they have no solver-mode effect.
Coupling Modes#
Inside one solver no coupling setup is needed. Rigid blocks and particle blocks are solved in the same iteration loop, rigid blocks first, and rigid-soft contacts and springs connect them as ordinary energy terms, and particle-body attachments as compliant augmented-Lagrangian rows (How it is solved). Cloth between the fingers of a VBD-integrated gripper is the plain case.
With integrate_with_external_rigid_solver=True the rigid bodies belong to
another solver. SolverVBD then skips the rigid solve, reads the body poses the
external solver wrote into the output state, and applies rigid-soft contact to
the particles only (one-way coupling). Attachments are not supported in this
mode. For two-way coupling through proxies, including contact forces reported
back to the source body, see Coupled Solvers.
Contact Handling#
The contact stack has three layers.
Detection. Self-contact pairs come from the triangle-mesh self-contact
detector: vertex-triangle and edge-edge candidates within
particle_self_contact_margin + particle_self_contact_gap, stored in fixed
buffers of particle_vertex_contact_buffer_size and
particle_edge_contact_buffer_size entries, and filtered by topological
distance (particle_topological_contact_filter_threshold) and rest-shape
distance (particle_rest_shape_contact_exclusion_radius). It runs on the
schedule of the SOFT_SELF_CONTACT collision slot, by default before and
after the inertial step of every substep. Rigid-soft records (particles, and
with enable_rigid_soft_full_surface_contact also soft edges and faces) and
body-body records come from the CollisionPipeline, which
detects out to its soft_contact_gap beyond particle_radius plus the
shape margin and carries the shape thicknesses of body-body pairs. With a
solver-owned pipeline the rigid slot follows collision_frequency_type as
well.
Forces. Each contact kind has its own law with friction and damping;
self-contact, rigid-soft contact and body-body contact are given in
The energy terms. Full-surface mesh records normally produce force only
when they pass the feature-region test of Offset Geometric Contact
(mesh_contact_valid); recovery records and rows stored after an overflow
bypass the test, and every stored record is kept for truncation.
Truncation. Divide and Truncate keeps self-contact, and with
rigid_soft_enable_dat also rigid-soft contact, free of new penetrations
between detections, through the planar and the radius truncation of
Penetration-free truncation (DAT). The guarantee holds for pairs that were separated at
detection and found within the query radius; the rigid trajectory check is
sampled unless rigid_soft_dat_use_interval_arithmetic is on.
How to set the distances, and what they cost, is on the tuning page: Contact Handling and Setting margin and gap. The contact methods are described in:
Anka H. Chen, Jerry Hsu, Youssef Ayman, and Miles Macklin. 2026. Divide and Truncate: A Penetration and Inversion Free Framework for Coupled Multi-physics System. In ACM SIGGRAPH 2026 Conference Papers. https://doi.org/10.1145/3799902.3811143
Anka He Chen, Jerry Hsu, Ziheng Liu, Miles Macklin, Yin Yang, and Cem Yuksel. 2025. Offset Geometric Contact. ACM Trans. Graph. 44, 4. https://doi.org/10.1145/3731205
Key Constructor Parameters#
The constructor has three groups. Common: iterations, friction_epsilon,
dat_conservative_bound_relaxation,
integrate_with_external_rigid_solver, and collision_pipeline with
collision_frequency_type and collision_frequency. particle_*:
self-contact on or off, margin and gap, the candidate buffer sizes, the
topological and rest-shape filters, and the tile solve. rigid_*: the
compliant ALM switch with its stabilization and retention factors, body-body
contact history, the per-body contact list capacities, the rigid-soft barrier
law and DAT, and the structural joint stiffness and damping. Defaults and exact
semantics are in SolverVBD; recommended values for
particle systems, and effects and symptoms for rigid bodies and cables, are
in the VBD Tuning parameter table.
How VBD Works#
VBD (Vertex Block Descent) treats every time step as one minimization and
solves it one small block at a time. Two kinds of blocks exist. A particle
\(i\) has 3 degrees of freedom, its position \(\mathbf{x}_i\), and a
mass \(m_i\). A rigid body \(b\) has 6 degrees of freedom, its pose
made of a centre of mass \(\mathbf{c}_b\) and a rotation
\(\mathbf{R}_b\), with a mass \(m_b\) and an inertia tensor
\(\mathbf{I}_b\). The notation follows the VBD paper cited above:
\(h\) is the time step and \(\mathbf{y}\) the inertial prediction.
The kernels named on this page live in newton/_src/solvers/vbd/
(solver_vbd.py, particle_vbd_kernels.py, rigid_vbd_kernels.py);
they are internal and can change without notice.
Step objective#
Every time step is one minimization, solved one block at a time. With the blocks and notation above, the step finds the positions and poses that minimize
The first sum is the particle inertia. Its target is the inertial prediction
computed by the forward_step kernel from the current velocity, gravity and
external forces. The second sum is the rigid inertia. Its target pose
\((\mathbf{c}_b^*, \mathbf{R}_b^*)\) comes from forward integration
(forward_step_rigid_bodies), and \(\boldsymbol{\theta}_b\) is the
rotation vector from the current rotation to the target rotation
(solve_rigid_body). The last sum is every energy term in the scene:
membrane, bending and volumetric elasticity, self-contact and rigid-soft
contact with friction and damping, springs, and the compliant
augmented-Lagrangian rows of joints, body-body contacts and particle-body
attachments.
How it is solved#
The solver never forms the full system. In each iteration it visits every block once, with the neighbours frozen at their current values. For a block it assembles the gradient \(\mathbf{g}_i\) and the Hessian \(\mathbf{H}_i\) of the terms that touch the block, and moves the block by one Newton step:
For a particle this is a 3x3 solve (solve_elasticity). For a body it is a
6x6 solve with the block Hessian
\(\operatorname{diag}(m_b/h^2\,\mathbf{1}_3,\ \mathbf{I}_b/h^2)\) plus the
energy terms, factored once per visit (solve_rigid_body, ldlt6_solve).
The energy blocks of \(\mathbf{H}_i\) are the usual stabilised second
derivatives, not always the exact ones: the membrane block is projected to be
positive semi-definite, the bending block keeps only the Gauss-Newton term
\(k\,\nabla\theta\,\nabla\theta^{\mathsf T}\), and the friction block
holds the contact normal fixed.
Blocks are grouped into colours so that no two blocks of one colour share a
mesh element: a triangle, a tetrahedron, a bending edge when
color() is called with include_bending=True, or
a joint for bodies. Contact pairs are not part of the colouring, so two blocks
of one colour can be in contact; they then update in Jacobi fashion. The
update is Gauss-Seidel across colours and Jacobi inside a colour. Rigid blocks
are solved before particle blocks in every iteration.
The real call tree, with the control flow written out:
SolverVBD.step(state_in, state_out, control, contacts, dt)
_initialize_rigid_bodies -> forward_step_rigid_bodies (predicted poses)
_initialize_particles -> forward_step (predicted positions y)
_initialize_particle_body_attachments -> step_particle_body_attachment_lambda_rho
(attachment rho from the endpoints' support; retained
multipliers scaled by rigid_avbd_gamma)
collision detection (default: self-contact before and after this initialization;
rigid contacts before it only, or before and after it with rigid-soft DAT)
for iteration in range(iterations):
if detection is due in this iteration (ITERATIONS schedule only):
detect again, refresh the DAT reference
_solve_rigid_body_iteration
for colour in body colours:
accumulate attachments, body-particle contacts, body-body contacts
solve_rigid_body (6x6 LDL^T per body)
update dual variables of the joint and contact rows
_solve_particle_iteration
for colour in particle colours:
accumulate attachment, body-particle contact, spring and self-contact terms
solve_elasticity (3x3 per particle: inertia + membrane + bending + volumetric)
DAT: planar truncation per pair, then the radius cap
update_duals_particle_body_attachments (dual step of the attachment rows,
after both primal blocks)
_finalize_rigid_bodies -> update_body_velocity
_finalize_particles -> update_velocity v = (x_new - x_old) / dt
Two facts follow directly. Information moves about one ring of neighbours per iteration: a block sees only the current values of its neighbours, so a change reaches a neighbour in the same iteration when that neighbour’s colour comes later in the sweep, and in the next iteration otherwise. And an unfinished solve is not discarded: whatever position the last iteration leaves becomes the new velocity through \(\mathbf{v} = (\mathbf{x}_{\text{new}} - \mathbf{x}_{\text{old}})/h\).
Rigid bodies and joints#
A rigid body \(b\) is one block with 6 degrees of freedom: the position \(\mathbf{c}_b\) of its centre of mass and a rotation vector about it. With the target pose \((\mathbf{c}_b^*, \mathbf{R}_b^*)\) from forward integration and \(\boldsymbol{\theta}_b\) the body-frame rotation vector from the current rotation \(\mathbf{R}_b\) to \(\mathbf{R}_b^*\), the inertia term of the step objective contributes
Every term that touches the body (joints, body-body contacts, rigid-soft
contacts, attachments, tendons) adds its force, its torque and the three
\(3 \times 3\) blocks of the Hessian. The block then solves one
\(6 \times 6\) system by \(LDL^{\mathsf T}\) (ldlt6_solve) and
applies the increments (solve_rigid_body; the rotation increment acts in
world space and is applied as a normalised small-angle quaternion, which
approximates the exponential map):
Joints and body-body contacts are constraint rows between two bodies, handled
as compliant augmented-Lagrangian (ALM) rows. A row has a residual \(C\),
a material stiffness \(k\), an internal metric \(\rho\) that the solver
sets when it initialises the row (once per step for joints, at every rigid
contact detection for body-body contacts), and a multiplier \(\lambda\).
With rigid_avbd_alpha above zero the residual is stabilised to
\(C - \alpha C_0\); the default under compliant ALM is zero, and the
formulas below use that default. Its force along the row and the dual update
after every iteration are
and \(k_{\text{eff}}\) enters the Hessian blocks through the row’s
Jacobian. A bilateral row keeps the full force. A unilateral row keeps only
the admissible sign: for contact \(f_n = \max(f, 0)\) and
\(\lambda_n \leftarrow \max(\,\cdot\,, 0)\). Joint multipliers persist
across steps, scaled by rigid_avbd_gamma (0.999); body-body contacts
restart theirs at zero unless rigid_contact_history=True.
Joints (evaluate_joint_force_hessian). The structural rows of a joint
between parent \(p\) and child \(c\) are the offset of the two joint
anchors and the rotation away from the rest relative orientation, each
projected onto the locked axes of the joint type
(compute_kappa_and_jacobian, build_joint_projectors). D6 joints with
two or three angular axes are the exception: their angular rows use the
absolute joint coordinates of eval_fk in the joint-axis basis
(compute_joint_kappa_and_jacobian):
so BALL keeps the three linear rows; REVOLUTE the three linear rows and two
angular rows; PRISMATIC two linear rows and the three angular rows; FIXED all
six; D6 whatever its free axes leave; FREE adds no row. The material stiffness
is rigid_joint_linear_ke or rigid_joint_angular_ke, and the metric is
\(\rho = \max(\text{support}, 9k)\), chosen so that \(k_{\text{eff}}
\ge 0.9\,k\) in every local solve;
rigid_joint_linear_kd and rigid_joint_angular_kd add absolute damping.
Drives are bilateral rows on the free axes of REVOLUTE, PRISMATIC and D6
joints, with the spring and the damper folded into one row:
with \(k_t = \mathtt{joint\_target\_ke}\), \(d_t =
\mathtt{joint\_target\_kd}\), and \(K_{\text{eff}}, s\) formed from
\(K\) and the axis metric as above. Limits are unilateral rows at the
nearer bound with joint_limit_ke and joint_limit_kd; only forces into
the bound are kept, and the multiplier is projected onto its half-line. ROD joints use four rows,
stretch, shear, bend and twist, with joint_target_ke and
joint_target_kd as their stiffness and damping. DISTANCE joints are not
supported.
Body-body contact is the unilateral row of the same family: \(C\) is the penetration of the two effective surfaces, \(k\) the pair average of the shape stiffness, \(\rho_n\) the pair’s inverse Delassus support over \(h^2\), \(f_n = \max(k_{\text{eff}}\,C + s\,\lambda_n,\ 0)\), and friction is a Coulomb projection carrying \(\boldsymbol{\lambda}_t\). The complete law, with the material mixing and the damping, is in The energy terms.
The deprecated legacy path (rigid_compliant_alm=False) replaces the
finite-material rows with ramped penalties and hard or soft slots.
The energy terms#
Each energy term depends on a few blocks and on model parameters:
Energy term |
Acts on |
Model parameters |
Kernel |
|---|---|---|---|
Membrane |
the 3 particles of a triangle |
|
|
Bending |
the 4 particles of an interior edge |
|
|
Volumetric |
the 4 particles of a tetrahedron |
|
|
Self-contact, with friction and damping |
a vertex and a triangle, or two edges, of one particle system |
|
|
Rigid-soft contact, with friction and damping |
a particle, or a point on a soft edge or face, and one shape |
|
|
Body-body contact, with friction and damping |
two shapes, on two rigid bodies or on a rigid body and a static shape |
shape |
|
Particle-body attachments |
a particle and a point fixed on a body |
attachment stiffness (finite: compliant; infinite: hard limit, damping ignored) and damping |
compliant augmented-Lagrangian row:
|
Springs |
two particles |
spring stiffness and damping |
|
Joints, drives and limits |
one or two bodies |
joint stiffness and damping |
compliant augmented-Lagrangian rows, |
Elasticity is three terms. The membrane term resists stretch and shear of a triangle. With \(\mathbf{F}\) the in-plane deformation gradient, \(I_C = \lVert \mathbf{F} \rVert^2\) and \(J_s\) the area ratio, the energy per unit rest area is the stable Neo-Hookean membrane energy
with \(\mu = \mathtt{tri\_ke}\) and \(\lambda_{\text{NH}} =
\mathtt{tri\_ka} + \mathtt{tri\_ke}\); the triangle energy is \(A_0\,\psi_{\text{m}}\).
tri_kd damps the rate of the metric \(\mathbf{F}^{\mathsf T}\mathbf{F}\),
so rigid rotations are not damped. The bending term resists a change of the
dihedral angle \(\theta\) across an interior edge of rest length
\(L_0\):
and edge_kd damps the angle rate with the force
\(\mathbf{f}_d = -\mathtt{edge\_kd}\,L_0\,\dot\theta\,\partial\theta/\partial\mathbf{x}\),
where \(\dot\theta\) is the change of \(\theta\) over the step divided
by \(h\). The volumetric term is the stable
Neo-Hookean energy of a tetrahedron of rest volume \(V_0\),
with \(\mu = \mathtt{k\_mu}\), \(\lambda_{\text{NH}} = \mathtt{k\_lambda} +
\mathtt{k\_mu}\), the same \(\alpha\), energy \(V_0\,\psi_{\text{v}}\),
and k_damp damping the rate of the metric \(\mathbf{F}^{\mathsf T}\mathbf{F}\)
like tri_kd. In both Neo-Hookean terms the user values are
Lamé parameters; the kernels convert them to the Neo-Hookean parameters, which
is why \(\lambda_{\text{NH}}\) carries the extra \(\mu\).
Contact is three different terms, one per contact kind. They share the friction smoothing length \(\epsilon_u = \mathtt{friction\_epsilon}\, h\) and nothing else: each kind has its own distance, its own stiffness and its own law.
Self-contact acts between a vertex and a triangle, or between two edges, of
one particle system (accumulate_self_contact_force_and_hessian). The
distance \(d\) is the closest-point distance of the pair, the interaction
distance is \(r = \mathtt{particle\_self\_contact\_margin}\) and the
stiffness is \(k = \mathtt{soft\_contact\_ke}\). The energy
(evaluate_self_contact_force_norm) is a quadratic penalty that turns into
a logarithmic barrier at half the interaction distance:
The gradient in the barrier region is \(-k\,(r/2)^2/d\), so the force
grows like \(1/d\) as the gap closes; below \(d_{\min}\) it continues
linearly, so the force stays finite. A pair whose closest points coincide
(\(d = 0\)) is skipped. The pair force is distributed to the
four vertices by the barycentric weights of the two closest points. Friction
follows the IPC form (compute_friction). With \(f_n = -\partial
E_{\text{sc}}/\partial d\) the normal force, \(\mathbf{u}\) the tangential
relative displacement of the two closest points over the step and
\(\mathbf{T}\) the tangent basis,
with \(\mu = \mathtt{soft\_contact\_mu}\). Unlike IPC, the contact normal
is held constant inside the block solve; the kernel notes that this improves
stability. Damping (damp_collision) acts on the gap rate \(\dot d =
\mathbf{n}\cdot\Delta\mathbf{x}_{\text{rel}}/h\), only while the gap closes:
\(\mathbf{f}_d = -\mathtt{soft\_contact\_kd}\;\dot d\,\mathbf{n}\). A
separating pair is not damped.
Rigid-soft contact acts between a soft contact point and one shape of a
rigid or static body. The point is a particle, or, with full-surface mesh
contact (enable_rigid_soft_full_surface_contact), a barycentric point on a
soft edge or face; the records come from the collision pipeline. The
distance is \(d = \mathbf{n}\cdot(\mathbf{x} - \mathbf{b})\), with
\(\mathbf{b}\) the shape’s contact point carried along with its body and
\(\mathbf{n}\) the stored contact normal. The interaction distance is
\(r = \mathtt{particle\_radius}\) plus the shape margin; an edge or
face point uses the largest radius among its corners. The material is mixed
per contact from the two sides (init_body_particle_contacts):
\(k = \tfrac{1}{2}(\mathtt{soft\_contact\_ke} + \mathtt{ke}_{\text{shape}})\),
\(k_d = \tfrac{1}{2}(\mathtt{soft\_contact\_kd} + \mathtt{kd}_{\text{shape}})\)
and \(\mu = \sqrt{\mathtt{soft\_contact\_mu}\;\mathtt{mu}_{\text{shape}}}\).
The normal law is the plain
quadratic penalty \(E = \tfrac{k}{2}(r - d)^2\) for \(d < r\);
rigid_soft_contact_use_log_barrier=True switches it to the self-contact
barrier law. Damping and friction have the same form as for self-contact,
written in world space (compute_projected_isotropic_friction): the damping force
\(-k_d\,\dot d\,\mathbf{n}\) acts only while the point approaches the
shape, and the friction force is \(-\mu\, f_n\,
f_1(\lVert\mathbf{u}_t\rVert)\,\hat{\mathbf{u}}_t\), with \(\mathbf{u}_t\)
the tangential slip of the point relative to the shape surface over the step.
The force is distributed to the particles by their barycentric weights, and
the equal and opposite force with its torque acts on the body
(accumulate_body_particle_contacts_per_body). With full-surface mesh
contact, every stored edge and face record is kept for DAT, but normally only
the records that pass the feature-region test of Offset Geometric Contact
produce force (mesh_contact_valid); penetration-recovery records and rows
stored after a buffer overflow bypass the test.
Body-body contact acts between two shapes, on two rigid bodies or on a
rigid body and a static shape. It is a constraint row of the rigid solve, not
a penalty energy. The gap is measured between the two effective surfaces:
with \(\mathbf{p}_0, \mathbf{p}_1\) the support points carried along with
their bodies and \(t_0, t_1\) the effective thicknesses (shape radius plus
margin, stored by the collision pipeline as rigid_contact_margin0 and
rigid_contact_margin1), the penetration is
The material is the pair average of the shape values:
\(k = \tfrac{1}{2}(\mathtt{ke}_0 + \mathtt{ke}_1)\),
\(k_d = \tfrac{1}{2}(\mathtt{kd}_0 + \mathtt{kd}_1)\) and
\(\mu = \sqrt{\mathtt{mu}_0\,\mathtt{mu}_1}\). On the default compliant
ALM path (rigid_compliant_alm=True) the normal force is
where \(\lambda_n\) is the multiplier of the row and \(\rho_n\) is an
internal metric that the solver sets every step from the pair’s inverse
Delassus support over \(h^2\) plus structural stiffness; the user does not set \(\rho_n\). The
force acts along \(\mathbf{n}\) at the support points with the matching
torques, and the row adds \(k_{\text{eff}}\,\mathbf{n}\mathbf{n}^{\mathsf T}\)
to the body Hessians (evaluate_rigid_contact_from_collision). After every
iteration the multiplier is advanced by the relaxed ascent
\(\lambda_n \leftarrow \max\big(s\,(\lambda_n + \rho_n C),\ 0\big)\)
(update_duals_body_body_contacts). At
the fixed point \(\lambda_n = k\,C\) and \(f_n = k\,C\): a converged
resting contact carries its load like a penalty of the material stiffness
\(k\), with the penetration \(f_n/k\), while each single iteration
uses the smaller \(k_{\text{eff}}\), which keeps the 6x6 systems well
conditioned. Friction is a Coulomb bound, not a smooth law: the tangential
force \(\rho_t\,\mathbf{u}_t + \boldsymbol{\lambda}_t\), with
\(\mathbf{u}_t\) the tangential offset between the two anchored surface
points and \(\rho_t\) the tangent metric, is projected onto the disk
\(\lVert\mathbf{f}_t\rVert \le \mu f_n\),
and \(\boldsymbol{\lambda}_t\) is advanced and projected onto the disk of
radius \(\mu\lambda_n\) in the same dual step. Damping is
\(-k_d\,(\mathbf{n}\cdot\mathbf{v}_{\text{rel}})\,\mathbf{n}\) while the
bodies approach and the row carries load. The multipliers start from zero in
every step unless rigid_contact_history=True warm-starts them from the
matched contacts of the previous step. The legacy AVBD path
(rigid_compliant_alm=False, deprecated) uses the pair stiffness as a fixed
penalty by default (rigid_avbd_beta is zero). With a positive beta
(rigid_avbd_linear_beta, falling back to rigid_avbd_beta) the penalty
starts at rigid_contact_k_start and grows by beta times the penetration
per iteration, capped at the pair stiffness. Its soft
mode uses the IPC friction law instead of the Coulomb projection.
The underlying methods are described in the references cited by
SolverVBD: VBD for the particle blocks, AVBD for the
rigid rows, Divide and Truncate (DAT) for penetration-free contact and Offset
Geometric Contact (OGC) for the feature regions of mesh contact. The two
contact papers are also listed under Contact Handling above.
Particle-body attachments are compliant augmented-Lagrangian rows of the
same form as the joint rows (Rigid bodies and joints). With \(C =
\mathbf{x} - \mathbf{a}\) the offset of the particle from its anchor point fixed
on the body, the row applies \(k_{\text{eff}}\,C + s\,\lambda\) to the
particle and the opposite force, with its torque, to the body;
\(k_{\text{eff}}\) and \(s\) follow from the authored stiffness
\(k\) and the metric \(\rho\) as above. The solver sets \(\rho\)
once per step to the row’s support, the endpoints’ reduced inverse mass over
\(h^2\) plus the particle’s incident material Hessian, and at the same time
scales the retained multiplier by rigid_avbd_gamma
(step_particle_body_attachment_lambda_rho). The multiplier advances once
per iteration after both primal blocks, \(\lambda \leftarrow s\,(\lambda +
\rho\,C)\) (update_duals_particle_body_attachments). Infinite stiffness is
the hard limit of the same row, \(k_{\text{eff}} = \rho\) and \(s =
1\), and ignores the attachment damping; a finite stiffness keeps its spring
offset under load and adds the damping force \(k_d\,\dot C\). Springs
between two particles stay plain penalty energies
(accumulate_spring_force_and_hessian).
Penetration-free truncation (DAT)#
Self-contact and rigid-soft contact are penalty energies, so a large step could jump through a surface between two detections. The Divide and Truncate scheme (DAT) prevents that for every pair that detection found separated and within the query radius. It works on the displacement each vertex has accumulated since the last detection and has two distinct parts.
Planar truncation acts per detected vertex-triangle or edge-edge pair
(apply_planar_truncation_parallel_by_collision). At the detection-time
positions a certified separating plane is found and placed between the two
primitives; it splits the gap in proportion to how far each side wants to move
toward the other, and each side keeps at least 5 percent of the gap
(place_dat_division_plane). Each vertex of the pair receives a factor
\(t \le 1\) so that its straight trajectory stays on its own side of the
plane, scaled by dat_conservative_bound_relaxation; motion away from the
plane or along it passes untouched (planar_truncation_t). A vertex in
several pairs keeps the smallest factor. A pair whose primitives already
intersect at detection time has no strict separator and is left unconstrained
in that pass.
Radius truncation acts on every particle (apply_truncation_ts). After
the planar factor is applied, the length of the accumulated displacement is
capped at one half of dat_conservative_bound_relaxation times the query
radius, isotropically. Two primitives farther apart than the query radius at
detection therefore cannot meet before the next detection.
Rigid-soft DAT (rigid_soft_enable_dat, off by default) applies the
same two parts to contacts between soft geometry and shapes. One thread per
contact record builds the division plane from the detection-time
configuration and truncates both sides against it
(apply_rigid_soft_truncation). The soft vertices of the
record move on straight rays, as in self-contact. The vertices of the rigid
primitive move on curved trajectories because the body rotates, so their
crossing times are found by sampling and bisection, with an interval bound on
the signed-distance derivative for rows that start inside the plane’s band;
rigid_soft_dat_use_interval_arithmetic additionally certifies the sampled
prefix with interval arithmetic (experimental). Without it, a crossing that
happens and reverses between two sample points can be missed. Full-surface
mesh contacts name the exact rigid triangle, edge or vertex through
soft_contact_rigid_indices; analytic SDF contacts and
penetration-recovery contacts truncate only their stored rigid surface point.
Each body then scales its whole pose update, translation and rotation about
the centre of mass, by the smallest factor of its contacts
(apply_body_truncation_ts). The radius truncation caps the motion of every
point of the body, bounded by the translation plus the rotation angle times
the body’s bounding radius, at one half of the relaxation factor times the
smallest rigid-soft query radius (pipeline gap plus soft feature radius plus
shape margin); that budget is fixed at construction. The pass runs after the
forward step and after every body colour of every iteration, and the particle
pass re-enforces the same planes after the particles move. Kinematic bodies
move outside the solver and are not truncated. A pair for which no separating
plane can be certified freezes. Contact Handling describes the
parts for tuning and illustrates the planes.
The guarantee has a price: a vertex cannot move faster than about one motion
budget per detection cycle. The motion that the truncation removes is called
DAT damping here and in VBD Tuning. It is not related to any kd
damping coefficient.
Examples#
Empirical parameter combinations on the tuning page lists one shipped example per scenario, from cloth drape to a gripper on a soft pad, and says what each one demonstrates; copy the values from the example file. For cables, see the cable example table.
Runnable examples using SolverVBD on GitHub:
VBD-specific examples (rigid and soft contacts, grippers, stiff materials)
Multiphysics examples (coupling with MuJoCo, XPBD, and MPM)