VDSim Plant — tutorial (closed-loop control plant)

A step-by-step guide to using VDSim as the vehicle-dynamics plant for a closed-loop controller (e.g. an MPC). You call step(); VDSim integrates one control period and returns ground-truth state + per-wheel tyre forces. Mid-fidelity (7-DOF + Pacejka MF tyres — MF96 by default, MF2002 .tir with load-dependent coefficients when you supply one — + load transfer), deterministic, sub-real-time — a drop-in upgrade from a single-track bicycle plant.

Status: shipped on main and published on PyPI since v0.7.0. The API below is the maintained public plant contract; changes are recorded in the project changelog.

---

0. What you get

• Control input each step: u = [delta_roadwheel_rad, Fx_total_N] — the road-wheel steer angle and a longitudinal force intent (no pedal/throttle calibration; +Fx drives, −Fx brakes).
• Observation each step: true vehicle state + per-wheel tyre forces in the contact frame, with the actual friction each wheel used. Everything ISO 8855, wheels FL=0,FR=1,RL=2,RR=3.
• Real tyre saturation: a too-aggressive command loses grip (the tyre slips), it is not silently clamped — this is the point of using VDSim over a toy plant.
• Spatially varying friction (low-μ patches) for "unseen low-μ" scenarios.

---

1. Build (one time)

cd /path/to/VDSim
cmake -B build -DCMAKE_BUILD_TYPE=Release -DVDSIM_BUILD_PYTHON=ON   # first time
cmake --build build -j                                              # core + python module

VDSIM_BUILD_PYTHON defaults OFF — you must pass -DVDSIM_BUILD_PYTHON=ON to get the vdsim pybind module under build/python/.

⚠️ Critical: Building for your Python environment (conda / venv)

The pybind module is ABI-locked to one Python interpreter. A shipped vdsim.cpython-38-*.so will not import under a different interpreter (e.g. a conda env on 3.11 — the typical acados / cvxpy MPC setup).

If your controller runs in its own environment (recommended), you MUST rebuild the module with that interpreter:

conda activate vla                       # your MPC env (e.g. python 3.11)
pip install pybind11                      # into that env
cmake -B build_vla -DCMAKE_BUILD_TYPE=Release -DVDSIM_BUILD_PYTHON=ON \
      -DPython3_EXECUTABLE="$(which python)" \
      -Dpybind11_DIR="$(python -m pybind11 --cmakedir)"
cmake --build build_vla -j

Then point sys.path at build_vla/python (not build/python).

Symptom of a mismatch: ImportError: ... undefined symbol or module compiled against a different Python. See the command above — do not skip the -DPython3_EXECUTABLE / -Dpybind11_DIR flags.

---

2. Import the plant

VDSimPlant lives in python/vdsim_plant.py; the compiled vdsim module in build/python/. Point Python at both:

import sys
sys.path += ["/path/to/VDSim/build/python", "/path/to/VDSim/python"]
from vdsim_plant import VDSimPlant

---

3. Quickstart (copy-paste, runs in seconds)

import sys
sys.path += ["/path/to/VDSim/build/python", "/path/to/VDSim/python"]
from vdsim_plant import VDSimPlant

plant = VDSimPlant(
    config="generic_sedan",  # bundled preset (configs/vehicles/generic_sedan.yaml)
    base_mu=0.9,             # uniform dry friction
    control_dt=0.05,         # one step() = 0.05 s of sim, held ZOH
    substep_dt=5e-4,         # internal fixed integration step (stiff wheel spin)
)

obs = plant.reset(state0=[0.0, 0.0, 0.0, 16.7, 0.0, 0.0])   # [X,Y,psi,vx,vy,r]; vx=16.7 m/s

for k in range(100):                     # 5 s
    delta = 0.03                         # rad, road-wheel
    Fx    = -1500.0                      # N, light brake (force intent)
    obs   = plant.step([delta, Fx])
    print(f"t={k*plant.control_dt:4.2f}  yaw_rate={obs['r']:+.3f}  ay={obs['ay']:+.2f}"
          f"  FL_kappa={obs['wheel'][0]['kappa']:+.3f}")

A ready-to-run version is examples/vla_plant_demo.py:

cd /path/to/VDSim && python3 examples/vla_plant_demo.py     # writes /tmp/vla_plant_demo.csv

---

4. The observation dict

reset() and step() both return one obs dict (true state — no sensor noise):

