quchip.backend

Backend selection and concrete backend implementations.

quchip ships two backends:

• QuTiPBackend — CPU, qutip.Qobj operators, process-parallel sweeps via loky. Default backend.

• DynamiqsBackend — JAX/dynamiqs, fully differentiable, native vmap batched solves. Optional extra (pip install quchip[dynamiqs]).

The public entry points are get_default_backend() / set_default_backend() / reset_default_backend(). Engine internals use the thread-safe _backend_context() to scope a backend for the duration of a single assembly pass (every backend call the engine makes is reentrant because the override is a ContextVar).

Example

>>> from quchip.backend import get_default_backend, set_default_backend
>>> backend = get_default_backend()  # QuTiPBackend by default
>>> set_default_backend("dynamiqs")

Functions

get_default_backend()	Return the active default backend.
reset_default_backend()	Clear the cached default backend; the next get_default_backend() rebuilds it.
set_default_backend(backend)	Set the active default backend by name ("qutip"/"dynamiqs") or instance.

class quchip.backend.Backend

Bases: ABC

Abstract contract implemented by every quchip operator backend.

A backend wraps a quantum-dynamics library (QuTiP, dynamiqs, …) and exposes three orthogonal surfaces:

• Operator algebra — destroy/create/number/identity, tensor, embed, dag, matmul, eigensystem helpers. Concrete operators support native Python arithmetic (+ - * @); quchip never introduces scale/add wrappers.

• State algebra — basis/coherent, overlap, ptrace, expect, ket ↔ density-matrix conversion.

• IR lowering & solve dispatch — prepare_hamiltonian / prepare_batch convert engine IR into native RHS form; sesolve / mesolve run the solver; solve_problem / solve_batch drive the full engine-to-result pipeline.

Defaults provided here are NumPy-based and correct for any backend, so concrete subclasses only override where they can do better (e.g. QuTiP reuses Qobj.eigenstates, dynamiqs reuses dq.expect).

abstract property array_module: Any

Return the array module used by backend-aware numeric code.

Returns numpy for CPU-only backends and jax.numpy for JAX-traceable backends. Engine code that must stay differentiable (signal evaluation, frame assembly) routes array ops through this module so the dynamiqs backend preserves JAX traceability end-to-end.

abstractmethod to_array(op)

Return a dense array (native to array_module) for op.

Parameters:

op (Any)

Return type:

Any

overlap(a, b)

Return the scalar inner product ⟨a|b⟩ for two kets.

Parameters:

• a (Any)

• b (Any)

Return type:

complex

norm(state_or_op)

Return the Frobenius / ℓ²-norm of a state or operator.

A concrete float for CPU-only backends; a JAX-traceable backend may return a native 0-d array instead so a traced amplitude stays differentiable — callers that need concreteness route through maybe_concrete_scalar().

Parameters:

state_or_op (Any)

Return type:

Any

trace(op)

Return the scalar trace Tr(op).

Parameters:

op (Any)

Return type:

complex

destroy(n)

Build the annihilation operator \(\hat a\) for an n-level Fock space.

Parameters:

n (int)

Return type:

Any

create(n)

Build the creation operator \(\hat a^\dagger\) for an n-level Fock space.

Parameters:

n (int)

Return type:

Any

number(n)

Build the number operator \(\hat n = \hat a^\dagger \hat a\).

Parameters:

n (int)

Return type:

Any

identity(n)

Build the identity operator for an n-level space.

Parameters:

n (int)

Return type:

Any

from_array(data, dims=None)

Construct a native operator from a dense matrix.

dims accepts quchip’s row/col layout ([[row_dims], [col_dims]]) or a flat list; None means a single subsystem of size data.shape[0].

Parameters:

• data (Any)

• dims (list[list[int]] | None)

Return type:

Any

diag(values, dims=None)

Construct a backend-native diagonal operator from main-diagonal values.

Backends that distinguish layouts (e.g. dynamiqs sparse-DIA) should return their sparse representation so element-wise composition with other diagonal operators stays sparse. Default falls through to from_array() after building a dense np.diag.

Parameters:

• values (Any)

• dims (list[list[int]] | None)

