First-principles calculations can accurately describe electron-phonon (e-ph) interactions and electronic transport in a wide range of materials, but are currently limited to unit cells with up to $\sim$100 atoms due to computational cost. Here, we develop an atomistic electronic potential with Holstein- and Peierls-like terms for modeling e-ph interactions and phonon-limited electronic transport that enables the study of moiré systems with thousands of atoms per unit cell. This method can accurately reproduce first-principles e-ph coupling and resistivity in graphene and large-angle twisted bilayer graphene (TBG). Using this approach, we study TBG over a range of twist angles down to 1.6$^\circ$ (5044-atom unit cell), and report the evolution of e-ph interactions and phonon-limited resistivity with twist angle. The predicted resistivity increases by two orders of magnitude between 13.2$^\circ$ and 1.6$^\circ$, driven by the progressive reduction of the electronic energy scale. Our calculations can predict key experimental trends in 2.0$^\circ$ and 1.6$^\circ$ TBG, including the resistivity and its dependence on temperature and band filling. Our work establishes a scalable approach for quantitative studies of e-ph interactions and transport in moiré materials and other systems with previously inaccessible length scales.
The Boltzmann transport equation (BTE) with electron-phonon (e-ph) interactions computed from first principles is widely used to study electronic transport and nonequilibrium dynamics in materials. Calculating the e-ph collision integral is the most important step in the BTE, but it remains computationally costly, even with current MPI+OpenMP parallelization. This challenge makes it difficult to study materials with large unit cells and to achieve high resolution in momentum space. Here, we show acceleration of BTE calculations of electronic transport and ultrafast dynamics using graphical processing units (GPUs). We implement a novel data structure and algorithm, optimized for GPU hardware and developed using OpenACC, to process scattering channels and efficiently compute the collision integral. This approach significantly reduces the overhead for data referencing, movement, and synchronization. Relative to the efficient CPU implementation in the open-source package Perturbo (v2.2.0), used as a baseline, this approach achieves a speed-up of 40 times for both transport and nonequilibrium dynamics on GPU hardware, and achieves nearly linear scaling up to 100 GPUs. The novel data structure can be generalized to other electron interactions and scattering processes. We released this GPU implementation in the latest public version (v3.0.0) of Perturbo. The new MPI+OpenMP+GPU parallelization enables sweeping studies of e-ph physics and electron dynamics in conventional and quantum materials, and prepares Perturbo for exascale supercomputing platforms.
The Migdal-Eliashberg (ME) formalism provides a reliable framework for describing phonon-mediated superconductivity in the adiabatic regime, where the electronic Fermi energy exceeds the characteristic phonon energy. In this work, we go beyond this limit by incorporating first-order vertex corrections to the electron-phonon (e-ph) interaction within the Eliashberg formalism and assess their impact on the superconducting properties of H3S and Pb using first-principles calculations. For H3S, where the adiabatic assumption breaks down, we find that vertex corrections to the e-ph coupling are substantial. When combined with phonon anharmonicity and the energy dependence of the electronic density of states, the predicted critical temperature (Tc) is in very good agreement with experimental observations. In contrast, for elemental Pb, where the adiabatic approximation remains valid, vertex corrections have a negligible effect, and the calculated Tc and superconducting gap closely match the predictions of the standard ME formalism. These findings demonstrate the importance of non-adiabatic corrections in strongly coupled high-Tc hydrides and establish a robust first-principles framework for accurately predicting superconducting properties across different regimes.
Matthias Ehrhardt, Michael Günther, Daniel Ševčovič
The port-Hamiltonian framework is a structure-preserving modeling approach that preserves key physical properties such as energy conservation and dissipation. When subsystems are modeled as port-Hamiltonian systems (pHS) with linearly related inputs and outputs, their interconnection remains port-Hamiltonian. This paper introduces a systematic method for transforming coupled port Hamiltonian ordinary differential equations systems (pHODE) into a single monolithic formulation, and for decomposing a monolithic system into weakly coupled subsystems. The monolithic representation ensures stability and structural integrity, whereas the decoupled form enables efficient distributed simulation via operator splitting or dynamic iteration.
The accuracy of phaseless auxiliary-field quantum Monte Carlo (ph-AFQMC) can be systematically improved with better trial states. Using multi-Slater determinant trial states, ph-AFQMC has the potential to faithfully treat strongly correlated systems, while balancing the static and dynamical correlations on an equal footing. This preprint presents an implementation and application of graphics processing unit-accelerated ph-AFQMC, for multi-Slater determinant trial wavefunctions (GPU-accelerated MSD-AFQMC), to enable efficient simulation of large-scale, strongly correlated systems. This approach allows for nearly-exact computation of ground state energies in multi-reference systems. Our GPU-accelerated MSD-AFQMC is implemented in the open-source code \texttt{ipie}, a Python-based AFQMC package [\textit{J. Chem. Theory Comput.}, 2022, 19(1): 109-121]. We benchmark the performance of the GPU code on transition-metal clusters like [Cu$_2$O$_2$]$^{2+}$ and [Fe$_2$S$_2$(SCH$_3$)]$^{2-}$. The GPU code achieves at least sixfold speedup in both cases, comparing the timings of a single A100 GPU to that of a 32-CPU node. For [Fe$_2$S$_2$(SCH$_3$)]$^{2-}$, we demonstrate that our GPU MSD-AFQMC can recover the dynamical correlation necessary for chemical accuracy with an MSD trial, despite the large number of determinants required ($>10^5$). Our work significantly enhances the efficiency of MSD-AFQMC calculations for large, strongly correlated molecules by utilizing GPUs, offering a promising path for exploring the electronic structure of transition metal complexes.
Benjamin Cohen-Stead, Shruti Agarwal, Sohan Malkaruge Costa
et al.
We introduce version 2.0 of the SmoQyDQMC.jl package, a Julia implementation of the determinant quantum Monte Carlo algorithm. SmoQyDQMC.jl supports generalized tight-binding Hamiltonians with local and extended Hubbard and generalized electron-phonon (e-ph) interactions, including non-linear e-ph coupling and anharmonic lattice potentials. Our implementation uses an optimized hybrid Monte Carlo method with exact forces to efficiently sample the phonon fields, enabling the 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, is found at https://smoqysuite.github.io/SmoQyDQMC.jl/stable/.
Fionn D. Malone, Ankit Mahajan, James S. Spencer
et al.
We report the development of a python-based auxiliary-field quantum Monte Carlo (AFQMC) program, ipie, with preliminary timing benchmarks and new AFQMC results on the isomerization of [Cu$_2$O$_2$$]^{2+}$. We demonstrate how implementations for both central and graphical processing units (CPUs and GPUs) are achieved in ipie. We show an interface of ipie with PySCF as well as a straightforward template for adding new estimators to ipie. Our timing benchmarks against other C++ codes, QMCPACK and Dice, suggest that ipie is faster or similarly performing for all chemical systems considered on both CPUs and GPUs. Our results on [Cu$_2$O$_2$$]^{2+}$ using selected configuration interaction trials show that it is possible to converge the ph-AFQMC isomerization energy between bis($μ$-oxo) and $μ$-$η^2$:$η^2$ peroxo configurations to the exact known results for small basis sets with $10^5$ to $10^6$ determinants. We also report the isomerization energy with a quadruple-zeta basis set with an estimated error less than a kcal/mol, which involved 52 electrons and 290 orbitals with $10^6$ determinants in the trial wavefunction. These results highlight the utility of ph-AFQMC and ipie for systems with modest strong correlation and large-scale dynamic correlation.
A semi-implicit fractional-step method that uses a staggered node layout and radial basis function-finite differences (RBF-FD) to solve the incompressible Navier-Stokes equations is developed. Polyharmonic splines (PHS) with polynomial augmentation (PHS+poly) are used to construct the global differentiation matrices. A systematic parameter study identifies a combination of stencil size, PHS exponent, and polynomial degree that minimizes the truncation error for a wave-like test function on scattered nodes. Classical modified wavenumber analysis is extended to RBF-FDs on heterogeneous node distributions and used to confirm that the accuracy of the selected 28-point stencil is comparable to that of spectral-like, 6th-order Padé-type finite differences. The Navier-Stokes solver is demonstrated on two benchmark problems, internal flow in a lid-driven cavity in the Reynolds number regime $10^2\leq$Re$\leq10^4$, and open flow around a cylinder at Re=100 and 200. The combination of grid staggering and careful parameter selection facilitates accurate and stable simulations at significantly lower resolutions than previously reported, using more compact RBF-FD stencils, without special treatment near solid walls, and without the need for hyperviscosity or other means of regularization.
Naman Bartwal, Shantanu Shahane, Somnath Roy
et al.
In recent years, a variety of meshless methods have been developed to solve partial differential equations in complex domains. Meshless methods discretize the partial differential equations over scattered points instead of grids. Radial basis functions (RBFs) have been popularly used as high accuracy interpolants of function values at scattered locations. In this paper, we apply the polyharmonic splines (PHS) as the RBF together with appended polynomial and solve the heat conduction equation in several geometries using a collocation procedure. We demonstrate the expected exponential convergence of the numerical solution as the degree of the appended polynomial is increased. The method holds promise to solve several different governing equations in thermal sciences.
We present an exponentially convergent semi-implicit meshless algorithm for the solution of Navier-Stokes equations in complex domains. The algorithm discretizes partial derivatives at scattered points using radial basis functions as interpolants. Higher-order polynomials are appended to polyharmonic splines (PHS-RBF) and a collocation method is used to derive the interpolation coefficients. The interpolating kernels are then differentiated and the partial-differential equations are satisfied by collocation at the scattered points. The PHS-RBF interpolation is shown to be exponentially convergent with discretization errors decreasing as a high power of a representative distance between points. We present here a semi-implicit algorithm for time-dependent and steady state fluid flows in complex domains. At each time step, several iterations are performed to converge the momentum and continuity equations. A Poisson equation for pressure corrections is formulated by imposing divergence free condition on the iterated velocity field. At each time step, the momentum and pressure correction equations are repeatedly solved until the velocities and pressure converge to a pre-specified tolerance. We have demonstrated the convergence and discretization accuracy of the algorithm for two model problems and simulated three other complex problems. In all cases, the algorithm is stable for Courant numbers in excess of ten. The algorithm has the potential to accurately and efficiently solve many fluid flow and heat transfer problems in complex domains. An open source code Meshless Multi-Physics Software (MeMPhyS) is available for interested users of the algorithm.
The study of chemical reactions in aqueous media is very important for its implications in several fields of science, from biology to industrial processes. Modelling these reactions is however difficult when water directly participates in the reaction. Since it requires a fully quantum mechanical description of the system, $\textit{ab-initio}$ molecular dynamics is the ideal candidate to shed light on these processes. However, its scope is limited by a high computational cost. A popular alternative is to perform molecular dynamics simulations powered by machine learning potentials, trained on an extensive set of quantum mechanical calculations. Doing so reliably for reactive processes is difficult because it requires including very many intermediate and transition state configurations. In this study, we used an active learning procedure accelerated by enhanced sampling to harvest such structures and to build a neural-network potential to study the urea decomposition process in water. This allowed us to obtain the free energy profiles of this important reaction in a wide range of temperatures, to discover a number of novel metastable states and to improve the accuracy of the kinetic rates calculations. Furthermore, we found that the formation of the zwitterionic intermediate has the same probability of occurring via an acidic or a basic pathway, which could be the cause of the insensitivity of reaction rates to the pH solution.
Meshless solution to differential equations using radial basis functions (RBF) is an alternative to grid based methods commonly used. Since the meshless method does not need an underlying connectivity in the form of control volumes or elements, issues such as grid skewness that adversely impact accuracy are eliminated. Gaussian, Multiquadrics and inverse Multiquadrics are some of the most popular RBFs used for the solutions of fluid flow and heat transfer problems. But they have additional shape parameters that have to be fine tuned for accuracy and stability. Moreover, they also face stagnation error when the point density is increased for accuracy. Recently, Polyharmonic splines (PHS) with appended polynomials have been shown to solve the above issues and give rapid convergence of discretization errors with the degree of appended polynomials. In this research, we extend the PHS-RBF method for the solution of incompressible Navier-Stokes equations. A fractional step method with explicit convection and explicit diffusion terms is combined with a pressure Poisson equation to satisfy the momentum and continuity equations. Systematic convergence tests have been performed for five model problems with two of them having analytical solutions. We demonstrate fast convergence both with refinement of number of points and degree of appended polynomials. The method is further applied to solve problems such as lid-driven cavity and vortex shedding over circular cylinder. We have also analyzed the performance of this approach for solution of Euler equations. The proposed method shows promise to solve fluid flow and heat transfer problems in complex domains with high accuracy.
We present an implementation of phaseless Auxiliary-Field Quantum Monte Carlo (ph-AFQMC) utilizing graphical processing units (GPUs). The AFQMC method is recast in terms of matrix operations which are spread across thousands of processing cores and are executed in batches using custom Compute Unified Device Architecture kernels and the hardware-optimized cuBLAS matrix library. Algorithmic advances include a batched Sherman-Morrison-Woodbury algorithm to quickly update matrix determinants and inverses, density-fitting of the two-electron integrals, an energy algorithm involving a high-dimensional precomputed tensor, and the use of single-precision floating point arithmetic. These strategies result in dramatic reductions in wall-times for both single- and multi-determinant trial wavefunctions. For typical calculations we find speed-ups of roughly two orders of magnitude using just a single GPU card. Furthermore, we achieve near-unity parallel efficiency using 8 GPU cards on a single node, and can reach moderate system sizes via a local memory-slicing approach. We illustrate the robustness of our implementation on hydrogen chains of increasing length, and through the calculation of all-electron ionization potentials of the first-row transition metal atoms. We compare long imaginary-time calculations utilizing a population control algorithm with our previously published correlated sampling approach, and show that the latter improves not only the efficiency but also the accuracy of the computed ionization potentials. Taken together, the GPU implementation combined with correlated sampling provides a compelling computational method that will broaden the application of ph-AFQMC to the description of realistic correlated electronic systems.
Mani Chandra, Gitansh Kataria, Deshdeep Sahdev
et al.
Electron transport in clean 2D systems with weak electron-phonon (e-ph) coupling can transition from an Ohmic to a ballistic or a hydrodynamic regime. The ballistic regime occurs when electron-electron (e-e) scattering is weak whereas the hydrodynamic regime arises when this scattering is strong. Despite this difference, we find that vortices and a negative nonlocal resistance believed to be quintessentially hydrodynamic are equally characteristic of the ballistic regime. These non-Ohmic regimes cannot be distinguished in DC transport without changing experimental conditions. Further, as our kinetic calculations show, the hydrodynamic regime in DC transport is highly fragile and is wiped out by even sparse disorder and e-ph scattering. We show that microwave-frequency AC sources by contrast readily excite hydrodynamic modes with current vortices that are robust to disorder and e-ph scattering. Indeed, current reversals in the non-Ohmic regimes occur via repeated vortex generation and mergers through reconnections, as in classical 2D fluids. Crucially, AC sources give rise to strong correlations across the entire device that unambiguously distinguish all regimes. These correlations in the form of nonlocal current-voltage and voltage-voltage phases directly check for the presence of a nonlocal current-voltage relation signifying the onset of non-Ohmic behavior as well as also for the dominance of bulk interactions, needed to confirm the presence of a hydrodynamic regime. We use these probes to demarcate all regimes in an experimentally realizable graphene device and find that the ballistic regime has a much larger extent in parameter space than the hydrodynamic regime.
We use determinant quantum Monte Carlo to study the single particle properties of quasiparticles and phonons in a variant of the two-dimensional Holstein model that includes an additional non-linear electron-phonon (e-ph) interaction. We find that a small positive non-linear interaction reduces the effective coupling between the electrons and the lattice, suppresses charge-density wave (CDW) correlations, and hardens the effective phonon frequency. Conversely, a small negative non-linear interaction can enhance the e-ph coupling resulting in heavier quasiparticles, an increased tendency towards a CDW phase at all fillings, and a softened phonon frequency. An effective linear model with a renormalized interaction strength and phonon frequency can qualitatively capture this physics; however, the quantitative effects of the non-linearity on both the electronic and phononic degrees of freedom cannot be captured by such a model. These results are significant for typical non-linear coupling strengths found in real materials, indicating that non-linearity can have a significant influence on the physics of many e-ph coupled systems.
DDSCAT 7.3 is an open-source Fortran-90 software package applying the discrete dipole approximation to calculate scattering and absorption of electromagnetic waves by targets with arbitrary geometries and complex refractive index. The targets may be isolated entities (e.g., dust particles), but may also be 1-d or 2-d periodic arrays of "target unit cells", allowing calculation of absorption, scattering, and electric fields around arrays of nanostructures. The theory of the DDA and its implementation in DDSCAT is presented in Draine (1988) and Draine & Flatau (1994), and its extension to periodic structures in Draine & Flatau (2008), and efficient near-field calculations in Flatau & Draine (2012). DDSCAT 7.3 includes support for MPI, OpenMP, and the Intel Math Kernel Library (MKL). DDSCAT supports calculations for a variety of target geometries. Target materials may be both inhomogeneous and anisotropic. It is straightforward for the user to "import" arbitrary target geometries into the code. DDSCAT automatically calculates total cross sections for absorption and scattering and selected elements of the Mueller scattering intensity matrix for user-specified scattering directions. DDSCAT 7.3 can efficiently calculate E and B throughout a user-specified volume containing the target. This User Guide explains how to use DDSCAT 7.3 to carry out electromagnetic scattering calculations, including use of DDPOSTPROCESS, a Fortran-90 code to perform calculations with E and B at user-selected locations near the target. A number of changes have been made since the last release, DDSCAT 7.2 .
Vladimír Chalupecký, Tasnim Fatima, Jens Kruschwitz
et al.
We consider a two-scale reaction diffusion system able to capture the corrosion of concrete with sulfates. Our aim here is to define and compute two macroscopic corrosion indicators: typical pH drop and gypsum profiles. Mathematically, the system is coupled, endowed with micro-macro transmission conditions, and posed on two different spatially-separated scales: one microscopic (pore scale) and one macroscopic (sewer pipe scale). We use a logarithmic expression to compute values of pH from the volume averaged concentration of sulfuric acid which is obtained by resolving numerically the two-scale system (microscopic equations with direct feedback with the macroscopic diffusion of one of the reactants). Furthermore, we also evaluate the content of the main sulfatation reaction (corrosion) product---the gypsum---and point out numerically a persistent kink in gypsum's concentration profile. Finally, we illustrate numerically the position of the free boundary separating corroded from not-yet-corroded regions.
The low energy scattering of heavy positively charged particles on hydrogen atoms (H) are investigated by solving the Faddeev equations in configuration space. A resonant value of the pH scattering length, $a=750\pm 5$ a.u., in the pp antisymmetric state was found. This large value indicates the existence of a first excited state with a binding energy B=1.14$\times10^{-9}$ a.u. below the H ground state. Several resonances for non zero angular momenta states are predicted.