# publications

List of all publications.

## 2024

- SmoQyDQMC.jl: A flexible implementation of determinant quantum Monte Carlo for Hubbard and electron-phonon interactionsBenjamin Cohen-Stead, Sohan Malkaruge Costa, James Neuhaus, and 5 more authors
*SciPost Phys. Codebases*, 2024We introduce the SmoQyDQMC.jl package, a Julia implementation of the determinant quantum Monte Carlo algorithm. SmoQyDQMC.jl supports generalized tight-binding Hamiltonians with on-site Hubbard and generalized electron-phonon (e-ph) interactions, including non-linear e-ph coupling and anharmonic lattice potentials. Our implementation uses hybrid Monte Carlo methods with exact forces for sampling the phonon fields, enabling efficient simulation of low-energy phonon branches, including acoustic phonons. The SmoQyDQMC.jl package also uses a flexible scripting interface, allowing users to adapt it to different workflows and interface with other software packages in the Julia ecosystem. The code for this package can be downloaded from our GitHub repository at https://github.com/SmoQySuite/SmoQyDQMC.jl or installed using the Julia package manager. The online documentation, including examples, can be obtained from our document page at https://smoqysuite.github.io/SmoQyDQMC.jl/stable/.

- Persistence of small polarons into the superconducting phase of Ba_(1-x) K_x Bi O_3Muntaser Naamneh, Eugenio Paris, Daniel McNally, and 8 more authors
*arXiv preprint arXiv:2408.00401*, 2024Bipolaronic superconductivity is an exotic pairing mechanism proposed for materials like Ba_(1-x) K_x Bi O_3 (BKBO); however, conclusive experimental evidence for a (bi)polaron metallic state in this material remains elusive. Here, we combine resonant inelastic x-ray and neutron total scattering techniques with advanced modelling to study the local lattice distortions, electronic structure, and electron-phonon coupling (e-ph) in BKBO as a function of doping. Data for the parent compound (x=0) indicates that the electronic gap opens in predominantly oxygen-derived states strongly coupled to a long-range ordered breathing distortion of the oxygen sublattice. Upon doping, short-range breathing distortions and sizable (e-ph) coupling persist into the superconducting regime (x=0.4). Comparisons with exact diagonalization and determinant quantum Monte Carlo calculations further support this conclusion. Our results provide compelling evidence that BKBO’s metallic phase hosts a liquid of small (bi)polarons derived from local breathing distortions of the lattice, with implications for understanding the low-temperature superconducting instability

- Kekulé valence bond order in the honeycomb lattice optical Su-Schrieffer-Heeger model and its relevance to grapheneSohan Malkaruge Costa, Benjamin Cohen-Stead, and Steven Johnston
*Phys. Rev. B*, Sep 2024We perform sign-problem-free determinant quantum Monte Carlo simulations of the optical Su-Schrieffer-Heeger (SSH) model on a half-filled honeycomb lattice. In particular, we investigate the model’s semi-metal (SM) to Kekulé Valence Bond Solid (KVBS) phase transition at zero and finite temperatures as a function of phonon energy and interaction strength. Using hybrid Monte Carlo sampling methods we can simulate the model near the adiabatic regime, allowing us to access regions of parameter space relevant to graphene. Our simulations suggest that the SM-KVBS transition is weakly first-order at all temperatures, with graphene situated close to the phase boundary in the SM region of the phase diagram. Our results highlight the important role bond-stretching phonon modes play in the formation of KVBS order in strained graphene-derived systems.

- SmoQyDEAC. jl: A differential evolution package for the analytic continuation of imaginary time correlation functionsJames Neuhaus, Nathan S Nichols, Debshikha Banerjee, and 4 more authors
*arXiv:2407.04568*, Sep 2024We introduce the SmoQyDEAC.jl package, a Julia implementation of the Differential Evolution Analytic Continuation (DEAC) algorithm [N. S. Nichols et al., Phys. Rev. E 106, 025312 (2022)] for analytically continuing noisy imaginary time correlation functions to the real frequency axis. Our implementation supports fermionic and bosonic correlation functions on either the imaginary time or Matsubara frequency axes, and treatment of the covariance error in the input data. This paper presents an overview of the DEAC algorithm and the features implemented in the SmoQyDEAC.jl. It also provides detailed benchmarks of the package’s output against the popular maximum entropy and stochastic analytic continuation methods.

