Table of Contents
Fetching ...

PyLIT: Reformulation and implementation of the analytic continuation problem using kernel representation methods

Alexander Benedix Robles, Phil-Alexander Hofmann, Thomas Chuna, Tobias Dornheim, Michael Hecht

TL;DR

This work addresses the ill-posed problem of analytic continuation from imaginary-time data to dynamic structure factors by reformulating the DSF as a linear combination of kernels with known Laplace transforms and solving a regularized inverse problem for the kernel coefficients. The PyLIT framework provides kernel construction, simulated-annealing-based hyperparameter tuning, and nonnegative gradient optimization, exploring entropic, Wasserstein, and $L^2$ regularizers within a Bayesian-prior setup. Key findings show that nonuniform kernel grids reduce the solution space, Wasserstein regularization performs on par with entropy while offering a linear gradient, and Gaussian kernels generally yield the best balance between stability and fidelity, with promising avenues to merge regularized and stochastic optimization. The approach is validated on synthetic UEG data and authentic PIMC data, revealing roton-like features in $S(q,\omega)$ and offering a practical, high-performance tool for quantum many-body dynamics and related scattering diagnostics, with potential extensions to multi-dimensional AC and joint stochastic-regularization frameworks.

Abstract

Path integral Monte Carlo (PIMC) simulations are a cornerstone for studying quantum many-body systems. The analytic continuation (AC) needed to estimate dynamic quantities from these simulations is an inverse Laplace transform, which is ill-conditioned. If this inversion were surmounted, then dynamical observables (e.g. dynamic structure factor (DSF) $S(q,ω)$) could be extracted from the imaginary-time correlation functions estimates. Although of important, the AC problem remains challenging due to its ill-posedness. To address this challenge, we express the DSF as a linear combination of kernel functions with known Laplace transforms that have been tailored to satisfy its physical constraints. We use least-squares optimization regularized with a Bayesian prior to determine the coefficients of this linear combination. We explore various regularization term, such as the commonly used entropic regularizer, as well as the Wasserstein distance and $L^2$-distance as well as techniques for setting the regularization weight. A key outcome is the open-source package PyLIT (\textbf{Py}thon \textbf{L}aplace \textbf{I}nverse \textbf{T}ransform), which leverages Numba and unifies the presented formulations. PyLIT's core functionality is kernel construction and optimization. In our applications, we find PyLIT's DSF estimates share qualitative features with other more established methods. We identify three key findings. Firstly, independent of the regularization choice, utilizing non-uniform grid point distributions reduced the number of unknowns and thus reduced our space of possible solutions. Secondly, the Wasserstein distance, a previously unexplored regularizer, performs as good as the entropic regularizer while benefiting from its linear gradient. Thirdly, future work can meaningfully combine regularized and stochastic optimization. (text cut for char. limit)

PyLIT: Reformulation and implementation of the analytic continuation problem using kernel representation methods

TL;DR

This work addresses the ill-posed problem of analytic continuation from imaginary-time data to dynamic structure factors by reformulating the DSF as a linear combination of kernels with known Laplace transforms and solving a regularized inverse problem for the kernel coefficients. The PyLIT framework provides kernel construction, simulated-annealing-based hyperparameter tuning, and nonnegative gradient optimization, exploring entropic, Wasserstein, and regularizers within a Bayesian-prior setup. Key findings show that nonuniform kernel grids reduce the solution space, Wasserstein regularization performs on par with entropy while offering a linear gradient, and Gaussian kernels generally yield the best balance between stability and fidelity, with promising avenues to merge regularized and stochastic optimization. The approach is validated on synthetic UEG data and authentic PIMC data, revealing roton-like features in and offering a practical, high-performance tool for quantum many-body dynamics and related scattering diagnostics, with potential extensions to multi-dimensional AC and joint stochastic-regularization frameworks.

Abstract

Path integral Monte Carlo (PIMC) simulations are a cornerstone for studying quantum many-body systems. The analytic continuation (AC) needed to estimate dynamic quantities from these simulations is an inverse Laplace transform, which is ill-conditioned. If this inversion were surmounted, then dynamical observables (e.g. dynamic structure factor (DSF) ) could be extracted from the imaginary-time correlation functions estimates. Although of important, the AC problem remains challenging due to its ill-posedness. To address this challenge, we express the DSF as a linear combination of kernel functions with known Laplace transforms that have been tailored to satisfy its physical constraints. We use least-squares optimization regularized with a Bayesian prior to determine the coefficients of this linear combination. We explore various regularization term, such as the commonly used entropic regularizer, as well as the Wasserstein distance and -distance as well as techniques for setting the regularization weight. A key outcome is the open-source package PyLIT (\textbf{Py}thon \textbf{L}aplace \textbf{I}nverse \textbf{T}ransform), which leverages Numba and unifies the presented formulations. PyLIT's core functionality is kernel construction and optimization. In our applications, we find PyLIT's DSF estimates share qualitative features with other more established methods. We identify three key findings. Firstly, independent of the regularization choice, utilizing non-uniform grid point distributions reduced the number of unknowns and thus reduced our space of possible solutions. Secondly, the Wasserstein distance, a previously unexplored regularizer, performs as good as the entropic regularizer while benefiting from its linear gradient. Thirdly, future work can meaningfully combine regularized and stochastic optimization. (text cut for char. limit)
Paper Structure (17 sections, 30 equations, 11 figures, 4 tables, 1 algorithm)

