Generating-Function Matrix Product States

• Generating-Function Matrix Product States is a framework that encodes quantum many-body states using generating functions and algebraic differentiation.
• The method efficiently extracts low-energy excited states and constructs finite-temperature Gibbs states through polynomial-time computations.
• GFMPS achieves exponential compression in free-fermion models and enables accurate simulation of many-body localized phases with constant bond dimension.

Generating-Function Matrix Product States (GFMPS) enable the efficient encoding and manipulation of large classes of quantum many-body states and their spectra within the matrix product state (MPS) and matrix product operator (MPO) frameworks by leveraging generating-function techniques and algebraic differentiation. GFMPS formalism has proved especially potent in the construction of low-temperature Gibbs states, the systematic extraction of excited-state manifolds, and the compression of Gaussian (free fermionic) states, as well as in the non-interacting and many-body localized (MBL) regimes.

1. Mathematical Foundation and Generating-Function Construction

The central principle of GFMPS is the replacement of explicit sums over exponentially many configurations (or tensor-network diagrams) with derivatives of a “generating-function” MPS or MPO. A key building block is the translationally invariant ground-state MPS,

$|\psi(A)\rangle = \prod_{j=0}^{N-1} A^{s_j}$

where $A$ is a rank-3 tensor of physical index $s_j$ (local Hilbert space dimension $d$) and bond dimension $\chi$.

To generate excited-state manifolds, a Bloch-type ansatz is constructed as

$$|G_\phi(\lambda, B)\rangle = \prod_{j=0}^{N-1} \left(A + \lambda e^{-ik r_j} B\right)_j$$

where $k$ is momentum, $\lambda$ a “bookkeeping” parameter, and $B$ an impurity tensor of the same shape as $A$. The first derivative at $\lambda=0$ recovers the single-mode excitation: $$|\phi_k(B)\rangle = \frac{\partial}{\partial \lambda} |G_\phi(\lambda, B)\rangle \big|_{\lambda=0}$$ This formalism generalizes: e.g., to compute norms, the generating function for the norm is $$G_{\|\phi\|}(\lambda) = \langle G_\phi(\lambda, \bar{B})|G_\phi(\lambda, B)\rangle$$, and

$$\|\phi_k(B)\|^2 = N\, \frac{\partial}{\partial \lambda} G_{\|\phi\|}(\lambda)\big|_{\lambda=0}$$

Similarly, operator insertions are encoded as product expansions with inserted sources, whose derivatives yield Fourier-mode observables. All generating functions are encoded as MPS or MPO networks, enabling polynomial-time computations in system size for suitably area-law states (Takahashi et al., 14 Jan 2026).

For non-interacting (Gaussian) fermionic states, GFMPS construction can be done at the level of correlation matrices. The generating function of the amplitudes is fully factorized: $$G(z_1, \dots, z_N) = [M^{[1]}(z_1) M^{[2]}(z_2) \cdots M^{[N]}(z_N)]_{1,1}$$ where $M^{[n]}(z_n)$ is the local generating matrix (Fishman et al., 2015, Schuch et al., 2019).

2. Excited States and Spectral Decomposition

Using the generating-function protocol, excited states are parametrized via impurity tensors $B$ in the single-mode ansatz. The optimal $B$ is found by solving the generalized eigenproblem: $H_{\mu\nu}\, B^\nu = E N_{\mu\nu} B^\nu$ where

$$H_{\mu\nu} = \frac{\partial^2}{\partial \bar{B}^\mu \partial B^\nu} \langle \phi_k(\bar{B})|H|\phi_k(B)\rangle, \quad N_{\mu\nu} = \frac{\partial^2}{\partial \bar{B}^\mu \partial B^\nu} \langle \phi_k(\bar{B})|\phi_k(B)\rangle$$

These derivatives are implemented efficiently by automatic differentiation of the generating functions, which are just translationally invariant diagrams, avoiding explicit $O(N)$ term expansions.

This procedure yields a set of eigenpairs $(E_{n,k}, B_n)$, producing an orthonormal set of Bloch-type states $|\Psi_k^n\rangle$. These states span a variationally controlled spectral window and are utilized in building truncated thermal ensembles or for dynamical observables (Takahashi et al., 14 Jan 2026).

An alternative construction for integrable or MBL models encodes all $2^n$ eigenstates via an $n$-variable polynomial MPS: $$M(\{\lambda_i\}) = \prod_{i=1}^n (G_i + \lambda_i E_i)$$ The extraction of any eigenstate is achieved via multivariate polynomial differentiation.

3. Finite-Temperature and Gibbs State Construction

The GFMPS approach enables explicit low-temperature Gibbs state construction via subspace expansion: $$\rho(\beta) \approx \frac{1}{Z} \sum_{n=0}^{M-1} e^{-\beta(E_n-E_0)} |\Psi^n\rangle\langle\Psi^n|$$ where $Z = \sum_{n=0}^{M-1} e^{-\beta(E_n-E_0)}$ and excited-state MPS are obtained as above. At inverse temperature $\beta \gg 1$, truncation error is exponentially suppressed since only low-lying states contribute to the thermal ensemble.

