Comparison with other packages

Several packages, across several languages, do symbolic second-quantized operator algebra. This page benchmarks SecondQuantizedAlgebra (SQA) head-to-head against four of them on a compact suite of practical tasks: building Hamiltonians, deriving Heisenberg equations of motion, nested commutators, operator powers, normal ordering, multi-mode chains, symbolic indexed sums, and mean-field expectation values.

The suite definition, canonical benchmark keys, and the runner scripts live in the repository under benchmark/comparison/; see Reproducing below.

Methodology

Further caveats, each disclosed rather than silently absorbed:

• Cross-language ratios are order-of-magnitude indicators. Every harness gives each scenario the same wall-clock budget (up to ~10 s) in a single process and reports the minimum per-call time (Julia via BenchmarkTools, Python via pyperf, Mathematica via a calibrated timing loop); fast scenarios stop early once ~50 samples are collected, and only scenarios slower than a few hundred ms per call use the full budget. Harness overhead and warm-up behavior still differ across runtimes, so read 90× as "roughly two orders of magnitude", not as a precise figure. The Julia-vs-Julia column is precise.
• OpenFermion uses numeric coefficients (ω = 1.0, J = 0.5, U = 0.25) because BosonOperator does not accept symbols. That is less work than the symbolic-coefficient arithmetic every other package performs, so OpenFermion's numbers are, if anything, flattered. It also cannot mix bosons with qubit operators in one expression, so the Jaynes–Cummings rows are not expressible.
• SymPy needs a small disclosed helper. normal_ordered_form does not know that boson and Pauli operators act on different Hilbert spaces (and hence commute), so mixed products fail to merge out of the box. The benchmark script includes a ~15-line helper that sorts boson factors before Pauli factors in each product; its cost is included in SymPy's timed path.
• sneg encodes the two-level atom as a single fermionic mode. For one site the fermion algebra ($f f† + f†f = 1$, $f² = 0$) is isomorphic to the σ algebra and bosons commute with fermions, so all Jaynes–Cummings results are identical to a spin encoding; this is also sneg's idiomatic usage in quantum impurity physics.
• Different canonical bases. QuantumAlgebra expands $σ⁺σ⁻$ into $σˣ/σʸ/σᶻ$ form and SymPy works in its Pauli basis, whereas SQA keeps the transition form. Results are physically equivalent but not term-by-term identical.

Results

Single Linux workstation, 2026-07-21: SQA 0.10.0 and QuantumAlgebra 1.6.0 on Julia 1.12.6, SymPy 1.14.0 and OpenFermion 1.8.1 on Python 3.14.6, sneg 2.0.21 on Mathematica 14.2. The SQA column is absolute time; every other cell shows that package's time with its ratio to SQA in parentheses, so (9.1×) means 9.1× slower than SQA. n/a means the scenario is not expressible in that package.

Reading the results

• SQA leads every operator-algebra row, against every package. The closest competitor throughout is QuantumAlgebra, which shares SQA's Julia substrate and a specialized bosonic normal-ordering engine. The general-purpose CAS approaches (SymPy, Mathematica-based sneg) sit well behind on the same physics.
• The lead widens exactly where real derivations spend their time. Nested commutators (the core of Heisenberg / Schrieffer–Wolff / cumulant derivations) are where the gap grows fastest with depth, as the top-left panel shows. Eager canonicalization keeps every intermediate compact, so the cost compounds more slowly.
• Operator powers separate width from depth. The top-right panel sweeps $Hⁿ$ (an eager fold, same contract as the table): a single wide product-and-collect rather than a deep recursion, and the place where general-CAS normal ordering degrades fastest while SQA and QuantumAlgebra stay near log-linear. It is the product-expansion counterpart of the nested-commutator panel and stresses like-term collection instead of cancellation.
• Chain length probes system size, not expression size. The bottom-right panel builds the Bose–Hubbard Hamiltonian on M = 2, 4, 8, 16, 32 modes; unlike the other three sweeps it grows the number of Hilbert-space factors rather than the depth or width of a single expression. The Hamiltonian has $O(M)$ terms, so every package scales near-linearly and the panel becomes a constant-factor race where the field bunches closest together. Building a Hamiltonian is the one axis with the least canonicalization work per term.
• OpenFermion is the strongest non-Julia contender on pure bosonic work, while computing with plain numeric coefficients, a strictly easier task than the symbolic arithmetic every other column performs.
• The one row SQA loses is measuring different objects. QuantumAlgebra's expval tags operators inside its own term type (a near-trivial field move), whereas SQA's average materializes a Symbolics.jl Num expression (splitting coefficients into real and imaginary parts and building a CAS term tree). That object is what feeds substitute, simplify, numeric conversion, and ModelingToolkit downstream, so the extra cost buys the bridge into the Julia symbolics ecosystem.
• Symbolic indexed sums are a two-package race. Only SQA and QuantumAlgebra can express $Σᵢ gᵢ(a†σ⁻ᵢ + aσ⁺ᵢ)$ with a symbolic atom number N, and the $[H, σ⁺ⱼσ⁻ⱼ]$ diagonal split is the operation that turns such sums into mesoscopic equations of motion.

Ergonomics: eager by default vs the lazy-product footgun

The contract above uses the eager workflow for lazy packages. Their default lazy workflow has a sharp footgun on powers, building the full product before canonicalizing: with QuantumAlgebra, normal_form(H⁴) for the Jaynes–Cummings H costs about 86 ms and 2.2 M allocations lazily versus 2.3 ms and 84 k allocations eagerly, a roughly 37× time difference from one workflow choice. SymPy's Bose–Hubbard H² row above is the same effect in a different package. SQA canonicalizes on every *, so this cliff does not exist for its users.

Reproducing

Each language has a runner script in benchmark/comparison/ that validates the known-answer identities, runs its scenarios, and writes one results/<package>.json per package (the Julia script covers both SQA and QuantumAlgebra in one run); make_table.jl merges whatever JSONs exist into the table and figure above. The scenario definitions and canonical keys live in benchmark/comparison/BENCHMARKS.md.

# Julia (SQA + QuantumAlgebra)
julia --project=benchmark benchmark/comparison/julia_bench.jl

# Python (each needs: pip install sympy openfermion pyperf)
python benchmark/comparison/sympy_bench.py
python benchmark/comparison/openfermion_bench.py

# Mathematica (needs sneg installed)
wolframscript -f benchmark/comparison/sneg_bench.wls

# Merge into the docs table + figure
julia --project=benchmark benchmark/comparison/make_table.jl