For quantum Monte Carlo simulations of molecular systems or supercells with thousands of electrons, matrix operations related to Slater determinants lead the computational cost. McDaniel et al. [ , 174107] proposed a delayed update algorithm to increase computational efficiency by using matrix-matrix multiplication when updating the inverse matrices of Slater determinants. However, preparing intermediate matrices for applying the Sherman-Morrison-Woodbury formula remained a bottleneck. In this work, we introduce an improved algorithm for CPUs and GPUs that (1) reduces this bottleneck by iteratively updating the intermediate matrices and (2) is efficient at any acceptance ratio, with no cost for rejected moves on CPUs and minimal cost on GPUs. We show the full scheme of integrating the delayed update algorithm into a single-electron move. The high efficiency of our algorithm is demonstrated on CPUs and GPUs for a 512 atom/6144 valence electron calculation, with 12× and 2× overall speed-up compared to traditional rank-1 update schemes in diffusion quantum Monte Carlo, respectively.

In the development of coarse-grained (CG) models, accurately capturing the electrostatic interactions between CG particles is crucial for enhancing the accuracy of these models. In this study, we developed a shell-space electrostatic potential (SS-ESP) method to derive the partial charges for CG particles by fitting the ESP of all-atom (AA) models. To validate the accuracy of the SS-ESP method, we systematically tested and compared the ESPs and electric dipole moments of CG charges obtained using four different methods, with AA models of amino acids serving as the reference. Based on the SS-ESP charges, we proposed a new Side-chain Reduced CG (SRCG) model for amino acids. The parametrization of the new SRCG force field follows a bottom-up strategy, where the bonding interaction parameters are inherited from the AA model, the CG charges are derived using the SS-ESP method, and the Lennard-Jones (LJ) parameters are fitted to reproduce the static potential for van der Waals (vdW) interactions. Molecular dynamics simulations of amino acids show that the test results for dipeptides and tripeptides based on this SRCG model are consistent with those obtained from the corresponding AA models. Furthermore, the SRCG force field can accurately describe the solvation free energy of amino acids and preserve the stability of protein secondary structures. Hence, the present SS-ESP charge fitting method provides a new framework for accurately describing the electrostatic interactions of CG particles for bottom-up development of protein CG models.

An overwhelming majority of molecules remain unexplored. This is mostly due to the sheer number of them, which prohibits any enumeration of chemical space, the set of all such molecules. In practice, only subsets of chemical space are considered, but those subsets exhibit substantial bias, prohibiting the data-driven characterization of chemical space itself. In this work, we provide a method to produce unbiased representative random samples of the chemical space without enumeration of constituent molecules and to estimate the number of molecules in any custom chemical space. The approach is applicable to molecules that can be represented as graphs and runs efficiently even for molecules of 30 atoms. We use it to estimate the representativeness of current databases with respect to their underlying chemical space and establish a necessary criterion for a lower bound of database sizes to be representative of an underlying chemical space.

We present GSCDB137, a rigorously curated benchmark library of 137 data sets (8377 entries) covering main-group and transition-metal reaction energies and barrier heights, (intra- and intermolecular) noncovalent interactions, dipole moments, polarizabilities, electric-field response energies, and vibrational frequencies. Legacy data from GMTKN55 and MGCDB84 have been updated to today's best reference values; redundant or low-quality points were removed, and many new, property-focused sets were added. Testing 29 popular density functional approximations (DFAs) confirms the expected Jacob's-ladder hierarchy overall but also reveals notable exceptions: functional performance for frequencies and electric-field properties correlates poorly with that for other ground-state energetics. ωB97M-V and ωB97X-V are the most balanced hybrid meta-GGA and hybrid GGA, respectively; B97M-V and revPBE-D4 lead the meta-GGA and GGA classes. Double hybrids lower mean errors by about 30% versus their hybrid analogues but demand careful frozen-core, basis set, and spin contamination treatment. GSCDB137 offers a comprehensive, openly documented platform for rigorous validation of DFA and universal machine learning potentials, and training of the next generation of exchange-correlation functionals.

