Skip to content

Trotterization

Choose the number of Trotter steps for a Pauli Hamiltonian, an evolution time and an operator-norm error budget, evaluate the product-formula error bound behind that choice, and build the product-formula circuit. Import the functions from nwqlib.subroutines.trotterization.

from qiskit.quantum_info import SparsePauliOp
from nwqlib.subroutines.trotterization import (
    build_trotter_evolution_circuit,
    select_trotter_step_count,
)

H = SparsePauliOp(["XI", "ZZ"], [0.5, 0.5])
selection = select_trotter_step_count(H, time=2.0, error_budget=0.01, order=2)
print(selection.step_count, round(selection.bound_value, 10))
circuit = build_trotter_evolution_circuit(H, selection)
print(selection.synthesis_method, circuit.num_qubits)
8 0.0078125
suzuki_trotter 2

The second-order coefficient of this H is 1/16, and 8 is the smallest r with (1/16) * 2**3 / r**2 <= 0.01.

The error bound

  • The first-order bound is Proposition 9, Eq. (120), and the second-order bound is Proposition 10, Eq. (121), of Childs, Su, Tran, Wiebe and Zhu, Phys. Rev. X 11, 011020 (2021), doi:10.1103/PhysRevX.11.011020. One application of the order-o formula obeys bound(t) = W_up * t**(o + 1).
  • For an evolution at one supplied time, the step count is the smallest one allowed by the rule of their Sec. V B: with r steps the total is at most r * bound(t / r).
  • The bound expands each commutator norm over Pauli strings by the triangle inequality, with ||[P, Q]|| = 2 when P and Q anticommute and ||[P, [Q, R]]|| = 4 when both commutators are nonzero. It upper-bounds the theorem's tail-sum norms and is not the evolution error itself. At second order the default "exact_census" uses the full Pauli-triangle expression. bound_variant="relaxed_prefix" uses pair tests and weighted degrees only, replacing each nested mass by min(S_i,D_i+D_j-a_i), where a_i is the coefficient magnitude of term i, S_i the later coefficient mass and D_v the mass anticommuting with term v. It can choose more steps than the full expression. Proposition 36 gives the proof. Its work is w*P+2*p+4*E for p Pauli terms, packed width w, P = p(p-1)/2 possible pairs and E anticommuting pairs. Its contraction workspace is 72*p+96*b+96*ceil(E/b)+65536 bytes at block size b, counted together with the complete input and the pair table against the byte limit.
  • The prefactor W_up is evaluated with scaled binary64 upper products and sums followed by rational rescaling, so it bounds the chosen expression from above and can exceed it. The bound and the step inversion use exact rationals formed from W_up and the binary64 time and budget, so no intermediate rounds to zero or overflows. The chosen count is minimal for W_up, and a returned bound is rounded upward to binary64, which keeps a nonzero bound positive and the bound at the chosen count at most the budget.

Functions and records

select_trotter_step_count

select_trotter_step_count(hamiltonian: SparsePauliOp, *, time: float, error_budget: float, order: int = 1, bound_variant: str = 'exact_census', max_work: int = DEFAULT_CENSUS_MAX_WORK, max_bytes: int = DEFAULT_MAX_BYTES) -> TrotterStepSelection

Return the smallest Trotter step count whose error bound meets an error budget.

Applies the rule of Childs et al., doi:10.1103/PhysRevX.11.011020, Sec. V B. With r steps the total additive error is at most r * bound(time / r) (triangle inequality over the unitary steps), and the smallest admissible r follows in closed form because bound(t) = W_up * t**(order + 1):

  • order 1: total = W_up * t^2 / r, so r = max(1, ceil(t^2 * W_up / epsilon)).
  • order 2: total = W_up * t^3 / r^2, so r = max(1, ceil(sqrt(ceil(t^3 * W_up / epsilon)))).

W_up is the outward upper bound of the bound_variant expression (full exact_census by default, relaxed_prefix on request), so the count is minimal for W_up, not for the exact triangle coefficient. At order 2, let W_2 and W_2,rel be the full and relaxed Pauli-triangle coefficients before outward numerical evaluation, and let r_exact and r_rel be their smallest positive step counts for the same time and positive error allowance. If W_2 > 0, then r_exact <= r_rel <= ceil(alpha * r_exact) for alpha = sqrt(W_2,rel / W_2) in exact arithmetic. Before converting labels or building its index tables, the function checks the work and bytes of counting anticommuting Pauli pairs and triples with a nonzero nested commutator against max_work and max_bytes, and a requested expression that does not fit raises before the stage that would exceed them.

