Source-linked AI summary
Observing ground-state properties of the Fermi-Hubbard model using a scalable algorithm on a quantum computer
Stasja Stanisic, Jan Lukas Bosse, Filippo Maria Gambetta, Raul A. Santos, Wojciech Mruczkiewicz, Thomas E. O'Brien, Eric Ostby, Ashley Montanaro
TL;DR
The paper addresses the lack of experimental demonstrations of Fermi-Hubbard ground-state preparation on quantum computers without compression. It uses Givens rotations for noninteracting-state preparation and relates BayesMGD parameter prediction to Kalman-filter state prediction; experiments show qualitative density and magnetic features, with limits in quantitative spin correlations and strong-interaction charge profiles.
Problem
To the authors’ knowledge, no experiment had demonstrated finding Fermi-Hubbard ground states on a quantum computer without introducing some form of compression.
Source: Appendix B · p. 11Method
The implementation prepares the noninteracting ground state with Givens rotations, and the paper relates BayesMGD’s prediction of optimal parameters to Kalman-filter state prediction.
Source: Appendix C · p. 12, §3 · p. 24Results
The VQE profiles show Friedel oscillations and antiferromagnetic correlations, although the total staggered spin correlation is not quantitatively reproduced.
Source: §IV · pp. 6–7Takeaways & Limitations
The 1D charge profiles suggest that a higher-depth variational ansatz is needed to capture strong-interaction effects fully.
Source: §IV · p. 6The VQE ground state follows the true ground state's staggered-spin-correlation trend but does not capture its value quantitatively.
Source: §IV · p. 7
Abstract
from arXiv · showhide
The famous, yet unsolved, Fermi-Hubbard model for strongly-correlated electronic systems is a prominent target for quantum computers. However, accurately representing the Fermi-Hubbard ground state for large instances may be beyond the reach of near-term quantum hardware. Here we show experimentally that an efficient, low-depth variational quantum algorithm with few parameters can reproduce important qualitative features of medium-size instances of the Fermi-Hubbard model. We address 1x8 and 2x4 instances on 16 qubits on a superconducting quantum processor, substantially larger than previous work based on less scalable compression techniques, and going beyond the family of 1D Fermi-Hubbard instances, which are solvable classically. Consistent with predictions for the ground state, we observe the onset of the metal-insulator transition and Friedel oscillations in 1D, and antiferromagnetic order in both 1D and 2D. We use a variety of error-mitigation techniques, including symmetries of the Fermi-Hubbard model and a recently developed technique tailored to simulating fermionic systems. We also introduce a new variational optimisation algorithm based on iterative Bayesian updates of a local surrogate model. Our scalable approach is a first step to using near-term quantum computers to determine low-energy properties of strongly-correlated electronic systems that cannot be solved exactly by classical computers.
I. VARIATIONAL ALGORITHM
The variational approach prepares a noninteracting ground state and applies low-parameter layers of Fermi-Hubbard time evolution. BayesMGD updates a local quadratic surrogate to guide noisy-measurement optimization.
I. VARIATIONAL ALGORITHM
- The Jordan-Wigner transform maps each fermionic mode to one qubit, but vertical hopping in 2D can produce multi-qubit terms.Snake ordering keeps horizontal terms on qubit pairs; a fermionic swap network rearranges modes so vertical interactions can use two-qubit gates.Source: §I · p. 2
- Each variational layer evolves sets of disjoint Hamiltonian terms for shared times, with 3, 4, or 5 parameters depending on lattice shape.The ansatz starts from the noninteracting ground state, prepared efficiently with Givens rotations.Source: §I · p. 2
- The ansatz’s small parameter count reduces optimization complexity, while its symmetries create opportunities for error mitigation.The approach uses the Hamiltonian variational ansatz within VQE.Source: §I · p. 2
- The Hamiltonian decomposition permits energy measurement with five distinct measurement settings, each requiring at most one additional two-qubit gate layer.The decomposition contains at most five sets of terms acting on disjoint modes.Source: §I · p. 3
- BayesMGD iteratively updates a Bayesian quadratic local surrogate from nearby noisy measurements, then uses its gradient for a descent step.The routine is introduced to optimize circuit parameters despite statistical and hardware measurement noise.Source: §I · p. 3
II. QUANTUM CIRCUIT IMPLEMENTATION
The experiments implement 1×8 and 2×4 Fermi-Hubbard instances on the Rainbow superconducting processor. Fermionic swap networks and a zig-zag layout accommodate the gate connectivity and reach up to eight sites.
II. QUANTUM CIRCUIT IMPLEMENTATION
- Each lattice site is encoded by two qubits, and the experiments use the Rainbow superconducting processor in Google Quantum AI’s Sycamore architecture.The processor configuration had 23 available qubits.Source: §II · p. 3
- The 1×8 zig-zag layout uses two nearby qubit lines to reach eight sites, with swaps needed to implement nonlocal onsite interactions.A more direct local mapping on this hardware would support only up to 1×4.Source: §II · p. 3
- For 2×Ly, a fermionic swap network converts otherwise four-qubit vertical hopping operations into two-qubit gates, at the cost of one extra swap layer per ansatz layer.The network constrains the order of vertical interactions.Source: §II · p. 3
- Near half-filling, one-layer circuits require at most 26 two-qubit gate depth and 140 gates for 1×8, or 32 depth and 176 gates for 2×4.The hardware configuration had 23 available qubits.Source: §II · p. 3
III. ERROR MITIGATION
The error-mitigation strategy combines hardware-oriented circuit changes, model symmetries, and a fermionic-system training method. These procedures use error detection, symmetry averaging, and classically tractable training data to improve estimates.
III. ERROR MITIGATION
- Alternating single- and two-qubit layers and spin-echo X gates substantially reduced experimental error.The authors attribute the reduction to possible inversion of accumulated single-qubit phase errors and modification of parasitic CPHASE errors.Source: §III · p. 4
- Number conservation detects errors by discarding runs whose final and initial occupations differ; time reversal, particle-hole, and lattice symmetries enable averaging with transformed partners.The methods exploit symmetries of the Fermi-Hubbard Hamiltonian.Source: §III · p. 4
- Training with Fermionic Linear Optics (TFLO) maps approximate quantum-computer energies to exact energies using paired training data from classically simulable fermionic linear-optics circuits.The training circuits should approximate the real algorithm circuits for the map to be accurate.Source: §III · p. 4
- For this model, TFLO training circuits can vary hopping parameters while setting onsite-term parameters to zero, because most VQE operations are fermionic linear optics.Onsite-term evolution is the circuit component that is not fermionic linear optics.Source: §III · p. 4
IV. RESULTS
The experiments reproduce several qualitative features expected of Fermi-Hubbard ground states, including energy behavior, charge-density patterns, and antiferromagnetic correlations. Agreement with the true ground state is not uniformly quantitative, and the metal-insulator-transition signal is unclear for 2×4.
IV. RESULTS
- BayesMGD converges to parameters for the VQE ground state within a small number of iterations.The resulting circuit is then used to measure properties of the approximate ground state.Source: §IV · p. 5
- Across occupations, VQE energies for 1×8, 2×4, and 1×4 agree quantitatively with the lowest energies achievable by the one-layer ansatz.For 1×4, a depth-2 ansatz reaches lower energy than is theoretically possible at depth 1.Source: §IV · p. 5
- In 1×8, the half-filling Mott-gap signal grows with U, while charge correlations decay most steeply at half-filling and more slowly away from it.The slower decay away from half-filling is described as a signature of increased conductivity; the 2×4 results are suggestive but do not clearly show the transition.Source: §IV · pp. 5–6
- The 1D charge profiles show Friedel oscillation peak counts, while Wigner oscillations are visible especially at Nocc≤4 and suggest that stronger interactions need greater ansatz depth.For even occupations the ground-state spin density is zero; the spin density is more substantial at odd occupations.Source: §IV · p. 6, Figure 5 · p. 6
- At half-filling with U=4, the 1×8 and 2×4 VQE states show antiferromagnetic correlations, stronger in 2D than 1D.The total staggered spin correlation is not captured quantitatively, though its trend matches the true ground state; Nocc=7 is less antiferromagnetic than half-filling for U=4 and 8.Source: §IV · p. 7
V. DISCUSSION
The authors extract qualitative features of medium-size Fermi-Hubbard instances with a low-depth variational ansatz, using four times the qubits of previously reported experiments. They report relatively high energy accuracy with one variational layer, while noting that quantitatively accurate true ground-state energies will require more layers as systems grow.
V. DISCUSSION
- A low-depth variational ansatz extracts fundamental qualitative features of medium-size Fermi-Hubbard instances using four times more qubits than previously reported experiments.The authors report relatively high energy accuracy for states produced with one variational layer.Source: §V · p. 7
- Higher-depth ansätze in larger systems are expected to demonstrate additional phenomena, while quantitatively accurate true ground-state energies are expected to require layers scaling with system size.The proposed future phenomena include Wigner oscillations, charge-density-wave ordering, and magnetic instabilities.Source: §V · p. 7
- The efficient algorithm and error-mitigation techniques provide a template that can be scaled as quantum-computing hardware improves.Source: §V · p. 7
SUPPLEMENTARY MATERIAL
The Fermi-Hubbard model describes interacting electrons in a narrow band and is used to study strongly correlated effects. Its phases depend on density and interaction strength, with higher-dimensional transition locations still unresolved.
Appendix A: Physics of the Fermi-Hubbard model
- The Fermi-Hubbard model captures inter-electron interactions in a narrow band and is a candidate for studying strongly correlated effects in solids.Its phenomenology has been linked to phases of cuprates, though whether the model supports high-temperature superconductivity remains open.Source: Appendix A · p. 11
- In the atomic limit of a one-band model without spin-orbit coupling, the tight-binding description reduces to the Fermi-Hubbard Hamiltonian with homogeneous nearest-neighbour hopping measured in units of t.The general description includes kinetic and ion-potential contributions in t and electron-electron interactions in U.Source: Appendix A · p. 11
- In one dimension, the model is Bethe-ansatz solvable and is a Mott insulator at half-filling for U > 0, while it is a Luttinger liquid away from half-filling.The model exhibits spin-charge separation in 1D.Source: Appendix A · p. 11
- For dimensions above one, transitions from metallic to Mott-insulating behavior and back to metallic behavior are expected as density changes, but their locations remain unresolved.The nature of the Mott transition is debated, with proposed contributions from quasiparticle localization and magnetic instabilities.Source: Appendix A · p. 11
Appendix B: Previous experimental implementations of quantum algorithms for the Fermi-Hubbard model and related problems
Prior Fermi-Hubbard experiments found ground states only in small, symmetry-compressed instances, while related work studied dynamics, low-depth ansätze, and other variational applications. The supplied comparisons distinguish those efforts by system size, compression, hardware, and task.
Appendix B: Previous experimental implementations of quantum algorithms for the Fermi-Hubbard model and related problems
- Earlier ground-state experiments reached up to four sites using up to four qubits, with symmetry-based compression limiting their near-term scalability.Examples include a 1×2 ion-trap experiment using five qubits for two ground-state copies and entanglement measurement [16], and compressed 1×2 and four-site VQE demonstrations [17, 18].Source: Appendix B · p. 11
- A separate 1×8 experiment simulated Fermi-Hubbard dynamics and demonstrated spin-charge separation, with meaningful results at two-qubit gate depth over 400.Unlike the present work, it did not optimise the ground state with VQE and did not study a 2×Ly system.Source: Appendix B · p. 11
- Related VQE experiments addressed small molecules and a three-qubit H3 metal-insulator transition, while QAOA targeted classical bit-string solutions and had reached 23 qubits in the cited work.Hartree-Fock preparation experiments on up to 12 qubits benchmarked hardware performance but addressed classically simulable states [40].Source: Appendix B · p. 12
- Independent low-depth ansatz work examined 1D eight-site chains by classical simulation, while the present paper identifies its own demonstration as the first uncompressed quantum-computer ground-state experiment for the model.The independent study examined spin correlations, double occupancies, energy, and ground-state fidelity [15].Source: Appendix B · pp. 11–12
Appendix C: Implementation of the Efficient Hamiltonian Variational ansatz
The implementation uses an Efficient Hamiltonian Variational ansatz built from Givens rotations and Hamiltonian time-evolution gates, with fermionic swaps supporting 2×Ly lattices. Its BayesMGD optimiser updates a local quadratic surrogate from noisy samples, and experiments separate optimisation from state preparation.
Appendix C: Implementation of the Efficient Hamiltonian Variational ansatz
- The Efficient Hamiltonian Variational ansatz prepares the noninteracting ground state with Givens rotations, then applies onsite and hopping time-evolution layers with variational times.For 2×Ly lattices, fermionic swap networks rearrange the Jordan-Wigner ordering so vertical interactions can use two-qubit gates; for the special 2×Ly case the swap-network overhead can be one layer.Source: Appendix C · pp. 12–13
- The circuit uses native-gate decompositions for its building blocks, and energy measurements use separate circuits for onsite and hopping terms, with some measurement basis changes merged into the final layer.Initial-state preparation has depth N − 1; circuit depth per ansatz layer can be as low as 2Lx + 1 for even Lx under the stated topology assumptions.Source: Appendix C · pp. 12–13
Appendix D: Variational optimiser
- BayesMGD updates a quadratic local surrogate using noisy samples and Bayesian parameter updates, then uses the surrogate's estimated gradient for a descent step.The method uses sample variance to estimate uncertainty in fitted parameters and model evaluations; its uncertainty growth after a parameter step is described as heuristic.Source: Appendix D · pp. 14–15
- The optimiser assumes finite variance for f(θ), enabling the central limit theorem for sample means, and a smooth mean objective µ(θ) locally described by a quadratic surrogate.VQE runs took under 30 minutes for 1 × 8 instances and under 70 minutes for 2 × 4 instances.Source: Appendix D · p. 14, §1 · p. 16
1. Details of implementation parameters
- The VQE optimisation and state-preparation stages can be performed at different times; the state-preparation stage uses all described error-mitigation techniques, while VQE omits particle-hole symmetry and TFLO.State-preparation energy estimates used 100,000 measurements per evaluation, and the procedure was repeated three times per instance.Source: §1 · pp. 15–16
Appendix E: Error mitigation
The appendix describes hardware-level, symmetry-based, and fermionic error-mitigation techniques. The symmetry methods use conserved occupations and transformations of the model to reduce experimental errors.
1. Low-level error mitigation
- The experiments used alternating single- and two-qubit gate layers, selected a qubit set using initial energy measurements, and chose the lowest-energy run from three trials.Averaging over qubit sets was not advantageous in these experiments.Source: §1 · p. 16
- Spin echo substantially reduced experimental errors, possibly by reversing phase errors or moving parasitic CPHASE errors to |00⟩.The technique sandwiches alternating two-qubit-gate layers with X gates.Source: §1 · p. 16, Figure 8 · p. 17
2. Symmetries of the Fermi-Hubbard model
- Time reversal, particle-hole, and lattice-reflection symmetries provide additional ways to average energies or observables and reduce worst-case errors.The particle-hole method averages suitably transformed observables at occupations Nocc and L − Nocc; reflection averaging uses two 1D or four 2D points.Source: §2 · pp. 17–18
- Occupation-number checks reject runs with incorrect spin-sector populations; assuming independent readout errors, the estimated undetected-error probability is below 6%.The experiments retained 7%–29% of runs after occupation-number postselection.Source: §2 · p. 17
- Particle-number conservation within each spin sector enables postselection because the variational ansatz preserves those occupations.The measured occupation number is compared with its known initial value.Source: §2 · p. 17
3. Training with fermionic linear optics
Training with fermionic linear optics uses efficiently simulable circuits to calibrate noisy VQE energies and observables. The method relies on the observed approximate linearity of noisy and exact energies and includes checks for cases where a fit is unsuitable.
3. Training with fermionic linear optics
- FLO circuits provide exact classical and approximate quantum expectations that serve as training pairs for mapping noisy VQE energies to approximations of exact energies.For the Hamiltonian variational ansatz, setting onsite parameters to zero yields FLO circuits that can be simulated efficiently.Source: §3 · p. 18
- Noisy and exact energies are approximately linearly related in practice, supporting a linear fit for energy correction.The authors use the Theil-Sen estimator to reduce the influence of outliers.Source: §3 · p. 19, Figure 9 · p. 19
- The training parameters are chosen to spread their exact energies, using a 16-point grid for one layer and 256 random choices for two or more layers.The grid-search cost grows exponentially with the number of layers.Source: §3 · p. 19
- A final correction subtracts residual error measured at the nearest FLO point, addressing coherent errors that depend on the circuit parameters.The nearby FLO point is obtained by setting the onsite parameter or parameters to zero.Source: §3 · p. 19
- For other observables, the method skips linear fitting when exact values are nearly constant and uses a fit only when it is sufficiently good; otherwise it returns the noisy value.The fit is accepted when R2 exceeds 0.7; observables invariant under FLO but changed by non-FLO gates are a failure case.Source: §3 · pp. 19–20
Appendix F: Classical comparators
The appendix describes classical methods used to construct simulated and ground-state comparisons, alongside Slater-determinant comparators. It also identifies the quantities shown in additional 2 × 4 results.
Appendix F: Classical comparators
- Ground-state energies, charges, spins, and correlations were computed by exact diagonalisation, while simulated VQE results used exact classical circuit simulation.The simulated VQE optimisation used BFGS, with results checked using Cobyla because BFGS can become trapped in local minima.Source: Appendix F · p. 20
- Slater determinants provide a classical comparator; the noninteracting U = 0 ground state is one trial state, while mean-field can yield a closer ground-state energy.The authors also classically optimised over Slater determinants to assess the best possible determinant comparator.Source: Appendix F · p. 20
- The Slater-determinant search parametrizes a single-particle Hamiltonian with an L × L Hermitian matrix, shared across spin-up and spin-down modes.The determinant is formed from the lowest-eigenvalue eigenvectors for the target occupation.Source: Appendix F · p. 20
- Figure 10 adds 2 × 4 chemical-potential results, an inset showing its derivative at even occupations, and total staggered spin correlation.The appendix also reports that the Slater-Condon rules and classical BFGS optimisation reliably converged to good solutions.Source: Appendix F · p. 21
Appendix G: Error analysis
The appendix estimates uncertainty for postselected observables and propagates it through TFLO using Monte Carlo sampling. It also distinguishes BayesMGD's internal certainty estimates from true measurement error bars.
Appendix G: Error analysis
- For postselected observables, the sample-mean variance accounts for the probability of retaining a run, assuming Gaussian measurements conditioned on correct spin-sector occupations.The stated approximation is σ^2/(pN)(1+(1−p)/(pN)+O(1/(pN)^2)).Source: Appendix G · p. 21
- TFLO error bars are estimated by sampling raw observables from Gaussian distributions 1,000 times and rerunning the TFLO procedure.The reported error bar is the sample standard deviation of those estimates.Source: Appendix G · p. 21
- For spin and charge correlations, the error-bar calculation assumes that the component observables are independent.Source: Appendix G · p. 21
- BayesMGD's reported energy error bars are internal certainty estimates and may not equal true energy-measurement error bars.Source: Appendix G · p. 21
Appendix H: Additional results for the 2 × 4 lattice
The 2 × 4 results show no onset of the metal-insulator transition, while charge and spin profiles retain features associated with one-dimensional physics. The appendix also details BayesMGD’s Bayesian updates, iterative optimization procedure, and mathematical equivalence to a Kalman filter.
Appendix H: Additional results for the 2 × 4 lattice
- The 2 × 4 results show no onset of the metal-insulator transition, with the authors citing processor errors and insufficient ansatz depth as possible causes.They suggest both may limit the circuit’s ability to capture correlations in this system.Source: Appendix H · p. 22
- At quarter filling, rung delocalisation explains the uniform charge profile and weakens the effect of onsite interactions.The explanation starts from two degenerate 1 × 4 bands at zero rung hopping; hopping splits them, and the lower band forms symmetric rung combinations.Source: Appendix H · p. 22
- For odd occupations, VQE and ground-state spin profiles are similar and near non-interacting; the authors infer that deeper circuits are needed to capture correlations.The ground-state spin is zero everywhere for even occupations.Source: Appendix H · p. 22, Figure 11 · p. 23
1. The Bayesian update step
- BayesMGD assumes a locally Gaussian model-parameter belief and uses new measurements to update the surrogate-model posterior.The appendix specifies a linear model in its parameters, a Gaussian prior, and an increase in parameter uncertainty after a gradient step.Source: §1 · pp. 22–23
2. Pseudocode for BayesMGD
- BayesMGD repeatedly samples near the current parameters, updates its surrogate model, and takes a gradient-descent step.Its pseudocode returns an estimate and uncertainty when the step is small or the evaluation limit is reached.Source: §2 · p. 24
3. Relation between BayesMGD and Kalman filters
- BayesMGD’s prediction and measurement-update equations are mathematically equivalent to Kalman-filter equations, with model parameters playing the role of the estimated state.In this application, the dynamics are trivial, there is no control input, and measurement outcomes, noise, and the design matrix correspond to the filter’s observation inputs.Source: §3 · p. 24
4. Experimental comparison of optimisation algorithms
On real hardware, BayesMGD and MGD outperform SPSA when the optimization budget is measured in shots, but their convergence rates are roughly similar under a wall-clock-time budget. In simulations, BayesMGD often outperforms MGD when quadratic fits are underdetermined.
4. Experimental comparison of optimisation algorithms
- The hardware comparison used well-understood problem instances and the same optimizer hyperparameter settings for BayesMGD and MGD.SPSA used a one-stage implementation with a reduced stability constant.Source: §4 · pp. 24–25
- With a shot budget, BayesMGD and MGD clearly outperform SPSA; with a wall-clock-time budget, their convergence rates are roughly similar.The comparison equalized shots per optimizer iteration; the figure shows energy difference from the best energy attainable at the given circuit depth.Source: §4 · p. 25, Figure 12 · p. 25
- For harder instances with six rather than four free parameters, (Bayes)MGD may stop improving after nearing the minimum.One SPSA run on the 2 × 4 lattice also reached parameters near the global minimum early and did not improve further.Source: §4 · p. 25
- In 12-qubit, 18-parameter kagome-Heisenberg simulations, BayesMGD often outperforms MGD for η < 1, while their performance is comparable for η ≥ 1.η is the ratio of evaluation points per iteration to the number needed for a fully determined quadratic fit; the figure compares simulated costs at selected η values.Source: §4 · p. 25, Figure 13 · p. 26