Accurately modeling the electronic structure of systems with many unpaired electrons remains a major challenge in quantum chemistry. Qualitatively correct electronic structures generally require large active space multireference wave functions, while dynamic correlation effects beyond the active space are crucial for quantitatively accurate descriptions of magnetic, catalytic and optical properties of such systems. Here, we present an uncontracted multireference perturbation theory based on the FCIQMC imaginary-time evolution of effective Hamiltonians, built upon the generalized active space concept and Löwdin's partitioning technique. The configurational interaction space is split into a reference space, consisting of the most important configurations, and a perturber space, containing the more numerous configurations responsible for dynamic correlation effects. The generalized active space algorithm allows the flexible partitioning of the configurational space. Löwdin's partitioning technique is then used to construct an effective Hamiltonian which is stochastically solved. This strategy allows us to apply perturbative corrections on large active space reference wave functions, without requiring high-order reduced density matrices, which have been found the bottleneck in other perturbation theory strategies. The capabilities of the resulting method, called Stochastic-SplitGAS, are demonstrated on the triplet-quintet spin gap of an Fe(II)-porphyrin model system and the spin ladder of a [Fe(III)S] complex.

Relative alchemical binding free energy calculations can be used to predict the effect of amino acid mutations on ligand binding affinities. However, these protocols are not well established for proteins containing intrinsically disordered regions (IDRs). In this work, we focus on the development of robust protein-free energy perturbation (FEP) protocols to reproduce experimental binding affinities that have been measured for a panel of mutants of the protein MDM2 against two ligands, AM-7209 and Nutlin-3a. We focus on mutations that occur in the N-terminal IDR lid of MDM2, which is known to undergo ligand-dependent folding upon binding. We systematically assess the effectiveness of both equilibrium and nonequilibrium alchemical protocols in reproducing these experimental binding affinities, in particular for mutations with slowly varying degrees of freedom. We show that the equilibrium protocol outperforms the nonequilibrium protocol in the precision of the free energy estimates obtained. In addition, we demonstrate the effect of the protein force field and the water model used to simulate the highly flexible IDR region. Overall, our findings demonstrate an accurate FEP protocol capable of reproducing these trends and further show the applicability of FEP protocols for elucidating the mutational effects on ligand binding affinity in highly dynamic intrinsically disordered protein regions.

We present an influence functional path integral framework for treating the coupled dynamics of solvated proton and electron transfer within a nonequilibrium open system. This method formulates a generalized Langevin equation describing dynamics in systems where proton is simultaneously coupled to fermionic (metal electrons) and bosonic (solvent phonons) reservoirs. It accounts for multiple dissipative channels without relying on phenomenological assumptions. With this scheme, we capture the relaxation of oscillations associated with large quantum zero-point fluctuations when protons are trapped in a harmonic potential. When the proton's translational motion is slow, the dynamics become effectively Markovian. In this regime, dissipation to the electronic reservoir is characterized by a position-dependent electronic friction. Using an effective electronic model Hamiltonian, we demonstrate that electronic friction introduces a sharp, localized resistance when the proton level crosses the Fermi level, effectively delaying the reaction. In contrast, solvent friction arising from assumed Caldeira-Leggett-type coupling, exerts a uniform, position-independent drag. Both mechanisms contribute comparable amounts to the overall energy dissipation. This framework offers a computationally efficient route to simulate complex electrochemical environments involving multiple dissipative baths.

Powered by dramatic advances in computer hardware, the advent of ultralarge make-on-demand virtual libraries, and a shift in small-molecule discovery toward more challenging targets with limited known actives, there has been a growing interest in the development of performant virtual screening methods that can reliably deliver novel hits. We report on a new method called Glide WS, that builds on our earlier efforts (WScore) to introduce an explicit representation of water structure and dynamics to an empirical scoring function suitable for high-throughput docking. This scoring function has been carefully tuned using absolute binding free energy perturbation calculations (ABFEP). Compared with Glide SP, Glide WS offers significant gains in the two primary tasks for molecular docking in drug discovery, pose prediction and virtual screening enrichment. For docking accuracy, Glide WS achieves a self-docking accuracy of 92% on a diverse set of 1477 protein ligand complexes as compared to 85% for Glide SP, using a criterion of 2.5 Å. We also demonstrate significantly improved virtual screening enrichment using a diverse data set covering of 38 targets together with three different computationally generated libraries of decoys, combined with standard known ChEMBL actives. We focus on ligands ranked in the top few percent of the database (the subset that is relevant to practical virtual screening efforts) and demonstrate that, along with improved enrichment of ChEMBL actives, Glide WS achieves a remarkable reduction in the number of poorly scoring decoys (as calibrated by ABFEP calculations), across a high percentage of targets, as compared to Glide SP. These results suggest that considerably higher hit rates will be observed, as compared to conventional rigid receptor docking, in practical virtual screening applications.