key	meaning	units / frame
X, Y	position (world)	m, ISO 8855 ENU
psi	yaw	rad
vx, vy	body-frame velocity	m/s (+vx fwd, +vy left)
r	yaw rate	rad/s (+ = left turn)
ax, ay	body specific force (accelerometer)	m/s²
beta	side-slip atan2(vy,vx)	rad
roll, pitch	quasi-static roll / pitch (L2–L5; L1 → 0)	rad
wheel[i]	per wheel i (FL0,FR1,RL2,RR3)	dict below
wheel[i].Fx, .Fy	tyre force, contact/wheel frame	N (+Fx drive, +Fy left)
wheel[i].Fz	vertical load	N
wheel[i].alpha	slip angle	rad
wheel[i].kappa	slip ratio	–
wheel[i].mu	road-surface friction at this wheel (contact μ)	–
wheel[i].mu_peak	realized tyre peak coefficient at this Fz (load-dependent)	–
wheel[i].alpha_peak	slip angle at the lateral force peak, this Fz/μ (inf=no peak)	rad
wheel[i].kappa_peak	slip ratio at the longitudinal force peak, this Fz/μ (inf=no peak)	–

The plant gives you the tyre's local operating-point characteristics as raw ground truth and leaves any metric to you. Per wheel you get:

• mu_peak — realized tyre peak coefficient at this Fz/μ (load-dependent). Rises above contact μ at low load (MF2002 PDY2 load sensitivity).
• alpha_peak, kappa_peak — slip position of the pure-axis force peak. This is the recommended core safety metric for limit-handling evaluation:

Recommended safety metric: slip-position classification

Compare current slip to peak slip per wheel. The force-slip curve peaks then falls, making a single "grip-usage" scalar ill-defined past the peak. Instead, expose the raw positions and classify directly:

|alpha| vs alpha_peak	interpretation
< alpha_peak	rising side — grip in reserve (understeer control room)
≈ alpha_peak	at peak — grip limit reached
> alpha_peak	past peak — drift / departure (sliding tail, no control recovery)

Same logic for |kappa| vs kappa_peak longitudinally. Why this works: - Unambiguous: Slip position directly tells you which side of the peak you are on (future grip available vs exhausted). - Load-dependent: Both alpha_peak and mu_peak account for vertical load transfer and MF2002 load sensitivity. - Control-relevant: Control laws can use reserve = (alpha_peak - |alpha|) / alpha_peak to detect margin loss before a slip event.

(For a friction-saturation ratio if you want one: sat = ‖[Fx,Fy]‖ / (mu_peak·Fz), bounded by 1. Use mu_peak, not mu: the realized peak rises above nominal at low load, so ‖F‖/(mu·Fz) can read >1 while the tyre is inside its own ellipse.)

Offline tyre query — combined-slip peak / friction ellipse

alpha_peak/kappa_peak in the obs are the pure-axis peaks (lateral peak at κ=0, long peak at α=0). Under combined slip the true peak shifts (lower slip, and the resultant grows toward the ellipse) — there is no single combined peak slip; it depends on the slip direction. Rather than bake that in, the plant exposes the tyre evaluator itself, so you reconstruct the exact combined force surface and find the peak in any direction yourself (same model the plant runs — no coefficient dump, nothing to re-implement):

Data access mirrors vehicle → part → physics (a read-only view for data delivery; the runtime itself stays flat):

plant.vehicle.params            # VehicleParams (mass/geometry/drivetrain/steering/brake/aero)
plant.vehicle.tire.params       # TireParams (backend, coefficients, lugre/belt)
plant.vehicle.tire.model        # live ITireModel  (== plant.tire_model shortcut)

import vdsim, math
# Preferred: the LIVE tyre the plant runs (single source of truth, no .tir re-load):
tire = plant.vehicle.tire.model            # or the shortcut plant.tire_model
# Or load a standalone instance from a .tir:
# tire = vdsim.create_magic_formula_tire_from_tir("external/chrono_parity/sample_pac02.tir")

def F(kappa, alpha, Fz=4000.0, mu=0.9):
    i = vdsim.TireInput(); i.Fz=Fz; i.kappa=kappa; i.alpha=alpha; i.mu_long=mu; i.mu_lat=mu
    o = tire.compute(i)            # also exposes o.mu_peak / o.alpha_peak / o.kappa_peak
    return o.Fx, o.Fy

