xad-rs

Exact automatic differentiation for Rust — forward-mode, reverse-mode, first- and second-order, with named variable support for ergonomic gradient readback.

Unofficial Rust port of the C++ XAD library. Not affiliated with the upstream project.

Choosing a mode

(1) Program your numerical code against the trait [Real]. (2) Instantiate it with the concrete mode that matches your problem shape.

use xad_rs::prelude::*;

fn quadratic<R: Real>(x: &R) -> R {
    x.clone() * x.clone() + R::from(2.0_f64) * x.clone() + R::from(1.0_f64)
}

// Same body, four call sites:
let v_passive = quadratic(&3.0_f64);                       // no AD
let v_jet1    = quadratic(&Jet1::new(3.0_f64, 1.0));       // forward 1st order
let v_jet2    = quadratic(&Jet2::variable(3.0_f64));       // forward 1st + 2nd
// Reverse mode: instantiate inside an active Tape.

The forward-mode types are named after k-jets — Taylor-coefficient bundles from differential geometry. The number is the order, and the Vec suffix means the bundle is propagated against a multi-variable tangent vector (full gradient / full Hessian) rather than a single seeded direction.

In v0.5.0, Real is implemented for f64, AReal<f64>, Jet1<f64>, and Jet2<f64>. Jet1Vec / Jet2Vec impls and the f32 flavours are tracked as follow-up changes (length-of-partial-vector and lossy-From<f64> issues respectively).

Every mode also has a named variant (NamedJet1Vec, NamedTape, etc.) that lets you read gradients by variable name instead of positional index.

Installation

[dependencies]
xad-rs = "0.6"

MSRV: 1.85 (Rust edition 2024).

Migration from 0.4 / 0.5 / 0.6

From 0.4 → 0.5

• Renamed trait: xad_rs::Scalar → xad_rs::Passive.
• New trait: xad_rs::Real — the unified active-scalar trait.
• New prelude: use xad_rs::prelude::*;.
• Bug fix: math::norm_pdf √2 magnitude error fixed.

From 0.5 → 0.6 (breaking changes)

• Removed: the deprecated alias xad_rs::Scalar. Use xad_rs::Passive directly.
• Trait shape: xad_rs::Real no longer requires Default as a supertrait. Generic code that wrote where R: Default should add the bound explicitly (the inherited bound is gone).
• Added: xad_rs::RealStats — extension trait providing erf, erfc, norm_cdf, inv_norm_cdf for f64, AReal<f64>, and Jet1<f64>. Use R: RealStats for option-style pricers.
• Added: Tape::activate_guard() returning an RAII TapeGuard<T> whose Drop deactivates the tape. The existing Tape::activate() / Tape::deactivate_all() APIs remain available.

From 0.6 → 0.7 (breaking changes)

This release trims the public surface to items with at least one in-repo caller and unblocks the build. Every removal below has a recommended replacement.

• Removed C++ XAD-compat alias modules at the crate root: xad_rs::adj, xad_rs::adjf, xad_rs::fwd, xad_rs::fwdf. Replace with the underlying types directly:
  • adj::TapeType → Tape<f64>
  • adj::ActiveType → AReal<f64>
  • adj::PassiveType → f64
  • adjf::* / fwdf::* — analogous f32 types
  • fwd::ActiveType → Jet1<f64> (note: xad_rs::math::fwd, the forward-mode AD math module, is unrelated and stays)
• Removed unused Tape methods:
  • Tape::push_statement — use Tape::push_unary / Tape::push_binary (the fixed-arity hot-path equivalents that the crate already uses everywhere internally).
  • Tape::compute_adjoints_to(pos) — no equivalent. Use the standard Tape::compute_adjoints to sweep the whole tape; the partial-sweep path had no callers and is dropped.
  • Tape::clear_derivatives_after(pos) — no equivalent (no callers).
  • Tape::increment_adjoint(slot, v) — no equivalent (no callers).
  • The pair Tape::get_position + Tape::reset_to(pos) is kept — that rollback flow is still supported.
