eaqsp: error-aware quantum signal processing

eaqsp finds QSP phases for a real polynomial and searches for a sequence with lower phase-error sensitivity or a lower Clifford+T gate count. It starts from pyqsp's symmetric solution for the same polynomial and keeps that baseline available as a fallback. Search is local: improvement and a global optimum are not guaranteed.

Version 0.2.1 adds an independently validated installed-wheel release, clearer benchmark domains, guarded numerical failure handling, and reproducibility tooling. See CHANGELOG.md.

Install

python -m pip install eaqsp

Python 3.12 or later is required; this release is validated on Python 3.12.14. Python 3.10 and 3.11 are not claimed as supported. NumPy, SciPy, pyqsp, and pygridsynth are dependencies. pyqsp provides baseline phases; pygridsynth synthesizes Clifford+T gates. Their code is imported rather than vendored; see NOTICE. The pyqsp version is pinned to 0.2.0 because the baseline bridge uses its symmetric-solver API. The tested pygridsynth version is 2.0.0. Exact dependency versions are recorded in validation/current_release/environment.json and requirements-reproduce.txt.

For development, clone the repository and install an editable copy with the test dependencies:

git clone https://github.com/kansari123/eaqsp.git
cd eaqsp
python -m pip install -e '.[test]'
python -m pytest tests

If you have downloaded the release wheel, install it directly with python -m pip install /path/to/eaqsp-0.2.1-py3-none-any.whl.

Start with a polynomial

Coefficients are in ascending Chebyshev order, so [0, 0.5, 0, 0.1] means 0.5 T₁(x) + 0.1 T₃(x). The polynomial must be real and have definite parity. Nonconstant targets must stay strictly below |p(x)| = 1 on [-1, 1]; constants may include the endpoints ±1.

import eaqsp

coefficients = [0, 0.5, 0, 0.1]
phases, report = eaqsp.solve(
    coefficients,
    objective="phase_error",
    tolerance=1e-3,
    error_model={"kind": "iid", "mode": "expected"},
    seed=0,
)
print(report.to_markdown())
assert report["total_verification"]["ok"]

The response is Re U₀₀(x), where

U(x) = exp(i φ₀ Z) W(x) exp(i φ₁ Z) … W(x) exp(i φd Z)
W(x) = [[x, i sqrt(1-x²)], [i sqrt(1-x²), x]]

A degree-d sequence has d + 1 phases. For custom polynomials, a larger explicit degree of the same parity pads the coefficients and returns the requested number of phases. A smaller degree or incompatible parity is rejected. Without an explicit degree, trailing zero coefficients are removed. The zero polynomial is supported. Degree-zero constants use an analytic baseline because the ordinary pyqsp symmetric solver does not handle that case.

Choose an objective

phases, report = eaqsp.solve(
    coefficients, objective="tcount", tolerance=1e-2, seed=0
)
print(report["comparison"]["T_base"], report["comparison"]["T_ours"])

phases, report = eaqsp.solve(
    coefficients,
    objective="mixed",
    weights={"phase_error": 0.7, "tcount": 0.3},
    tolerance=1e-2,
    seed=0,
)

rows, target = eaqsp.tradeoff(
    coefficients, "phase_error", "tcount",
    tolerance=1e-2, n_points=5, seed=0,
)

A trade-off call selects from a finite pool of locally discovered candidates and returns one feasible selection per requested weight. Duplicate selections are allowed: different weights may prefer the same candidate. The rows are a sampled trade-off, not a certified Pareto frontier.

The T-count covers phase rotations only. Signal operations W(x), state preparation, controls, and measurement are outside that count. solve returns ideal phase angles; the compiled implementation is in report["synthesis"]["ours"]. Use the full synthesis record, including phase corrections, when reconstructing it:

import numpy as np
from eaqsp.objectives.tcount import compiled_response

x = np.linspace(-1, 1, 101)
y_compiled = compiled_response(report["synthesis"]["ours"], x)

Named targets and tolerance

# Choose the degree automatically from the approximation budget.
phases, report = eaqsp.solve(
    "inverse", kappa=2.0, objective="phase_error", tolerance=6e-2
)

# Or specify a degree that can meet that budget.
phases, report = eaqsp.solve(
    "sign", degree=21, delta=0.7,
    objective="phase_error", tolerance=5e-2
)
phases, report = eaqsp.solve(
    "cos", degree=12, tau=3.0,
    objective="phase_error", tolerance=1e-3
)

Named targets are scaled functions. sign approximates scale × sign(x) outside the excluded central interval. inverse approximates s/x on |x| ≥ 1/kappa; its scale s is reported. cos and sin approximate scale × cos(tau x) and scale × sin(tau x). Read report["target"] for the exact target parameters and approximation error.

