The ALF (Algorithms for Lattice Fermions) project release 2.4. Documentation for the auxiliary-field quantum Monte Carlo code

TL;DR

The paper presents ALF 2.4, a versatile open-source code base implementing finite-temperature and projective auxiliary-field quantum Monte Carlo for a broad class of lattice fermion models. It formalizes a general Hamiltonian built from one-body terms, squares of one-body terms, and bosonic couplings, and supports five predefined model families (e.g., SU(N) Hubbard, Kondo, t-V, long-range Coulomb, and Z2 lattice gauge theories), together with an integrated stochastic Maximum Entropy solver for analytic continuation. It details multiple updating schemes (local and global), parallel tempering, Langevin dynamics, and symmetric Trotter decompositions, along with stabilization strategies for the fermionic determinant, enabling scalable simulations on modern HPC platforms. The package provides both finite-temperature and projective (ground-state) workflows, a detailed data-structure library, and a Python interface, making it a powerful framework for benchmarking new algorithms and exploring correlated fermion phenomena across diverse lattice geometries. The work emphasizes careful error analysis, autocorrelation handling, and symmetry-based optimizations to mitigate the sign problem and ensure reliable results for complex models.

Abstract

The Algorithms for Lattice Fermions package provides a general code for the finite-temperature and projective auxiliary-field quantum Monte Carlo algorithm. The code is engineered to be able to simulate any model that can be written in terms of sums of single-body operators, of squares of single-body operators and single-body operators coupled to a bosonic field with given dynamics. The package includes five pre-defined model classes: SU(N) Kondo, SU(N) Hubbard, SU(N) t-V and SU(N) models with long range Coulomb repulsion on honeycomb, square and N-leg lattices, as well as $Z_2$ unconstrained lattice gauge theories coupled to fermionic and $Z_2$ matter. An implementation of the stochastic Maximum Entropy method is also provided. One can download the code from our Git instance at https://git.physik.uni-wuerzburg.de/ALF/ALF/-/tree/ALF-2.4 and sign in to file issues.

Figures (5)

• Figure 1: Total energy for the 6-site Hubbard chain at $U/t=4$, $\beta t = 4$ and with open boundary conditions. For this system it can be shown that the determinant is always positive, so that no singularities occur in the action and, consequently, the Langevin dynamics works very well. The reference data point at $\delta t_{l} =0$ comes from the discrete field code for the field coupled to the z-component of the magnetization and reads $-2.8169 \pm 0.0013$, while the extrapolated value is $-2.8176 \pm 0.0010$. Throughout the runs the maximal force remained bellow the threshold of 1.5. The displayed data has been produced by the pyALF script https://git.physik.uni-wuerzburg.de/ALF/pyALF/-/blob/\pyALFbranch/Scripts/Langevin.py.
• Figure 2: Analysis of Trotter systematic error. Left: We consider a $6 \times 6$ Hubbard model on the Honeycomb lattice, $U/t=2$, half-band filling, inverse temperature $\beta t =40$, and we have used an HS transformation that couples to the density. The figure plots the local-time displaced Green function. Right: Here we consider the $6\times 6$ Hubbard model at $U/t=4$, half-band filling, inverse temperature $\beta t =5$, and we have used the HS transformation that couples to the $z$-component of spin. We provide data for the four combinations of the logical variables Symm and Checkerboard, where Symm=.true. (.false.) indicates a symmetric (asymmetric) Trotter decomposition has been used, and Checkerboard=.true. (.false.) that the checkerboard decomposition for the hopping matrix has (not) been used. The large deviations between different choices of Symm are here $\sim [T, [T,H]]$ as detailed in goth2020.
• Figure 3: Comparison between the finite-temperature and projective codes for the Hubbard model on a $6 \times 6$ Honeycomb lattice at $U/t=2$ and with periodic boundary conditions. For the projective code (blue and black symbols) $\beta t = 1$ is fixed, while $\theta$ is varied. In all cases we have $\Delta \tau t = 0.1$, no checkerboard decomposition, and a symmetric Trotter decomposition. For this lattice size and choice of boundary conditions, the non-interacting ground state is degenerate, since the Dirac points belong to the discrete set of crystal momenta. In order to generate the trial wave function we have lifted this degeneracy by either including a Kékulé mass term Lang13 that breaks translation symmetry (blue symbols), or by adding a next-next nearest neighbor hopping (black symbols) that breaks the symmetry nematically and shifts the Dirac points away from the zone boundary Ixert14. As apparent, both choices of trial wave functions yield the same answer, which compares very well with the finite temperature code at temperature scales below the finite-size charge gap.
• Figure 4: The autocorrelation function $S_{\hat{O}}(t_{\textrm{Auto}})$ (a) and the scaling of the error with effective bin size (b) of three equal-time, spin-spin correlation functions $\hat{O}$ of the Hubbard model in the $M_z$ decoupling (see Sec. \ref{['sec:hubbard']}). Simulations were done on a $6 \times 6$ square lattice, with $U/t=4$ and $\beta t = 6$. We used $\textrm{N\_auto}=500$ (see Sec. \ref{['sec:running']}) and a total of approximately one million bins. The original bin contained only one sweep and we calculated around one million bins on a single core. The different autocorrelation times for the $xy$-plane compared to the $z$-direction can be detected from the decay rate of the autocorrelation function (a) and from the point where saturation of the error sets in (b), which defines the required effective bin size for independent measurements. The improved estimator $(S_{\hat{S}^{x}} + S_{\hat{S}^{y}}+ S_{\hat{S}^{z}})/3$ appears to have the smallest autocorrelation time, as argued in the text.
• Figure 5: Structure of HDF5 output file data.h5. In parameters all $n$ namelists connected with the simulated Hamiltonian can be found.