A first-pass Julia port of autoLie focused on the core Lie derivative / Lie bracket workflow and a few representative control-affine dynamics.
Attribution — The original
autoLiePython library was written by Vineet Tiruvadi and is released under the MIT License. This Julia port inherits that license; see LICENSE for details.
Most physical systems of interest can be written in control-affine form:
ẋ = f(x) + Σ gᵢ(x) uᵢ
where f is the autonomous drift vector field, gᵢ are control vector fields, uᵢ are scalar inputs, and x ∈ ℝⁿ is the state. Understanding whether the system can be steered between arbitrary states is the controllability problem.
The Lie derivative of a scalar function h along a vector field f is defined as:
Lf h(x) = ∇h(x) · f(x)
For a vector field g, the iterated Lie derivative Lf^k g describes how g evolves under the flow of f. Numerically, this is simply the directional derivative of one field with respect to another.
The Lie bracket of two vector fields f and g is:
[f, g](x) = (∂g/∂x) f(x) − (∂f/∂x) g(x)
It measures the infinitesimal rotation between the two flows. Physically, if you follow f, then g, then −f, then −g for equal small times, the net displacement is approximately [f, g] scaled by the time squared.
A nonlinear control-affine system is locally accessible at x₀ if the involutive closure (Lie algebra generated by iterated brackets of f and all gᵢ) spans ℝⁿ at x₀. This is a generalisation of the classical Kalman rank condition to nonlinear settings:
rank( g₁, [f,g₁], [f,[f,g₁]], …, [gᵢ,gⱼ], … ) = n
The condition is satisfied generically (at almost every point) for most well-posed systems, but fails when control authority is structurally restricted — exactly the kind of situation that motivates the analysis in autoLie.
The dual question is observability: can the state be inferred from a measurement output y = h(x)? The Lie derivative Lf h encodes how much information the drift dynamics contribute to the observed output, while Lg h encodes the measurement–control coupling. When Lg h = 0 the control input has no immediate effect on the measurement, which directly informs sensor/actuator placement strategies.
| Domain | What Lie analysis tells you |
|---|---|
| Adaptive DBS (deep brain stimulation) | Optimal placement of recording and stimulation probes for maximal disease-state coverage; whether a given electrode configuration can observe/control a target symptom subspace. |
| Network neuroscience | Which nodes in a network are structurally controllable; rank of the controllability Gramian for reduced brain-network models. |
| Robotics | Kinematic controllability of underactuated manipulators and mobile robots (e.g., the car-parking problem is solved via the Lie bracket of steering and forward motion). |
| Chemical processes | Accessibility analysis of reaction networks with few actuated species. |
| Power systems | Nonlinear reachability for grid stabilisation controllers. |
L_dandL_bracketequivalents for numerical Lie analysis- A lightweight
Operablewrapper mirroring the original operator-composition idea - Representative 3D dynamics from the Python codebase (
f1,f4,f8,f9,g1,h1,f_main, etc.) - A simplified DBS-style control system example
- Makie-based plotting helpers that save PNG figures from 2D slices of the 3D fields
- Dense-grid field evaluation via
Tullio.jl(slice_field_dense,slice_scalar_dense) - Numeric controllability rank checks via Chow's theorem (
controllability_distribution,controllability_rank,chow_rank_test) over random sample points Graphs.jl/GraphMakie.jlnetwork visualisation for richer DBS layout plots (seeexamples/dbs_network_graphmakie.jl)Symbolics.jlsymbolic codepath —sym_L_d,sym_L_bracket,sym_controllability_matrix, andsym_to_numericfor exact closed-form Lie derivatives, brackets, and generated numeric functions (seesrc/symbolic.jl)- Richer trajectory simulation (
simulate_trajectory) and parameterized DBS presets (dbs_scenario,DBSScenario) - Expanded control-affine catalog including canonical nonholonomic and SE(2) systems (
f_brockett_drift,g_brockett_1,g_brockett_2,f_dubins_drift,g_dubins_drive,g_dubins_turn) plus a controlled Lorenz drift (f_lorenz,g_lorenz_z) - A placement-focused Makie figure that maps
L_g hacross control/readout node pairs (examples/sensor_actuator_placement.jl)
# generate basic Lie derivative figures
julia --project=. examples/basic.jl
# generate DBS network figures (Makie version)
julia --project=. examples/dbs_network_ctrl.jl
# run expanded control-affine examples + Chow rank checks
julia --project=. examples/control_affine_examples.jl
# generate richer DBS network figure (Graphs.jl + GraphMakie.jl)
julia --project=. examples/dbs_network_graphmakie.jl
# generate sensor/actuator placement coupling map from L_g h
julia --project=. examples/sensor_actuator_placement.jl
# benchmark dense-grid backends + compare parameterized DBS trajectories
julia --project=. examples/improvements_demo.jl
# run tests
julia --project=. -e 'using Pkg; Pkg.test()'This dashboard is the primary placement workflow and has four panels:
- A (strict toy model) confirms the structural rule: off-diagonal entries are zero (
L_g h = 0), so control at nodeihas no immediate measurement effect at nodej != i. - B (kernel-spread model, nominal state) relaxes that strict assumption by allowing spatially spread control and readout kernels, producing graded pairwise coupling.
- C (robustness map) averages
|L_g h|over drift-trajectory state samples to favor placements that remain informative away from a single operating point. - D (ranked candidates) scores each pair with
mean(|L_g h|)/(1 + std)and highlights the top practical sensor/actuator choices.
Together, these panels move from a binary structural check to a ranking that directly supports placement decisions.
The graph view provides structure (who is connected), while the L_g h matrix above provides immediate measurement-coupling analytics (who is directly observable from control action at first order).
The left panel shows the in-plane (x-y) drift projection and is near-zero on z = 0 for this field. The right panel shows the out-of-plane component (w = z-dot), which explains why a pure in-plane quiver can look blank even when drift is active.
The scalar heatmap shows Lf h(x) = ∇h · f where h(x) = 2x₁ + 3x₃. This is the drift-only contribution to measurement change; together with L_g h, it separates autonomous measurement motion from control-induced measurement motion.
The top panel compares dense field-slice evaluation runtime on a 181 x 181 grid using the legacy nested-loop backend and the new Tullio.jl backend. The bottom panel shows readout trajectories from three parameterized DBS presets (:baseline, :early_drive, :high_gain) integrated with the new simulate_trajectory API.
The Python repository is intentionally exploratory and numerically direct. This port keeps that spirit by using explicit function arguments, a small catalog of sample fields, and grid-slice plotting rather than a new abstraction-heavy framework.
- Formalize a reproducible
BenchmarkTools.jlbenchmarking plan.- Scope and benchmark groups: split into
Lie derivative,Lie bracket,controllability_distribution,slice_field_dense, andsimulate_trajectorygroups so each algorithmic path is measured independently. - Deterministic fixtures: freeze RNG seeds, fixed state dimensions (
n = 3, 10, 50), and representative field sets (legacyf8, DBSf_trivial/g_mono, symbolic-generated numeric functions) to avoid drift in timing comparisons. - Parameter sweeps: benchmark grid sizes (
n = 41, 81, 161, 321) and Lie depth (1..4) to expose scaling behavior, not just single-point latency. - Metrics: capture median time, memory, allocations, and variance (
Trialstatistics + confidence intervals) and export machine-readable summaries (CSV/JSON) for regression tracking. - CI integration strategy: keep heavy benchmarks optional (nightly/manual), while adding lightweight guardrail benchmarks in PR checks with conservative thresholds to catch major regressions.
- Reporting: auto-generate a benchmark report figure/table in
figures/and include hardware/Julia-version metadata for reproducibility.
- Scope and benchmark groups: split into
- Current scope constraints (kept intentionally narrow for this first pass). No GUI frontend, no Mayavi/full 3D interactive stack, no broad graph abstraction layer, no large API-shape refactor, and no attempt to generalize far beyond the curated control-affine example catalog.