We present a gridless framework for computing high-dimensional conformational free energy surfaces (FES) of flexible molecules using enhanced sampling trajectories. By combining concurrent well-tempered metadynamics with Density Peaks Advanced (DPA) clustering, our approach bypasses the dimensionality limitations of conventional grid-based FES reconstruction. Free energies are assigned on a per-configuration basis via local density estimation and Zwanzig reweighting, allowing for a direct, resolution-independent mapping of the conformational ensemble. Conformers are identified as density peaks in torsional angle space, and convergence is assessed via systematic consistency metrics. We validate this approach by reproducing the paradigmatic FES of alanine dipeptide and extend it to explore molecules with 4-, 7-, and 11-dimensional torsional angle spaces. As a key application, we investigate the solvent-dependent conformational preferences of bicalutamide in vacuum, chloroform, and DMSO. The predicted global minima reflect the known solvent-induced conformational shift between open and closed forms, in agreement with NMR and crystallographic data. These results demonstrate that our workflow provides a scalable route to high-dimensional conformational free energy landscapes, with direct relevance for polymorphism, solvation, and drug design.

Characterizing electronic thermal states at low temperatures is an important but challenging task in quantum chemistry and condensed matter physics, making it a prime candidate for a useful application in quantum computing. One of the most successful methods for state preparation on quantum computers is the Adaptive, Problem-Tailored (ADAPT) Variational Quantum Eigensolver (VQE), which has recently been generalized to treat excited states within a state-averaged framework as well as Gibbs states. In this work, we introduce Helmholtz-Optimized Thermal (HOT) ADAPT-VQE, an ancilla-free strategy for preparing Gibbs states that directly minimizes the Helmholtz free energy by targeting the dominant eigenstates of the thermal ensemble. We demonstrate the usefulness of HOT-ADAPT-VQE by predicting the free energy of two model systems with strongly correlated ground states: (1) the Fe cation in a magnetic field and (2) a [CuO] fragment of the Mott insulator LaCuO. Our results demonstrate that HOT-ADAPT-VQE significantly improves upon Gibbs-state estimates from multistate variants of ADAPT-VQE, often with substantially shallower quantum circuits, making it a promising candidate for thermal-state calculations.

Extensive benchmarks and reviews of time-dependent density functional theory (TDDFT) have been published, covering at least 50 functionals. Here, we do not use one of the (meta)GGA or hybrid functionals but highlight the particular TDDFT method that does not use the Kohn-Sham (KS) potential of some exchange-correlation functional, but uses the exact Kohn-Sham potential or a close approximation to it in the SCF calculations. Such KS-potential-based TDDFT results have been proven to yield excellent results. For routine application, it is required that a computationally simple approximation (as a density functional) to the exact KS potential is available. In this paper, we benchmark TDDFT calculations with a recently developed model KS potential and compare them to advanced quantum chemical methods and experimental data. The target systems here are medium-sized molecules. These TDDFT calculations based on good KS potentials prove to be competitive in accuracy with sophisticated methods. An advantage is the possibility to use large basis sets (also for large molecules), enabling a description of valence excitations and Rydberg (or mixed valence-Rydberg) excitations on the same footing. This is necessary for high accuracy; a sophisticated but expensive method that cannot handle large basis sets cannot achieve high accuracy. An advantage of the use of a (close to) exact Kohn-Sham potential is the realistic nature (shape and energy) of both the occupied and the virtual KS orbitals, affording an interpretation of excitations in terms of one or a few single orbital-to-orbital transitions. This obviates the need for extensive optimization of the virtual orbitals for the purpose of interpreting excitations. To highlight the importance of the realistic nature of the KS virtuals, including the Rydberg orbitals, we discuss and reassess the nature of the 2 state of s-trans-1,3-butadiene, which has been widely considered as a prototype excitation with large double excitation character. Charge-transfer and true double excitations cannot be handled by the simple ALDA kernel used here.

