Py-DiSMech: Python Soft Robotics Simulator

• Py-DiSMech is a Python-based simulation framework that leverages discrete differential geometry to model soft robotic rods, shells, and hybrid structures.
• It integrates fully vectorized NumPy operations, an implicit penalty-energy contact model, and natural strain-based feedback control for reliable simulation and design.
• Building on prior DisMech frameworks, Py-DiSMech offers enhanced modularity, computational efficiency, and a unified approach to sim-to-real research in soft robotics.

Searching arXiv for Py-DiSMech and related DisMech papers. Py-DiSMech is a Python-based, open-source simulation framework for modeling and control of soft robotic structures grounded in the principles of Discrete Differential Geometry (DDG). By discretizing geometric quantities such as curvature and strain directly on meshes, it represents rods through the Discrete Elastic Rods (DER) formulation, shells through hinge-based or mid-edge-normal-based bending models, and hybrid rod–shell structures through shared nodes and coupled stencils. The framework combines a fully vectorized NumPy implementation, a fully implicit penalty-energy contact model, and a natural-strain-based feedback-control module intended for simulation-driven design, control validation, and sim-to-real research in soft robotics (Lahoti et al., 10 Dec 2025).

1. Position within the DiSMech family

Py-DiSMech belongs to a sequence of DDG-based soft-robotics simulators that progressively expanded modeling scope, implementation environment, and control support. DisMech introduced a simulator for highly dynamic motions of rod-like soft continuum robots and structures, using a fully implicit discrete differential geometry-based solver, fast and accurate contact handling, and natural-curvature actuation for elastic rods with arbitrary interconnections (Choi et al., 2023). MAT-DiSMech extended the DDG approach to rods, shells, and rod–shell systems in an open-source MATLAB-based framework with implicit integration, modular external forces, and examples such as a PneuNet-like actuator, manta ray, and parachute (Lahoti et al., 24 Apr 2025). Py-DiSMech retains the DDG foundation while shifting to a Python-native, vectorized, object-oriented implementation and broadening the contact model, software extensibility, and closed-loop control interface (Lahoti et al., 10 Dec 2025).

The progression is most visible in the representational scope. DisMech is centered on rod-like soft continuum robots and structures; MAT-DiSMech formalizes rods, shells, and their combinations; Py-DiSMech presents a unified framework for rods, shells, and hybrid structures with direct interoperability with Python scientific and machine-learning ecosystems (Choi et al., 2023, Lahoti et al., 24 Apr 2025, Lahoti et al., 10 Dec 2025).

Framework	Modeled structures	Salient traits
DisMech	Rod-like soft continuum robots and structures	Fully implicit DDG solver; arbitrary interconnections; natural curvature actuation
MAT-DiSMech	Rods, shells, and rod–shell hybrids	MATLAB-based; implicit integration; modular external forces and contact
Py-DiSMech	Rods, shells, and hybrid rod–shell structures	Fully vectorized NumPy; fully implicit penalty-energy contact; natural-strain-based feedback control

This lineage also clarifies Py-DiSMech’s technical orientation. It is not a volumetric FEM or FEA package; it is a reduced-order, energy-based mechanics framework whose state variables, residuals, and Jacobians are defined directly on discrete geometric primitives. That positioning is explicit throughout the DiSMech literature and is central to the framework’s balance between physical accuracy and computational tractability (Choi et al., 2023, Lahoti et al., 24 Apr 2025, Lahoti et al., 10 Dec 2025).

2. Discrete geometric mechanics

Py-DiSMech discretizes a structure using $N$ nodes in 3D, $E$ edges, and $T$ triangles. The user supplies nodal positions as an $N \times 3$ array, edge connectivity as an $E \times 2$ array, and triangle connectivity as a $T \times 3$ array. For rods, the generalized coordinates are

$$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$

so the state dimension is $3N + E$. For shells using the mid-edge-normal bending model, additional shell-edge rotational DOFs $\xi$ are introduced: $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E, \xi^1, \dots, \xi^Z]^\top,$$ where $E$0 is the number of shell edges (Lahoti et al., 10 Dec 2025).

The governing equations of motion are written as

$E$1