- Fluctuating charge-density-wave correlations in the three-band Hubbard modelP. Mai, B. Cohen-Stead, T. A. Maier, and 1 more author
*arXiv:2405.13164*, Sep 2024The high-temperature superconducting cuprates host unidirectional spin- and charge-density-wave orders that can intertwine with superconductivity in non-trivial ways. While the charge components of these stripes have now been observed in nearly all cuprate families, their detailed evolution with doping varies across different materials and at high and low temperatures. We address this problem using non-perturbative determinant quantum Monte Carlo calculations for the three-band Hubbard model. Using an efficient implementation, we can resolve the model’s fluctuating spin and charge modulations and map their evolution as a function of the charge transfer energy and doping. We find that the incommensurability of the charge modulations is decoupled from the spin modulations and decreases with hole doping, consistent with experimental measurements at high temperatures. These findings support the proposal that the high-temperature charge correlations are distinct from the intertwined stripe order observed at low-temperature and in the single-band Hubbard model.

## 2023

- Comparative study of the superconductivity in the Holstein and optical Su-Schrieffer-Heeger modelsAndy Tanjaroon Ly, Benjamin Cohen-Stead, Sohan Malkaruge Costa, and 1 more author
*Phys. Rev. B*, Nov 2023Theoretical studies suggest that Su-Schrieffer-Heeger-like electron-phonon (e-ph) interactions can mediate high-temperature bipolaronic superconductivity that is robust against repulsive electronelectron interactions. Here we present a comparative analysis of the pairing and competing charge/bond correlations in the two-dimensional Holstein and optical Su-Schrieffer-Heeger (SSH) models using numerically exact determinant quantum Monte Carlo. We find that the SSH interactions support light bipolarons and strong superconducting correlations out to relatively large values of the e-ph coupling λand densities near half-filling, while the Holstein interaction does not due to the formation of heavy bipolarons and competing charge-density-wave order. We further find that the Holstein and SSH models have comparable pairing correlations in the weak coupling limit for dilute carrier concentrations, where competing orders and polaronic effects are absent. These results support the proposal that SSH (bi)polarons can support superconductivity to larger values of λin comparison to the Holstein polaron, but that the resulting T_c gains are small in the weak coupling limit. We also find that the SSH model’s pairing correlations are suppressed after including a weak on-site Hubbard repulsion. These results have important implications for identifying and engineering bipolaronic superconductivity.

- Comparative determinant quantum Monte Carlo study of the acoustic and optical variants of the Su-Schrieffer-Heeger modelSohan Malkaruge Costa, Benjamin Cohen-Stead, Andy Tanjaroon Ly, and 2 more authors
*Phys. Rev. B*, Oct 2023We compare the acoustic Su-Schrieffer-Heeger (SSH) model with two of its optical variants where the phonons are defined on either on the sites or bonds of the system. First, we discuss how to make fair comparisons between these models in any dimension by ensuring their dimensionless coupling λand relevant phonon energies are the same. We then use determinant quantum Monte Carlo to perform non-perturbative and sign-problem-free simulations of all three models on one-dimensional chains at and away from half-filling. By comparing the results obtained from each model, we demonstrate that the optical and acoustic models produce near identical results within error bars for suitably chosen phonon energies and λat half-filling. In contrast, the bond model has quantitatively different behavior due to its coupling to the q = 0 phonon mode. These differences also manifest in the total length of the chain, which shrinks for the bond model but not for the acoustic and optical models when λ\ne 0. Our results have important implications for quantum Monte Carlo modeling of SSH-like interactions, where these models are sometimes regarded as being interchangeable.