# combined peak along a fixed slip direction (kappa:alpha ratio): 1-D search
def combined_peak(ratio, Fz=4000.0, mu=0.9):
    best, bestm = (0.0, 0.0), 0.0
    s = 0.0
    while s <= 0.5:
        fx, fy = F(ratio*s, s, Fz, mu); m = math.hypot(fx, fy)
        if m > bestm: bestm, best = m, (ratio*s, s)
        s += 0.001
    return best, bestm                # ((kappa*, alpha*), |F|*)

For an MF2002 .tir, the combined peak at κ:α=1 sits at a smaller α than the pure-lateral peak, with a larger resultant — the shift the pure alpha_peak cannot show.

---

5. Control input u = [delta, Fx_total]

• delta — road-wheel steer angle [rad] (what an MPC outputs; no hand-wheel ratio).
• Fx_total — total longitudinal force intent [N], +drive / −brake. The plant converts it to per-wheel drive/brake torque (50/50 axle split, drive_split_front) and applies it to the wheel-spin dynamics, so the tyre's combined-slip limit is real:
• moderate command → realised ΣFx ≈ Fx_total;
• over-command → wheels slip/lock, realised ΣFx < Fx_total = physical grip loss.

---

6. Low-μ patches (spatially varying friction)

6a. 1-D bands along the road (friction_map)

Pass friction_map as a list of (x_start, x_end, mu) along the road (world-x; straight road ⇒ s=x). Outside the patches, base_mu applies. At each patch boundary the grip transitions linearly over 1.0 m (entry and exit); x < x0 − 1 or x > x1 + 1 uses base_mu. Overlapping patches take the minimum μ. Friction is per-wheel by contact x, so the front axle enters a patch before the rear (real time lag).

plant = VDSimPlant(
    config="generic_sedan",
    base_mu=0.9,
    friction_map=[(40.0, 60.0, 0.5)],   # μ=0.5 ice patch from x=40 m to x=60 m
    control_dt=0.05, substep_dt=5e-4,
)

6b. 2-D polygon patches (friction_map_2d)

For arbitrary regions in the ground plane (e.g. a localized ice spot offset from the lane centre), pass friction_map_2d — a list of {"polygon": [(x0,y0), ...], "mu": 0.5} entries (world frame, ISO 8855). Each polygon needs ≥3 vertices; edges close implicitly (last → first).

Per-wheel μ is evaluated at the wheel contact (x, y). Inside a polygon: patch μ. Outside: linear blend to base_mu over 1.0 m from the polygon boundary (global blend_distance). Overlapping patches take the minimum μ (most restrictive grip). Mutually exclusive with friction_map (pass one or the other).

plant = VDSimPlant(
    config="generic_sedan",
    base_mu=0.9,
    friction_map_2d=[{
        "polygon": [(80.0, -2.0), (120.0, -2.0), (120.0, 2.0), (80.0, 2.0)],
        "mu": 0.5,
    }],
    control_dt=0.05, substep_dt=5e-4,
)

Low-level access: vdsim.create_polygon_friction_ground(z, base_mu, patches, blend_distance=1.0).

---

7. Steering actuator lag (optional)

By default the plant applies your road-wheel steer command instantaneously (ZOH at control_dt) — the right baseline for MPC. For robustness studies you can add a first-order lag on the realized steer angle using VehicleParams.steer_deadtime_s [s] as the time constant τ:

\[ \tau \dot\delta_\text{real} + \delta_\text{real} = \delta_\text{cmd} \]

(discretized with an exact exponential step each substep_dt). Set steer_deadtime_s = 0 (default) to recover the instantaneous path.

VDSimPlant does not expose steer_deadtime_s in its constructor yet. Use one of these manual workflows before the session is created:

A — add to the vehicle YAML (simplest):

steer_deadtime_s: 0.05    # [s] 1st-order steer lag (τ)

B — load, patch, write a temp config:

import tempfile
import vdsim
from vdsim_plant import VDSimPlant, resolve_vehicle_config

vp_path = resolve_vehicle_config("generic_sedan")
vp = vdsim.VehicleParams.from_yaml(str(vp_path))
vp.steer_deadtime_s = 0.05
with tempfile.NamedTemporaryFile(suffix=".yaml", delete=False) as f:
    tmp = f.name
vp.to_yaml(tmp)
plant = VDSimPlant(config=tmp, control_dt=0.05, substep_dt=5e-4)

Observe the delay — step response after a steer step:

