multicalc

Scientific computing that fits on a microcontroller, built and tested from scratch in one integrated package. Estimation, control, kinematics, Lie groups, calculus, autodiff and linear algebra in stable no_std Rust with no heap, no panics and no unsafe.

combined_demo_reel.mp4

A reel of the live showcase demos: a quadcopter performing state estimation, planning and control to fly waypoints, 2D robot running particle filter localization + EKF sensor fusion + obstacle avoidance over a 1kHz loop rate; then a Franka Panda arm, loaded from its MuJoCo model file, tracking a moving 3D pose. Every number on screen is measured live, inside a 1 ms tick.

Highlights

• 1 kHz loop rates: No heap, fixed-size types, bounded work per call. Results in a full robotics control loop at 1 kHz.
• Tested on six embedded targets: Every commit is built and tested on six targets: the x86_64 and aarch64 Linux hosts and on four bare-metal ABIs (thumbv7em soft-float, thumbv7em hardware-FPU, thumbv6m, and riscv32imc), running the real math under QEMU. no_std, no-alloc, and no-panic rules hold on each target.
• Measured against external references: Each module's results are verified against established libraries like numpy, scipy, and filterpy fixtures within ~1 ulp, thus validating the rust implementation. See the benchmarks.
• Pure safe and panic-free.#![forbid(unsafe_code)], no C dependencies, and unwrap/ panic denied on library paths; every fallible call returns a typed error. Types are fixed-size and stack-allocated, and iteration counts are bounded.

What it does

Robotics and control

• Estimation: linear, extended, and unscented KalmanFilters (autodiff Jacobians, no hand-derived ones; the unscented one needs no derivatives at all), an ErrorStateKalmanFilter that fuses an IMU with position and heading fixes, MahonyFilter and MadgwickFilter for attitude estimation, and a ParticleFilter for nonlinear, non-Gaussian problems (alloc only), with a Monte Carlo Localization built on top of it.
• Control: Pid control, infinite horizon Lqr, GeometricAttitudeController for drones, the pure pursuit path-following law, and FollowTheGap reactive obstacle avoidance. Model-based torque control on any ArticulatedBody: ComputedTorqueController (feedback linearization through H(q)), JointImpedanceController and CartesianImpedanceController (spring-damper in joint and tool space), and JointPdController (PD with optional gravity compensation).
• Spatial math: Quaternion, the SO2/SE2/SO3/SE3 Lie groups with left/right Jacobians and inverses, and Twist/Wrench spatial algebra in [v; ω] — Plücker transforms, 6×6 adjoints, motion and force cross products, and SpatialInertia momentum, bias wrench, energy and composition.
• Rigid-body dynamics: RigidBody for a single body, ArticulatedBody for a jointed robot — inverse dynamics, the joint-space inertia matrix and forward dynamics (RNEA, CRBA, ABA) with armature, viscous damping and Coulomb friction, checked against Pinocchio. Both read from an MJCF or URDF file with multicalc-robot-model.
• Plant: Actuator models between a command and the wrench a body feels — MultirotorMixer maps a wrench to rotor thrusts and back, RotorLag a rotor's first-order thrust lag, and PositionServo a position-commanded joint's second-order servo, discretized exactly.
• Kinematics: KinematicTree for revolute/prismatic/continuous/fixed/floating chains and forward and inverse kinematics, generic over the scalar with autodiff, built by hand or read from an MJCF or URDF file; a damped-least-squares SE(3) pose solver with joint limits and null-space redundancy resolution.
• Collision checking: CollisionQuery for sphere/capsule proximity — primitives on tree frames against each other and against world-fixed obstacles, with pair exclusions and fixed capacities.
• Motion: PolylinePath for waypoint paths with arc-length, closest-point, and lookahead queries, MinimumSnapPlanner for the smoothest trajectory through them, and MotionProfilePlanner for jerk-limited point-to-point moves with multi-axis synchronization.
• Mapping: 2D OccupancyGrid and ScanGeometry