Return type:

Any

abstractmethod to_canonical_operator(op)

Serialize a native operator into the backend-agnostic canonical IR.

Parameters:

op (Operator)

Return type:

CanonicalOperator

abstractmethod from_canonical_operator(canonical)

Reconstruct a native operator from the canonical IR payload.

Parameters:

canonical (CanonicalOperator)

Return type:

Operator

coerce_operator(op)

Coerce an array-like local operator into backend-native form.

The operator-side mirror of coerce_state(). Device-side physical-coupling accessors (e.g. charge_coupling_operator()) return dense, trace-safe arrays rather than backend-native operators; composition entry points (tensor(), dag()) route their operands through this hook so both forms compose freely. Backends whose native operators are already arrays override with a trace-safe passthrough.

Parameters:

op (Any)

Return type:

Any

dag(op)

Return the Hermitian conjugate op†.

Parameters:

op (Any)

Return type:

Any

matmul(a, b)

Return the matrix product a @ b.

Parameters:

• a (Any)

• b (Any)

Return type:

Any

eigenenergies(op)

Return the ascending eigenvalues of a Hermitian operator.

Parameters:

op (Any)

Return type:

Any

eigenstates(op)

Return (eigenvalues, eigenstates) of a Hermitian operator, ascending.

Parameters:

op (Any)

Return type:

tuple[Any, Any]

eigensystem_data(op)

Return ascending eigenvalues, stacked eigenvector matrix, and lazy eigenstates.

Parameters:

op (Any)

Return type:

EigensystemData

expect(op, state)

Return the expectation value ⟨state|op|state⟩ — accepts ket or density matrix.

Parameters:

• op (Any)

• state (Any)

Return type:

complex

ptrace(state, keep, dims)

Reduce a composite state onto subsystem(s) keep via partial trace.

Parameters:

• state (Any) – Composite ket or density matrix over the subsystems in dims.

• keep (int | list[int]) – Subsystem index (or indices) to retain; all others are traced out.

• dims (list[int]) – Subsystem dimensions of the full composite space, in order.

Return type:

Any

permute_state(state, dims, order)

Reorder a composite state’s subsystems.

dims are the subsystem levels in state’s current tensor order. order follows numpy.transpose convention: order[i] is the current-order subsystem index that becomes position i after permuting (so order == list(range(len(dims))) is a no-op).

Default: densify, reshape to one axis per subsystem, transpose, reshape back — correct for any backend. A ket reshapes/transposes once; a density matrix applies the same subsystem permutation to both the row and column index groups. Concrete backends may override this with a native reorder (e.g. QuTiP’s Qobj.permute) that stays sparse-friendly and keeps dims metadata.

Every array op routes through array_module (never bare numpy), so a traced state stays traced end-to-end under the dynamiqs backend. dims/order are always concrete Python ints (subsystem structure, never a physics parameter), so their product is taken with math.prod() rather than array_module.prod — routing a shape value through the array module would turn it into a tracer under jit and make the subsequent reshape fail.

Parameters:

• state (Any)

• dims (Sequence[int])

• order (Sequence[int])

Return type:

Any

stack_states(states)

Stack a trajectory’s saved states into a single (T, …) array.

Default densifies each state and ``np.stack``s along a new leading time axis. Backends whose solver already returns a stacked native state (dynamiqs) override this to a no-op pass-through.

Parameters:

states (Any)

Return type:

Any

expect_over_time(op, stacked_states)

Return ⟨op⟩(t) for every save point, as one (T,) array in array_module.

Accepts a ket stack (T, n, 1) (returns ⟨ψ|op|ψ⟩) or a density-matrix stack (T, n, n) (returns Tr(op·ρ)). One einsum replaces the per-point expect loop. Every intermediate stays in array_module, so the JAX backend keeps the trace differentiable.

Parameters:

• op (Any)

• stacked_states (Any)

Return type:

Any

overlap_over_time(target, stacked_states)

Return ⟨target|ψ(t)⟩ for every save point, as one (T,) array.

