qarp.algorithms

Algorithms. Public depth: flat — every primitive and composite algorithm is qarp.algorithms.<Name>, together with the eigenspectrum helpers.

Submodules are private implementation and may be reorganised without notice.

class qarp.algorithms.AdaptVQD(reference_block: Block, hamiltonian: QubitOperator, excitation_pool: List[QubitOperator], orthogonal_states: List[Block], betas: List[float], gradient: bool | str = False, optimizer: Optimizer | None = None, diminishing: bool = False, exc_per_iter: int = 1, verbose: bool = True, term_gradient: str = 'ADAPT-VQD', convergence_thresh=1e-07, gradient_thresh=1e-07, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: AdaptAnsatzMixin, CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block()

Use the current object reference, ansatz excitations and symbols to construct a wavefunction object.

Returns:

A wavefunction object constructed with the current AdaptVQE attributes.

iterate()

property optimal_parameters: dict

Optimized parameters keyed by symbol — the order-proof surface (alias of ansatz_parameters_dict, populated by run()).

pool_scan()

run(*, max_iter: int = 25) → Tuple[float, List[float] | None]

Run the AdaptVQD algorithm for a maximum of max_iter iterations.

Returns:

The final energy and final parameters as a float and list of floats.

class qarp.algorithms.AdaptVQE(reference_block: Block, system_hamiltonian: QubitOperator, excitation_pool: List[QubitOperator], optimizer: Optimizer | None = None, gradient: bool | str = False, verbose: bool = True, diminishing: bool = False, exc_per_iter: int = 1, qubit_adapt: bool = False, gradient_thresh: float = 1e-05, convergence_thresh: float = 1e-06, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: AdaptAnsatzMixin, CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block()

Use the current object reference, ansatz excitations and symbols to construct a wavefunction object.

Returns:

A wavefunction object constructed with the current AdaptVQE attributes.

iterate()

property optimal_parameters: dict

Optimized parameters keyed by symbol — the order-proof surface (alias of ansatz_parameters_dict).

pool_scan()

run(*, max_iter: int = 25) → Tuple[float, List[float] | None]

Run the AdaptVQE algorithm for a maximum of max_iter iterations.

Returns:

The final energy and final parameters as a float and list of floats.

class qarp.algorithms.AmplitudeAmplification(state_preparation: Block, oracle: Block, n_iterations: int, *, good_states: Sequence[int] | None = None, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

Prepare a state and apply a fixed number of amplification iterates.

The built circuit is A followed by n_iterations applications of AmplitudeAmplificationBlock. The supplied oracle must implement exactly I - 2 Pi_good; see the block documentation for the phase convention.

good_states is optional classical reporting metadata. It contains full-register integer labels and is used only to sum the returned sampling probabilities. It neither defines nor modifies the quantum oracle.

The iterate and the probability law sin((2 * n_iterations + 1) * theta)**2 follow Brassard, Hoyer, Mosca, and Tapp, Quantum Amplitude Amplification and Estimation, arXiv:quant-ph/0005055, Eqs. (1), (5), and (8).

Parameters:

• state_preparation – Unitary A preparing the initial state.

• oracle – Good-state phase oracle on the same positive-width register.

• n_iterations – Explicit non-negative number of amplification iterates.

• good_states – Optional unique full-register integer labels whose returned probabilities are summed into success_probability.

• primitive – Sampling primitive. Defaults to a private Sampler. A sampler carrying initial_state is rejected: the amplification law assumes the circuit starts in |0...0>, so A must act on that state and not on a seeded one.

• engine – Execution engine. Defaults to QarpEngine through the composite-algorithm base class.

build() → Self

Build and compile A followed by the requested iterates.

run() → dict[tuple[int, ...], float]

Execute once and return the LSB-first sampling distribution.

class qarp.algorithms.AmplitudeEstimation(state_preparation: Block, oracle: Block, n_ancilla: int, *, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

Estimate the initial good-state probability with canonical QAE.

The supplied oracle must implement exactly I - 2 Pi_good. With M = 2**n_ancilla, a sampled QPE label y maps to sin(pi*y/M)**2. The conjugate phase branches are folded together at z = min(y, M-y) before selecting the most probable estimate.

This is the canonical estimator of Brassard, Hoyer, Mosca, and Tapp, arXiv:quant-ph/0005055.

Parameters:

• state_preparation – Unitary preparing the state whose amplitude is estimated.

• oracle – Phase-exact good-state oracle on the same register.

• n_ancilla – Positive number of estimation qubits.

• primitive – Sampling primitive. Defaults to a private Sampler.

• engine – Execution engine.

build() → Self

Build and compile the canonical QAE circuit.

run() → float

Execute QAE once and return the modal folded-bin estimate.

class qarp.algorithms.BasisRotationAveraging(ket: Block | None = None, integrals: tuple | None = None, constant: float = 0.0, tolerance: float = 1e-08, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

gradient_kind: str = 'expectation'

property n_groups: int

run(results: list) → float

constant + Σ_ℓ c_ℓ · (const_ℓ + Σ_masks coeff · ⟨Z-parity⟩).

In group ℓ’s rotated basis every mask expectation is (1/n_shots) Σ_outcome counts[outcome] · (−1)^popcount(outcome & mask).

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE})

class qarp.algorithms.CompositeAlgorithm(primitive: PrimitiveAlgorithm, engine: Engine | None = None)

Bases: ABC

abstractmethod build() → Self

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

abstractmethod run() → Any

Execute the algorithm after build().

The return type is algorithm-specific (documented per class); subclasses may add optional keyword-only arguments such as max_iter but take no positional arguments.

class qarp.algorithms.CuttingPrimitive(ket: Block | None = None, operator: QubitOperator | None = None, max_subcircuit_qubits: int | None = None, n_shots: int = 5000, force_max_number_cuts: bool = False, experiment_fraction: float | None = None, rng_seed: int | None = None)

Bases: PrimitiveAlgorithm

PrimitiveAlgorithm that evaluates expectation values via QPD circuit cutting.

Works with any OpenQARP engine without modification.

build() → Self

Partition the circuit and build all measurement sub-blocks.

Uses streaming experiment generation and QWC grouping so that sub_blocks contains n_experiments × n_subcircuits × n_qwc_groups blocks (instead of × n_obs_terms).

gradient_kind: str = 'expectation'

run(results: list) → float

Reconstruct the expectation value from engine sampling results.

Parameters:

results – list[qx.SamplingResult] in the same order as sub_blocks.

Returns:

Reconstructed expectation value (float).

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE})

supports_exact: bool = False

class qarp.algorithms.DOSQPE(unitary: Block, n_ancilla: int, hamming_weight: int | None = None, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

build()

Build the DOS Phase Estimation circuit.

First offers the problem to the engine’s structured fast path (Engine.prepare_structured_qpe — matrix exponentiation, the controlled-U ladder is never compiled). Engines without the path, or refusing it (EXACT readout, noise, routing, parametric U — see QarpEngine.prepare_structured_qpe), return None and the full DOSQPEBlock circuit is built instead.

Returns:

The instance of the class.

Return type:

self

property freqs

Get the frequency grid for the ancilla register.

plot(figsize=None, return_fig=False, max_xticks=11)

Plot the results of the DOS Phase Estimation algorithm.

Parameters:

• figsize (tuple) – Size of the figure. If None, scales with n_ancilla.

• return_fig (bool) – If True, return (fig, ax) for external saving/customization.

• max_xticks (int) – Upper bound on the number of x-axis ticks.

Returns:

(fig, ax) if return_fig is True, otherwise None.

plot_against_spectrum(unique_eigs, normalized_degeneracy, unique_occ_numbers, figsize=None, return_fig=False, max_xticks=11)

Plot the results of the DOS Phase Estimation algorithm against the spectrum.

Parameters:

• unique_eigs (list) – Unique eigenvalues.

• normalized_degeneracy (list) – Normalized degeneracy.

• unique_occ_numbers (list) – Unique occupation numbers.

• figsize (tuple) – Size of the figure. If None, scales with n_ancilla.

• return_fig (bool) – If True, return (fig, ax) for external saving/customization.

• max_xticks (int) – Upper bound on the number of x-axis ticks.

Returns:

(fig, ax) if return_fig is True, otherwise None.

run()

Run the DOS Phase Estimation algorithm.

Returns:

The distribution of the measurement results.

Return type:

distribution

class qarp.algorithms.GroupedTransitionHadamardTest(bra: Block | None = None, operator: QubitOperator | None = None, ket: Block | None = None, real: bool = True, imaginary: bool = True, n_shots: int | Shots | None = None, grouping: GroupingStrategy | None = None)

Bases: PrimitiveAlgorithm

Estimate a full-operator transition amplitude by QWC grouping.

The circuit for a group samples the correlation between the ancilla quadrature and every Pauli parity in that group. The primitive returns <bra|H|ket> directly; it does not expose or run one sub-primitive per Pauli term.

Parameters:

• bra – State-preparation block for <bra|.

• operator – Real-coefficient QubitOperator to measure.

• ket – State-preparation block for |ket>.

• real – Include real-quadrature circuits.

• imaginary – Include imaginary-quadrature circuits.

• n_shots – Shots per group/quadrature; None defers to the engine.

• grouping – Term grouping strategy. The initial circuit primitive supports strategies with qubit_wise=True; default is QubitWiseCommuting.

General (non-QWC) commuting groups require an entangling diagonalisation and are intentionally rejected until a transition-aware Clifford path is added. This keeps the first implementation exact and easy to validate.

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property expectation_type: str

gradient_kind: str = 'expectation'

property n_groups: int

property n_terms: int

run(results: list) → float | complex

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.TRANSITION_AMPLITUDE})

