Research
Lindblad Multiproduct Formulas
Overview Research area: Quantum error mitigation and classical simulation of quantum circuits, specifically combining two-dimensional tensor network methods with noisy quantum hardware to estimate noi

- arXiv
- 2609.05024
- Published
- 2026-09-04
- Authors
- Niall F. Robertson, Andrea D'Urbano, Anton Dekusar, Max Rossmannek, Eric D. Switzer, James R. Garrison, Nicolas Lorente, Igor Pasichnyk, Paolo D'Alberto, Muhammad Osama, Zachary Streeter, Yasuko Eckert, Constantinos Evangelinos, Sergiy Zhuk
AI summary
Overview
Research area: Quantum error mitigation and classical simulation of quantum circuits, specifically combining two-dimensional tensor network methods with noisy quantum hardware to estimate noiseless expectation values.
Technical level: Advanced. The paper assumes familiarity with Lindblad master equations, Pauli-Lindblad noise models, probabilistic error cancellation/amplification, zero noise extrapolation, tensor network states and operators, belief propagation, and discrete time crystal physics.
Scope (one sentence): The paper introduces Lindblad Multiproduct Formulas (L-MPF), an error mitigation workflow in which coefficients for a linear combination of noisy expectation values are computed with two-dimensional tensor networks contracted via loop-corrected belief propagation, and demonstrates it on a 65-qubit two-dimensional discrete time crystal experiment on ibm_basquecountry.
What This Paper Is About
Estimating noiseless expectation values from noisy quantum hardware usually relies on methods such as probabilistic error cancellation (PEC) or zero noise extrapolation (ZNE), but PEC's sampling cost grows quickly with circuit size and depth, and ZNE suffers from bias in fitting expectation values against the noise rescaling parameter. This paper asks whether state-of-the-art two-dimensional tensor network techniques can be embedded directly into an error mitigation workflow, so that advances in classical simulation translate into more accurate and reliable mitigated quantum results. The goal is to produce an estimate of the noiseless observable that is more accurate than either a purely classical tensor network simulation or the noisy quantum computer alone.
Key Contributions
-
The L-MPF algorithm itself. A scheme that estimates a noiseless observable as a weighted sum of noisy expectation values, ∑ᵢ cᵢ⟨O⟩{λᵢ}, where the coefficients cᵢ are computed classically by minimizing a cost function built from overlaps Fᵢⱼ = ⟨⟨O|Λᵀ{λᵢ}Λ_{λⱼ}|O⟩⟩, subject to the constraint ∑ᵢ cᵢ = d fixed by Clifford rescaling.
-
A contraction-path-free use of two-dimensional tensor networks. Rather than storing Λᵀ_{λᵢ}Λ_{λⱼ} as a matrix product operator via a "middle-out" Heisenberg-picture contraction path (as in prior tensor-network multiproduct work), the authors implement Λ_{λⱼ} as tensor network operators, store the resulting states as two-dimensional tensor network states, and evaluate the overlaps with loop-corrected belief propagation, using bond dimension extrapolation to estimate the infinite-bond-dimension limit.
-
Ergodic amplification. A new amplification method that perturbs parameters of the Hamiltonian H (for example, the spin-flip coupling ε of the time crystal model) rather than the dissipative noise rates λₐ, producing progressively dampened expectation values with significantly reduced sampling cost compared to probabilistic error amplification (PEA). For k = 3 the authors typically choose γ₁ = 1, γ₂ = 1.4, γ₃ = 1.8, or slight variants.
-
An error-propagation analysis and error bars. An argument, formalized in an expansion of the L-MPF estimate, showing why tensor network truncation error in the coefficients need not degrade the mitigated observable as much as a purely classical workflow would, plus incorporation of Clifford rescaling and a heuristic error bar into the workflow.
Main Findings
-
L-MPF outperforms both ingredients alone on the tested model: Applied to the two-dimensional discrete time crystal studied previously, the workflow produces more accurate estimates of the noiseless observable than either tensor networks or the noisy quantum computer alone.
-
Resolution of the apparent paradox: The quantities needed to implement error mitigation (the overlaps Fᵢⱼ) can be less computationally expensive to evaluate with tensor networks than the expectation values themselves, so a problem that is hard to simulate classically does not automatically inherit an equally large truncation error in the mitigated result. This statement is claimed to be independent of whether a middle-out contraction path is used.
-
Errors enter as products: In the reported error expansion, the tensor network errors σ_{cᵢ} on the coefficients appear only multiplied by the hardware noise errors ε^O_{λᵢ} on the expectation values, so for sufficiently small values of each the product sum can be smaller than either quantity in isolation.
-
Poor conditioning is largely harmless here: Although tensor network errors on the coefficients can be magnified when the matrix of Fᵢⱼ entries for i, j ≥ 1 is approximately singular, the direction of magnification is a vector with components εᵢ satisfying ∑εᵢO_{λᵢ} ≈ 0, so the effect on the final weighted sum is negligible.
-
Projection error is not the bottleneck: The residual projection error ε_proj (the error even with an exponentially large bond dimension) is argued in an appendix to be smaller than the average ZNE error.
-
Ergodic amplification is well suited to this model: The authors report that it works particularly well for the time crystal studied, and that the classical cost of computing the coefficients does not depend on whether noise or ergodicity is amplified. Results obtained with PEA are reported in an appendix rather than the main text.
-
Hardware scale: The quantum experiments were run on 65 qubits arranged in a 3×3 heavy-hexagonal topology on
ibm_basquecountry. -
GPU acceleration: A GPU implementation of the classical part of the workflow achieves a speedup of up to 5.6×.
-
Not reported in the available content: Specific numerical accuracy values for the mitigated observables (for example, absolute errors or improvements over baselines), the bond dimensions χ used, the number of Floquet cycles, and the exact values of ε studied are not given in the text available here.
Methodology in Plain English
The workflow has three stages.
-
Learn the noise and measure. The noise on
ibm_basquecountryis characterized first, following established Pauli-Lindblad noise learning, producing a set of noise rates λₐ and corresponding Pauli strings. The device is then used to measure noisy expectation values of the observable at several settings. Those settings come either from probabilistic error amplification, which deliberately scales up the noise, or from the new ergodic amplification, which instead nudges the circuit's Hamiltonian parameters (such as the spin-flip coupling ε) so that the observable is progressively dampened — qualitatively similar to adding noise but without extra sampling overhead. -
Solve for the coefficients classically. On CPUs/GPUs, the overlaps between the noisy and noiseless evolved observables are computed as two-dimensional tensor networks, using belief propagation with loop corrections, and a classical optimizer finds the coefficients that minimize the mismatch to the noiseless observable, subject to a normalization fixed by a Clifford-circuit rescaling benchmark. Truncation error from finite bond dimension is handled by repeating at increasing bond dimension and extrapolating.
-
Combine. The hardware measurements are multiplied by the classically optimized coefficients and summed to give the mitigated estimate, along with a heuristic error bar.
The physics testbed is a periodically driven, many-body-localized spin model — a discrete time crystal — defined by a Floquet operator composed of three XXZ coupling layers plus a global X-rotation "kick" of angle φ, with disordered couplings Jᵢⱼ = 1 + δᵢⱼ where δᵢⱼ is sampled uniformly from [−0.5, 0.5]. Accuracy is judged on the order parameter Δ(t), the average over sites of the two-time Zᵢ correlation, which is a sum of single-qubit observables.
Why This Matters
Impact on research. The paper reframes error mitigation as a co-design problem between classical tensor network methods and quantum hardware, rather than a purely sampling-based procedure. It shows that the classical objects needed for mitigation need not be as hard to compute as the observable itself, which suggests a route to applying mitigation on problems where classical methods struggle to give the answer directly. It also extends tensor-network-enhanced multiproduct ideas from one-dimensional circuits and Trotter-error reduction to two-dimensional circuits and device noise, and it comes with an error bar — something that practical rescaling approaches typically lack.
Real-world applications (as motivated in the paper):
- Many-body physics and mapping phase diagrams of strongly correlated and periodically driven systems, including discrete time crystals.
- Simulation of photonic systems.
- Solving partial differential equations via quantum algorithms.
- Benchmarking and validating noisy quantum hardware on problems whose phase structure is already understood, such as the two-dimensional time crystal used here.
Industry relevance. The work is a collaboration between IBM Research and AMD, spanning quantum hardware (ibm_basquecountry), error mitigation software, and GPU acceleration, with a reported GPU speedup of up to 5.6× on the classical portion. That coupling of quantum devices with accelerated classical tensor network kernels is directly relevant to vendors building hybrid quantum-HPC stacks, since error mitigation accuracy and classical runtime are both practical bottlenecks for useful quantum simulation.
Future Directions
- Generalizing ergodic amplification beyond this model. The authors explicitly leave the general study of how to apply ergodic amplification to other classes of quantum circuits for future work, noting that any parameter-tuning heuristic will change the resulting error bars, and that poorly separated amplified circuits make the overlap matrix near-singular and the error bars impractical.
- Better conditioning heuristics. Because the choice of amplification strategy strongly affects how well-conditioned the coefficient problem is, selecting or designing heuristics per problem class — or falling back to PEA when no good heuristic exists — remains open.
- Revisiting structured contraction paths. The authors flag exploring contraction paths such as storing Λᵀ_{λᵢ}Λ_{λⱼ} as a matrix product operator or general tensor network operator, noting these can stay close to the identity even for long evolution times but cost a bond dimension quadratic in that of the state, implying large memory overhead in two dimensions and beyond.
- Scaling and comparison. Larger lattices, longer Floquet evolution, and systematic comparison of the ergodic-amplification results in the main text against the PEA results collected in the appendix would clarify when each amplification strategy is preferable.
Target Audience
Quantum error mitigation researchers and practitioners who already work with PEC, PEA, or ZNE and want a mitigation scheme that trades sampling cost for classical tensor network computation; tensor network developers interested in two-dimensional loop-corrected belief propagation applied to noisy quantum data; quantum simulation researchers studying discrete time crystals and many-body localization; and hybrid quantum-HPC engineers at hardware and accelerator vendors who care about GPU-accelerated classical components in quantum workflows. Readers without a background in open quantum dynamics and tensor networks will find the paper difficult, though the high-level workflow and the error-propagation argument are accessible with effort.
Authors’ abstract
We introduce Lindblad Multiproduct Formulas: a quantum error mitigation technique that uses two-dimensional tensor networks contracted with loop-corrected belief propagation. The quantities required to implement the error mitigation scheme that are evaluated with tensor networks can be less computationally expensive to calculate than the expectation values themselves, thus allowing for the possibility of applying our method to certain systems for which tensor network methods may struggle to calculate the observable quantities of interest. The workflow incorporates Clifford rescaling techniques and outputs an estimated error bar. We apply our method to a model of two-dimensional discrete time crystals studied previously and implement it on $65$ qubits arranged in a $3\!\times\!3$ heavy-hexagonal topology on the quantum computer ibm_basquecountry. We show that a GPU implementation of the classical part of our workflow achieves a speedup of up to $5.6\times$.