Coarse-grained (CG) molecular dynamics (MD) simulations have emerged as a powerful and cost-effective approach for modeling materials by simplifying atomic structures into CG beads. However, accurately parametrizing interatomic potential models (force fields, FFs) that can reliably reproduce material properties and quantifying the uncertainties associated with both the model parameters and their predictions remains a major challenge. In this study, we developed coarse-grained embedded atom method (CG EAM) potentials to model interatomic interactions in face-centered cubic (FCC) metals, including palladium (Pd), gold (Au), silver (Ag), copper (Cu), and platinum (Pt). The CG EAM potentials combine the physical interpretability of a traditional EAM with the computational efficiency of coarse-graining. We first employed a Particle Swarm Optimization (PSO) framework integrated with CG MD simulations to explore a 14-dimensional parameter space and identify CG EAM parameters that reproduce key physical, mechanical, and thermodynamic properties, such as cohesive energy, lattice constants, and elastic moduli. These parameters were subsequently refined using a Bayesian uncertainty quantification (BUQ) approach, which allowed the systematic assessment of uncertainties in both the FF parameters and the predicted properties. For all five metals, this framework yielded robust parameter ranges within which the predicted properties generally remained within their 95% confidence intervals. Overall, this integrated parameter optimization and BUQ approach provides an effective strategy for developing accurate and reliable interatomic potentials while offering a generalizable framework for designing both hard and soft materials with targeted properties.

Modern spectroscopic techniques enable the determination of the spacing between rovibrational levels of H with a relative accuracy of approximately 10. At this extreme level of precision, subtle quantum electrodynamic (QED) effects, such as electron self-interaction and vacuum polarization, are probed. A theoretical model aiming to achieve similar accuracy must precisely describe not only these relatively small QED effects but also the more significant contributions related to electron correlation, coupling between electronic and nuclear motions, and relativistic effects. Although the hydrogen molecule exhibits most of the phenomena found in larger molecules, it is simple enough to meet the requirements mentioned above. In this article, we report on enhancements to the current capabilities of quantum mechanical calculations for the hydrogen molecule. We present a method based on exponential functions that fully captures electron correlation or, more broadly, interparticle correlation, enabling a comprehensive description of effects related to nuclear motion. Specifically, we solve the four-particle Schrödinger equation without invoking commonly used approximations such as the one-electron or the Born-Oppenheimer approximation. The only source of nonrelativistic energy error comes from the finite size of the basis set. The explicitly correlated nonadiabatic wave function used here is then employed to determine the relativistic and QED effects. As a result, the dissociation energy for the lowest rovibrational levels in the electronic ground state of H has been obtained with a relative accuracy of 7 × 10, while the frequencies of intervals between these levels have been determined with sub-MHz accuracy, corresponding to a relative accuracy of 3 × 10. In consequence, the discrepancies between the highest precision measurements and earlier theoretical predictions have been resolved.

Markov state models (MSMs) are widely employed to analyze the kinetics of complex systems. But despite their effectiveness in many applications, MSMs are prone to systematic or statistical errors, often exacerbated by suboptimal hyperparameter choice. In this article, we attempt to understand how these choices affect the error of estimates of mean first-passage times and committors, key quantities in chemical rate theory. We first evaluate the performance of the recently introduced "stopped-process estimator" [Strahan, J. Long-time-scale predictions from short-trajectory data: A benchmark analysis of the trp-cage miniprotein. 2021, 17, 2948-2963. 10.1021/acs.jctc.0c00933.] that attempts to reduce error caused by choosing a too-large lag time. We then study the effect of statistical errors on Markov state model construction using the condition number, which measures an MSM's sensitivity to perturbation. This analysis helps give an insight into which factors cause an MSM to be more or less sensitive to statistical error. Our work highlights the importance of choosing a good sampling measure, the measure from which the initial points are drawn, and has implications for recent work applying a variational principle for evaluating the committor.