class qarp.algorithms.Grover(oracle: Block, n_marked: int = 1, *, good_states: Sequence[int] | None = None, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

Search a uniform basis-state space with a known marked-state count.

good_states is optional reporting metadata. It is never used to build or modify the oracle; when supplied, its returned probability mass is exposed as success_probability. It is not cross-checked against the oracle (an opaque unitary cannot be inspected without an exponential dense matrix), so declaring labels the oracle does not actually mark yields a misleading success_probability — the labels are the caller’s promise.

most_likely_states is only meaningful once amplification has concentrated the distribution. When the optimal iteration count is zero (n_marked >= N/2, so the uniform state is already at or past the amplification optimum) the distribution stays (near-)uniform and most_likely_states is the entire register; read predicted_success_probability to see that no concentration occurred.

The algorithm follows Grover, arXiv:quant-ph/9605043, with the known-count analysis of Boyer, Brassard, Hoyer, and Tapp, arXiv:quant-ph/9605034.

Parameters:

• oracle – Phase oracle implementing exactly I - 2 Pi_good.

• n_marked – Number of basis states marked by the oracle.

• good_states – Optional explicit marked integer labels for reporting only.

• primitive – Sampling primitive. Defaults to a private Sampler.

• engine – Execution engine.

build() → Self

Build and compile uniform preparation plus amplification.

run() → dict[tuple[int, ...], float]

Execute once and return the full LSB-first search distribution.

class qarp.algorithms.HadamardTest(bra: Block | None = None, operator: Block | None = None, ket: Block | None = None, real: bool = True, imaginary: bool = True, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property expectation_type: str

gradient_kind: str = 'expectation'

run(results: list) → complex

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE, Target.OVERLAP, Target.TRANSITION_AMPLITUDE})

class qarp.algorithms.InterferometricTest(bra: ComputationalBasisStateBlock | None = None, operator: Block | None = None, ket: ComputationalBasisStateBlock | None = None, real: bool = True, imaginary: bool = True, sampling_algorithm: PrimitiveAlgorithm | None = None, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property expectation_type: str

gradient_kind: str = 'expectation'

run(results: list) → complex

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE})

class qarp.algorithms.MMQCELS(operator: QubitOperator | ndarray | TrotterBlock, state: ndarray | Block, *, execution_mode: Literal['classical', 'statevector', 'hadamard'], T0: float, parameter_mode: Literal['standard', 'error_rate'] = 'standard', N0: int | None = None, Nj: int | Sequence[int] | None = None, n_dominant_eigenvalues: int = 1, error_rate: float = 0.001, n_levels: int | None = None, q: float | None = None, gamma: float = 1.0, initial_eigenvalues: Sequence[float] | None = None, n_initial_guesses: int = 10, n_shots: int | Shots | None = None, lam_min: float = -3.141592653589793, lam_max: float = 3.141592653589793, optimizer: Optimizer | None = None, seed: int | None = None, verbose: bool = True, engine: Engine | None = None)

Bases: CompositeAlgorithm

Estimate several dominant eigenvalues with MM-QCELS.

The implementation follows Algorithm 2 and Eq. (13) of Ding and Lin, Quantum 7, 1136 (2023): times are independent truncated-Gaussian draws, observations are complex, and the time scale doubles between levels.

parameter_mode="standard" is the recommended default. It requires T0, N0, Nj, and exactly one of q or n_levels. parameter_mode="error_rate" retains OpenQARP’s historical epsilon-derived level and sample-count helpers as an explicitly opt-in convenience. It is not a theorem-backed MM-QCELS parameter calibration.

Parameters:

• operator – Hamiltonian in classical mode, or a symbolic-time TrotterBlock in circuit modes.

• state – Statevector/block in classical mode, or a built block in circuit modes.

• execution_mode – Required data-generation strategy: "classical", "statevector", or "hadamard".

• T0 – Positive initial Gaussian time scale.

• parameter_mode – "standard" or the opt-in "error_rate" compatibility preset.

• N0 – Number of samples at level zero. Required in standard mode.

• Nj – Positive sample count reused after level zero, or one count per subsequent level. Required in standard mode, even though a scalar value is unused when n_levels=1.

• n_dominant_eigenvalues – Number of dominant modes K.

• error_rate – Accuracy epsilon used by parameter_mode="error_rate" and by standard mode’s q schedule. It does not affect a standard schedule with an explicit n_levels.

• n_levels – Explicit number of evaluated levels. Mutually exclusive with q in standard mode.

• q – Theorem-1 schedule parameter. Standard mode derives l=max(ceil(log2(q/(epsilon*T0))), 1) and evaluates l+1 levels.

• gamma – Truncation radius in standard deviations. Sampled times lie in [-gamma*T_j, gamma*T_j].

• initial_eigenvalues – Optional length-K first-level phase guess.

• n_initial_guesses – Total number of first-level optimization starts, including initial_eigenvalues when supplied.

• n_shots – Hadamard shots per quadrature. None means one shot, as in the paper; qarp.EXACT evaluates the exact protocol.

• lam_min – Initial lower eigenvalue bound.

• lam_max – Initial upper eigenvalue bound.

• optimizer – SciPy-compatible bounded optimizer.

• seed – Seed for time sampling and generated optimizer starts. When no engine is supplied, it also seeds the default QarpEngine.

• verbose – Print per-level progress.

• engine – Optional configured execution engine.

Variables:

• result – Eigenvalues returned by the latest successful run().

• eigenvalues – Alias of result, following the QPE output style.

• amplitudes – Complex amplitudes fitted alongside the eigenvalues.