• Removed AD-wrapped admin functions that propagated no derivative information and only shadowed f64::* methods: xad_rs::math::ad::is_nan, is_infinite, is_finite, is_normal, signum, floor, ceil, round, trunc, fract. Replace with x.value().is_nan() (etc.).
• Removed Jet2Vec::symmetrize — Jet2Vec ops produce structurally symmetric Hessians; if a downstream user accumulates a Hessian outside the Jet2Vec hot path, average H and H.t() directly with ndarray.
• Removed Jet2Vec::powd — for an active Jet2Vec exponent, route through exp(y * ln(x)) explicitly (y.clone() * x.clone().ln()).exp()).
• Removed impl From<i32> for AReal<f32> — Real is implemented only for AReal<f64> so the f32 conversion was unconstrained. Use AReal::new(value as f32).
• Build hygiene: Cargo.toml no longer declares the 5 dead [[bench]] targets and the criterion dev-dependency is gone. The README's benchmarking section is also removed.
• Internal refactor (no public behavior change): the four ref-permutation impls per binary +/-/*// on AReal<T> are now stamped via impl_areal_binop! macros, mirroring the existing __named_areal_binop! and impl_unary_math! patterns. All operator paths produce identical numerical results.
• Examples: examples/hessian.rs now uses the exact compute_full_hessian helper (Jet2Vec based, machine-precision) instead of the finite-difference compute_hessian. The finite-difference helper is still in the public API for consumers who specifically want it.

From 0.7 → 0.8 (breaking changes)

This release consolidates the named-forward AD path: one lifecycle for every named-forward type, type-generic vec types, and a uniform declare_* naming. Every removal/rename below has a recommended replacement.

• Unified named-forward lifecycle. The eager-construction half of NamedForwardTape is gone; everything now goes through declare → into_scope:

  • Removed NamedForwardTape::input_jet1<T>, NamedForwardTape::constant_jet1<T>, NamedForwardTape::freeze(), NamedForwardTape::is_frozen(), the tape-side NamedForwardTape::registry() accessor.
  • Replaced by: NamedForwardTape::<T>::declare_jet1(name, value) returning a Jet1Handle<T>, then tape.into_scope() to obtain a NamedForwardScope<T>, then scope.jet1(handle) for the value. scope.constant_jet1(value) produces a constant; scope.registry() is the new accessor.

  Before:

  let mut ft = NamedForwardTape::new();
  let x: NamedJet1<f64> = ft.input_jet1("x", 3.0);
  let _registry = ft.freeze();
  let f = &x * &x + &x;
  

  After:

  let mut ft = NamedForwardTape::<f64>::new();
  let x_h = ft.declare_jet1("x", 3.0);
  let scope = ft.into_scope();
  let x = scope.jet1(x_h);
  let f = x * x + x;
  

• NamedForwardTape and NamedForwardScope are now generic over T: Passive (with T = f64 default). A single tape/scope picks one underlying scalar type for the whole problem. Jet1Handle<T>, Jet1VecHandle<T>, and Jet2Handle<T> are typed so the wrong scope cannot be indexed with a foreign-T handle. Existing call sites that used the f64-only path keep working without turbofish thanks to the default; explicit <f64> (or other T) only needed if a literal is ambiguous.

• Jet1Vec and Jet2Vec are now generic over T: Passive (with T = f64 default). NamedJet1Vec<T> similarly. The genericity opens the f32-mode AD path for multi-variable forward mode that already existed for the single-direction Jet1<T> / Jet2<T> siblings. Most call sites are unchanged; ambiguous-literal cases (e.g. Jet1Vec::variable(2.0, 0, 2) with no other constraint) need an explicit _f64 literal suffix or : Jet1Vec<f64> annotation.

• declare_jet2_f64 renamed to declare_jet2 (and now generic over T). Same for constant_jet2_f64 → constant_jet2. The _f64 suffix existed because the previous tape stored only f64; with the type-generic tape it's a misnomer.

• Internal: forward_tape.rs was rewritten as a mono-T builder. The TLS generation guard (check_gen in debug builds) is unchanged; the cross-registry-op debug panic still catches mixing values from different tape scopes.

Quick start

Reverse mode

use xad_rs::{AReal, Tape, math};

let mut tape = Tape::<f64>::new(true);
tape.activate();

let mut x = AReal::new(3.0);
let mut y = AReal::new(4.0);
AReal::register_input(std::slice::from_mut(&mut x), &mut tape);
AReal::register_input(std::slice::from_mut(&mut y), &mut tape);

// f(x, y) = x^2 * y + sin(x)
let mut f = &(&x * &x) * &y + math::ad::sin(&x);
AReal::register_output(std::slice::from_mut(&mut f), &mut tape);
f.set_adjoint(&mut tape, 1.0);
tape.compute_adjoints();

println!("df/dx = {}", x.adjoint(&tape));  // 2xy + cos(x)
println!("df/dy = {}", y.adjoint(&tape));  // x^2

Forward mode (full gradient)

use xad_rs::Jet1Vec;

let (x, y) = (Jet1Vec::variable(3.0, 0, 2), Jet1Vec::variable(4.0, 1, 2));
let f = &(&x * &x) * &y;  // x^2 * y

assert_eq!(f.partial(0), 24.0);  // df/dx = 2xy
assert_eq!(f.partial(1),  9.0);  // df/dy = x^2

Second-order derivatives

use xad_rs::Jet2;

let x = Jet2::variable(2.0_f64);
let y = x * x * x;  // x^3
assert_eq!(y.first_derivative(),  12.0);  // 3x^2
assert_eq!(y.second_derivative(), 12.0);  // 6x

Named variables

Access derivatives by name — useful in financial models with many risk factors:

use xad_rs::{NamedForwardTape, NamedForwardScope};

let mut ft = NamedForwardTape::<f64>::new();
let spot_h   = ft.declare_jet1vec("spot",   100.0);
let strike_h = ft.declare_jet1vec("strike", 105.0);
let scope: NamedForwardScope<f64> = ft.into_scope();

let spot   = scope.jet1vec(spot_h);
let strike = scope.jet1vec(strike_h);
let ratio  = spot / strike;

assert!((ratio.partial("spot") - 1.0 / 105.0).abs() < 1e-14);

Named reverse mode returns gradients as IndexMap<String, f64>:

use xad_rs::NamedTape;

let mut tape = NamedTape::new();
let x = tape.input("x", 3.0);
let y = tape.input("y", 4.0);
let _registry = tape.freeze();

let f = &(&x * &x) * &y + x.sin();
let grad = tape.gradient(&f);

assert!((grad["x"] - (2.0 * 3.0 * 4.0 + 3.0_f64.cos())).abs() < 1e-12);
assert!((grad["y"] - 9.0).abs() < 1e-12);

Jacobian and Hessian

use xad_rs::{compute_jacobian_rev, compute_hessian};

// f: R^2 -> R^2, f(x, y) = [x*y, x + y]
let jac = compute_jacobian_rev(&[3.0, 5.0], |v| {
    vec![&v[0] * &v[1], &v[0] + &v[1]]
});

// g: R^2 -> R, g(x, y) = x^2 * y + y^3
let hess = compute_hessian(&[2.0, 3.0], |v| {
    let x2 = &v[0] * &v[0];
    let y3 = &v[1] * &v[1] * &v[1];
    x2 * &v[1] + y3
});

Dense full Hessian (Jet2Vec)

use xad_rs::Jet2Vec;

let x = Jet2Vec::variable(1.0, 0, 2);
let y = Jet2Vec::variable(2.0, 1, 2);
let f = &(&(&x * &x) * &y) + &(&(&y * &y) * &y);

assert_eq!(f.hessian()[[0, 0]], 4.0);   // d2f/dx2 = 2y
assert_eq!(f.hessian()[[0, 1]], 2.0);   // d2f/dxdy = 2x
assert_eq!(f.hessian()[[1, 1]], 12.0);  // d2f/dy2 = 6y

Per-op cost is O(n^2). For n > ~50, prefer seeded Jet2<T> with n passes.

Crate structure

src/
  real.rs           The unified active-scalar trait Real
  passive.rs        The passive-scalar bound Passive (f32, f64) — was scalar.rs
  prelude.rs        Real, Passive, AReal, Jet1, Jet2, Tape, TapeStorage
  forward/          Jet1, Jet1Vec, Jet2, Jet2Vec + Named wrappers
  reverse/          AReal, NamedAReal, NamedTape
  ops/              compute_jacobian_*, compute_hessian, compute_full_hessian
  math.rs           AD-aware transcendentals (sin, exp, erf, norm_cdf, ...)
  tape.rs           Reverse-mode tape and thread-local active-tape slot
  registry.rs       VarRegistry — ordered name-to-index map
  forward_tape.rs   NamedForwardTape / NamedForwardScope setup

Examples

cargo run --release --example swap_pricer

Design notes

• Tape storage is thread-local. One Tape<T> per thread; NamedTape is !Send.
• Forward mode is allocation-light. Jet1Vec keeps tangents in a single Vec<f64> with fused, autovectorizable loops.
• Zero-alloc operator fast paths. Every AReal binary op uses fixed-arity Tape::push_binary / push_unary — no intermediate Vec per op.

Tests

cargo test

165 tests covering operator correctness, transcendentals, second-order derivatives, Jacobian/Hessian helpers, named variable readback, and cross-mode consistency.

License

AGPL-3.0-or-later, matching the upstream XAD project. See LICENSE.md.

Acknowledgements

• The C++ XAD library — architectural inspiration and source of the financial examples.
• QuantLibAAD — XAD-instrumented QuantLib build; reference for the AAD-on-quant-finance patterns the financial examples in this crate are modelled after.
• num-traits for generic scalar plumbing.