For a ket stack this is the phase-sensitive complex amplitude ⟨target|ψ(t)⟩ (never |·|² — phase-dependent gradients flow through it). For a density-matrix stack there is no single phase, so it returns the target-state population ⟨target|ρ(t)|target⟩. Stays in array_module for JAX traceability.

Parameters:

• target (Any)

• stacked_states (Any)

Return type:

Any

populations_over_time(stacked_states)

Return full-chip diagonal populations (T, ∏dims) (real) for every save point.

For a ket stack this is |ψ(t)|² along the level axis (the density matrix is never built); for a DM stack it is the real diagonal. Stays in array_module for JAX traceability.

Parameters:

stacked_states (Any)

Return type:

Any

ptrace_over_time(stacked_states, keep, dims)

Reduce onto subsystem(s) keep at every save point, returning a (T, k, k) DM stack.

Reduces a ket or DM trajectory onto subsystem(s) keep without a per-point Python loop.

Parameters:

• stacked_states (Any)

• keep (int | list[int])

• dims (list[int])

Return type:

Any

embed(op, device_index, dims)

Embed a single-device operator at device_index into the full tensor space.

Validation is shared; the identity factors and tensor product dispatch through the backend’s own identity() / tensor(), so each backend keeps its native-optimal layout.

Parameters:

• op (Any) – Single-device operator whose dimension matches dims[device_index].

• device_index (int) – Position of op within the tensor product.

• dims (Sequence[int]) – Subsystem dimensions of the full space, in order.

Return type:

Any

abstractmethod embed_two_body(op_ab, index_a, index_b, dims)

Embed a two-body operator on devices index_a ⊗ index_b into the full space.

Handles subsystem SWAP when index_a > index_b and identity-padding for non-adjacent devices, without materializing the dense full matrix when the backend supports sparse layouts.

Parameters:

• op_ab (Any)

• index_a (int)

• index_b (int)

• dims (Sequence[int])

Return type:

Any

basis(n, k)

Build the Fock basis ket \(|k\rangle\) in an n-level space.

Parameters:

• n (int)

• k (int)

Return type:

Any

tensor_states(*states)

Return the tensor product of states (defaults to tensor()).

Parameters:

states (Any)

Return type:

Any

coherent(n, alpha)

Build the coherent state \(|\alpha\rangle\) truncated to n Fock levels.

Built from the analytic series \(|\alpha\rangle = e^{-|\alpha|^2/2} \sum_k \alpha^k/\sqrt{\Gamma(k+1)}\,|k\rangle\).

Parameters:

• n (int)

• alpha (complex)

Return type:

Any

state_to_dm(state)

Return a density matrix; pass through if state is already one.

Parameters:

state (Any)

Return type:

Any

is_ket(state)

Return whether state is a column vector (ket) rather than a density matrix.

Parameters:

state (Any)

Return type:

bool

as_density_matrix(state)

Promote a ket to its density matrix; pass a density matrix through unchanged.

Parameters:

state (Any)

Return type:

Any

abstractmethod tensor(*operators)

Return the tensor product of operators, preserving subsystem dims metadata.

Parameters:

operators (Any)

Return type:

Any

abstractmethod sesolve(H, psi0, tlist, e_ops=None, options=None)

Solve the Schrödinger equation \(i\hbar \partial_t |\psi\rangle = H |\psi\rangle\).

Parameters:

• H (Any) – Native Hamiltonian (Qobj/QobjEvo or dynamiqs TimeQArray).

• psi0 (Any) – Initial ket.

• tlist (Any) – 1D array of save times in ns.

• e_ops (list[Any] | None) – Optional observables to accumulate ⟨ψ|Oᵢ|ψ⟩ along the trajectory.

• options (dict[str, Any] | None) – Backend-specific options dict. See each backend’s resolve_solver_options for supported keys.

Return type:

SolverResult

abstractmethod mesolve(H, rho0, tlist, c_ops=None, e_ops=None, options=None)

Solve the Lindblad master equation (Lindblad 1976; Breuer & Petruccione 2002).

Integrates \(\dot\rho = -i[H,\rho] + \sum_k \mathcal D[L_k]\rho\) with \(\mathcal D[L]\rho = L\rho L^\dagger - \tfrac12\{L^\dagger L, \rho\}\).

