Efficient Pseudomode Representation and Complexity of Quantum Impurity Models

TL;DR

This work develops a rigorous and practical framework to efficiently simulate non-equilibrium quantum impurity models by representing the bath with a finite set of pseudomodes. It proves that, under mild analytic assumptions on the bath spectrum, the required pseudomode count grows polylogarithmically with the evolution time $T$ and the inverse error $1/\varepsilon$, and shows that interpolative decomposition (ID) can reduce this count to $N_{\text{ID}} \sim \log(T)\log(1/\varepsilon)$. To broaden applicability, the authors augment the analytic construction with the adaptive AAA algorithm to obtain rational approximants of general spectral densities, which, when compressed by ID, yield comparable scaling. They also connect the pseudomode representation to a Liouvillian form that governs the impurity+pseudomode dynamics, enabling straightforward integration with tensor-network techniques. Across various baths (flat, linear, Gaussian-gapped, semicircular), the results indicate near-optimal, scalable bath representations with strong potential for DMFT and HEOM-type approaches in non-equilibrium settings.

Abstract

Out-of-equilibrium fermionic quantum impurity models (QIM), describing a small interacting system coupled to a continuous fermionic bath, play an important role in condensed matter physics. Solving such models is a computationally demanding task, and a variety of computational approaches are based on finding approximate representations of the bath by a finite number of modes. In this paper, we formulate the problem of finding efficient bath representations as that of approximating a kernel of the bath's Feynman-Vernon influence functional by a sum of complex exponentials, with each term defining a fermionic pseudomode. Under mild assumptions on the analytic properties of the bath spectral density, we provide an analytic construction of pseudomodes, and prove that their number scales polylogarithmically with the maximum evolution time $T$ and the approximation error $\varepsilon$. We then demonstrate that the number of pseudomodes can be significantly reduced by an interpolative matrix decomposition (ID). Furthermore, we present a complementary approach, based on constructing rational approximations of the bath's spectral density using the ``AAA'' algorithm, followed by compression with ID. The combination of two approaches yields a pseudomode count scaling as $N_\text{ID} \sim \log(T)\log(1/\varepsilon)$, and the agreement between the two approches suggests that the result is close to optimal. Finally, to relate our findings to QIM, we derive an explicit Liouvillian that describes the time evolution of the combined impurity-pseudomodes system. These results establish bounds on the computational resources required for solving out-of-equilibrium QIMs, providing an efficient starting point for tensor-network methods for QIMs.

Figures (17)

• Figure 1: Schematic illustration of the approach used in this work. Starting from a continuous spectral density (a), we determine a set of complex modes through either analytic construction (b, top) or the AAA algorithm (b, bottom). The yellow dotted line indicates the integration contour used for the Fourier integrals that define the hybridization funtion. In the analytic construction, pseudomodes are obtained via exponential frequency parametrization along a rotated contour in the complex plane (black dots). Blue stars represent poles from the spectral density and the Fermi-Dirac distribution, which must be considered when rotating the contour. The AAA algorithm identifies pseudomodes as sets of poles and residues in the complex plane (red stars). c) In order to compress the set of modes, we consider the hybridization function which encodes the sum of temporal correlations induced by each mode. By employing interpolative matrix decomposition, we extract a subset of pseudomodes with renormalized couplings which approximates the hybridization function with a controlled error. d) Illustration of the Keldysh contour.
• Figure 2: Integration contour and different groups of poles.
• Figure 3: Schematic illustration of the mode compression: Performing an interpolative decomposition of the kernel matrix $\mathcal{K}$ selects a subset of modes $\mathcal{J}$ with renormalized couplings encoded in the weights $\alpha_j.$
• Figure 4: Before compression: Scaling of $N_\text{bath}$ for a wide flat band, $\Gamma_\text{flat}^\Lambda(\omega),$ with width $\Lambda/\Gamma = 10^5$ and sharpness $\nu = 20 /\Lambda$ at inverse temperatures $\beta = 0$ (top) and $\beta = 10^4$ (bottom). This figure refers to the positive-frequency branch of the particle component, $\Delta^p_+(t)$ in Eq. (\ref{['eq:cont_integral']}). a) Scaling of $N_\text{bath}$ with the error $\epsilon.$ Inset: Squareroot of same data, illustrating the scaling law $N_\text{bath}\sim \log^2(1/\epsilon).$ The dotted line serves as guide to the eye. b) Left column: Scaling of number of modes with the time for fixed error (interpolated). The curves are fitted with a function $a \log(T) + \text{const.}$ (dotted). Right column: Dependence of fit parameter $a$ on the error.
• Figure 5: After compression: Scaling of $N_\text{ID}$ for a wide flat band. All parameters are as specified in Fig. \ref{['fig:scaling_wideband']}. a) Scaling of the number of modes retained after compression, $N_\text{ID},$ with the error $\epsilon$. b) Left column: Scaling of number of modes with time for fixed error (interpolated). The curves are fitted with a function $a \log(T) + \text{const.}$ (dotted). Right column: Dependence of the fit parameter $a$ on the error.

Theorems & Definitions (6)

• Theorem 3.1
• proof
• Theorem B.1
• proof
• Lemma B.2
• proof