Parameters:

  • hamiltonian (SparsePauliOp) –

    H = sum_j c_j P_j in the term order of the product formula. Coefficients must be finite and real (an imaginary part of at most 1e-12 in absolute value is accepted), and identity terms are rejected, so remove them and keep their global phase exp(-i t c) yourself.

  • time (float) –

    Positive finite evolution time t.

  • error_budget (float) –

    Positive finite allowed operator-norm error epsilon of the whole evolution.

  • order (int, default: 1 ) –

    Default 1. Product-formula order, 1 (Lie-Trotter) or 2 (second-order Suzuki).

  • bound_variant (str, default: 'exact_census' ) –

    Default "exact_census", the full Pauli-triangle expression. At order 2, "relaxed_prefix" selects the degree-capped suffix relaxation, which needs pair tests only and can give a larger bound.

  • max_work (int, default: DEFAULT_CENSUS_MAX_WORK ) –

    Default 1_000_000_000. Work limit of label preparation, pair/triple tests and coefficient contraction.

  • max_bytes (int, default: DEFAULT_MAX_BYTES ) –

    Default 10 GB (decimal, 10_000_000_000 bytes). Byte limit of label conversion, index tables and bounded coefficient-reduction workspace.

Returns:

  • selection ( TrotterStepSelection ) –

    selection.step_count is the step count, and selection.bound_value is the certified total bound at that count (<= error_budget). selection.synthesis_method names the Qiskit product formula ("lie_trotter" or "suzuki_trotter").

Raises:

  • ValueError –

    If time or error_budget is not positive and finite, order is not 1 or 2, bound_variant is unknown, a coefficient is not finite or has an imaginary part above 1e-12, hamiltonian has an identity term, or as for trotter_bound_coefficient.

Examples:

For H = X + Z the second-order coefficient is ||[Z, [Z, X]]||/12 + ||[X, [X, Z]]||/24 = 1/2, so at t = 1 and epsilon = 0.01 the smallest r with (1/2)/r**2 <= 0.01 is 8:

>>> from qiskit.quantum_info import SparsePauliOp
>>> from nwqlib.subroutines.trotterization import (
...     select_trotter_step_count)
>>> H = SparsePauliOp(["X", "Z"], [1.0, 1.0])
>>> selection = select_trotter_step_count(
...     H, time=1.0, error_budget=0.01, order=2)
>>> print(selection.step_count, round(selection.bound_value, 10))
8 0.0078125

build_trotter_evolution_circuit

build_trotter_evolution_circuit(hamiltonian: SparsePauliOp, selection: TrotterStepSelection) -> QuantumCircuit

Build the product-formula circuit S_order(time / step_count)**step_count of a step-count selection.

Checks that hamiltonian has finite real coefficients and no identity term, as select_trotter_step_count requires, and builds the circuit with Qiskit's product-formula synthesis (LieTrotter or SuzukiTrotter, as selection.synthesis_method names) repeated step_count times. It does not check that hamiltonian is the operator selection was computed for, and selection.bound_value holds only for that operator in the same term order. The summands are applied in term order, matching the ordering the recorded bounds are stated for.

Parameters:

  • hamiltonian (SparsePauliOp) –

    The Hamiltonian the selection was computed for.

  • selection (TrotterStepSelection) –

    Output of select_trotter_step_count.

Returns:

  • circuit ( QuantumCircuit ) –

    The evolution circuit on hamiltonian.num_qubits qubits.

Raises:

  • ValueError –

    If a coefficient is not finite, has an imaginary part above 1e-12 in absolute value, or hamiltonian has an identity term.

trotter_bound_coefficient

trotter_bound_coefficient(hamiltonian: SparsePauliOp, *, order: int = 1, bound_variant: str = 'exact_census', max_work: int = DEFAULT_CENSUS_MAX_WORK, max_bytes: int = DEFAULT_MAX_BYTES) -> float

Return the time-independent prefactor W_up of the Trotter error bound W_up * t**(order + 1).

W_up is an outward upper bound on the Pauli-triangle expression chosen by bound_variant, which bounds the commutator sum of Childs et al., doi:10.1103/PhysRevX.11.011020, Prop. 9 (order=1) or Prop. 10 (order=2), so one product-formula application obeys bound(t) = W_up * t**(order + 1) and the Sec. V B r-step total is W_up * t**(order + 1) / r**order. The function evaluates the coefficient without dense matrices. W_up can exceed the exact triangle coefficient, and it is returned rounded upward to binary64. Before converting labels or building its index tables, the function checks the work and bytes of counting anticommuting Pauli pairs and triples with a nonzero nested commutator against max_work and max_bytes, and a requested expression that does not fit raises before the stage that would exceed them.