Parameters:

• H (Any)

• rho0 (Any)

• tlist (Any)

• c_ops (list[Any] | None)

• e_ops (list[Any] | None)

• options (dict[str, Any] | None)

Return type:

SolverResult

batched_sesolve(problems, *, n_jobs=-1, progress=True)

Solve multiple sesolve problems. Default: sequential dispatch.

Backends that can vectorize (dynamiqs vmap) or parallelize (QuTiP via loky) override this method. Each problem dict carries the same kwargs as sesolve().

Parameters:

• problems (list[dict[str, Any]])

• n_jobs (int)

• progress (bool)

Return type:

list[SolverResult]

batched_mesolve(problems, *, n_jobs=-1, progress=True)

Solve multiple mesolve problems. Default: sequential dispatch.

Parameters:

• problems (list[dict[str, Any]])

• n_jobs (int)

• progress (bool)

Return type:

list[SolverResult]

prepare_hamiltonian(description, tlist)

Convert a HamiltonianDescription into a native solver RHS.

The engine passes frequencies already converted by 2π (angular rad/ns); backends must not rescale. Concrete backends override this method to choose their native time-dependence representation.

Parameters:

• description (HamiltonianDescription)

• tlist (Any)

Return type:

PreparedHamiltonian

resolve_solver_options(options, *, metadata, tlist)

Merge user options with backend-side defaults and metadata heuristics.

Default is a shallow copy; concrete backends use metadata (e.g. spectral_bound_ghz) to fill in integrator step budgets when the user did not specify one.

Parameters:

• options (dict[str, Any])

• metadata (dict[str, Any])

• tlist (Any)

Return type:

dict[str, Any]

coerce_state(state, dims=None)

Convert a foreign-native state into this backend’s native form.

Called at the solve boundary so a per-call backend= override (see QuantumSequence.simulate()) accepts initial states that were built under another backend — e.g. a QuTiP Qobj from chip.state(...) handed to a dynamiqs gradient solve. Default: pass through unchanged. dims are the chip’s subsystem levels, for backends whose state type carries tensor structure.

Parameters:

• state (Any)

• dims (tuple[int, ...] | None)

Return type:

Any

solve_problem(problem)

Lower and solve a SolveProblem — the single-element entry point.

Picks mesolve when collapse operators are present (open system) or sesolve otherwise, unless problem.solver forces a choice.

Parameters:

problem (SolveProblem)

Return type:

SolverResult

prepare_batch(description, tlist)

Lower a BatchedHamiltonianDescription into a prepared batch.

The return type declares the batching strategy: EagerBatch (one RHS per element), VmappedBatch (one natively batched RHS — dynamiqs), or DeferredBatch (backend-private payload consumed by an overridden solve_batch() — QuTiP). Default: lowers each element independently via prepare_hamiltonian() into an EagerBatch.

Parameters:

• description (BatchedHamiltonianDescription)

• tlist (Any)

Return type:

PreparedBatch

solve_batch(batch, *, progress=True)

Lower and solve a SolveBatch.

Default: prepares the batch then dispatches each element’s RHS through batched_sesolve() / batched_mesolve(). Native backends override to avoid the per-element unpack.

Parameters:

• batch (SolveBatch)

• progress (bool)

Return type:

list[SolverResult]

quchip.backend.Operator

alias of Any

quchip.backend.State

alias of Any

class quchip.backend.SolverResult(times, states=None, expect=None, final_state=None, stats=<factory>, solver='')

Bases: object

Backend-agnostic container for a single time-evolution solve.

Parameters:

• times (Any)

• states (list[Any] | None)

• expect (list[list[float]] | dict[str | tuple[str, str], Any] | None)

• final_state (Any | None)

• stats (dict[str, Any])

• solver (str)

times

1D array of save times (ns) matching the solver’s tlist.

Type:

Any

states

Saved states, one per save time, in the backend’s native state type (Qobj / QArray). None if store_states was disabled.

Type:

list[Any] | None

expect

Expectation-value traces. list[list[float]] indexed by e_ops-then-time, or a dict when the engine supplied labelled observables (keys are either strings or (drive_label, op) pairs).