tolerance is the total nominal error budget. For named targets, one third is reserved for function-to-polynomial approximation. The rest goes to phase recovery, and is split between phase recovery and gate synthesis when synthesis is included. Each applicable stage and the total must pass. An insufficient fixed degree raises eaqsp.VerificationError; it cannot produce a successful result merely because phases reproduce the wrong polynomial accurately.

Verification includes conservative analytic or coefficient-based bounds, sampled response checks, and matrix norm checks for gate synthesis. Bounds are evaluated in floating point; they are not interval-arithmetic certificates. The nominal tolerance does not guarantee accuracy under arbitrary noise.

Reading robustness results

All robustness costs and percentages are evaluated on the target's declared domain: |x| >= delta for sign, |x| >= 1/kappa for inverse, and [-1, 1] for custom and trigonometric targets. They use equally weighted grid points, including interval endpoints. A gain measured on this domain need not hold on the entire interval.

The default objective minimizes the first-order squared RMS coefficient under independent phase jitter. Noise tables report RMS itself, jointly over sample draws and input points. Ratios of these two quantities need not be equal.

The table also evaluates a common phase offset and bounded phase errors. Small-error rows use derivatives; larger-error rows evaluate the perturbed QSP product directly. A direct bounded-error search supplies a lower bound on the true worst case. Comparing two lower bounds cannot establish a robustness improvement, so such rows do not claim that the selected sequence helps.

error_model={"size": 0.03} requests a noise table at that single size. Pass sizes=(1e-4, 1e-3, 1e-2, 3e-2, 1e-1) explicitly to compare a range. Size changes reporting, not the minimizer of the homogeneous first-order objective.

Seeds determine the randomized calculations reproducibly in a fixed dependency environment. Runtime fields naturally vary. The report distinguishes the optimized first-order model from finite-noise measurements, and numerical evidence from worst-case guarantees.

Examples and limitations

python examples/case_custom.py
python examples/case_sign.py
python examples/case_inverse.py
python examples/case_hamiltonian.py
python examples/case_tradeoff.py
python examples/summarize_results.py

Examples write complete reports and concise summaries to examples/results/. Rerun evidence is in examples/results/combined_table.md. Those results describe their specific parameters; they do not promise improvement on all polynomials. See docs/walkthrough.md for the API and docs/plan.md for implementation methods and boundaries.

Geometry and limitations

Fiber implements the coefficient motion dψ/dt = B J(ψ).T α, numerical projection, and constrained local descent. The application uses these directions and projected starts to search the space of implementations. It does not implement holonomy steering, enumerate topological components, or prove access to a global optimum. Fiber.tangent returns the numerical Jacobian kernel; at a singular point those vectors need not integrate to actual fiber curves, so proposed steps require retraction and target verification.

The ideal phase-error search aims at the same polynomial to numerical precision. Gate-cost search may use a nonzero recovery budget, so its returned polynomial can differ slightly within the declared total tolerance. Reported gate savings therefore apply to verified approximate implementations with a common target and tolerance. This distinction is explicit in each report's recovery bound.

Runtime depends on degree, target, search settings, compilation precision and processor. The fresh degree-21 sign example used 12.6 seconds for phase-error search and 10.7 seconds for T-count search with one BLAS thread; the degree-25 inverse used 16.1 and 17.9 seconds. These are example timings, not a scaling guarantee: Jacobian-kernel factorizations can require cubic work, while finite-difference objectives, ODE steps and gate synthesis add variable costs. Full settings and timings accompany the reports and validation/current_release/.

Reproduce and cite

python -m pip install -r requirements-reproduce.txt
python -m pip install -e . --no-deps
python -m pytest tests -q
python validation/independent/verify_release.py --output validation/current_release/independent_results.json

The independent script imports the installed package, uses its own QSP product and elementary H/S/T/X/W matrices, and records the import path. The CI workflow builds and installs the wheel before running the tests. Validation covers floating-point software results; it includes no QPU experiment or interval-arithmetic proof.

For the geometry and coefficient-flow method, cite Kamran Ansari, Geometry, topology, and navigation of quantum signal processing phase fibers (accompanying manuscript, 2026). Kamran Ansari, Selecting quantum signal processing implementations for robustness and compilation cost (accompanying manuscript, 2026), describes the constrained search, verification contract and numerical experiments. These companion manuscripts are not included in the package; no publication venue or DOI is implied. Package citation metadata is in CITATION.cff.

Please also cite the software and methods used by its dependencies:

• pyqsp, version 0.2.0; Y. Dong, X. Meng, K. B. Whaley and L. Lin, Efficient phase-factor evaluation in quantum signal processing, Physical Review A 103, 042419 (2021), doi:10.1103/PhysRevA.103.042419. Further references are in pyqsp's CITATION file.
• pygridsynth, version 2.0.0; N. J. Ross and P. Selinger, Optimal ancilla-free Clifford+T approximation of z-rotations, Quantum Information and Computation 16, 901–953 (2016), arXiv:1403.2975.