Core math

• Automatic differentiation: Exact autodiff of any order (total and partial), plus Jacobian and Hessian matrices.
• Linear algebra: fixed-size, stack-allocated Matrix and Vector with LU, Cholesky, column-pivoted QR, SVD, symmetric eigendecomposition, and the matrix exponential expm. General N×N determinant and inverse, pseudo-inverse, eigenvalue clamping, zero-copy MatrixView / VectorView, solve_discrete_riccati and solve_discrete_lyapunov.
• Least-squares optimization: LevenbergMarquardt and GaussNewton solvers for nonlinear curve fitting.
• Root finding: bracketed bisection and Newton solvers for scalar equations and square systems, with an optional damped line search.
• Polynomials: Polynomial for evaluation with any number of derivatives in one pass, arithmetic, calculus, fitting and real roots; PiecewisePolynomial for curves made of pieces; and MultivariatePolynomial for several variables with symbolic partial derivatives.
• Integration: iterative Newton-Cotes rules (Boole, Simpson, Trapezoidal) and Gaussian quadrature (Legendre, Hermite, Laguerre) over finite, semi-infinite, and infinite limits.
• ODE integrators: fixed-step Rk4 and adaptive Rk45 (Dormand-Prince 5(4)) with PI step control and dense output, plus ExponentialMap, which is a purely orientation integrator.
• Discretization: zero-order hold, Van Loan, and discrete white-noise models for continuous-time linear systems.
• Signal processing: Biquad low-pass, high-pass, band-pass, and notch filters; with cascades, motor-harmonic notches, and per-channel filtering. Plus MovingAverage, RunningMedian, SavitzkyGolay smoothing, Deadband, Hysteresis and SlewRateLimiter conditioning.
• Vector calculus: curl, divergence, and line and flux integrals.
• Approximation: linear and quadratic Taylor models with goodness-of-fit metrics.
• Random: Pcg32 and the RandomSource trait, a seedable no_std generator for the particle filter and for stochastic models.

Quick start

Two formulas, written once, carried through six modules, each step feeding the next:

use multicalc::prelude::*;use multicalc::{Hessian,Jacobian,KalmanFilter,KalmanModel,Matrix,Newton,SE3,SO3,Vector, constant };use multicalc::{scalar_fn, scalar_fn_vec};fnmain() -> Result<(),CalcError>{// Written once, evaluated at f64 here and at an autodiff number wherever a derivative is asked// for — the formula text never changes.let f = scalar_fn!(|x| x * x * x - constant(2.0)* x);// f(x) = x³ - 2xlet g = scalar_fn!(|v:&[f64;2]| v[0]* v[0]* v[1] + v[0].sin());// g(x, y) = x²y + sin x// Derivatives — exact, by forward-mode autodiff. No step size, no truncation error.let single_point = 2.0_f64;let slope = derivative(&f, single_point);// f'(2) = 10let bend = second_derivative(&f, single_point);// f''(2) = 12let point = [1.0_f64,2.0];let x_index = 0;let dg_dx = partial(&g, x_index,&point)?;// The derivative matrices of those same two formulas.let hessian = Hessian::new().evaluate(&g,&point)?;// 2x2 second derivativeslet both = scalar_fn_vec!(|v:&[f64;2]| [
v[0]* v[0]* v[1] + v[0].sin(),
v[0]* v[0]* v[0] - constant(2.0)* v[0],]);let jacobian = Jacobian::new().evaluate(&both,&point)?;// 2x2 first derivatives// Integration — f again, this time over an interval.let limits = [0.0,2.0];let area = integral(&|x:f64| f.eval(x), limits)?;// ∫₀² f = 0// Linear algebra — solve H·x = b with the Hessian computed three lines up.let b = Vector::new([1.0,2.0]);let x = hessian.solve(b)?;// Root finding — Newton on the same f, its derivative supplied by autodiff.let initial_guess = 2.0;let root = Newton::new().solve(&f, initial_guess)?.root;// √2 ≈ 1.41421356// Rigid-body motion — SO(3)/SE(3), generic over the scalar like everything above.let quarter_turn_about_z = Vector::new([0.0,0.0, core::f64::consts::FRAC_PI_2]);let translation = Vector::new([1.0,2.0,3.0]);let start = Vector::new([1.0,0.0,0.0]);let pose = SE3::from_parts(SO3::exp(quarter_turn_about_z), translation);let moved = pose.act(start);// rotate, then translate → (1, 3, 3)// Estimation — a Kalman filter recovering the velocity it never measures.let initial_state = Vector::new([0.0,0.0]);// [position, velocity]let initial_covariance = Matrix::new([[1.0,0.0],[0.0,1.0]]);let model = KalmanModel{state_transition:Matrix::new([[1.0,1.0],[0.0,1.0]]),measurement_model:Matrix::new([[1.0,0.0]]),// position onlyprocess_noise:Matrix::new([[0.01,0.0],[0.0,0.01]]),measurement_noise:Matrix::new([[0.1]]),};letmut filter = KalmanFilter::new(initial_state, initial_covariance, model);
filter.predict();let measurement = Vector::new([1.0]);// the target moved about 1 m
filter.update(measurement)?;let velocity = filter.state()[1];// recovered, though never measuredOk(())}