Type:

list[list[float]] | dict[str | tuple[str, str], Any] | None

final_state

Final state (states[-1] when states are stored). Useful for multi-segment protocols without retaining full trajectories.

Type:

Any | None

stats

Solver-specific diagnostics (integrator step count, batch index, etc.). Purely informational; no physics depends on it.

Type:

dict[str, Any]

solver

The dispatched solver name ("sesolve" or "mesolve").

Type:

str

times: Any

states: list[Any] | None = None

expect: list[list[float]] | dict[str | tuple[str, str], Any] | None = None

final_state: Any | None = None

stats: dict[str, Any]

solver: str = ''

class quchip.backend.PreparedHamiltonian(rhs, metadata=<factory>)

Bases: object

Backend-native Hamiltonian produced by Backend.prepare_hamiltonian().

rhs is whatever the backend’s solver accepts directly — a Qobj / QobjEvo for QuTiP, a dynamiqs TimeQArray / sum of them for dynamiqs. metadata passes engine-level hints (e.g. spectral_bound_ghz for integrator step heuristics) through opaquely.

Parameters:

• rhs (Any)

• metadata (dict[str, Any])

rhs: Any

metadata: dict[str, Any]

class quchip.backend.EagerBatch(rhs_list, batch_size, metadata=<factory>, tlist=None)

Bases: object

Batched-solve payload with one native RHS per element.

The default Backend.solve_batch() dispatches each element’s RHS through Backend.batched_sesolve() / Backend.batched_mesolve().

Parameters:

• rhs_list (list[Any])

• batch_size (int)

• metadata (dict[str, Any])

• tlist (Any | None)

rhs_list: list[Any]

batch_size: int

metadata: dict[str, Any]

tlist: Any | None = None

class quchip.backend.VmappedBatch(rhs, batch_size, metadata=<factory>, tlist=None)

Bases: object

Batched-solve payload with a single natively batched RHS.

rhs covers every element at once — the dynamiqs path, where the per-element signals are stacked along a leading batch axis and the solver runs one vmapped call.

Parameters:

• rhs (Any)

• batch_size (int)

• metadata (dict[str, Any])

• tlist (Any | None)

rhs: Any

batch_size: int

metadata: dict[str, Any]

tlist: Any | None = None

class quchip.backend.DeferredBatch(shared, batch_size, metadata=<factory>, tlist=None)

Bases: object

Batched-solve payload whose RHS construction is deferred.

shared carries backend-private state; the producing backend must override Backend.solve_batch() to consume it (the QuTiP path — final QobjEvo assembly happens inside loky workers so main-process overhead stays O(1) in batch size).

Parameters:

• shared (Any)

• batch_size (int)

• metadata (dict[str, Any])

• tlist (Any | None)

shared: Any

batch_size: int

metadata: dict[str, Any]

tlist: Any | None = None

class quchip.backend.EigensystemData(eigenvalues, eigenvector_matrix, _states_builder=None, _states_cache=None)

Bases: object

Hermitian eigensystem returned in a single diagonalization call.

eigenvalues is ascending. eigenvector_matrix stacks the eigenvectors as columns in the bare-product basis used for the diagonalization. The per-column backend-native kets are exposed lazily via the eigenstates property so the dressing / sweep hot path (which reads only eigenvalues + eigenvector_matrix + labeling) never pays for D backend-ket allocations and a second O(D**2) densification it does not use.

Backends populate either _states_builder (a callable that materializes the ket list on demand) or prime _states_cache directly (when the diagonalizer already produced the kets, e.g. QuTiP’s Qobj.eigenstates).

Parameters:

• eigenvalues (Any)

• eigenvector_matrix (Any)

• _states_builder (Callable[[], list[Any]] | None)

• _states_cache (list[Any] | None)

eigenvalues: Any

eigenvector_matrix: Any

property eigenstates: list[Any]

Return per-column backend-native eigenstate kets (built on first access).

Memoizes only when the materialized kets are tracer-free: under jit/grad/vmap the states carry tracers bound to the current trace, so caching them would let a stale tracer escape into a later trace. Concrete states are cached.