Nanoscale systems often contain weakly coupled components, as exemplified by layered materials. Time-domain atomistic modeling of excited state processes in such systems with nonadiabatic (NA) molecular dynamics (MD) runs into severe challenges due to the divergence of the NA coupling. At the same time, standard NAMD methods work well within each component. We develop an efficient ab initio NAMD methodology using a mixed diabatic-adiabatic representation (dNAMD), implement decoherence-induced surface hopping (DISH) within the dNAMD framework, and demonstrate its utility with long-range charge transfer in 2D perovskites taking place on nano- to microsecond time scales. The dNAMD method bypasses the trivial state crossing issue of traditional NAMD by using a diabatization technique to derive diabatic electronic coupling integrals between weakly coupled components, while employing adiabatic representation within each component. We demonstrate the approach by application to 2D perovskites, which are promising materials for optoelectronic applications, but show limited efficiencies because of the insulating nature of organic spacer cations and slow interlayer charge transport. The interlayer charge transfer time scales predicted by DISH-dNAMD are consistent with experimental data and Marcus rate constants. The simulations show that phenethylammonium spacers enhance inorganic lattice rigidity via strong hydrogen bonding and π-π stacking interactions, and reduce electron-vibrational coupling while increasing interlayer spacing and charge localization. These effects significantly reduce the electronic couplings, yielding charge transfer rates that are 1-2 orders of magnitude lower than those for the more structurally flexible butylammonium spacers. The DISH-dNAMD simulations highlight the critical role of the spacer rigidity in the interlayer charge transport of 2D perovskites. The developed dNAMD framework provides an efficient and versatile tool for simulating and elucidating excited state dynamics in weakly coupled nanoscale and condensed phase systems at the atomistic level and in the time domain as it occurs in nature and experiments, advancing the design of next-generation optoelectronic devices.

The calculation of rate constants for proton-coupled electron transfer (PCET) reactions is a challenging task in quantum chemistry. This task involves identifying the mechanism of the process, that can take place either adiabatically or nonadiabatically, and calculating the necessary quantities, such as vibronic couplings, according to the mechanism. Due to different electronic configurations involved, it almost becomes inevitable to use wave function-based multireference methods. However, the accurate prediction of rate constants for large molecular systems is limited by the high computational cost of common choices such as complete active space self-consistent field (CASSCF). Since PCET reactions occur in a wide range of biological processes, the development of alternatives with large scalability is of particular interest. A promising alternative is the multistate density-functional theory method based on block-localized Kohn-Sham (BLKS) orbitals. This gives access to the diabatic donor and acceptor states, adiabatic ground and first excited states, and the electronic coupling. In this work, different operators for the construction of BLKS orbitals are considered. A comparison with CASSCF and -electron valence state second-order perturbation theory shows that accurate vibronic couplings can be obtained using a non-Hermitian operator. As the method relies on a fragmentation of the system, spectator fragments can be explicitly included with a convenient computational cost. This is demonstrated by the example of a DNA-acrylamide complex.

We formulate and implement the core-valence separated multireference equation-of-motion-driven similarity renormalization group method (CVS-IP-EOM-DSRG) for simulating X-ray photoelectron spectra (XPS) of strongly correlated molecular systems. This method is numerically robust and computationally efficient, delivering accurate core-ionization energies with () scaling relative to basis set size in the EOM step. To ensure rigorous core intensivity, we propose a simple modification of the ground-state MR-DSRG formalism. We develop and compare three variants of the theory based on different approximations of the effective Hamiltonian: two derived from low-order perturbative methods (DSRG-MRPT2 and DSRG-MRPT3) and one from a nonperturbative scheme truncated to 1- and 2-body operators [MR-LDSRG(2)]. We benchmark the CVS-IP-EOM-DSRG methods by computing vertical core-ionization energies for a representative molecular test set and comparing the results against the established single-reference and multireference methods. To demonstrate the applicability of CVS-IP-EOM-DSRG to strongly correlated systems, we compute the potential energy curves and vibrationally resolved XPS of N and CO and the XPS of ozone. Comparison with experimental data and other high-level theoretical results shows that all three CVS-IP-EOM-DSRG variants accurately predict vertical ionization energies but only those based on the DSRG-MRPT3 and MR-LDSRG(2) levels of theory reliably capture the full dissociation behavior and reproduce the experimental vibrational structure.