Parameters:

  • hamiltonian (SparsePauliOp) –

    H = sum_j c_j P_j in the term order of the product formula. Coefficients must be finite and real (an imaginary part of at most 1e-12 in absolute value is accepted), and identity terms are rejected, so remove them and keep their global phase exp(-i t c) yourself.

  • order (int, default: 1 ) –

    Default 1. Product-formula order, 1 (Lie-Trotter) or 2 (second-order Suzuki).

  • bound_variant (str, default: 'exact_census' ) –

    Default "exact_census", the full Pauli-triangle expression. At order 2, "relaxed_prefix" selects the degree-capped suffix relaxation, which needs pair tests only and can give a larger bound.

  • max_work (int, default: DEFAULT_CENSUS_MAX_WORK ) –

    Default 1_000_000_000. Work limit of label preparation, pair/triple tests and coefficient contraction.

  • max_bytes (int, default: DEFAULT_MAX_BYTES ) –

    Default 10 GB (decimal, 10_000_000_000 bytes). Byte limit of label conversion, index tables and bounded coefficient-reduction workspace.

Returns:

  • coefficient ( float ) –

    W_up, rounded upward to binary64.

Raises:

  • ValueError –

    W_up exceeds the largest binary64 number, a limit is not a positive integer, or the requested pair and triple tests do not fit max_work and max_bytes.

evaluate_trotter_bound

evaluate_trotter_bound(hamiltonian: SparsePauliOp, *, time: float, steps: int, order: int = 1, bound_variant: str = 'exact_census', max_work: int = DEFAULT_CENSUS_MAX_WORK, max_bytes: int = DEFAULT_MAX_BYTES) -> float

Return the error bound W_up * time**(order + 1) / steps**order of a product formula with a given step count.

For r=steps, the triangle-inequality composition of Childs et al. (doi:10.1103/PhysRevX.11.011020) gives r * bound(time / r) = W_up * time**(order + 1) / r**order with the outward coefficient W_up of the bound_variant expression. With the step count fixed, the relaxed expression changes only the certified bound, which can be larger. Before converting labels or building its index tables, the function checks the work and bytes of counting anticommuting Pauli pairs and triples with a nonzero nested commutator against max_work and max_bytes, and a requested expression that does not fit raises before the stage that would exceed them.

Parameters:

  • hamiltonian (SparsePauliOp) –

    H = sum_j c_j P_j with finite real coefficients and no identity term, in the term order of the product formula, as for select_trotter_step_count.

  • time (float) –

    Total evolution time.

  • steps (int) –

    Positive number of product-formula steps.

  • order (int, default: 1 ) –

    Default 1. Product-formula order, 1 (Lie-Trotter) or 2 (second-order Suzuki).

  • bound_variant (str, default: 'exact_census' ) –

    Default "exact_census", or "relaxed_prefix" for the second-order degree-capped suffix relaxation.

  • max_work (int, default: DEFAULT_CENSUS_MAX_WORK ) –

    Default 1_000_000_000. Work limit of label preparation, pair/triple tests and coefficient contraction.

  • max_bytes (int, default: DEFAULT_MAX_BYTES ) –

    Default 10 GB (decimal, 10_000_000_000 bytes). Byte limit of label conversion, index tables and bounded coefficient-reduction workspace.

Returns:

  • bound ( float ) –

    The total bound, rounded upward to binary64.

Raises:

  • ValueError –

    If time is not positive and finite, steps is not a positive integer, order is not 1 or 2, bound_variant is unknown, a coefficient is not finite or has an imaginary part above 1e-12, hamiltonian has an identity term, or as for trotter_bound_coefficient.

dense_trotter_bound_coefficient

dense_trotter_bound_coefficient(hamiltonian: SparsePauliOp, *, order: int = 1) -> float

Return the dense tail-sum coefficient of Props. 9-10 for validation on small instances.

This evaluates the spectral norms in Eqs. (120)-(121) from dense matrices, so its cost is exponential in the qubit count. The other functions of this module use the Pauli-triangle upper bound, which is never smaller.

TrotterStepSelection