quchip.backend.compute_two_body_permutation(idx_a, idx_b, dims)

Compute permutations moving two subsystems to the front of a tensor product.

Returns (order, inverse) where order lists the original subsystem positions in their reordered layout (the two target subsystems first, in ascending index order, followed by spectators) and inverse maps each original position to its new slot so order[inverse[i]] == i.

Used by backends to embed a two-body operator acting on arbitrary device indices into the full tensor product without densifying spectators.

Parameters:

• idx_a (int) – Subsystem positions the two-body operator acts on, in any order.

• idx_b (int) – Subsystem positions the two-body operator acts on, in any order.

• dims (Sequence[int]) – Subsystem dimensions of the full tensor product; only len(dims) is read.

Return type:

tuple[list[int], list[int]]

Examples

>>> from quchip.backend._dims import compute_two_body_permutation
>>> compute_two_body_permutation(2, 0, [2, 2, 2, 2])
([0, 2, 1, 3], [0, 2, 1, 3])

quchip.backend.default_solver_steps(metadata, tlist)

Return a heuristic step budget: ~100 steps per period of the fastest scale, min 200,000.

nsteps is an abort ceiling for adaptive integrators, not a step-size choice — the integrator picks its own step and the ceiling only decides when to give up, so a generous value costs nothing on a converging solve. The previous tight budget (10/period, min 500) sat at the same order as the true step count of high-order methods at tight tolerance over long spans, which forced user code to carry {"nsteps": 200_000} overrides; the floor now bakes that in.

The fastest oscillation is the larger of the static spectral span (spectral_bound_ghz) and the largest time-dependent carrier (max_carrier_freq_ghz). Both are ordinary GHz. Considering the carrier is essential in a rotating frame, where the static span is near-zero but inter-mode / drive carriers still drive fast dynamics.

Parameters:

• metadata (dict[str, Any]) – Hamiltonian metadata; reads spectral_bound_ghz and max_carrier_freq_ghz (both ordinary GHz).

• tlist (Any) – Save-time grid; only its span tlist[-1] - tlist[0] (ns) is used.

• available (Returns None when neither scale is)

• is (the tlist span)

• non-positive

• be (or any input is a JAX tracer (which must not)

• default). (concretized — the caller falls back to the solver library's)

Return type:

int | None

quchip.backend.normalize_dims_from_list(dims, fallback=None)

Normalize quchip-style dims metadata to a flat subsystem-dim tuple.

Accepts [[rows], [cols]], [[rows]], or a flat list. None means a single subsystem of size fallback.

Parameters:

• dims (list[list[int]] | list[int] | None) – [[rows], [cols]], [[rows]], a flat list of subsystem dimensions, or None.

• fallback (int | None) – Single-subsystem dimension used when dims is None.

Return type:

tuple[int, …]

Examples

>>> from quchip.backend._dims import normalize_dims_from_list
>>> normalize_dims_from_list([[2, 3], [2, 3]])
(2, 3)
>>> normalize_dims_from_list(None, fallback=4)
(4,)

quchip.backend.validate_two_body_indices(index_a, index_b, dims)

Raise ValueError when two-body device indices are equal or out of range.

Parameters:

• index_a (int)

• index_b (int)

• dims (Sequence[int])

Return type:

None

quchip.backend.get_default_backend()

Return the active default backend.

Returns (in order) the nearest _backend_context() override, the user-set default from set_default_backend(), or a freshly instantiated QuTiPBackend (cached for subsequent calls).

Return type:

Backend

quchip.backend.set_default_backend(backend)

Set the active default backend by name ("qutip"/"dynamiqs") or instance.

Parameters:

backend (str | Backend)

Return type:

None

quchip.backend.reset_default_backend()

Clear the cached default backend; the next get_default_backend() rebuilds it.

Return type:

None

Modules

containers	Backend-agnostic solver-result and IR-lowering containers.
dynamiqs	Dynamiqs backend — JAX-native solver for differentiable, vmappable simulation.
protocol	Backend protocol and backend-agnostic solver-result containers.
qutip	QuTiP backend implementation of the quchip.backend.protocol contract.