with $E$2 a diagonal lumped mass matrix, $E$3 the internal elastic force, $E$4 the external force, and $E$5 the self-contact and friction force. Py-DiSMech computes deformation energy through local discrete stencils, or “springs,” associated with stretching, rod bending, rod twisting, hinge-based shell bending, and mid-edge shell bending (Lahoti et al., 10 Dec 2025).

For rods, the DDG structure follows DER. Each rod edge stores a reference frame $E$6, a material frame $E$7, and a tangent

$E$8

A key DER feature is that the reference frame is updated by parallel transport. For consecutive edges $E$9 and $T$0, the curvature binormal is

$T$1

Projected onto the material directors, this yields the scalar curvatures $T$2 and $T$3, and the bending strain is defined relative to the natural curvatures $T$4. Twist strain is represented as

$T$5

which gives a discrete torsional response compatible with the same energy-based assembly (Lahoti et al., 10 Dec 2025).

For shells, Py-DiSMech includes two bending formulations. The hinge-based model uses the dihedral angle between adjacent triangle normals, with strain

$T$6

The mid-edge-normal-based model introduces a mid-edge normal $T$7, shell-edge DOFs $T$8, and a discrete shape operator $T$9, with the bending strain expressed through

$N \times 3$0

The mid-edge formulation is presented as more geometrically sophisticated and more robust to mesh quality, though more expensive because of the additional edge DOFs (Lahoti et al., 10 Dec 2025).

Hybrid rod–shell coupling is realized through joint nodes shared by rod edges and shell triangles. In these regions, shell edges near the joint are treated as rod edges, material frames and twist DOFs are assigned, and rod-like bending–twisting stencils are constructed across valid local edge combinations. For the mid-edge shell model, a joint edge may carry both $N \times 3$1 and $N \times 3$2 DOFs. This extends the rod-only and MATLAB-based DiSMech formulations into a single hybrid representation (Lahoti et al., 24 Apr 2025, Lahoti et al., 10 Dec 2025).

3. Implicit dynamics, Newton solves, and contact

Py-DiSMech is primarily integrated with implicit Euler, while also providing Implicit Midpoint and Newmark-$N \times 3$3 steppers. Under implicit Euler, with velocity $N \times 3$4, the residual is

$N \times 3$5

with

$N \times 3$6

At each step, $N \times 3$7 is solved by Newton–Raphson,

$N \times 3$8

where $N \times 3$9 is usually $E \times 2$0, but may be reduced by line search. The Jacobian is

$E \times 2$1

This fully implicit formulation is presented as the key to stable integration for stiff systems, large elastic moduli, and contact-rich motion (Lahoti et al., 10 Dec 2025).

The contact model generalizes the earlier rod-centric Implicit Contact Model used in DisMech and MAT-DiSMech. DisMech already used an implicit penalty-energy method for rod–rod self-contact, rod-on-floor contact, and frictional interactions in rod-based robots (Choi et al., 2023). MAT-DiSMech retained IMC for rod self-contact and indicated that shell contact could be approximated by edge–edge interactions on triangle edges (Lahoti et al., 24 Apr 2025). Py-DiSMech explicitly extends fully implicit penalty-energy contact to rod–rod, rod–shell, and shell–shell interactions (Lahoti et al., 10 Dec 2025).

The simulator defines contact pairs analogously to elastic springs and classifies closest-feature interactions as Point-to-Point, Point-to-Edge, Edge-to-Edge, or Point-to-Triangle. Distances are computed analytically in batch; for edge–edge interactions the paper uses Lumelsky’s algorithm, and for triangle–triangle proximity it uses a batched barycentric method inspired by NVIDIA PhysX. At the beginning of each time step, during the first Newton iteration, a candidate contact set is built to restrict collision checks to likely interacting pairs and reduce cost (Lahoti et al., 10 Dec 2025).

Given separation $E \times 2$2, the contact energy is

$E \times 2$3

The contact force on node $E \times 2$4 is

$E \times 2$5

Friction is modeled as a regularized Coulomb-like force,

$E \times 2$6

with $E \times 2$7, $E \times 2$8 a slipping tolerance, and $E \times 2$9 a smooth activation function of tangential relative velocity. The paper states that this regularization avoids a discontinuous static/kinetic switching law (Lahoti et al., 10 Dec 2025).