• level_losses – Final complex least-squares loss at every level.

• sample_counts – Number of observations used at every level.

• sampled_max_times – Largest absolute sampled time at every level.

• sampled_total_times – Sum of absolute sampled times at every level.

• optimizer_failed_starts – Number of discarded optimizer starts at every level. A level for which every start fails raises instead.

Every evaluated level must contain more observations than dominant modes; otherwise the variable-projection residual cannot identify the frequencies.

build() → MMQCELS

Build and cache the selected dataset evaluator.

generate_dataset(T: float, n: int | None = None) → list[tuple[float, complex]]

Generate one independent MM-QCELS dataset at scale T.

generate_times(T: float, n: int) → NDArray[float64]

Draw n independent samples from Eq. (3)’s distribution.

objective(eigenvalues: ArrayLike, dataset: Sequence[tuple[float, complex]]) → float

Evaluate the full complex Eq. (13) loss after eliminating amplitudes.

optimal_dataset_length() → int

Return OpenQARP’s historical epsilon-only sample-count heuristic.

optimal_number_of_iterations() → int

Return the predecessor-QCELS epsilon-only level heuristic.

run() → NDArray[float64]

Run all levels and return the sorted eigenvalues.

The returned array is also stored in result and eigenvalues. Fitted amplitudes and per-level diagnostics are available on the corresponding instance attributes.

time_scale(level: int) → float

Return the MM-QCELS scale T_j = 2**j * T0.