Every fallible call propagates with ?: each module has its own error enum, and all of them convert into the CalcError umbrella, so one return type covers a program that mixes modules.

Full tutorial

Refer to the tutorials for a comprehensive tutorial for each module. They show the full imports, expected outputs in comments, error-path notes, and pointers to runnable demos. Start there when you need the complete picture of a feature.

Accuracy

Verified against external-library fixtures (mpmath, numpy, scipy, filterpy) in the multicalc-qa crate, with per-module tables generated from those fixtures. See benchmarks/README.md for the index, or go straight to calculus, linear_algebra, optimization, ode, estimation, kinematics, dynamics, or root_finding.

Runnable demos

Runnable, self-contained programs for each module live in the demos/ crate. See demos/README.md. Run one with:

cargo run -p multicalc-demos --example <name>

Documentation

• Tutorials: A worked page for every module: imports, a runnable snippet, error paths, and a demo pointer.
• Crate README / API docs: the crates.io page and full API reference, with notes on no_std, error handling, and heap allocation.
• Examples: Self-contained, self-checking programs for each module in the demos/ crate. Run one with cargo run -p multicalc-demos --example <name>.
• Benchmarks: Per-module accuracy tables and latency measurements, generated from the QA fixtures and checked in CI.
• Live showcases: Six animated Rerun demos, including a quadcopter using particle filter localization then using state estimation to fly a planned loop, a 2D robot using monte carlo localization, then fusing wheel odometry+IMU+GPS for state estimation to lap a course of obstacles. The others are a Franka Panda arm, loaded from its MuJoCo model file tracking a moving 3D pose, a Newton fractal, Fourier epicycles drawing Ferris, and gradient-driven marbles, each streaming live-measured speed and accuracy.
• QA crate: multicalc-qa holds the CI-enforced accuracy fixtures and generates the benchmarks tables from them.

Repository layout

The published library crate lives in crates/multicalc; the repository root is a Cargo workspace. Runnable demos live in the dev-only demos/ crate (basics and live Rerun showcases), and tools/embedded-smoke runs multicalc on the four bare-metal targets (three Cortex-M targets + riscv32imc) under QEMU every PR. crates/multicalc-robot-model reads MuJoCo MJCF and URDF model files into multicalc's robot types, and ships a model_viewer binary that draws any of them in a Rerun viewer. See README for more details.

Contributing

Contributions are welcome. See CONTRIBUTING.md.

Acknowledgements

The least-squares solvers and QR factorization port the public-domain MINPACK routines (Moré, Garbow, Hillstrom; netlib); the full citation is in the crate README.

License

Licensed under the MIT License.

Contact

anmolkathail@gmail.com