Several limitations are explicit. Because the method is a penalty method, it does not guarantee no penetration. Edge-to-triangle contact is not expected in well-tuned scenarios but can occur under interpenetration and is handled if it arises. The candidate contact set is only updated at the start of a time step or first Newton iteration, which improves speed but may miss some extreme within-step changes if motion is too aggressive. Robustness therefore depends on proper tuning of penalty stiffness and tolerance parameters (Lahoti et al., 10 Dec 2025).

4. Actuation, feedback control, and software design

Py-DiSMech includes a control module based on natural strain fields. Instead of directly commanding nodal positions or generalized forces, actuation is expressed as changes in the structure’s intrinsic or natural strains. The paper discusses 2D rod examples in which the control inputs are natural longitudinal strain $T \times 3$0 and natural curvature $T \times 3$1, while the controller acts on the residual between desired and measured strain. The built-in feedback law is a proportional–integral controller for shape regulation and trajectory tracking, with stabilizing features that include spatial smoothing of the actuation field, rate limits on curvature change, and anti-windup in the PI integrator (Lahoti et al., 10 Dec 2025).

Two closed-loop tasks are presented. In shape regulation, a rod simply supported at both ends is driven to an S-shaped target under gravity and damping; the paper reports that strain residuals converge to zero and that the shape RMSE between actual and target nodal positions converges to zero. In trajectory tracking, a rod performs serpentine locomotion in viscous fluid using curvature-based actuation; after about $T \times 3$2 s the tracked and reference configurations overlap closely, and the strain residual converges to within $T \times 3$3 of zero (Lahoti et al., 10 Dec 2025).

The control interface generalizes an earlier actuation philosophy in DisMech. DisMech actuated rods by changing natural curvature and introduced a gradient descent approach that mapped the motions of hardware robot prototypes to control inputs in simulation. In a shape-memory-alloy-actuated dual-limb soft manipulator, that real2sim procedure reported an average tip position error of $T \times 3$4 mm, corresponding to $T \times 3$5 of robot length (Choi et al., 2023). Py-DiSMech shifts the emphasis from offline curvature fitting to a built-in feedback-control module, but the common abstraction remains natural strain rather than low-level force prescription (Choi et al., 2023, Lahoti et al., 10 Dec 2025).

The software architecture is explicitly modular and object-oriented. The framework centers on a SoftRobot object configured with Geometry, Mesh, Material, SimParams, and Environment. Runtime variables are stored in a RobotState frozen dataclass containing the DOF vector $T \times 3$6, velocity $T \times 3$7, acceleration $T \times 3$8, reference frames $T \times 3$9, and material frames $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$0. Local mechanics are represented by stretch_springs, bend_springs, twist_springs, hinge_springs, and triangle_springs, each carrying nat_strain and inc_strain for strain-based actuation. The code also defines ElasticEnergy, ContactEnergy, FrictionForce, ExternalForce, TimeStepper, and Control abstractions, as well as a modular linearSolve layer supporting a dense NumPy pseudoinverse and sparse solve via PyPardiso (Lahoti et al., 10 Dec 2025).

The implementation strategy is heavily vectorized. Rather than evaluating spring and contact stencils in Python loops, Py-DiSMech computes force and Jacobian contributions in batched array operations. The authors present this as a major engineering contribution because it reduces Python overhead, improves cache-friendly computation, and maps naturally to large stencil sets. The source code is reported as open source at https://github.com/StructuresComp/dismech-python, and the examples are implemented as Jupyter notebooks (Lahoti et al., 10 Dec 2025).

5. Benchmarking, validation, and example systems

The main quantitative comparison in Py-DiSMech is against PyElastica on cantilever, helix, and snake-locomotion benchmarks. All wall-clock times are averaged over five runs on an AMD Ryzen 9 9950X CPU with 64 GB DRAM. The runtime table includes, for PyElastica, a helix case with $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$1 Pa, Verlet, $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$2, taking $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$3 s; for Py-DiSMech, the corresponding helix case with Implicit Euler and $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$4 takes $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$5 s, while Newmark-$$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$6 takes $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$7 s. Across the reported cases, the paper summarizes the overall gain as $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$8 to $$\mathbf{q} = [\mathbf{x}_1, \dots, \mathbf{x}_N, \theta^1, \dots, \theta^E]^\top,$$9 faster than PyElastica while maintaining comparable physical accuracy (Lahoti et al., 10 Dec 2025).

