Summing Real Time Feynman Paths of Lattice Polaron with Matrix Product States

TL;DR

The paper tackles the challenging problem of real-time polaron dynamics by merging a real-time Feynman path integral with matrix product state techniques. It shows that the phonon influence functional can be expressed as an exponential of an MPS, $I({s}) = e^{\lambda^2 C({s})}$, where $C({s})$ has bond dimension $N+2$, and propagates the world-line amplitude $W({s})$ via a flow equation $dW/d(\lambda^2) = C({s}) W$, solved with TDVP. This approach enables sign-tolerant, scalable summation of high-dimensional path integrals and yields accurate 1D spectral functions that agree with diagrammatic Monte Carlo and variational results, while extending to 2D without phonon truncation and enabling diffusion studies. Overall, the method broadens the computational toolkit for real-time polaron dynamics and offers a pathway to tackling higher-dimensional and more complex polaron problems with controlled approximations.

Abstract

We study numerically the real time dynamics of lattice polarons by combining the Feynman path integral and the matrix product state (MPS) approach. By constructing and solving a flow equation, we show that the integrand, viewed as a multivariable function of polaron world line parameters, can be compressed as a low bond dimension MPS, thereby allowing for efficient evaluation of various dynamical observables. We establish the effectiveness of our method by benchmarking the calculated polaron spectral function in one dimension against available results. We further demonstrate its potential by presenting the polaron spectral function in two dimensions and simulating polaron diffusion in both one and two dimensions.

Figures (11)

• Figure 1: MPS summation of real time Feynman paths. (a) World line of a polaron is parametrized by its displacements $\{s\}$ between successive times. The amplitude $W$ as a function of $\{s\}$ may be viewed as a high dimensional tensor. (b) $W$ can be efficiently summed provided that it is well approximated by a low bond dimension MPS. (c) $W$ is given by the free amplitude $W_0$ dressed by phonon influence functional $I$, which, in turn, is the exponential of a low bond dimension matrix product operator (d). (e) $W$ satisfies an ordinary differential equation with respect to electron-phonon coupling $\lambda^2$. Integrating it from $\lambda^2=0$ yields the MPS approximation to $W$.
• Figure 2: Polaron spectral function in one dimension. (a)(b) Spectral function $A(k,\omega)$ for the Holstein model $\omega_q = \omega_0 = 1$, $f_q = 1$, and $\lambda = \sqrt{2}$ calculated with system size $N=16$ and bond dimension $\chi =120$, compared with the result from MPS Chebyshev expansion Zhao2023. Green crosses mark the dispersion of lowest polaron band from a variational method calculation Bonca1999. (c)(d)(e) Constant-$k$ cuts of $A(k,\omega)$ at representative momenta. (f)(g) Polaron ground state energy $E_\mathrm{GS}$ and spectral weight $Z$ for the parameter set $\omega_q = \omega_0 = 1$ and $f_q = 1$ extracted from our results, benchmarked against diagrammatic Monte Carlo data Ragni2020. (h)(i) Spectral function for a model with dispersive phonons. $\omega_q = 1+0.4\cos(q)$, $f_q = 1$, and $\lambda = \sqrt[4]{3.36}$, compared with the variational result Bonca2021. (j)(k)(l) Representative constant-$k$ cuts of $A(k,\omega)$.
• Figure 3: Polaron spectral function in two-dimensional square lattice. (a) Spectral function $A(k,\omega)$ along high symmetry directions of the first Brillouin zone for the Holstein model $\omega_q = \omega_0 = 2$, $f_q = 1$, $\lambda = 2$. The system size is $6\times6$. Maximal bond dimension is $\chi = 120$. Dashed line shows the tight binding dispersion. (b) Similar to (a) but for $\lambda = 2\sqrt{2}$. (c)(d) Convergence of the spectral function with respect to the system size and the bond dimension $\chi$ at the X point, which we find empirically exhibits the slowest convergence rate.
• Figure 4: Top row: the zeroth (a) and first (b) order moments of the spectral function along the high symmetry directions of the first Brillouin zone for the Holstein model $\omega_q = \omega_0 = 2$, $f_q = 1$, $\lambda = 2$, corresponding to Fig. \ref{['fig:spectra_2d']}a. Blue open circles and red dots are the numerically evaluated value and the exact sum rule, respectively. Insets show the absolute error. Bottom row: similar to the top row but for $\lambda=2\sqrt{2}$, which is to be cross-referenced with Fig. \ref{['fig:spectra_2d']}b.
• Figure 5: Polaron diffusion. (a) Electron density $\rho(x,t)$ as a function of coordinate $x$ and time $t$ for the one-dimensional Holstein model. $\omega_q = \omega_0 = 1$, $f_q = 1$, $\lambda=\sqrt{2}$. The system size $N = 16$. Bond dimension $\chi = 400$. (b) Evolution of $\rho(x,t)$ for selected sites. (c) Snapshots of $\rho(x,t)$ at different times for the two-dimensional square lattice Holstein model. $\omega_q = \omega_0 = 2$, $f_q=1$, $\lambda = 2\sqrt{2}$. We use a $7\times 7$ system with $\chi = 200$. Insets of (b) and (d) show the mean displacement squared (MSD) as a function of time (blue). The ballistic (red) and diffusive (yellow) behaviors are displayed as reference.