- Flexible class of exact Hubbard-Stratonovich transformationsSeher Karakuzu, Benjamin Cohen-Stead, Cristian D. Batista, and 2 more authors
*Phys. Rev. E*, May 2023We consider a class of Hubbard-Stratonovich transformations suitable for treating Hubbard interactions in the context of quantum Monte Carlo simulations. A tunable parameter p allows usto continuously vary from a discrete Ising auxiliary field (p = ∞) to a compact auxiliary field that couples to electrons sinusoidally (p = 0). In tests on the single-band square and triangular Hubbard models, we find that the severity of the sign problem decreases systematically with increasing p. Selecting p finite, however, enables continuous sampling methods like the Langevin or Hamiltonian Monte Carlo methods. We explore the tradeoffs between various simulation methods through numerical benchmarks.

- A hybrid Monte Carlo study of bond-stretching electron–phonon interactions and charge order in BaBiO_3Benjamin Cohen-Stead, Kipton Barros, Richard Scalettar, and 1 more author
*npj Computational Materials*, May 2023The relationship between electron-phonon (e-ph) interactions and charge-density-wave (CDW) order in the bismuthate family of high-temperature superconductors remains unresolved. We address this question using nonperturbative hybrid Monte Carlo calculations for the parent compound BaBiO_3. Our model includes the Bi 6s and O 2p_σorbitals and coupling to the Bi-O bond-stretching branch of optical phonons via modulations of the Bi-O hopping integral. We simulate three-dimensional clusters of up to 4000 orbitals, with input model parameters taken from ab initio electronic structure calculations and a phonon energy \hbar\Omega_0 = 60 meV. Our results demonstrate that the coupling to the bond-stretching modes is sufficient to reproduce the CDW transition in this system, despite a relatively small dimensionless coupling. We also find that the transition deviates from the weak-coupling Peierls’ picture. This work demonstrates that off-diagonal e-ph interactions in orbital space are vital in establishing the bismuthate phase diagram.

- Charge order in the kagome lattice Holstein model: a hybrid Monte Carlo studyOwen Bradley, Benjamin Cohen-Stead, Steven Johnston, and 2 more authors
*npj Quantum Materials*, May 2023The Holstein model is a paradigmatic description of the electron-phonon interaction, in which electrons couple to local dispersionless phonon modes, independent of momentum. The model has been shown to host a variety of ordered ground states such as charge density wave (CDW) order and superconductivity on several geometries, including the square, honeycomb, and Lieb lattices. In this work, we study CDW formation in the Holstein model on the kagome lattice, using a recently developed hybrid Monte Carlo simulation method. We present evidence for CDW order at an average electron filling of ⟨n ⟩= 2/3 per site, with an ordering wavevector at the K-points of the Brillouin zone. We estimate a phase transition occurring at T_c ≈t/18, where t is the nearest-neighbor hopping parameter. Our simulations find no signature of CDW order at other electron fillings or ordering momenta for temperatures T ≥t/20.

- Charge correlations suppress unconventional pairing in the Holstein modelPhilip M Dee, Benjamin Cohen-Stead, Steven Johnston, and 1 more author
*Physical Review B*, May 2023In a recent work by Schrodi et al. [Phys. Rev. B 104, L140506 (2021)], the authors find an unconventional superconducting state with a sign-changing order parameter using the Migdal-Eliashberg theory, including the first vertex correction. This unconventional solution arises despite using an isotropic bare electron-phonon coupling in the Hamiltonian. We examine this claim using hybrid quantum Monte Carlo for a single-band Holstein model with a cuprate-like noninteracting band structure and identical parameters to Schrodi et al. Our Monte Carlo results for these parameters suggest that unconventional pairing correlations do not exceed their noninteracting values at any carrier concentration we have checked. Instead, strong charge-density-wave correlations persist at the lowest accessible temperatures for dilute and nearly half-filled bands. Lastly, we present arguments for how vertex-corrected Migdal-Eliashberg calculation schemes can lead to uncontrolled results in the presence of Fermi surface nesting.

## 2022