class qarp.algorithms.MirrorTest(bra: Block | None = None, operator: Block | None = None, ket: Block | None = None, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

gradient_kind: str = 'expectation'

returns_probability: bool = True

run(results: list) → float

P(all zeros) = the squared overlap estimate.

The single-circuit engine pipeline returns one qx.SamplingResult per sub_blocks entry; we read its counts and divide by n_shots. All-zeros corresponds to outcome integer 0 regardless of register size.

supported_targets: frozenset[Target] = frozenset({Target.OVERLAP})

class qarp.algorithms.MonteCarlo(hamiltonian: QubitOperator | ndarray, approx_ground_state_energy: float | List[float], total_time: float, time_step: float, reference_walker_label: int | List[int], unitary_block: Block, num_target_states: int = 1, initial_walker_count: int | List[int] = 200, walker_basis: List[WalkerState] | None = None, shift_damping: float | List[float] = 0.1, population_threshold: int = 350, mode: str = 'Semiclassical', primitive: PrimitiveAlgorithm | None = None, n_shots: int = 1000, num_trajectories: int = 10, ham_cache: Dict[Tuple[int, int], Tuple[complex, float]] | None = None, save_walker_history: bool = False, history_save_interval: int = 1, qdrift: bool = False, qdrift_samples: int = 100, qdrift_ratio: float | None = None, verbose: bool = False, engine: Engine | None = None, seed: int | None = None)

Bases: CompositeAlgorithm

build()

Build the simulator by setting up walker states and caches.

Returns:

Self for method chaining

estimate_ground_state_energy(walkers_by_label: Dict[str, List[WalkerState]], initial_energy: float, visited_states: Set[int]) → float

Calculate energy using pre-grouped walkers (semiclassical mode).

estimate_ground_state_energy_quantum(walkers: List[WalkerState], initial_energy: float, visited_states: Set[int]) → float

Calculate ground state energy from walker list (quantum mode).

iterate() → None

Run single Monte Carlo simulation iteration.

prepare_vector_and_projector(walkers) → Tuple[ndarray, ndarray]

From list of walkers, return vector in walker basis and projector of that vector

project_for_orthogonalization(current_excitation, projectors) → Tuple[List[WalkerState], float]

Apply projectors onto current ES walker vector to ensure orthogonalisation with lower energy states

remove_opposite_sign_pairs(walkers: List[WalkerState]) → List[WalkerState]

Remove walkers with opposite signs in the same state (annihilation).

Parameters:

walkers – List of walker states

Returns:

Walkers after annihilation

run() → List[float]

Run multiple Monte Carlo trajectories.

Returns:

Final energy estimate from last trajectory

class qarp.algorithms.PCE(graph: Graph, order: int, ket: Block, merging: bool = False, classical_function: Callable = <function classical_function_max_cut>, quantum_function: Callable = <function quantum_function_max_cut>, initial_parameters: Iterable[float] | None = None, optimizer: Optimizer | None = None, verbose: bool = True, gradient: bool | str = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

build()

Build the PCE object with UCC ansatz

If not all attributes are set, some defaults are loaded. These are: - COBYLA as optimizer - No gradients used - StateVector for all measurement strategies - QarpEngine() for all backends

Returns:

self

get_final_state_block() → Block

The ket bound at the optimal parameters (set by run()).

get_gradient_function()

Builds a gradient function as an option for certain optimizers to minimize

Returns:

The gradient function assisting the optimizer

get_objective_function()

Builds the objective function needed for the optimizer to minimize

Returns:

The objective function to be minimized

property optimal_parameters: Dict

Optimized parameters keyed by symbol — the order-proof surface.

run()

Run the PCE algorithm by calling the optimizer with the built objective function and optionally gradients.

Returns:

The final energy, final parameters and solution to graph problem as a float, list of floats and list of node names.

class qarp.algorithms.PauliAveraging(bra: Block | None = None, operator: QubitOperator | None = None, ket: Block | None = None, n_shots: int | Shots | None = None, grouping: GroupingStrategy | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

gradient_kind: str = 'expectation'

property n_groups: int

property n_terms: int

run(results: list) → float

Per-group expectation, weighted sum, plus the constant term.

After the group’s diagonalising Clifford, term i equals s_i · Z^{mask_i} (s_i = ±1), so from the measurement counts:

E[term_i] = s_i · (1/n_shots) Σ_outcome counts[outcome] · (-1)^popcount(outcome & mask_i)

⟨H⟩ = constant + Σ_groups Σ_terms c_i · E[term_i].

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE})

class qarp.algorithms.PauliKernel(n_qubits: int)

Bases: ShadowKernel

Random-Pauli (local-Clifford) inverse channel.

A setting is one axis per qubit (0/1/2 = X/Y/Z). The inverse channel factorizes per qubit: for a Pauli term P a snapshot contributes prod_{q in supp(P)} 3 * (-1)^{b_q} iff the measured axis on every qubit of supp(P) matches P’s axis there, and 0 otherwise.

capabilities: frozenset[str] = frozenset({'expval'})

ensemble: str = 'pauli'

snapshot_estimate(setting: ndarray, outcome: int, term) → float

Single-snapshot estimate of one Pauli term under one (setting, outcome). term is a tuple of (qubit, axis_char) pairs (empty for identity — never passed here; the estimator handles the constant).

class qarp.algorithms.PauliShadow(operator: str | QubitOperator, ket: Block | None = None, *, dataset: ShadowDataset | None = None, n_settings: int = 1000, n_shots: int = 1, seed: int | None = None)

Bases: ShadowProtocol

Random-Pauli classical shadows (Huang–Kueng–Preskill, 2020).

class qarp.algorithms.PrimitiveAlgorithm(ket: Any = None, bra: Any = None, operator: Any = None, n_shots: int | Shots | None = None, target: Target = Target.SAMPLING)

Bases: ABC

Pure Python base class for qarpx-backed primitive algorithms.

Uses the C++ kernel (qarpx) for compilation and simulation. All algorithm logic — circuit construction and result post-processing — lives in Python, keeping the boundary simple.

Subclasses must implement:

• build() — populate self.sub_blocks with qx.Block instances

• run(results) — post-process the engine’s per-sub-block results (qx.SamplingResult, or ExactResult under n_shots=qarp.EXACT) → scalar

Variables:

• sub_blocks – Populated by build(); list of qx.Block to compile.

• compiled_circuits – Set by QarpEngine.build(); list of command lists.

• _n_qubits_list – Set by QarpEngine.build(); n_qubits per circuit.

Class capability metadata (override on the subclass; class-level so it survives the deepcopy-then-mutate pattern QSE/MonteCarlo/QMEGS use):

• supported_targets: The Target values this primitive can estimate.

• consumes: Consumes.COUNTS — estimator reads sampled measurement statistics (the default); Consumes.AMPLITUDES — estimator contracts simulator statevectors directly, and engines dispatch to run_from_amplitudes(compiled_circuits) instead of sampling.

• supports_exact: True iff the protocol has a meaningful ∞-shot limit. Consulted only for COUNTS consumers (an AMPLITUDES primitive is always exact); False for inherently stochastic protocols (CuttingPrimitive).

• supports_backprop_gradient: True iff the primitive is eligible for the engine’s adjoint backprop gradient path (Engine.run_gradient). Default False — sampling primitives use the parameter-shift fallback instead.

• gradient_kind: what a shift rule may assume about run() as a function of each compiled circuit’s state, separately: "expectation" (bilinear — every sampled estimator), "amplitude" (linear or antilinear in one circuit’s amplitudes), "squared_overlap" (run() returns ⟨bra|ket⟩, the gradient is of |·|²), or "none" (non-linear, e.g. a median-of-means or a ratio: parameter shift refuses, finite differences still apply). Default "none" — a new primitive loses speed, never correctness, until it declares its class.

• requires_noiseless: True iff the estimator’s post-processing is only valid under ideal (noiseless) measurement. Enforced in Engine._validate_primitive: a primitive with this flag is rejected on an engine whose noise model is active. Default False (noise is fine). PauliShadow sets it — the shadow inverse channel assumes ideal readout, so a noisy campaign fed through it is silently biased.

• returns_probability: True iff an OVERLAP estimator’s run() already returns |⟨bra|ket⟩|² rather than the amplitude ⟨bra|ket⟩. Consumers that need the squared overlap (VQD/ADAPT-VQD deflation) must not square such a value again. Default False (amplitude-returning).

accepts_initial_state: bool = False

property bra: Any

abstractmethod build() → PrimitiveAlgorithm

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

compiled_circuits: list[list]

consumes: Consumes = 0

gradient_kind: str = 'none'

infer_target() → Target

Detect which quantity to compute from the {ket, bra, operator} inputs.

Single source of truth for input-driven target detection, shared by the primitives that vary their operation with their inputs (HadamardTest, StateVector). Sets and returns self.target, and defaults bra to ket for the expectation-value case. Fixed-purpose primitives (Sampler, SWAPTest, MirrorTest, …) pass an explicit target to __init__ and never call this.

Rules:

• ket only → SAMPLING

• ket + operator → EXPECTATION_VALUE (bra ≔ ket)

• bra ≠ ket + operator → TRANSITION_AMPLITUDE

• bra ≠ ket, no operator → OVERLAP

The bra ≔ ket default is live: while bra has never been explicitly assigned, every call re-derives it from the current {ket, operator}, so rebinding only ket (e.g. on a deepcopy) keeps an EXPECTATION_VALUE primitive an expectation value. Any explicit assignment takes ownership; bra = None hands it back.

Raises:

RuntimeError – if ket is missing, or bra was explicitly set to the same object as ket with no operator (a degenerate ⟨ψ|ψ⟩ request).

initial_state: Any = None

requires_noiseless: bool = False

returns_probability: bool = False

abstractmethod run(results: list) → float | complex | dict[tuple[int, ...], float]

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

run_from_amplitudes(compiled_circuits: list, simulator=None) → float | complex

Evaluate the target directly from simulator amplitudes.

Only called by engines when self.consumes is Consumes.AMPLITUDES. Default raises — sampling primitives don’t override this.

Parameters:

• compiled_circuits – One command stream per sub_blocks entry, already parameter-substituted.

• simulator – Optional statevector backend (anything implementing statevector(commands, n_qubits)). Engines pass their own simulator so amplitude primitives run on the same backend — e.g. CudaqEngine passes its GPU qx.CudaqSimulator. Defaults to a CPU qx.QarpSimulator when None.

sub_blocks: list

supported_targets: frozenset[Target] = frozenset({})

supports_backprop_gradient: bool = False

supports_exact: bool = True

class qarp.algorithms.ProjectedVQE(operator: QubitOperator, ket: Block, projector: Block | Sequence[Block] | None = None, projector_matrix: ndarray | None = None, initial_parameters: NDArray[float64] | None = None, optimizer: Optimizer | None = None, verbose: bool = False, engine: Engine | None = None, postselection_tol: float = 1e-14, save_energy_history: bool = False)

Bases: VQE

Variation-after-projection VQE using a fast StateVector objective.

Parameters:

• operator – Hamiltonian as a qarp QubitOperator on the system register, with no ancilla shift.

• ket – Parameterized ansatz block on the system register.

• projector – One symmetry-projector block, or a sequence of them (applied as a product). Known projector blocks are converted to dense system-register projector matrices and are not simulated as circuits during optimization.

• projector_matrix – Dense projector matrix in qarpx LSB ordering, for custom projectors. Mutually exclusive with projector.

• postselection_tol – Squared-norm threshold below which the projected state counts as annihilated; the objective returns +inf there.

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

run()

Run the VQA algorithm by calling the optimizer with the built objective function and optionally gradients.

Returns:

The final energy and final parameters as a float and list of floats.

class qarp.algorithms.QAOA(problem: Graph | QubitOperator, n_layers: int = 1, use_rzz: bool = True, initial_parameters: NDArray[float64] | None = None, gradient: bool | str = False, optimizer: Optimizer | None = None, verbose: bool = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None, save_energy_history: bool = False)

Bases: VQA

class qarp.algorithms.QITE(hamiltonian: QubitOperator, initial_block: Block, dtau: float, n_steps: int, *, pool: List[QubitOperator] | None = None, regularization: float = 1e-06, trotter_steps: int = 1, trotter_order: int = 2, engine: Engine | None = None, verbose: bool = False)

Bases: CompositeAlgorithm

build() → QITE

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block() → Block

The accumulated initial_block ∘ layer₁ ∘ … ∘ layerₙ (after run).

run() → Tuple[float, ndarray]

Sweep n_steps imaginary-time steps; return (final energy, ψ).

class qarp.algorithms.QMEGS(unitary: Block, state: Block, n_shots: int | Shots | None, target_indices: List[int], sigma: float = 1.0, eta: float = 0.01, T: int = 100, filtering_function: Callable | None = None, mode_dataset: Literal['analytical', 'sampling'] = 'analytical', mode_time: Literal['rvs', 'rejection_sampling'] = 'rvs', overlaps: Literal['classical'] | Tuple[float, float] = 'classical', verbose: bool = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

Quantum Multiple Eigenvalue Gaussian filtered Search (QMEGS) algorithm.

Implementation based on arxiv2402.01013. This algorithm finds multiple eigenphases of a unitary operator using Hadamard tests with Gaussian filtered time sampling.

Parameters:

• unitary – Block representing the time-evolution unitary.

• state – Block representing the trial state.

• n_shots – Number of measurement shots — same meaning in both modes: an int adds shot noise (synthetic in “analytical” mode, sampled in “sampling” mode with a sampling primitive); qarp.EXACT means no shot noise; None means the default shot count (10_000 / engine default). The default exact primitive ignores it in “sampling” mode.

• target_indices – List of indices of target eigenvalues to find (Dset in the paper).

• sigma – Truncation level for time sampling.

• eta – Error tolerance for the algorithm.

• T – Time window for time sampling.

• filtering_function – Custom filtering function for time sampling. Defaults to truncated_gaussian_density.

• mode_dataset – Dataset generation mode. - “analytical”: Generate data using classical simulation (faster, for testing). - “sampling”: Generate data using actual sampling execution.

• mode_time – Time sampling mode. - “rvs”: Use scipy’s truncnorm.rvs for fast sampling. - “rejection_sampling”: Use rejection sampling (enables arbitrary filtering functions).

• overlaps –

  Where the algorithm’s a-priori inputs p_min / p_tail (Theorem 1 of arXiv:2402.01013) come from. - "classical" (default): build() dense-diagonalises the

  Hamiltonian (O(4^n)) and reads the trial statevector through the engine’s simulator — a validation harness, refused on a noisy or routed engine. Required by mode_dataset="analytical".

  • (pmin, ptail): the two overlaps supplied by the user; no diagonalisation, no statevector. Only len(target_indices) is consulted (the number of eigenvalues to extract) — the index values themselves are not used.

• verbose – Whether to print progress information.

• primitive – Primitive used to evaluate ⟨ψ|U(t)|ψ⟩ in “sampling” mode. Defaults to StateVector() — exact, noise-free amplitudes. Pass HadamardTest() for finite-shot protocol realism.

• engine – Quantum engine for circuit execution.

build()

Build the algorithm by computing necessary parameters.

Calculates overlaps, probabilities, and algorithm parameters based on the input unitary and trial state.

Returns:

Self for method chaining.

Raises:

ValueError – If the algorithm’s assumptions are not satisfied (pmin must be greater than ptail).

filtered_density_function(dataset: List[tuple], theta_js: ndarray, n_samples: int) → ndarray

Compute the filtered density function G_j for given theta values.

Parameters:

• dataset – List of (time, measurement) tuples from quantum measurements.

• theta_js – Array of theta values to evaluate.

• n_samples – Number of samples for normalization.

Returns:

filtered density values for each theta.

Return type:

G_js

generate_data() → List[tuple]

Generate the measurement dataset according to selected mode.

Returns:

List of (time, measurement) tuples.

run() → List[float]

Run the QMEGS algorithm to find eigenphases.

Follows Algorithm 2 in arxiv2402.01013.

Returns:

List of found eigenphases.

validate_inputs()

Validate input parameters for correctness.

Raises:

ValueError – If any input parameter is invalid.

class qarp.algorithms.QPE(state: Block, unitary: Block, n_ancilla: int, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

build()

Build the Canonical Phase Estimation circuit.

First offers the problem to the engine’s structured fast path (Engine.prepare_structured_qpe — matrix exponentiation, the controlled-U ladder is never compiled). Engines without the path, or refusing it (EXACT readout, noise, routing, parametric U — see QarpEngine.prepare_structured_qpe), return None and the full QPE circuit is built instead.

Returns:

The instance of the class.

Return type:

self

estimate_phase(fit_range=None, verbose=False)

Estimation of the phase done by fitting the Dirichlet kernel squared to the distribution obtained from the run.

Parameters:

• fit_range (tuple or None) – Optional (min, max) range between 0 and 1 for fitting. If None, use on the entire histogram.

• verbose (bool) – If True, print the fitted phase.

Returns:

The estimated phase.

Return type:

phi_fit

plot(figsize=None, return_fig=False, max_xticks=11)

Plot the results of the Canonical Phase Estimation algorithm.

Parameters:

• figsize (tuple) – Size of the figure. If None, scales with n_ancilla.

• return_fig (bool) – If True, return (fig, ax) for external saving/customization.

• max_xticks (int) – Upper bound on the number of x-axis ticks.

Returns:

(fig, ax) if return_fig is True, otherwise None.

run()

Run the Canonical Phase Estimation algorithm.

Returns:

The estimated eigenvalue.

Return type:

result

class qarp.algorithms.QSE(hamiltonian: QubitOperator | Block, ground_state: Block, primitive: PrimitiveAlgorithm, overlap_primitive: PrimitiveAlgorithm, excitation_operators: List[QubitOperator], engine: Engine | None = None, real_symmetric: bool = False, verbose: bool = False)

Bases: CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

compute_hamiltonian()

Constructs the Hamiltonian matrix H in HC = SCe from the primitives.

Note

It is assumed that the QSE Hamiltonian is symmetric or Hermitian.

Returns:

A numpy array representing the Hamiltonian matrix H.

compute_overlap()

Constructs the overlap matrix S from the primitives.

Note

Assumes that the QSE overlap matrix is real valued and symmetric or Hermitian.

Returns:

A numpy array representing the overlap matrix S.

property hamiltonian_primitives: List[PrimitiveAlgorithm]

property overlap_primitives: List[PrimitiveAlgorithm]

run() → Tuple[ndarray, ndarray]

Runs the QSE algorithm.

Returns:

A tuple containing the eigenvalues and eigenvectors of the QSE Hamiltonian.

solve() → Tuple[NDArray, NDArray]

Solves the generalized eigenvalue problem HC = SCe.

Uses the QSE Hamiltonian stored in self.hamiltonian_matrix and overlap matrix stored in self.overlap_matrix.

Returns:

A tuple containing the eigenvalues and eigenvectors.

class qarp.algorithms.SSVQE(operator: QubitOperator | Block, basis_state_blocks: List[Block], ansatz_block: Block, weights: List[float], initial_parameters: Mapping | NDArray[float64] | None = None, optimizer: Optimizer | None = None, verbose: bool = True, gradient: bool | str = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block(index: int) → Block

Built block for optimized state index: basis state + ansatz bound at the optimal parameters. No symbol handling required — evaluate ⟨N⟩, ⟨S²⟩, etc. directly on the returned block.

objective(x)

objective_gradient(x)

property optimal_parameters: Dict[Symbol, float]

Optimized parameters keyed by symbol — the order-proof surface.

Use this (not manual zip against a symbol list) to evaluate properties of the optimized states.

run()

Run the SSVQE algorithm to find multiple eigenstates.

Returns:

• List of computed energies for each eigenstate.

• Optimized parameters as a NumPy array.

Return type:

A tuple containing

class qarp.algorithms.SWAPTest(bra: Block | None = None, operator: Block | None = None, ket: Block | None = None, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

gradient_kind: str = 'expectation'

returns_probability: bool = True

run(results: list) → float

2·P(ancilla=0) − 1 = |⟨bra|ket⟩|².

SamplingResult.counts is keyed on the full all-qubit outcome integer; the ancilla is qubit 0 so its measured value is bit 0.

supported_targets: frozenset[Target] = frozenset({Target.OVERLAP})

class qarp.algorithms.Sampler(ket: Block | None = None, n_shots: int | Shots | None = None, measured_qubits: list[int] | None = None, initial_state=None)

Bases: PrimitiveAlgorithm

qarpx-backed sampling primitive.

Parameters:

• ket – Block to execute.

• n_shots – Number of measurement shots (default: use engine default). qarp.EXACT returns the exact Born distribution |ψ|² — note: probabilities, not amplitudes (phases are discarded; raw amplitudes are simulator internals, sim.statevector).

• measured_qubits – Qubit indices whose outcomes appear in the output tuple. Defaults to range(ket.n_qubits). Non-measured qubits are marginalised out by summing probabilities over their bit values.

• initial_state – Optional LSB-indexed amplitudes seeding the register instead of |0…0⟩ (length 2**n_qubits, unit norm within 1e-10 — ValueError at run otherwise). QarpEngine only; other engines raise CapabilityError at build. Mutable between runs for step → snapshot → re-seed loops.

accepts_initial_state: bool = True

build() → Self

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

run(results: list) → dict[tuple[int, ...], float]

Convert the first result (qx.SamplingResult, or the duck-typed ExactResult under n_shots=qarp.EXACT) to a {bitstring-tuple: probability} dict.

supported_targets: frozenset[Target] = frozenset({Target.SAMPLING})

class qarp.algorithms.ShadowDataset(kernel: ShadowKernel, n_qubits: int, records: list[tuple[ndarray, dict[int, float]]], shot_exact: bool = False)

Bases: object

A collected campaign: settings + outcomes + the ensemble kernel.

estimator()

Return a ShadowEstimator over this dataset.

classmethod from_dict(d: dict) → ShadowDataset

classmethod load(path: str | PathLike[str]) → ShadowDataset

merge(other: ShadowDataset) → ShadowDataset

Concatenate two compatible campaigns (same ensemble, n_qubits, shot_exact). Batching is round-robin by setting index, so the merged batches stay evenly mixed between the two campaigns.

property n_settings: int

save(path: str | PathLike[str]) → None

Write to an NPZ file (schema-versioned). Accepts any os.PathLike.

to_dict() → dict

Schema-versioned logical form: metadata + flat integer arrays.

class qarp.algorithms.ShadowEstimate(value: float, error: float, delta: float, n_batches: int, ensemble_only: bool)

Bases: NamedTuple

One estimate with the HKP Theorem-1 accuracy half-width.

error is the eps = sqrt(34 * V_hat / N) of arXiv:2002.08953 Theorem 1 (V_hat the empirical single-setting variance, N the per-batch size) — the half-width for which Pr[|value - true| >= error] <= delta at the delta that fixed n_batches. That promise holds for the true variance sigma^2; V_hat is a plug-in for it, so the reported bar inherits whatever error that estimate carries (and HKP’s own Remark that the constant 34 is a loose worst case).

delta: float

Alias for field number 2

ensemble_only: bool

Alias for field number 4

error: float

Alias for field number 1

n_batches: int

Alias for field number 3

value: float

Alias for field number 0

class qarp.algorithms.ShadowEstimator(dataset)

Bases: object

Estimate observables from a collected ShadowDataset.

expval(observable: str | QubitOperator, *, delta: float = 0.05, n_batches: int | None = None) → ShadowEstimate

Median-of-means estimate of <observable> (sum-inside for a multi-term operator).

expval_many(observables: Sequence[str | QubitOperator], *, delta: float = 0.05, n_batches: int | None = None, marginal: bool = False) → list[ShadowEstimate]

Estimate a family of observables from the one dataset.

By default this is exactly Theorem 1 of HKP applied to the family: the batch count uses the union bound K = ~2 ln(2M/delta) over the M supplied, so each error is the half-width for which all M estimates lie within their band jointly with probability at least 1 - delta. marginal=True opts out to a per-observable delta (K = ~2 ln(2/delta), no union correction).

Looping expval() yourself does not give the joint guarantee — use this method for a family. And choosing observables after seeing results from the same dataset voids it: collect fresh or split.

An explicit n_batches overrides the derived count for every estimate (power-user knob; it then defeats the union bound this method exists to apply, so the joint guarantee no longer holds at delta).

fidelity(pure_state)

Fidelity to a pure state (global-Clifford ensemble; future).

one_rdm()

1-RDM (matchgate ensemble; future).

purity()

Tr(ρ²) (global-Clifford ensemble; future).

renyi2_entropy()

Rényi-2 entropy (global-Clifford ensemble; future).

two_rdm()

2-RDM (matchgate ensemble; future).

class qarp.algorithms.ShadowKernel(n_qubits: int)

Bases: ABC

Stateless inversion kernel for one measurement ensemble.

Variables:

• ensemble (str) – registry tag used by ShadowDataset.from_dict().

• n_qubits – register width the kernel was built for.

• capabilities (frozenset[str]) – estimator features this ensemble supports (e.g. {"expval"}; matchgate would add "rdm"). The estimator raises CapabilityError for anything outside this set.

capabilities: frozenset[str]

ensemble: str

classmethod from_params(n_qubits: int, params: dict) → ShadowKernel

params() → dict

Ensemble parameters for the serialized descriptor (besides n_qubits/ensemble). Default: none.

abstractmethod snapshot_estimate(setting: ndarray, outcome: int, term) → float

Single-snapshot estimate of one Pauli term under one (setting, outcome). term is a tuple of (qubit, axis_char) pairs (empty for identity — never passed here; the estimator handles the constant).

class qarp.algorithms.ShadowProtocol(operator: str | QubitOperator, ket: Block | None = None, *, dataset: ShadowDataset | None = None, n_settings: int = 1000, n_shots: int = 1, seed: int | None = None)

Bases: PrimitiveAlgorithm

Base collector for a randomized-measurement shadow ensemble.

build()

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property dataset: ShadowDataset

The collected (or passed) ShadowDataset. Raises if neither a dataset was passed nor a campaign collected yet.

gradient_kind: str = 'none'

release() → None

Free the compiled circuits in place (the detached dataset stays valid).

run() never does this automatically — the engine reuses the compiled circuits under rebuild=False and across batch_run parameter sets. After release(), a re-run needs a rebuild.

requires_noiseless: bool = True

run(results: list) → float

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE})

supports_exact: bool = True

class qarp.algorithms.Shor(number: int, *, base: int | None = None, n_counting_qubits: int | None = None, max_attempts: int = 8, base_seed: int | None = None, force_quantum: bool = False, primitive: Sampler | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

Exact small-integer/reference implementation of Shor factoring.

This implementation uses an exponentially synthesized basis permutation for modular arithmetic and supports at most ModularMultiplicationBlock.MAX_REFERENCE_WORK_QUBITS work qubits (number <= 64 for the shipped value of six). It establishes a correct reference workflow for small examples; it does not claim cryptographic-scale performance or an asymptotic quantum speedup.

build() prepares the order-finding samplers and run() produces the factor pair, classical shortcuts included. run() executes the batch once; a repeated call returns the stored outcome, None included, so a finite-shot retry is a new instance. Bases are processed in attempt order, as in Shor’s sequential algorithm: every coprime base gets an order-finding circuit, and the first non-coprime base ends the attempt list with gcd(base, number) as a classical fallback that run() returns only if every quantum attempt was inconclusive. Even and perfect-power inputs are factored classically without any circuit.

Parameters:

• number – Composite integer greater than one to factor.

• base – Optional first modular-order-finding base. Must satisfy 1 < base < number. A non-coprime base is a valid classical gcd shortcut unless force_quantum is set.

• n_counting_qubits – Counting-register width. Defaults to twice the work width and must be at least that large.

• max_attempts – Maximum number of distinct bases to try.

• base_seed – Seed for the algorithm-local base-selection RNG.

• force_quantum – Guarantee the quantum path. Disables the gcd shortcut (random bases are drawn coprime; a supplied non-coprime base is rejected) and rejects even or perfect-power inputs, which lie outside the preconditions of Shor’s order-finding theorem: for N = 2p and odd N = p**k every even order gives a**(r/2) == -1 (mod N), so no base can succeed.

• primitive – Sampling primitive. A private deep copy is used, and Shor owns its ket and measured_qubits (the counting register); the caller’s sampler contributes shot settings only.

• engine – Execution engine. Defaults to QarpEngine.

build() → Self

Choose bases and compile the order-finding samplers; never factors.

run() → tuple[int, int] | None

Execute the batch once and return the first recovered factor pair.

class qarp.algorithms.StateVector(bra: Block | None = None, operator: QubitOperator | Block | None = None, ket: Block | None = None, initial_state: ndarray | None = None)

Bases: PrimitiveAlgorithm

accepts_initial_state: bool = True

build() → Self

Populate sub_blocks with [ket] (or [bra, ket] for overlap-type).

Re-disambiguates self.target from the current bra / operator / ket attributes so callers (VQA-family composites) that mutate the primitive after construction get the right target without having to re-instantiate.

consumes: Consumes = 1

run(results: list) → complex | float

Convenience entry point: evaluate the build-time compiled circuits.

Engines dispatch directly to run_from_amplitudes() (because consumes is Consumes.AMPLITUDES), so this method exists mainly for callers driving the primitive without an engine. Symbolic parameters must be substituted into self.compiled_circuits before calling.

run_from_amplitudes(compiled_circuits: list, simulator=None) → complex | float

Simulate the substituted commands and contract the statevectors.

Parameters:

• compiled_circuits – One command stream per entry in self.sub_blocks, already parameter-substituted by the engine (or by the caller).

• simulator – Optional statevector backend implementing statevector(commands, n_qubits). Defaults to a CPU qx.QarpSimulator; CudaqEngine injects its GPU simulator so the statevectors are computed on-device. A simulator that also provides transition(bra, ket, n_qubits, observable) (QarpSimulator does) contracts the operator in C++; otherwise the numpy sweep pauli_expectation() is used.

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE, Target.OVERLAP, Target.TRANSITION_AMPLITUDE})

supports_backprop_gradient: bool = True

class qarp.algorithms.Target(*values)

Bases: Enum

Enumeration of measurement targets for primitive algorithms.

Target specifies the type of quantum measurement being performed by a primitive algorithm. Different targets determine how measurement circuits are constructed and how results are post-processed to compute the desired quantum mechanical quantity.

Variables:

• SAMPLING – Return raw measurement distribution as a dictionary of bitstrings to probabilities.

• EXPECTATION_VALUE – Compute the expectation value <ket|operator|ket⟩.

• OVERLAP – Compute the overlap (inner product) ⟨bra|ket⟩.

• TRANSITION_AMPLITUDE – Compute the transition amplitude ⟨bra|operator|ket⟩.

EXPECTATION_VALUE = 1

OVERLAP = 2

SAMPLING = 0

TRANSITION_AMPLITUDE = 3

class qarp.algorithms.TermwiseHadamardTest(bra: Block | None = None, operator: Block | List[Block] | QubitOperator | None = None, ket: Block | None = None, coefficients: List[float] | None = None, real: bool = True, imaginary: bool = True, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → TermwiseHadamardTest

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property expectation_type: str

get_hadamard_test(index: int) → HadamardTest

gradient_kind: str = 'expectation'

property n_operators: int

property n_sub_algorithms: int

run(results: list) → complex

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.EXPECTATION_VALUE, Target.OVERLAP, Target.TRANSITION_AMPLITUDE})

class qarp.algorithms.TermwiseSWAPTest(bra: Block | List[Block], ket: Block | None = None, coefficients: List[float] | None = None, n_shots: int | Shots | None = None)

Bases: PrimitiveAlgorithm

build() → TermwiseSWAPTest

Construct circuits: populate self.sub_blocks with qx.Block instances.

Returns self so callers can chain primitive.build().

property expectation_type: str

get_swap_test(index: int) → SWAPTest

gradient_kind: str = 'expectation'

property n_bra: int

property n_sub_algorithms: int

returns_probability: bool = True

run(results: list) → float

Post-process the engine’s sampling output into a result.

Contract: a pure function of results — no simulator or engine access. Keeps every estimator unit-testable with synthetic results (anything exposing counts / n_shots / n_qubits).

Parameters:

results – One qx.SamplingResult per entry in self.sub_blocks, as produced by the engine’s sampler.

Returns:

A scalar (float / complex) for expectation/overlap-style primitives, or a SamplingDictionary ({bitstring-tuple: probability}) for Sampler.

supported_targets: frozenset[Target] = frozenset({Target.OVERLAP})

class qarp.algorithms.VFF(operator: QubitOperator | Block, ansatz_block: Block, use_trotter: bool = False, t_time: float = 1.0, t_steps: int = 1, t_order: int = 1, initial_parameters: NDArray[float64] | None = None, optimizer: Optimizer | None = None, verbose: bool = True, primitive: PrimitiveAlgorithm | None = None, engine=None)

Bases: CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

objective(x)

property optimal_parameters

Optimized parameters keyed by symbol — the order-proof surface.

run()

Execute the algorithm after build().

The return type is algorithm-specific (documented per class); subclasses may add optional keyword-only arguments such as max_iter but take no positional arguments.

class qarp.algorithms.VQA(operator: QubitOperator | Block, ket: Block, name: str, initial_parameters: Mapping | NDArray[float64] | None = None, gradient: bool | str = False, optimizer: Optimizer | None = None, verbose: bool = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None, save_energy_history: bool = False)

Bases: CompositeAlgorithm

build()

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block() → Block

The ket bound at the optimal parameters (set by run()).

property optimal_parameters: Dict[Symbol, float]

Optimized parameters keyed by symbol — the order-proof surface.

run()

Run the VQA algorithm by calling the optimizer with the built objective function and optionally gradients.

Returns:

The final energy and final parameters as a float and list of floats.

class qarp.algorithms.VQD(operator: QubitOperator | Block, kets: List[Block], weights: List[float], initial_parameters: Iterable[Iterable[float]] | None = None, verbose: bool = False, gradient: bool | str = False, optimizer: Optimizer | None = None, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: CompositeAlgorithm

build() → Self

Build the composite algorithm by constructing all sub-algorithms.

This method should: 1. Create and configure all necessary sub-algorithms 2. Build each sub-algorithm

Returns:

Self for method chaining

get_final_state_block(index: int) → Block

Ket index bound at its optimal parameters.

iterate()

objective(theta: Iterable[float])

objective_gradient(theta: Iterable[float])

property optimal_parameters: List[Dict[Symbol, float]]

Per-state optimized parameters keyed by symbol (alias of state_parameters — already the order-proof surface).

run()

Run the VQD algorithm to find multiple eigenstates.

Returns:

• List of computed energies for each eigenstate.

• List of optimized parameter dictionaries for each eigenstate.

Return type:

A tuple containing

class qarp.algorithms.VQE(operator: QubitOperator | Block, ket: Block, initial_parameters: NDArray[float64] | None = None, gradient: bool | str = False, optimizer: Optimizer | None = None, verbose: bool = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None, save_energy_history: bool = False)

Bases: VQA

class qarp.algorithms.WalkerState(state_data: ndarray | Block, sign: float, label: str)

Bases: NamedTuple

Represents a quantum walker state.

Variables:

• state_data (numpy.ndarray | qarpx.Block) – State vector (np.ndarray) for semiclassical or Block for quantum mode

• sign (float) – Sign of the walker (+1.0 or -1.0)

• label (str) – String label identifying which basis state this walker represents

label: str

Alias for field number 2

sign: float

Alias for field number 1

state_data: ndarray | Block

Alias for field number 0

class qarp.algorithms.cPCE(graph: Graph, order: int, ket: Block, merging: bool = False, classical_function: Callable = <function classical_function>, quantum_function: Callable = <function quantum_function>, coefficients: Iterable[float] | None = None, initial_parameters: Iterable[float] | None = None, optimizer: Optimizer | None = None, verbose: bool = True, gradient: bool | str = False, primitive: PrimitiveAlgorithm | None = None, engine: Engine | None = None)

Bases: PCE

Continuous Pauli Correlation Encoding (cPCE).

cPCE reuses the same measurement and optimization machinery as PCE, but drops the tanh + sign thresholding step used to turn Pauli expectation values into a discrete bitstring. Instead, the expectation values are used directly (optionally rescaled by coefficients) as the solution to a real-valued optimization problem.

This implementation follows the approach described in https://arxiv.org/abs/2604.05637.

qarp.algorithms.calculate_qubits(n_nodes: int, order: int, merging: bool = False) → int

Compute the needed number of qubits for a given number of graph nodes and order

Parameters:

• n_nodes – Number of nodes to map into the qubits.

• order – Order to consider in the Pauli correlators.

• merging – If True, considers merging of same Pauli terms of different correlators.

Returns:

Number of qubits required.

qarp.algorithms.classical_function_max_cut(graph: Graph, binary_string: list[int]) → float

Calculate the cut size of the given partitions of the graph.

Parameters:

• graph – Original graph over which to compute the cut size.

• binary_string – Binary representation of solution, where each position represents a node, and it value (0 or 1) represents the subset location.

Returns:

Sum of weighted edges cut due to the partitions.

qarp.algorithms.dirichlet_kernel_squared(x, phi, N)

Dirichlet kernel squared function. Used to reconstruct the phase information in the QPE algorithm.

Parameters:

• x (array-like) – Input values (e.g. eigenvalues).

• phi (float) – Phase to fit.

• N (int) – Number of qubits.

Returns:

Squared Dirichlet kernel values.

Return type:

array-like

qarp.algorithms.find_eigenspectrum_degeneracy(eigs, tolerance=1e-15, verbose=False)

Finds the degeneracy of the eigenvalues in the eigenspectrum.

Parameters:

• eigs (np.ndarray) – The eigenvalues of the Hamiltonian.

• tolerance (float) – The tolerance for determining degeneracy.

• verbose (bool) – If True, prints the eigenvalues and their degeneracy.

Returns:

A dictionary with eigenvalues as keys and their degeneracy as values.

Return type:

dict

qarp.algorithms.find_occupation_numbers(hamiltonian, n_qubits, tolerance=1e-15, verbose=False)

Finds the occupation numbers of the eigenstates of a Hamiltonian.

Parameters:

• hamiltonian (FermionOperator) – The Hamiltonian operator.

• n_qubits (int) – The number of qubits.

• tolerance (float) – The tolerance for determining occupation numbers.

• verbose (bool) – If True, prints the eigenvalues and occupation numbers.

Returns:

An array of occupation numbers for the eigenstates.

Return type:

np.ndarray

qarp.algorithms.find_unique_eigs_and_occupation_numbers(eigs, occ_numbers, select_occ=None, tolerance=1e-15, verbose=False)

Finds unique eigenvalues and their corresponding occupation numbers.

Parameters:

• eigs (np.ndarray) – The eigenvalues of the Hamiltonian.

• occ_numbers (np.ndarray) – The occupation numbers corresponding to the eigenvalues.

• select_occ (int, optional) – If provided, filters the unique eigenvalues by this occupation number.

• tolerance (float) – The tolerance for determining degeneracy.

• verbose (bool) – If True, prints the unique eigenvalues and their occupation numbers.

Returns:

A tuple containing:

• dict: A dictionary with unique eigenvalues as keys and their degeneracy as values.

• list: A list of unique occupation numbers corresponding to the unique eigenvalues.

Return type:

tuple

qarp.algorithms.generate_states_new_basis(U: Block, hamming_weight: int | List[int] | None = None, get_statevector: bool = True) → Tuple[List[ndarray], List[Block], List[int]]

Apply a unitary block U to every computational basis state and (optionally) compute the resulting statevectors.

The new basis is \(\{ U |i\rangle : i \in \mathrm{indices} \}\).

WARNING: exponential scaling — the indices set is 2^n_qubits unless hamming_weight filters it down.

Parameters:

• U – A built qarp.blocks.AnyBlock whose unitary is the change of basis.

• hamming_weight – If int, restrict to computational basis states with exactly that Hamming weight. If list[int], restrict to the union of those weights. If None, use every basis state.

• get_statevector – If True, compute each new basis state’s statevector via the qarpx simulator (still exponential). If False, return an empty statevector list.

Returns:

• new_basis_states — list of length-2^n_qubits ndarrays (only populated when get_statevector=True).

• new_basis_blocks — Block per basis state: ComputationalBasisStateBlock(bitstring) + U.

• basis_states_indices — the integer indices the basis was built from.

Return type:

(new_basis_states, new_basis_blocks, basis_states_indices) where

qarp.algorithms.get_overlaps(trial_state: Block, hamiltonian: QubitOperator, return_eigenvalues: bool = False) → ndarray | Tuple[ndarray, ndarray]

Compute overlaps of trial state with eigenstates of the Hamiltonian.

Parameters:

• trial_state – Block representing the trial state.

• hamiltonian – QubitOperator representing the Hamiltonian.

• return_eigenvalues – If True, also return the eigenvalues.

Returns:

Array of overlap probabilities with each eigenstate, sorted by ascending eigenvalue. If return_eigenvalues is True, returns a (overlaps, eigenvalues) tuple instead.

Note

Uses np.linalg.eigh which returns eigenvectors as COLUMNS of the matrix. The overlaps are computed correctly by iterating over columns.

class qarp.algorithms.iterativePCE(graph: Graph, order: int, ket: Block, pce_cls: type = <class 'qarp.algorithms._composite.pce.PCE'>, alpha0: float = 1.0, threshold: float = 0.9, max_outer_iterations: int = 50, stabilized_update: bool = False, verbose: bool = True, **pce_kwargs)

Bases: CompositeAlgorithm

Iterative-alpha Pauli Correlation Encoding.

Wraps a PCE-like algorithm (by default PCE, but compatible with any class sharing its constructor/build/run/get_objective_function surface) and repeatedly re-runs it, progressively increasing the sharpness parameter alpha used to binarize Pauli-correlator expectation values via tanh(alpha * <Pi>). At each round, only the least-binarized correlator is pushed just past the binarization threshold M, and the ansatz parameters are warm-started from the previous round’s optimum.

This avoids the tradeoff of a single fixed alpha: too small and the loss/constraint is evaluated on an under-binarized (and thus misleading) solution; too large and tanh’s vanishing derivative stalls the optimizer everywhere at once.

Implements Algorithm 1 (“Iterative-alpha PCE heuristic”) from https://arxiv.org/abs/2602.17479.

After run(), alpha_history, raw_expectations_history and solution_history hold, per round, the alpha used, the raw <Pi_i> expectation values obtained at convergence, and the rounded solution. There is exactly ONE inner algorithm — built once and re-run — so sub_algorithms has a single entry, not one per round.

build()

No-op: each round’s inner algorithm is built when it is constructed in run.

Returns:

self

run()

Run the Iterative-alpha PCE heuristic.

At every round, runs pce_cls to convergence with the current alpha, then rescales alpha to push the least-binarized correlator just past threshold, warm-starting the next round from the converged ansatz parameters. Stops once every correlator is binarized or max_outer_iterations is reached.

Returns:

The final energy, final parameters and solution, as returned by the last round’s pce_cls.run().

qarp.algorithms.map_binary_to_integer_keys(probs)

Maps binary keys to integer keys in a dictionary.

Parameters:

probs (dict) – A dictionary with binary keys and float values.

Returns:

A new dictionary with integer keys and float values.

Return type:

dict