Double excitations in organic molecules have garnered significant interest as a result of their importance in singlet fission and photophysics. These excitations play a crucial role in understanding the photoexcitation processes in polyenes. To describe photoexcited states with both single and double excitation character, we use a first-principles many-body theory that combines the GW/Bethe-Salpeter equation and the configuration interaction (CI) methods. Specifically, we develop and employ two CI-based methods: screened configuration interaction singles and doubles (scrCISD) and screened configuration interaction singles with perturbative doubles (scrCIS(D)), applied to an effective many-body Hamiltonian that incorporates screening. We apply these methods to Thiel's set of molecules, which exhibit excited states predominantly characterized by single excitations with a partial double excitation character. Our results indicate that the scrCISD method systematically underestimates the excitation energies compared to the best theoretical estimates, while the scrCIS(D) method shows good agreement with these estimates. Furthermore, we used the scrCISD method to calculate the binding energies of the dominantly doubly excited correlated triplet pair states, TT, in pentacene dimers, finding that the TT binding energies agree well with empirical calculations.

An economic modeling approach for cavity quantum electrodynamics is provided by mean-field dynamics, wherein the optical field is described classically while a self-consistent interaction with quantum emitters is incorporated through the Ehrenfest theorem. However, conventional implementations of mean-field dynamics are known to suffer from a catastrophic leakage of zero-point energy, to lose accuracy in the short-cavity limit, and to require large numbers of trajectories to be sampled. Here, we address these three shortcomings within a single integrated approach. This approach builds on our recently proposed modification of the Ehrenfest theorem, referred to as decoupled mean-field (DC-MF) dynamics, in combination with a focused sampling scheme that enforces zero-point energy at the single-trajectory level. The approach is shown to yield high accuracy in both short and long-cavity limits while reaching convergence within a minimal amount of trajectories.

We present a new nonadiabatic event sampling method that combines the fewest switches surface hopping (FSSH) method with enhanced sampling technology. In this new method, a biasing function is introduced to drive dynamic trajectories to cross high-energy barriers and toward product regions on potential energy surfaces, enabling this method to effectively sample rare nonadiabatic transition events and different decay channels between reactants and products. The biasing function also provides additional kinetic energy to the system, thereby accelerating the simulations and reducing the required computational cost. Two types of biasing functions, based on bond length differences and root-mean-square deviation, are designed to accommodate different reaction systems. The test results indicate that the new method can effectively sample rare nonadiabatic transition events and decay channels that are difficult to capture using the standard FSSH method.

Molecular platforms for optically addressable spin states are emerging as fascinating alternatives to solid-state spin centers, offering scalable synthesis, structural tunability, and chemical versatility. Here, we present a molecular design strategy for achieving photoinduced spin polarization in organic diradicals bridged by systems featuring an inverted singlet-triplet (InveST) energy gap. These InveST units possess HOMO and LUMO orbitals localized on complementary atomic sites. By covalently linking the non-SOMO-bearing positions of alternant hydrocarbon radicals to the LUMO-localized atoms of the InveST bridge, we construct diradicals in which the radical centers remain electronically decoupled in the ground state, yielding degenerate singlet and triplet configurations. Upon photoexcitation, the population of the InveST LUMO activates an excited-state exchange interaction between the radicals, generating a finite singlet-triplet gap and enabling spin-selective intersystem crossing to polarized triplet states. Using a combination of model Hamiltonians and multireference ab initio calculations, we establish design principles for tuning exchange interactions and spin-orbit coupling to achieve molecular-level control over optical-spin interfaces. The resulting InveST-bridged diradicals have emerged as promising scaffolds for molecular quantum technologies.