obs = plant.reset(state0=[0, 0, 0, 16.7, 0, 0])
r0 = obs["r"]
for _ in range(20):
    obs = plant.step([0.05, 0.0])   # δ_cmd = 0.05 rad
print(obs["r"] - r0)              # yaw rate builds more slowly than with τ=0

Compare against steer_deadtime_s: 0 (or default) on the same manoeuvre; larger τ ⇒ slower obs["r"] / lateral response for the same delta command.

---

8. Closed-loop with your controller (drop-in)

The plant replaces a Python bicycle plant with no controller change. Pattern:

plant = VDSimPlant(config="generic_sedan", base_mu=0.9,
                   friction_map=[(40.0, 60.0, 0.5)],
                   control_dt=0.1, substep_dt=1e-3)   # control_dt = your MPC sample
obs = plant.reset(state0=[0, 0, 0, 16.7, 0, 0])
for k in range(N):
    # build the MPC state from obs (ISO 8855, same convention)
    x = [obs["X"], obs["Y"], obs["psi"], obs["vx"], obs["vy"], obs["r"]]
    delta, Fx = mpc.solve(x)             # your acados MPC — UNCHANGED
    obs = plant.step([delta, Fx])

To replace ~/vla_design/sim/closed_loop_sim.py's plant: keep the MPC; swap the plant object for VDSimPlant and feed step([delta, Fx]) instead of the bicycle plant_deriv/rk4.

Note: control_dt (the MPC sample, held ZOH) must be ≥ substep_dt, and substep_dt must divide control_dt evenly. The MPC at dt=0.1 with control_dt=0.1, substep_dt=1e-3 runs 100 internal substeps per control step.

---

9. Parameter file — generic_sedan.yaml

The default vehicle preset (configs/vehicles/generic_sedan.yaml) holds the solver defaults of vdsim.VehicleParams() — a generic mid-size sedan, not a measured car:

mass: 1500.0
mass_sprung: 1350.0
inertia_diag: [500.0, 2000.0, 2500.0]   # [Ixx, Iyy, Izz]; Izz is the yaw inertia
wheelbase: 2.7
cg_to_front: 1.2                          # lf
cg_to_rear: 1.5                           # lr
track_front: 1.55
track_rear: 1.55
cg_height: 0.55
wheel_radius_nominal: 0.32
drive_type: RWD
steering_ratio: 15.0
aero_drag_coeff: 0.3
frontal_area: 2.2
tire_yaml: parts/tire/generic_pacejka.yaml  # the tyre sidecar

(excerpt; the file also carries suspension, drivetrain and brake fields). Its tyre configs/parts/tire/generic_pacejka.yaml is the MF96 default tyre (constant B/C/D/E).

MF2002 (.tir) tyres

A tyre YAML with backend: mf2002 and tir_path: <file>.tir evaluates the full Pacejka-2002 model from a coefficient file. No .tir ships in the package: .tir files usually hold measured data, so they stay with their owner (see Private vehicles below). The repository's public synthetic external/chrono_parity/sample_pac02.tir (hand-authored, no measured data) is the test input. The load-shaping coefficients:

coeff	role
PKY1, PKY2	cornering-stiffness magnitude + its saturation with load (Ca concave in Fz)
PDY1, PDY2	lateral peak factor + peak-μ fall with load (μ_peak ↓ as Fz↑)
PCY1, PEY1	lateral shape / curvature
PKX1, PDX1, PDX2, PCX1	longitudinal stiffness, peak, load sensitivity, shape

Making your own vehicle

Copy generic_sedan.yaml, edit mass/geometry. The MF2002 cornering stiffness per wheel is Kya(Fz) = PKY1 · FNOMIN · sin( 2·atan( Fz / (PKY2·FNOMIN) ) ) · LKY — load-dependent, so a single .tir reproduces different front/rear axle stiffness automatically from each axle's static Fz (no per-axle file needed). To target known axle stiffnesses, pick PKY1 (overall magnitude) and PKY2 (where it saturates) so that 2·Kya(Fz_front)=Ca_f and 2·Kya(Fz_rear)=Ca_r.

Private vehicles

A car whose parameters must not be published lives outside the package, in a root with the configs/ layout: <root>/vehicles/<name>.yaml, <root>/parts/tire/<tyre>.yaml (+ its .tir). Set VDSIM_PRIVATE_CONFIGS=<root>; VDSimPlant(config="<name>.yaml") and vdsim_rl then search that root before the bundled presets, and the vehicle's tire_yaml is resolved in the same root first.

