Quantum Krylov Algorithm

Estimate ground-state energies by building a small Krylov subspace from time-evolved reference states and diagonalizing the projected Hamiltonian classically. It trades a much shallower quantum circuit for a single small classical generalized eigensolve.

The problem it solves

One of the central tasks in quantum chemistry and materials science is finding the ground-state energy of a Hamiltonian H — its smallest eigenvalue. The variational quantum eigensolver (VQE) approaches this by parameterizing a deep ansatz circuit and repeatedly optimizing it to push the energy down; but deep circuits are easily swallowed by decoherence on noisy hardware, and the optimization often stalls on barren plateaus. The quantum Krylov method takes a different road: it does no variational optimization at all, and instead leans on the classical power of linear algebra, using only shallow circuits to collect a little information.

The intuition — spanning a good subspace with time evolution

Start from a reference state |ψ⟩ that is easy to prepare and has non-zero overlap with the true ground state. Apply the time-evolution operator e^−iHt for several different durations to obtain a family of states

{ |ψ⟩, e−iHt|ψ⟩, e−2iHt|ψ⟩, … }

These states span a Krylov subspace. The trick is this: the time-evolution operator shares its eigenstates with H, but stamps each eigencomponent with a different phase rotation set by its eigenenergy. By forming linear combinations of this family, we can constructively amplify the low-energy components and destructively suppress the high-energy ones. If the reference state has sufficient nonzero overlap with the ground state, and the evolution times and subspace dimension are chosen well, this Krylov subspace can effectively approximate the ground-state direction — even though every individual generating state may still be far from it. All that remains is to find the lowest-energy combination inside the subspace.

How it works — a generalized eigenproblem on the subspace

Minimizing the energy within the subspace reduces to solving a generalized eigenproblem

H c = E S c

Here there are two small matrices, with dimension equal to the subspace size (far smaller than the full Hilbert space): the Hamiltonian matrix H_ij = ⟨ψ_i|H|ψ_j⟩ and the overlap matrix S_ij = ⟨ψ_i|ψ_j⟩. The matrix S appears because the Krylov basis states are not orthogonal — they overlap one another, so this is not an ordinary eigenproblem but a generalized one carrying the metric S. The quantum computer is responsible for one thing only: estimating the entries of these two matrices, ⟨ψ_i|ψ_j⟩ and ⟨ψ_i|H|ψ_j⟩, with the Hadamard test. For a quantum-chemistry Hamiltonian, H is usually decomposed into a sum of Pauli strings H = Σℓ hℓ Pℓ. What is actually measured are the terms ⟨ψi|Pℓ|ψj⟩, then weighted by hℓ and summed to give Hij. Solving H c = E S c is then entirely a classical job, and its smallest eigenvalue E is the ground-state energy estimate.

An honest look at ill-conditioning and noise

• The overlap matrix S is often badly ill-conditioned. The Krylov basis states are highly similar to one another, which drives S toward singular with a large condition number. This means small errors in the entries of H and S get amplified into large errors in the energy estimate — the thorniest difficulty of the method in practice.
• Noisy matrix elements break the structure of the problem. The entries from the Hadamard test carry sampling noise and hardware error; when S is ill-conditioned, that noise can rob the generalized eigenproblem of physical meaning (even yielding an "energy" below the true ground state). A common remedy is thresholding: discarding the directions of S with eigenvalues that are too small, to stabilize the solve.
• Growing the subspace is a double-edged sword. Adding more time-evolved states can in principle improve accuracy, but it also worsens the ill-conditioning of S. The usable subspace dimension is often modest in practice.
• It still needs reliable time evolution. The whole method presupposes a faithful hardware implementation of e^−iHt; errors from its Trotter–Suzuki decomposition seep directly into the matrix elements and thus into the final energy.

The algorithm, step by step

1. Choose a reference state |ψ⟩ that is easy to prepare and overlaps the ground state.
2. Generate the Krylov basis {|ψ⟩, e^−iHt|ψ⟩, …} from a set of evolution times.
3. Estimate the matrix elements S_ij and H_ij with the Hadamard test.
4. Classically solve the generalized eigenproblem H c = E S c (thresholding S if needed to stabilize it).
5. Take the smallest eigenvalue E as the ground-state energy estimate.

Seeing it in code

Read the snippet below against the theory above: the double loop fills the overlap matrix S and Hamiltonian matrix H element by element, where qalgora_overlap estimates ⟨ψ_i|ψ_j⟩ with the Hadamard test and qalgora_matrix_element estimates ⟨ψ_i|H|ψ_j⟩ with H inserted; finally scipy.linalg.eigh(H, S, …) solves exactly the generalized eigenproblem H c = E S c, and energies[0] is the ground-state energy estimate.

import qalgora
import numpy as np

# overlap S_ij = <psi_i|psi_j>  and  H_ij = <psi_i|H|psi_j>
# estimated on hardware with the Hadamard test
dim = 4
S = np.zeros((dim, dim), dtype=complex)
H = np.zeros((dim, dim), dtype=complex)
for i in range(dim):
    for j in range(dim):
        S[i, j] = qalgora_overlap(i, j)        # Hadamard test
        H[i, j] = qalgora_matrix_element(i, j)  # Hadamard test with H

from scipy.linalg import eigh
energies = eigh(H, S, eigvals_only=True)        # generalized eigenproblem
print("ground-state estimate:", energies[0])