High-temperature series expansion of the dynamic Matsubara spin correlator

TL;DR

We extend the high-temperature expansion to the Matsubara spin-spin correlator (Dyn-HTE) for Heisenberg spins with a single exchange $J$ and spins $S\in\{1/2,1\}$, obtaining exact rational coefficients up to order $n_{max}=12$. The coefficients are precomputed on ~10^6 embeddable graphs and released with open-source tools, enabling static susceptibilities and momentum-resolved results for arbitrary site-pairs; a kernel-trick framework makes the imaginary-time integrals tractable and links high-frequency Dyn-HTE coefficients to spectral moments for real-frequency dynamics. Benchmarking against exact solutions, QMC, and diagrammatic Monte Carlo data for chain and triangular lattices demonstrates reliable accuracy, and the static susceptibility can be described by a renormalized mean-field form with small beyond-rMF corrections. The work provides a scalable, high-precision route to analyze dynamical properties in frustrated and high-dimensional spin systems, with clear pathways to extensions to multiple couplings, anisotropies, fields, and higher-order correlators.

Abstract

The high-temperature series expansion for quantum spin models is a well-established tool to compute thermodynamic quantities and equal-time spin correlations, in particular for frustrated interactions. We extend the scope of this expansion to the dynamic Matsubara spin-spin correlator and develop an algorithm that yields exact expansion coefficients in the form of rational numbers. We focus on Heisenberg models with a single coupling constant J and spin lengths S=1/2,1. The expansion coefficients up to 12th order in J/T are precomputed on all possible $\sim 10^6$ graphs embeddable in arbitrary lattices and are provided in a repository. This enables calculation of static momentum-resolved susceptibilities for arbitrary site-pairs or wavevectors. We test our results for the antiferromagnetic S=1/2 chain and triangular lattice model. An important application that we discuss in a companion letter is the calculation of real-frequency dynamic structure factors. This is achieved by identifying the high-frequency expansion coefficients of the Matsubara correlator with frequency moments of the spectral function.

Figures (5)

• Figure 1: All required graphs $g^{(n)}_t$ with $n=0,1,2,3$ edges (black lines) and arbitrary numbered vertices (blue circles). Terminal vertices $j,j^\prime$ (also blue circles) are indicated by their red-square terminal flags attached with thin gray lines. The symmetry factor is $s[g^{(n)}_t]=1$ for all graphs shown except for $(n,t)=(3,5)$ where it is two (exchange of vertices $2 \leftrightarrow 3$).
• Figure 2: Subtractions from a graph $g^{(6)}$ shown in the upper left. All connected sub-graphs $g^{(k)}$ with $k<6$ and non-zero graph evaluation are shown in the top row, second to last column. The vacuum graphs $g^{(6)}\backslash g^{(k)}\equiv v$ are shown in the bottom row. These vacuum graphs are possibly disconnected, see third column. The multiplicity factors $f$ for the particular subtraction are also given.
• Figure 3: Matsubara correlators of the Heisenberg $S=1/2$ AFM chain for real-space distances $i-i^\prime=0,1,2$ (left) and wavevectors $k/\pi=0.4,1$ (right) as obtained from Dyn-HTE (lines). The frequencies are $\nu_m=2\pi m T$ with $m=0,1,2$ (top to bottom). Markers show benchmark results from QMC simulations of a 256-site ring with error bars smaller than the symbol size. The thin gray line denotes the evaluation of the bare Dyn-HTE series truncated at order $n_{max}=12$. Symmetric Padé approximants ([6,6] and [5,5]) of the x-series are shown by blue lines (full and dashed), for the transformed series in $u=\mathrm{tanh}(fx)$ with $f=0.205$ they are indicated in purple.
• Figure 4: Static susceptibility $\chi_\mathbf{k}$ for the Heisenberg $S=1/2$ AFM on the triangular lattice. Results from Dyn-HTE (lines) are compared against the bold-line diagrammatic Monte Carlo (BDMC) data obtained by Kulagin et al. in Ref. kulagin_bold_2013 (dots). Top: $\chi_{\mathbf{k}=K}$ at the $K$-point in the corner of the hexagonal BZ versus $T$. The convergence of various Padé approximants denoted by blue linestyles is decent and still improves by using the $u$-series (for $f=0.25$, purple lines, same Padé approximants as for $x$-series). Bottom: $\chi_\mathbf{k}$ for a path through the BZ at various $T$. The $\Gamma$-point is the center of the BZ and $M$ denotes the center of the BZ-edge. Away from the $K$-point there is good agreement between BDMC (dots) and Dyn-HTE (full lines). Here the [6,6] $u$-Padé with $f=0.25$ is shown. The dashed lines represent the best fit of the renormalized mean-field form in Eq. \ref{['eq:QCC']} to the Dyn-HTE results. The associated fit parameters $\{f_x,g_x\}$ are given in the legend.
• Figure 5: Amplitudes of parameter-ratios of the rMF form \ref{['eq:QCC']}: $x/f_x$, $g_x/f_x$ and the first few corrections $\epsilon_{2,3,4}/f_x$ [c.f. Eq. \ref{['eq:susc_exact']}] for the nearest-neighbor $S=1/2$ Heisenberg AFM model on the chain, triangular, kagome and pyrochlore lattice from Dyn-HTE. For resummation, we used $[6,6]$ and $[5,5]$ u-Padé approximants. The data is multiplied by the maximum (eigenvalue) of $-\gamma_1(\mathbf{k})$ given in the panels alongside the u-Padé parameter $f$. For the chain geometry QMC data (dots) lodepolletLodePolletWorm2024 matches well with the Dyn-HTE data (lines). For triangular, kagome and pyrochlore lattice $\epsilon_{2,3,4}/f_x$ are of order $10^{-3}$ rendering the rMF an excellent approximation. Whenever the [6,6] or [5,5] Padé approximant is anomalous (with obvious poles), we use a stable lower-order Padé approximant.