This paper contains 17 sections, 30 equations, 11 figures, 4 tables, 1 algorithm.

Figures (11)

  • Figure 1: Each plot shows various PyLIT DSF estimates obtained with: 50 uniform kernels, the RPA default model, denoted (prior), and Gull's Bayesian averaging over $\lambda$gull1989MEMBayesianWeighting. The static approximation, denoted with (input), was used to construct the synthetic data and is the true solution. We compute the DSF at different values of the wavenumber normalized by the Fermi wavenumber $q/q_F$, as indicated by dash shadow line underneath each curve. The $\omega$ axis is normalized by the uniform electron gas's plasma frequency $\omega_{p,e}$. The $r_s$ and $\Theta$ values are indicated by the plot title. The regularization term is indicated by plot legend.; from left to right the cross entropy regularizer (Table \ref{['tab:regularizers_analytic']}-entropy), the Wasserstein distance regularizer (Table \ref{['tab:regularizers_analytic']}-WD), the $L^2$-distance regularizer (Table \ref{['tab:regularizers_analytic']}-$L^2$-distance). Each plot contains different reconstructions of $S(q,\omega)$ obtained varying noise levels $\sigma_0 = 10^{-1}, 10^{-2}, 10^{-3}$.
  • Figure 2: Each plot shows various PyLIT DSF estimates obtained with: 50 uniform kernels, the RPA default model, denoted (prior), and the $\chi^2$-kink regularization weight $\lambda$KAUFMANN2023ana_cont. The static approximation, denoted with (input), was used to construct the synthetic data and is the true solution. We compute the DSF at different values of the wavenumber normalized by the Fermi wavenumber $q/q_F$, as indicated by dash shadow line underneath each curve. The $\omega$ axis is normalized by the uniform electron gas's plasma frequency $\omega_{p,e}$. The $r_s$ and $\Theta$ values are indicated by the plot title. The regularization term is indicated by plot legend.; from left to right the entropy regularizer (Table \ref{['tab:regularizers_analytic']}-entropy), the Wasserstein distance regularizer (Table \ref{['tab:regularizers_analytic']}-WD), the $L^2$-distance regularizer (Table \ref{['tab:regularizers_analytic']}-$L^2$-distance). Each plot contains different reconstructions of $S(q,\omega)$ obtained varying noise levels $\sigma_0 = 10^{-1}, 10^{-2}, 10^{-3}$.
  • Figure 3: Heatmaps of the DSF $S(q,\omega)$ estimated by PyLIT. The $\omega$ axis is normalized by the plasma frequency $\omega_{p,e}$ and the $q$ values are normalized by the Fermi wavenumber $q_F$. The $r_s$ and $\Theta$ values indicated by the plot title, and the regularization weight procedure, regularization term, and default model are indicated by the plot subtitle. The dispersion relation for the RPA is plotted in all PyLIT estimates as a dashed cyan line and the dispersion relation for the static approximation is plotted as a solid cyan line.
  • Figure 4: Cross sections from the DSF heat maps plotted in Fig. \ref{['fig:heatmaps-authentic-Bayesian']}. The cross sections are offset along the $y$-axis and the associated $q$ values are indicated by dash shadow line underneath each curve. The $r_s$ and $\Theta$ values indicated by the plot title, and regularization term indicated by plot legend. From left to right the entropy regularizer (Table \ref{['tab:regularizers_analytic']}-CE)), the Wasserstein distance regularizer (Table \ref{['tab:regularizers_analytic']}-WD), the $L^2$ distance regularizer (Table \ref{['tab:regularizers_analytic']}-$L^2$)
  • Figure 5: Cross sections from the DSF heat maps, $q$ values indicated by dash shadow line underneath each curve, $r_s$ and $\Theta$ values indicated by the plot title, regularization term indicated by plot legend. From left to right the entropy regularizer (Table \ref{['tab:regularizers_analytic']}-CE)), the Wasserstein distance regularizer (Table \ref{['tab:regularizers_analytic']}-WD). We do not include $q=0.63 q_F$ in this plot because the width of the peak is Dirac-delta-like and cannot be resolved by our $\omega$ resolution.
  • ...and 6 more figures