Observable expectation values are evaluated using precomputed matrix elements via further derivatives of generating functions. For dynamical structure factors or frequency-resolved quantities, the necessary correlators are assembled from matrix elements $\langle \Psi^m|O|\Psi^n\rangle$, efficiently matched by momentum conservation rules (Takahashi et al., 14 Jan 2026).

4. Algorithmic Complexity, Scalability, and Benchmarks

The methodological pipeline for GFMPS at finite temperature is as follows:

1. Variationally optimize the translationally invariant ground state MPS $A$.
2. Form $|G_\phi(\lambda, B)\rangle$ and compute the required generalized eigenproblems for $B$.
3. Diagonalize to obtain $M$ lowest excited states for relevant momenta $k$.
4. Assemble the truncated Gibbs density matrix using those eigenstates.
5. Compute observables using generating-function derivatives.

The cost of the most demanding steps scales as $O(M^2 N d \chi^5)$ for observable evaluation and $O(M N d \chi^5)$ for preparing the spectrum, with typical $\chi \approx 20-40$ and $M$ up to several thousand feasible for $N \lesssim 100$ (Takahashi et al., 14 Jan 2026).

In the free fermion (Gaussian) setting, the GFMPS protocol realizes exponential compression: the conventional MPS scaling $O(N D^3)$ with $D \sim e^S$ (entanglement entropy $S$) is replaced with $O(N \chi^3)$ with $\chi \sim S$ (Schuch et al., 2019, Fishman et al., 2015). In the many-body localized phase, the entire spectrum can be encoded with a bounded bond dimension, while delocalized (thermal) phases feature growth with system size (Pekker et al., 2014).

Benchmark studies for the S=½ Heisenberg chain and S=1 chain at low temperature demonstrate that retaining sufficiently many states reproduces exact diagonalization to high fidelity ($\geq 0.998$), rendering GFMPS highly accurate in the low-energy window (Takahashi et al., 14 Jan 2026).

5. Comparison with Alternative Finite-Temperature Tensor Network Methods

Standard tensor-network finite temperature methods include purification (which encodes the Gibbs state as a pure state in an enlarged Hilbert space evolved via imaginary time) and the minimally entangled typical thermal states (METTS) protocol, which samples random product states and applies imaginary-time evolution. Both approaches suffer from entanglement growth at low temperatures, necessitating large bond dimensions and/or extensive sampling.

In contrast, GFMPS circumvents real or imaginary time evolution of mixed states entirely, instead handling entanglement by expanding the eigenstate subspace. The area-law constraint only enters at the level of individual excited states, not the full density matrix, allowing highly efficient and accurate simulations at low temperature for one-dimensional models (Takahashi et al., 14 Jan 2026).

6. Extensions, Limitations, and Outlook

The current GFMPS implementation is based on the single-mode excitation ansatz, which robustly captures the essential low-temperature spectrum in gapless 1D chains (e.g., the deconfined spinon continuum in S=½). At higher temperatures, capturing observables accurately requires systematic inclusion of multi-mode excitations; strategies based on two- or three-particle ansätze or excitation-space algorithms have been proposed.

A direct extension to higher dimensions is in principle straightforward: the MPS tensors are replaced by their projected entangled pair state (PEPS) analogs, and the generating-function protocol is adapted with contraction cost determined by the PEPS contraction scheme.

Error control and convergence are regulated by two tunable parameters: the MPS bond dimension $\chi$ and the number of retained states $M$. At low temperatures, moderate $M$ suffices due to exponential decay in $w_n(\beta)$, but for temperatures above the spectral gap, $M$ must grow rapidly and costs eventually become prohibitive.

For interacting disordered systems in the MBL phase, GFMPS encapsulates the entire many-body spectrum with only two rank-3 tensors per site, rendering many observables and l-bit operators analytically tractable at constant bond dimension. Approaching the thermal phase, bond dimension increases rapidly, and full-spectrum representation loses its efficiency advantages (Pekker et al., 2014).

GFMPS, as a “third paradigm” for finite-temperature tensor network simulations alongside purification and METTS, is particularly well-adapted for extremely low temperatures and applications where the spectral window is dominated by few-body excitations (Takahashi et al., 14 Jan 2026). For Gaussian fermionic systems, the exponential compression continues to make it the method of choice for large-scale simulations (Schuch et al., 2019, Fishman et al., 2015).

7. Representative Applications

• Low-temperature dynamical and thermodynamic response of quantum spin chains, with benchmark-level accuracy versus exact diagonalization (Takahashi et al., 14 Jan 2026).
• Time-dependent and finite-temperature correlation functions in free fermionic models with system sizes up to $10^6$ sites; application to quantum impurity models and quasi-1D geometries (Schuch et al., 2019).
• Encoding the entire many-body spectrum in MBL phases, enabling extraction of all eigenstates and their observables at constant bond dimension (Pekker et al., 2014).
• Efficient MPS-based representation and manipulation of ground states of quadratic Hamiltonians, allowing scalable studies of free-fermion topological and symmetry-broken phases (Fishman et al., 2015).

The GFMPS formalism therefore integrates seamlessly with established tensor-network algorithms, and its algebraic, differentiable structure directly supports analytical computation, scalable numerics, and generalization to excited-state and finite-temperature properties across a wide class of quantum lattice models.