TrotterStepSelection(*, formula_order: int, evolution_time: float, error_budget: float, bound_value: float, step_count: int, pauli_term_count: int = 0, pair_commutation_checks: int = 0, nested_commutation_checks: int = 0, bound_variant: str, common_tau: float | None = None, common_g: int | None = None, common_m: int | None = None, common_step_time: float | None = None, common_prefix_steps: tuple[tuple[int, int], ...] | None = None, common_power_allowances: tuple[tuple[int, float], ...] | None = None, common_recheck_outcome: str | None = None, common_prefix_bounds: tuple[tuple[int, float], ...] | None = None)

The step count chosen for one product-formula error budget, with its error bound.

select_trotter_step_count returns it. Its step_count field holds the result, and bound_value is its certified error bound. The fields below are read-only.

The fields follow the quantities of Childs, Su, Tran, Wiebe and Zhu (doi:10.1103/PhysRevX.11.011020). bound_value is the certified total additive Trotter error bound step_count * bound(evolution_time / step_count) at the chosen step count (at most error_budget), and the properties bound_formula and step_formula give the per-application bound and the smallest-r arithmetic as text with their sources.

The Pauli-triangle coefficient used for the bound is evaluated with scaled binary64 products and sums rounded upward, followed by exact rational rescaling. bound_variant identifies the expression that this coefficient bounds.

A common selection, made by select_common_step for QPE powers on a common time grid, uses one ordered product-formula step throughout the trajectory. At each requested nonnegative integer power p, the target time is p*val(common_tau), where val denotes the exact real value of a stored binary64 number. Its circuit prefix consists of the first r(p) steps, with r(p) recorded in common_prefix_steps. select_common_step accepts the selection only after the product-formula, pruning, time-displacement, identity-phase and rotation-angle bounds for every prefix fit that power's allowance, evaluated with the binary64 step and angles that the circuit uses. The common_* fields are None ("not recorded") on a selection for one evolution time. The maximum target time p_max*val(common_tau), where p_max is the largest requested power, and each represented prefix time r(p)*val(common_step_time) are computed from the fields and not recorded. On construction, the record computes g = common_g as the greatest common divisor of the positive powers, requires every power to be divisible by it and checks r(p) = (p/g)*m with m = common_m, the positive subdivision count. The all-zero schedule has no g, m or step and no common selection.

Attributes:

  • formula_order (int) –

    Lie-Trotter order 1 or symmetric Suzuki order 2. On a common selection: order of the shared ordered product-formula step. Common QPE trajectories use the symmetric second-order Suzuki formula.

  • evolution_time (float) –

    Total evolution time t > 0 of the product formula. On a common selection: binary64 display of the largest common-grid target time, whose exact value is defined by the maximum power and stored tau.

  • error_budget (float) –

    Requested operator-error allowance for this evolution. On a common selection: effective product-formula allowance for the full common trajectory after the other structural contributions are subtracted at every requested power. At each positive requested power p, let epsilon_p be its allowance in common_power_allowances, and let P_p, T_p, I_p and A_p be its pruning, time-displacement, identity-phase and rotation-angle error bounds, respectively. Let r_p be its cumulative step count in common_prefix_steps and r_max the largest of those counts. The allowance is min_(p>0) (epsilon_p-P_p-T_p-I_p-A_p)*r_max/r_p rounded downward to binary64.

  • bound_value (float) –

    Upward-rounded total operator-error bound at step_count, computed from the supplied upper coefficient and evolution_time. On a common selection: upward-rounded product-formula bound W_up*r_max*abs(val(h_hat))**3 for the common step that the circuit uses.

  • step_count (int) –

    Smallest positive integer allowed by the supplied upper coefficient and remaining error_budget. On a common selection: number r(p_max) of shared steps through the last requested power. This is a sufficient common-grid count and need not be the independently smallest count for that power.

  • pauli_term_count (int) –

    Number of nonidentity terms kept, after pruning, in this Hamiltonian or in the generator of this LCHS quadrature node.

  • pair_commutation_checks (int) –

    Number of new pair commutation tests that this selection performed, zero when the coefficient or the Pauli commutation table (anticommuting pairs and, for the full second-order expression, triples with a nonzero nested commutator) was reused.

  • nested_commutation_checks (int) –

    Number of new nested commutation tests that this selection performed, zero for relaxed_prefix or a reused Pauli commutation table.

  • bound_variant (str) –

    Pauli-triangle expression used, "exact_census" for the Pauli commutation table of anticommuting pairs and, at second order, triples with a nonzero nested commutator, or "relaxed_prefix" for the second-order degree-capped suffix relaxation. Order 1 uses "exact_census".

  • common_tau (float | None) –

    Finite positive binary64 time unit. The exact target time at integer power p is p times its represented value.

  • common_g (int | None) –

    Greatest common divisor of the positive requested integer powers. It identifies the exact target-time grid.

  • common_m (int | None) –

    Positive integer subdivisions of g times the represented time unit. The ideal common step is g*val(tau)/m.

  • common_step_time (float | None) –

    Binary64 step h_hat that the circuit uses for every shared step. Its represented value can differ from the ideal rational step.

  • common_prefix_steps (tuple[tuple[int, int], ...] | None) –

    Integer cumulative counts r(p)=(p/g)*m at each requested power, with r(0)=0, as (power, count) pairs. Execution between adjacent points uses the difference of these counts.

  • common_power_allowances (tuple[tuple[int, float], ...] | None) –

    Per-power structural operator-error allowances against which the bounds of the circuit's step were rechecked, as (power, allowance) pairs. The binary64 allowances are compared through their exact represented values.

  • common_recheck_outcome (str | None) –

    Outcome of comparing the complete subtotal of every circuit prefix with its allowance. An accepted common selection passed at every power. None means the recheck was not recorded.

  • common_prefix_bounds (tuple[tuple[int, float], ...] | None) –

    Upward-rounded values of the exact pruning, product-formula, time-displacement, identity-phase and rotation-angle subtotals for the circuit prefixes, as (power, bound) pairs.