---

10. API reference

VDSimPlant(config="generic_sedan",   # vehicle preset name or path
           friction_map=None,     # [(x0, x1, mu), ...] 1-D bands along x
           friction_map_2d=None,  # [{"polygon": [(x,y),...], "mu": ...}, ...]
           base_mu=0.9,           # (0, 1.2]
           control_dt=0.05,       # s, one step() period (ZOH); >= substep_dt
           substep_dt=5e-4)       # s, internal fixed integration step; divides control_dt

plant.reset(state0=[X, Y, psi, vx, vy, r]) -> obs   # seeds wheel spin = vx / R_eff
plant.step([delta_rad, Fx_total_N])         -> obs   # advance one control period

Inputs are validated; bad ones raise a clear ValueError/TypeError (see Troubleshooting).

---

11. What's validated (so you can trust it)

Gated by ctest -R VlaPlant; the current full-suite result and measurement configuration are tracked in Validation: - Dry handling = linear bicycle within 15 % (steady yaw-rate gain vs the single-track model with the same Caf/Car) — your bicycle results carry over. Checked for both the analytic-B tyre and an MF2002 .tir (the public synthetic sample_pac02.tir, axle stiffness taken from the same file at the static loads). - Load-dependent tyre (MF2002): cornering stiffness concave in Fz (Kya(2Fz0) < 2·Kya(Fz0)) and peak μ falls with load — the load-transfer limit behaviour MF96 could not represent. - Friction circle held: ‖[Fx,Fy]‖ ≤ mu·Fz every step; over-command ⇒ real slip + the vehicle cannot deliver the commanded Fx (grip loss), and a brake never spins a wheel backwards. - Force ↔ accel consistency: ΣFy ≈ m·ay. - Throttle bypass: Fx is applied as torque to the wheel spin, not a pedal map. - Speed: a 5 s trajectory runs in ~25 ms (≫ real time → cheap sweeps).

This is an independent-of-the-controller mid-fidelity plant; the MF2002 load-dependent tyre ≠ your bicycle internal (linear-Ca) model is an intended plant–model mismatch — a stronger robustness story than a matched plant, and exactly the Fz-driven nonlinearity a σ_Fz / chance-constraint study wants to stress.

---

12. Tips & gotchas

• Determinism: no RNG; identical inputs ⇒ identical outputs. Keep substep_dt fixed for reproducible figures.
• Two-tier dt: control_dt = your control rate; substep_dt ≤ 1 ms (wheel spin is stiff). Coarsening substep_dt can make hard-brake locking inaccurate.
• Frames: wheel.Fx/Fy are in the contact (wheel) frame — correct for the friction circle; do not re-rotate by steer.
• Low-speed: the plant uses a tyre VLOW floor (no kinematic-blend fade), so the saturating tyre stays active even near a stop.
• Beta: report rough edges — API and defaults may be tuned from your usage.

---

13. Troubleshooting

symptom	cause / fix
ValueError: friction_map[i]: x0>=x1	patch start must be < end
ValueError: ... mu=... outside (0, 1.2]	friction must be in (0, 1.2]
ValueError: substep_dt must divide control_dt evenly	pick e.g. control_dt=0.05, substep_dt=5e-4 (×100)
ValueError: step(u): u must be [delta_rad, Fx_total_N]	pass a 2-element [delta, Fx]
FileNotFoundError: ... config	use a bundled name (generic_sedan) or a car under $VDSIM_PRIVATE_CONFIGS/vehicles/
import vdsim fails	build first (cmake --build build -j) and add build/python to sys.path
car won't slow under braking	Fx is a force intent; on low μ the realised brake force is capped at grip — expected

---

14. Planned features (P2 — next release)

Based on BETA #1 customer feedback (docs/PLANT_CUSTOMER_FEEDBACK.md):

Split-μ friction map (left/right different grip)

friction_map_left=[(x0, x1, mu), ...], friction_map_right=[...] - C++ backend exists (create_split_mu_ground); not yet exposed to Python plant API. - Use case: ABS/ESC algorithm testing, split-μ braking scenarios.

---

15. Related documents

• Customer feedback and validation: docs/PLANT_CUSTOMER_FEEDBACK.md
• Full spec and design decisions: docs/design/VLA_THESIS_PLANT.md
• Computational cost per dynamics level: docs/CATALOG_AND_PHYSICS.md
• Bundled example: examples/vla_plant_demo.py