The step pipeline

A step has two halves. The first runs once per articulation and is classical Featherstone dynamics, with a Cholesky factor of the mass matrix. The second runs once per world and turns contacts, limits, drives, mimic couplings and loop closures into rows of a velocity-level projected Gauss-Seidel (PGS) solve. Positions then follow from the solved velocities by symplectic Euler.

Figure 1One step. The front half is grouped by articulation size so robots of one size share kernels. The green box is the constraint solve, where most of the options in the reference act.

With joint-space inertia $H(q)$, bias forces $C(q,\dot q)$, applied generalized force $\tau$ and step $h$, the unconstrained velocity is

$$\tilde H = H + \operatorname{diag}(\text{armature}) + \operatorname{diag}(h\,k_d + h^2 k_e),\qquad \hat v = \dot q + h\,\tilde H^{-1}\big(\tau + u_0 - C(q,\dot q)\big),$$

where the extra diagonal folds the default implicit drives into the mass matrix and $u_0 = -k_e(q - q^* + h\dot q) - k_d(\dot q - \dot q^*)$ is the drive force at the start of the step, clamped to the effort limit. Constraint impulses $\lambda$ act through row Jacobians $J$, so $v^+ = \hat v + \tilde H^{-1}J^\top\lambda$. With $u = Jv^+$ the rows’ physical velocity and $b$ their bias (position correction, speculative gaps and bounce targets), normal rows satisfy the complementarity condition $0 \le \lambda_n \perp u_n + b_n \ge 0$ and friction rows stay in the Coulomb disk $\|\lambda_t\| \le \mu\lambda_n$.

The bias forces $C(q,\dot q)$ and the passive joint springs and damping in $\tau$ are taken at the start of the step, so they are explicit, as in SolverFeatherstone. That sets two limits on $h$:

• Passive joint damping $c$ is stable only while $h\,c\,\lambda_{\max}(\tilde H^{-1}) < 2$. Light links make $\lambda_{\max}$ large: for a six-link chain of 2 cm capsules it is about $4\times10^4$, so at 240 Hz damping above about 0.012 N·m·s/rad diverges. Drive gains do not have this limit, because they are folded into $\tilde H$. A change that folds passive damping into $\tilde H$ the same way is in review in skild-ai/newton#18.
• Velocity-product forces of fast, light chains add energy at a rate that is first order in $h$; Featherstone and MuJoCo’s euler integrator behave the same. Iterations do not help, since no rows are involved. A smaller $h$ reduces the gain; armature large against the links’ inertia also does, but it changes the model. skild-ai/newton#19 documents this in the solver. Tutorial: Armature · Joint springs and damping