bound_formula

bound_formula: str

Coefficient expression used, with its source theorem and inequality.

For a common selection, it names the composition of the step that the circuit uses.

step_formula

step_formula: str

Integer inversion of W_up*abs(t)**(order+1)/r**order against the remaining error budget.

synthesis_method

synthesis_method: str

Qiskit product formula used for the circuit: "lie_trotter" at order 1, "suzuki_trotter" at order 2.

bound_method

bound_method: str

Name of the bound, "pauli_triangle": commutator norms bounded by the Pauli triangle inequality.

bound_value_status

bound_value_status: str

Meaning of bound_value, "structural_upper_bound": an upper bound on the product-formula operator error.

to_dict

to_dict() -> dict[str, Any]

Return a JSON-like record of the selection.

Cost and limits

  • A coefficient or fixed-step bound above the largest binary64 number raises ValueError.
  • The keyword arguments max_work (default 1,000,000,000) and max_bytes (default 10,000,000,000) limit the tests for anticommuting Pauli pairs and for triples with a nonzero nested commutator. The pair stage is checked before label conversion. The requested expression is checked again with the counted number of anticommuting pairs before the nested tests or the contraction.
  • The full second-order expression stores three NumPy intp indices per triple with a nonzero nested commutator and can require cubic storage in the number of terms. The contraction block limits the coefficient-reduction workspace and does not limit the stored index count. If the requested expression does not fit, the function raises before the stage that would exceed the limits.
  • A work- or byte-limit refusal names the failed stage and reports its required work and the current max_work and max_bytes. A byte-fit failure is reported only when no contraction block size tried fits the current byte limit.
  • An order-two refusal of the full expression also names bound_variant="relaxed_prefix" as an explicit alternative, which is checked against the same limits.

Common time grid of QPE powers

A QPE trajectory on a common time grid uses one step for all power positions. At each requested power it checks the circuit prefix, meaning all steps up to that power position, against that power's allowance for the sum of the product-formula, pruning, time-displacement, identity-phase and rotation-angle error bounds. The QPE implementation map lists where controlled product-formula powers use these bounds.

  • val(tau) denotes the exact real value of the stored binary64 time unit, and a nonnegative integer power p targets the time p * val(tau).
  • If h is the exact real value of the binary64 step that the circuit uses and r is its cumulative count, the second-order product-formula bound is r * W * abs(h)**3. The discrepancy between r * h and the target contributes a separate operator bound.
  • A complete controlled-power bound also includes pruning and the applicable identity-phase and rotation-angle formation terms. The common_* fields of TrotterStepSelection record this choice.

Proposition numbering

Unless marked as arXiv:1912.08854v3, every proposition, equation and section number on this page and in the docstrings follows the Phys. Rev. X version, doi:10.1103/PhysRevX.11.011020. The arXiv:1912.08854v3 preprint, titled "A Theory of Trotter Error", numbers the same results Proposition 15/Eq. (145) (p. 38), Proposition 16/Eq. (152) (p. 39) and Sec. 5.2 (p. 40).