The cantilever benchmarks validate static and dynamic accuracy against Euler–Bernoulli beam theory. A damped cantilever rod with 101 nodes, $3N + E$0, radius $3N + E$1 m, and length $3N + E$2 m is tested for several Young’s moduli. For the cantilever with $3N + E$3 Pa, PyElastica is reported at $3N + E$4 s with Verlet and $3N + E$5, whereas Py-DiSMech is reported at $3N + E$6 s with Implicit Euler and $3N + E$7. The helix under gravity exercises large geometric nonlinearity, and the snake benchmark exercises dynamic actuation in viscous fluid with hydrodynamic drag (Lahoti et al., 10 Dec 2025).

The case studies broaden the modeled class beyond rod-only systems. Py-DiSMech presents a jellyfish example in which a hemispherical shell forms the bell and slender rods act as tentacles, with periodic upward propulsion and negligible lateral motion. It also presents rectangular shell folding with self-contact, circular cloth or shell folding with self-contact, and the closed-loop control demonstrations already noted. These examples emphasize shell mechanics, hybrid structures, and the extended contact model rather than only rod dynamics (Lahoti et al., 10 Dec 2025).

The earlier DiSMech papers provide the immediate precursors for these use cases. DisMech demonstrated a cantilever beam, a helical rod under gravity, a frictional rod-on-floor test, a four-legged spider-like soft robot built from interconnected rods, an active entanglement gripper, and a shape-memory-alloy-actuated dual-limb soft manipulator (Choi et al., 2023). MAT-DiSMech added a PneuNet-like actuator represented as a rod with natural-curvature actuation, a manta ray modeled as a shell, a parachute modeled as a rod–shell hybrid, a rod dropped on the ground, and cantilever validation using rod, hinge-shell, and mid-edge-shell models (Lahoti et al., 24 Apr 2025). Py-DiSMech can therefore be read as the Python-native consolidation of a benchmark and example suite that had already been distributed across rod-only and MATLAB-based predecessors (Choi et al., 2023, Lahoti et al., 24 Apr 2025, Lahoti et al., 10 Dec 2025).

6. Scope, limitations, and name ambiguity

Py-DiSMech is designed for soft robots dominated by large geometric deformation, stiff elastic response, and complex contact between slender or thin structures. Its strengths lie in reduced-order but high-fidelity DDG mechanics, stable large-step implicit integration, unified implicit contact for rod–rod, rod–shell, and shell–shell interactions, and a control interface based on natural strains rather than ad hoc actuation forces (Lahoti et al., 10 Dec 2025).

The limitations stated in the literature are equally important. In the original DisMech formulation, the rod model is Kirchhoff rather than full Cosserat: it models stretching, bending, and twisting but does not model shear, which the authors identify as a major scope condition for sufficiently slender rods and a disqualifying issue where shear is important (Choi et al., 2023). In Py-DiSMech, the authors explicitly note that penalty contact allows penetration, that no differentiable simulation or gradient-based control or optimization pipeline is demonstrated, that no direct sim-to-real experimental validation is presented, that no GPU backend is presented, and that shell-heavy and hybrid-heavy benchmark scaling is not tabulated against baselines. They also note that the more advanced mid-edge shell formulation adds DOFs and computational cost (Lahoti et al., 10 Dec 2025).

A recurrent source of confusion is the label “Py.” In Py-DiSMech, the prefix denotes the Python-based simulator. In unrelated magnetic literature, however, “Py” denotes Permalloy. The arXiv literature includes a Co/Py bilayer microstructure for element-specific visualization of dynamic magnetic coupling and a Py/Gd/Py/Gd/Py/SiN$3N + E$8 artificial ferrimagnet studied through anisotropic magnetoresistance; in both cases, “Py” refers to Permalloy rather than a Python implementation (Feggeler et al., 2019, Zhang et al., 2024). For soft-robotics and DDG research, Py-DiSMech therefore designates a software framework, not a Py-based magnetic material system.

Taken together, the DiSMech papers define Py-DiSMech as a DDG-based, fully implicit, contact-capable, control-ready Python platform for rods, shells, and hybrid soft structures. Its central contribution is not a single constitutive law, but a unified software realization in which geometric discretization, implicit mechanics, contact handling, and strain-based control are organized into a scalable research workflow for soft robotics (Lahoti et al., 10 Dec 2025).