- Fast and scalable quantum Monte Carlo simulations of electron-phonon modelsBenjamin Cohen-Stead, Owen Bradley, Cole Miles, and 3 more authors
*Phys. Rev. E*, Jun 2022We introduce methodologies for highly scalable quantum Monte Carlo simulations of electron-phonon models, and report benchmark results for the Holstein model on the square lattice. The determinant quantum Monte Carlo (DQMC) method is a widely used tool for simulating simple electron-phonon models at finite temperatures, but incurs a computational cost that scales cubically with system size. Alternatively, near-linear scaling with system size can be achieved with the hybrid Monte Carlo (HMC) method and an integral representation of the Fermion determinant. Here, we introduce a collection of methodologies that make such simulations even faster. To combat "stiffness" arising from the bosonic action, we review how Fourier acceleration can be combined with time-step splitting. To overcome phonon sampling barriers associated with strongly-bound bipolaron formation, we design global Monte Carlo updates that approximately respect particle-hole symmetry. To accelerate the iterative linear solver, we introduce a preconditioner that becomes exact in the adiabatic limit of infinite atomic mass. Finally, we demonstrate how stochastic measurements can be accelerated using fast Fourier transforms. These methods are all complementary and, combined, may produce multiple orders of magnitude speedup, depending on model details.

- Dynamical tuning of the chemical potential to achieve a target particle number in grand canonical Monte Carlo simulationsCole Miles, Benjamin Cohen-Stead, Owen Bradley, and 3 more authors
*Phys. Rev. E*, Apr 2022We present a method to facilitate Monte Carlo simulations in the grand canonical ensemble given a target mean particle number. The method imposes a fictitious dynamics on the chemical potential, to be run concurrently with the Monte Carlo sampling of the physical system. Corrections to the chemical potential are made according to time-averaged estimates of the mean and variance of the particle number, with the latter being proportional to thermodynamic compressibility. We perform a variety of tests, and in all cases find rapid convergence of the chemical potential – inexactness of the tuning algorithm contributes only a minor part of the total measurement error for realistic simulations.

## 2021

- Quantum Monte Carlo study of an anharmonic Holstein modelG Paleari, F Hébert, B Cohen-Stead, and 3 more authors
*Physical Review B*, Apr 2021We study the effects of anharmonicity on the physics of the Holstein model, which describes the coupling of itinerant fermions and localized quantum phonons, by introducing a quartic term in the phonon potential energy. We find that the presence of this anharmonic term reduces the extent of the charge density wave phase (CDW) at half-filling as well as the transition temperature to this phase. Doping away from half-filling, we observe a first order phase transition between the CDW and a homogeneous phase which is also present in the harmonic model. In addition, we study the evolution of the superconducting susceptibility in the doped region and show that anharmonicity can enhance the superconducting response.

## 2020

- Langevin simulations of the half-filled cubic Holstein modelBenjamin Cohen-Stead, Kipton Barros, Zi Yang Meng, and 3 more authors
*Physical Review B*, Apr 2020Over the past several years, reliable Quantum Monte Carlo results for the charge density wave transition temperature T_cdw of the half-filled two dimensional Holstein model in square and honeycomb lattices have become available for the first time. Exploiting the further development of numerical methodology, here we present results in three dimensions, which are made possible through the use of Langevin evolution of the quantum phonon degrees of freedom. In addition to determining T_cdw from the scaling of the charge correlations, we also examine the nature of charge order at general wave vectors for different temperatures, couplings, and phonon frequencies, and the behavior of the spectral function and specific heat.

## 2019

- Effect of strain on charge density wave order in the Holstein modelBenjami Cohen-Stead, NC Costa, Ehsan Khatami, and 1 more author
*Physical Review B*, Apr 2019We investigate charge ordering in the Holstein model in the presence of anisotropic hopping, t_x, t_y = 1−δ, 1+δ, as a model of the effect of strain on charge density wave (CDW) materials. Using Quantum Monte Carlo simulations, we show that the CDW transition temperature is relatively insensitive to moderate anisotropy δ≤0.3, but begins to decrease more rapidly at δ≥0.4. However, the density correlations, as well as the kinetic energies parallel and perpendicular to the compressional axis, change significantly for moderate δ. Accompanying mean-field theory calculations show a similar qualitative structure, with the transition temperature relatively constant at small δ and a more rapid decrease for larger strains. We also obtain the density of states N(ω), which provides clear signal of the charge ordering transition at large strain, where finite size scaling of the charge structure factor is extremely difficult because of the small value of the order parameter.