Quantum-Inspired Algorithm for Classical Spin Hamiltonians Based on Matrix Product Operators

We propose a tensor-network (TN) approach for solving classical optimization problems that is inspired by spectral filtering and sampling on quantum states. We first shift and scale an Ising Hamiltonian of the cost function so that all eigenvalues become non-negative and the ground states correspond to the the largest eigenvalues, which are then amplified by power iteration. We represent the transformed Hamiltonian as a matrix product operator (MPO) and form an immense power of this object via truncated MPO-MPO contractions, embedding the resulting operator into a matrix product state for sampling in the computational basis. In contrast to the density-matrix renormalization group, our approach provides a straightforward route to systematic improvement by increasing the bond dimension and is better at avoiding local minima. We also study the performance of this power method in the context of a higher-order Ising Hamiltonian on a heavy-hexagonal lattice, making a comparison with simulated annealing. These results highlight the potential of quantum-inspired algorithms for solving optimization problems and provide a baseline for assessing and developing quantum algorithms.

💡 Research Summary

The paper introduces a quantum‑inspired tensor‑network algorithm for solving classical spin‑glass optimization problems, particularly Ising‑type Hamiltonians. The authors start by encoding a combinatorial cost function C(z) as a diagonal Hamiltonian Ĥ acting on a real Hilbert space spanned by computational basis states |z⟩. Because the ground‑state configurations of interest correspond to the smallest eigenvalues of Ĥ, they apply a linear shift and scaling, defining a new operator Ĝ = aĤ + b I with a = −1 and b = Λ, where Λ is chosen larger than the spectral radius of Ĥ. This transformation makes Ĝ positive semidefinite and flips the spectrum so that the original ground states become the largest eigenvalues of Ĝ.

Ĝ is then represented as a matrix‑product operator (MPO). Each site i contributes a local tensor W_i with virtual indices (α_{i‑1},α_i); the MPO bond dimension R is bounded by the number of interaction terms and the operator‑Schmidt rank, but in practice it is much smaller because of locality and sparsity (e.g., for QUBO problems R ≤ min(⌊N/2⌋+2, r+2), where r is the rank of the coupling matrix).

The core of the method is a power‑iteration (spectral filtering) performed directly on the MPO. The authors construct Ĝ^K by repeatedly multiplying the MPO with itself. Two contraction schedules are discussed: a “linear” schedule that applies Ĝ K times sequentially (cost O(N d² R χ³ K)) and a “doubling” schedule that squares the MPO repeatedly, achieving K = 2^m with only m multiplications (higher χ scaling but logarithmic multiplication count). After each multiplication, a full singular‑value decomposition (SVD) truncates the bond dimension to a user‑defined maximum χ, discarding singular values below a fixed threshold ε = 10⁻¹⁵. The resulting powered MPO remains diagonal in the computational basis.

Because the powered MPO is diagonal, it can be embedded into a matrix‑product state (MPS) simply by acting on the uniform product state |+⟩⊗N, where |+⟩ = (1/√d)∑_{z_i}|z_i⟩. The resulting state |Ψ⟩ ∝ Ĝ^K|+⟩⊗N has amplitudes proportional to (Λ − C(z))^K. Sampling this MPS in the computational basis yields configurations z with probabilities proportional to the square of those amplitudes, i.e., p(z) ∝