Efficient quantum analytic nuclear gradients with double factorization.

Efficient representations of the Hamiltonian, such as double factorization, drastically reduce the circuit depth or the number of repetitions in error corrected and noisy intermediate-scale quantum (NISQ) algorithms for chemistry. We report a Lagrangian-based approach for evaluating relaxed one- and two-particle reduced density matrices from double factorized Hamiltonians, unlocking efficiency improvements in computing the nuclear gradient and related derivative properties. We demonstrate the accuracy and feasibility of our Lagrangian-based approach to recover all off-diagonal density matrix elements in classically simulated examples with up to 327 quantum and 18 470 total atoms in QM/MM simulations with modest-sized quantum active spaces. We show this in the context of the variational quantum eigensolver in case studies, such as transition state optimization, ab initio molecular dynamics simulation, and energy minimization of large molecular systems.

Rank-reduced coupled-cluster. III. Tensor hypercontraction of the doubles amplitudes.

We develop a quartic-scaling implementation of coupled-cluster singles and doubles (CCSD) based on low-rank tensor hypercontraction (THC) factorizations of both the electron repulsion integrals (ERIs) and the doubles amplitudes. This extends our rank-reduced (RR) coupled-cluster method to incorporate higher-order tensor factorizations. The THC factorization of the doubles amplitudes accounts for most of the gain in computational efficiency as it is sufficient, in conjunction with a Cholesky decomposition of the ERIs, to reduce the computational complexity of most contributions to the CCSD amplitude equations. Further THC factorization of the ERIs reduces the complexity of certain terms arising from nested commutators between the doubles excitation operator and the two-electron operator. We implement this new algorithm using graphical processing units and demonstrate that it enables CCSD calculations for molecules with 250 atoms and 2500 basis functions using a single computer node. Furthermore, we show that the new method computes correlation energies with comparable accuracy to the underlying RR-CCSD method.

An ab initio exciton model for singlet fission.

We present an ab initio exciton model that extends the Frenkel exciton model and includes valence, charge-transfer, and multiexcitonic excited states. It serves as a general, parameter-free, yet computationally efficient and scalable approach for simulation of singlet fission processes in multichromophoric systems. A comparison with multiconfigurational methods confirms that our exciton model predicts consistent energies and couplings for the pentacene dimer and captures the correct physics. Calculations of larger pentacene clusters demonstrate the computational scalability of the exciton model and suggest that the mixing between local and charge-transfer excitations narrows the gap between singlet and multiexcitonic states. Local vibrations of pentacene molecules are found to facilitate singlet-multiexcitonic state-crossing and hence are important for understanding singlet fission. The exciton model developed in this work also sets the stage for further implementation of the nuclear gradients and nonadiabatic couplings needed for first principles nonadiabatic quantum molecular dynamics simulations of singlet fission.

Psi4 1.4: Open-source software for high-throughput quantum chemistry.

PSI4 is a free and open-source ab initio electronic structure program providing implementations of Hartree-Fock, density functional theory, many-body perturbation theory, configuration interaction, density cumulant theory, symmetry-adapted perturbation theory, and coupled-cluster theory. Most of the methods are quite efficient, thanks to density fitting and multi-core parallelism. The program is a hybrid of C++ and Python, and calculations may be run with very simple text files or using the Python API, facilitating post-processing and complex workflows; method developers also have access to most of PSI4's core functionalities via Python. Job specification may be passed using The Molecular Sciences Software Institute (MolSSI) QCSCHEMA data format, facilitating interoperability. A rewrite of our top-level computation driver, and concomitant adoption of the MolSSI QCARCHIVE INFRASTRUCTURE project, makes the latest version of PSI4 well suited to distributed computation of large numbers of independent tasks. The project has fostered the development of independent software components that may be reused in other quantum chemistry programs.

Quantum Computation of Electronic Transitions Using a Variational Quantum Eigensolver.

We develop an extension of the variational quantum eigensolver (VQE) algorithm-multistate contracted VQE (MC-VQE)-that allows for the efficient computation of the transition energies between the ground state and several low-lying excited states of a molecule, as well as the oscillator strengths associated with these transitions. We numerically simulate MC-VQE by computing the absorption spectrum of an ab initio exciton model of an 18-chromophore light-harvesting complex from purple photosynthetic bacteria.

Observation of Ultrafast Intersystem Crossing in Thymine by Extreme Ultraviolet Time-Resolved Photoelectron Spectroscopy.

We studied the photoinduced ultrafast relaxation dynamics of the nucleobase thymine using gas-phase time-resolved photoelectron spectroscopy. By employing extreme ultraviolet pulses from high harmonic generation for photoionization, we substantially extend our spectral observation window with respect to previous studies. This enables us to follow relaxation of the excited state population all the way to low-lying electronic states including the ground state. In thymine, we observe relaxation from the optically bright 1ππ* state of thymine to a dark 1nπ* state within 80 ± 30 fs. The 1nπ* state relaxes further within 3.5 ± 0.3 ps to a low-lying electronic state. By comparison with quantum chemical simulations, we can unambiguously assign its spectroscopic signature to the 3ππ* state. Hence, our study draws a comprehensive picture of the relaxation mechanism of thymine including ultrafast intersystem crossing to the triplet manifold.

Rank reduced coupled cluster theory. II. Equation-of-motion coupled-cluster singles and doubles.

Equation-of-motion coupled-cluster singles and doubles (EOM-CCSD) is a reliable and popular approach to the determination of electronic excitation energies. Recently, we have developed a rank-reduced CCSD (RR-CCSD) method that allows the ground-state coupled-cluster energy to be determined with low-rank cluster amplitudes. Here, we extend this approach to excited-state energies through a RR-EOM-CCSD method. We start from the EOM-CCSD energy functional and insert low-rank approximations to the doubles amplitudes. The result is an approximate EOM-CCSD method with only a quadratic number (in the molecular size) of free parameters in the wavefunction. Importantly, our formulation of RR-EOM-CCSD preserves the size intensivity of the excitation energy and size extensivity of the total energy. Numerical tests of the method suggest that accuracy on the order of 0.05-0.01 eV in the excitation energy is possible with 1% or less of the original number of wavefunction coefficients; accuracy of better than 0.01 eV can be achieved with about 4% or less of the free parameters. The amount of compression at a given accuracy level is expected to increase with the size of the molecule. The RR-EOM-CCSD method is a new path toward the efficient determination of accurate electronic excitation energies.

Rank reduced coupled cluster theory. I. Ground state energies and wavefunctions.

We propose a compression of the opposite-spin coupled cluster doubles amplitudes of the form τij ab≡Uia VTVWUjb W, where Uia V are the nV-highest magnitude eigenvectors of the MP2 or MP3 doubles amplitudes. Together with a corresponding parameterization of the opposite-spin coupled cluster Lagrange multipliers of the form λab ij≡Uia VLVWUjb W, this yields a fully self-consistent parameterization of reduced-rank coupled cluster equations in terms of the Lagrangian L0TVW,LVW. Making this Lagrangian stationary with respect to the LVW parameters yields a perfectly determined set of equations for the TVW equations and coupled cluster energy. These equations can be solved using a Lyapunov equation for the first-order amplitude updates. We test this "rank-reduced coupled cluster" method for coupled cluster singles and doubles in medium sized molecules and find that substantial compression of the T^2 amplitudes is possible with acceptable accuracy.

The Surprising Importance of Peptide Bond Contacts in Drug-Protein Interactions.

The study of noncovalent interactions, notably including drug-protein binding, relies heavily on the language of localized functional group contacts: hydrogen bonding, π-π interactions, CH-π contacts, halogen bonding, etc. Applying the state-of-the-art functional group symmetry-adapted perturbation theory (F-SAPT) to an important question of chloro versus methyl aryl substitution in factor Xa inhibitor drugs, we find that a localized contact model provides an incorrect picture for the origin of the enhancement of chloro-containing ligands. Instead, the enhancement is found to originate from many intermediate-range contacts distributed throughout the binding pocket, particularly including the peptide bonds in the protein backbone. The contributions from these contacts are primarily electrostatic in nature, but require ab initio computations involving nearly the full drug-protein pocket system to be accurately quantified.

Communication: A difference density picture for the self-consistent field ansatz.

We formulate self-consistent field (SCF) theory in terms of an interaction picture where the working variable is the difference density matrix between the true system and a corresponding superposition of atomic densities. As the difference density matrix directly represents the electronic deformations inherent in chemical bonding, this "difference self-consistent field (dSCF)" picture provides a number of significant conceptual and computational advantages. We show that this allows for a stable and efficient dSCF iterative procedure with wholly single-precision Coulomb and exchange matrix builds. We also show that the dSCF iterative procedure can be performed with aggressive screening of the pair space. These approximations are tested and found to be accurate for systems with up to 1860 atoms and >10 000 basis functions, providing for immediate overall speedups of up to 70% in the heavily optimized TeraChem SCF implementation.

Communication: Practical intramolecular symmetry adapted perturbation theory via Hartree-Fock embedding.

We develop a simple methodology for the computation of symmetry-adapted perturbation theory (SAPT) interaction energy contributions for intramolecular noncovalent interactions. In this approach, the local occupied orbitals of the total Hartree-Fock (HF) wavefunction are used to partition the fully interacting system into three chemically identifiable units: the noncovalent fragments A and B and a covalent linker C. Once these units are identified, the noninteracting HF wavefunctions of the fragments A and B are separately optimized while embedded in the HF wavefunction of C, providing the dressed zeroth order wavefunctions for A and B in the presence of C. Standard two-body SAPT (particularly SAPT0) is then applied between the relaxed wavefunctions for A and B. This intramolecular SAPT procedure is found to be remarkably straightforward and efficient, as evidenced by example applications ranging from diols to hexaphenyl-ethane derivatives.

Quantum-mechanical evaluation of π-π versus substituent-π interactions in π stacking: direct evidence for the Wheeler-Houk picture.

The influence of substituents on π-stacking interactions has previously been explained by two competing hypotheses: a nonlocal effect in which tuning of the π density by the substituent alters the interaction (the Hunter-Sanders picture) or a local effect in which the direct interaction of the added substituent and changed polarity of the phenyl-substituent σ bond alter the interaction (the Wheeler-Houk picture). In this work, we applied the recently developed functional-group partition of symmetry-adapted perturbation theory (F-SAPT) to directly quantify these two effects in situ. The results show that both pictures contribute to the change in interaction energy but that the Wheeler-Houk picture is usually dominant.

Spatial assignment of symmetry adapted perturbation theory interaction energy components: The atomic SAPT partition.

We develop a physically-motivated assignment of symmetry adapted perturbation theory for intermolecular interactions (SAPT) into atom-pairwise contributions (the A-SAPT partition). The basic precept of A-SAPT is that the many-body interaction energy components are computed normally under the formalism of SAPT, following which a spatially-localized two-body quasiparticle interaction is extracted from the many-body interaction terms. For electrostatics and induction source terms, the relevant quasiparticles are atoms, which are obtained in this work through the iterative stockholder analysis (ISA) procedure. For the exchange, induction response, and dispersion terms, the relevant quasiparticles are local occupied orbitals, which are obtained in this work through the Pipek-Mezey procedure. The local orbital atomic charges obtained from ISA additionally allow the terms involving local orbitals to be assigned in an atom-pairwise manner. Further summation over the atoms of one or the other monomer allows for a chemically intuitive visualization of the contribution of each atom and interaction component to the overall noncovalent interaction strength. Herein, we present the intuitive development and mathematical form for A-SAPT applied in the SAPT0 approximation (the A-SAPT0 partition). We also provide an efficient series of algorithms for the computation of the A-SAPT0 partition with essentially the same computational cost as the corresponding SAPT0 decomposition. We probe the sensitivity of the A-SAPT0 partition to the ISA grid and convergence parameter, orbital localization metric, and induction coupling treatment, and recommend a set of practical choices which closes the definition of the A-SAPT0 partition. We demonstrate the utility and computational tractability of the A-SAPT0 partition in the context of side-on cation-π interactions and the intercalation of DNA by proflavine. A-SAPT0 clearly shows the key processes in these complicated noncovalent interactions, in systems with up to 220 atoms and 2845 basis functions.

Levels of symmetry adapted perturbation theory (SAPT). I. Efficiency and performance for interaction energies.

A systematic examination of the computational expense and accuracy of Symmetry-Adapted Perturbation Theory (SAPT) for the prediction of non-covalent interaction energies is provided with respect to both method [SAPT0, DFT-SAPT, SAPT2, SAPT2+, SAPT2+(3), and SAPT2+3; with and without CCD dispersion for the last three] and basis set [Dunning cc-pVDZ through aug-cc-pV5Z wherever computationally tractable, including truncations of diffuse basis functions]. To improve accuracy for hydrogen-bonded systems, we also include two corrections based on exchange-scaling (sSAPT0) and the supermolecular MP2 interaction energy (δMP2). When considering the best error performance relative to computational effort, we recommend as the gold, silver, and bronze standard of SAPT: SAPT2+(3)δMP2/aug-cc-pVTZ, SAPT2+/aug-cc-pVDZ, and sSAPT0/jun-cc-pVDZ. Their respective mean absolute errors in interaction energy across the S22, HBC6, NBC10, and HSG databases are 0.15 (62.9), 0.30 (4.4), and 0.49 kcal mol(-1) (0.03 h for adenine·thymine complex).

Communication: Acceleration of coupled cluster singles and doubles via orbital-weighted least-squares tensor hypercontraction.

We apply orbital-weighted least-squares tensor hypercontraction decomposition of the electron repulsion integrals to accelerate the coupled cluster singles and doubles (CCSD) method. Using accurate and flexible low-rank factorizations of the electron repulsion integral tensor, we are able to reduce the scaling of the most vexing particle-particle ladder term in CCSD from O(N6) to O(N5), with remarkably low error. Combined with a T1-transformed Hamiltonian, this leads to substantial practical accelerations against an optimized density-fitted CCSD implementation.

Accurate description of torsion potentials in conjugated polymers using density functionals with reduced self-interaction error.

We investigate the torsion potentials in two prototypical π-conjugated polymers, polyacetylene and polydiacetylene, as a function of chain length using different flavors of density functional theory. Our study provides a quantitative analysis of the delocalization error in standard semilocal and hybrid density functionals and demonstrates how it can influence structural and thermodynamic properties. The delocalization error is quantified by evaluating the many-electron self-interaction error (MESIE) for fractional electron numbers, which allows us to establish a direct connection between the MESIE and the error in the torsion barriers. The use of non-empirically tuned long-range corrected hybrid functionals results in a very significant reduction of the MESIE and leads to an improved description of torsion barrier heights. In addition, we demonstrate how our analysis allows the determination of the effective conjugation length in polyacetylene and polydiacetylene chains.

Reducing the Runtime of Fault-Tolerant Quantum Simulations in Chemistry through Symmetry-Compressed Double Factorization.

Quantum phase estimation based on qubitization is the state-of-the-art fault-tolerant quantum algorithm for computing ground-state energies in chemical applications. In this context, the 1-norm of the Hamiltonian plays a fundamental role in determining the total number of required iterations and also the overall computational cost. In this work, we introduce the symmetry-compressed double factorization (SCDF) approach, which combines a CDF of the Hamiltonian with the symmetry shift technique, significantly reducing the 1-norm value. The effectiveness of this approach is demonstrated numerically by considering various benchmark systems, including the FeMoco molecule, cytochrome P450, and hydrogen chains of different sizes. To compare the efficiency of SCDF to other methods in absolute terms, we estimate Toffoli gate requirements, which dominate the execution time on fault-tolerant quantum computers. For the systems considered here, SCDF leads to a sizable reduction of the Toffoli gate count in comparison to other variants of DF or even tensor hypercontraction, which is usually regarded as the most efficient approach for qubitization.

Estimation of Electrostatic Interaction Energies on a Trapped-Ion Quantum Computer.

We present the first hardware implementation of electrostatic interaction energies by using a trapped-ion quantum computer. As test system for our computation, we focus on the reduction of NO to N2O catalyzed by a nitric oxide reductase (NOR). The quantum computer is used to generate an approximate ground state within the NOR active space. To efficiently measure the necessary one-particle density matrices, we incorporate fermionic basis rotations into the quantum circuit without extending the circuit length, laying the groundwork for further efficient measurement routines using factorizations. Measurements in the computational basis are then used as inputs for computing the electrostatic interaction energies on a classical computer. Our experimental results strongly agree with classical noise-less simulations of the same circuits, finding electrostatic interaction energies within chemical accuracy despite hardware noise. This work shows that algorithms tailored to specific observables of interest, such as interaction energies, may require significantly fewer quantum resources than individual ground state energies would require in the straightforward supermolecular approach.

Quantum-mechanical analysis of the energetic contributions to π stacking in nucleic acids versus rise, twist, and slide.

Symmetry-adapted perturbation theory (SAPT) is applied to pairs of hydrogen-bonded nucleobases to obtain the energetic components of base stacking (electrostatic, exchange-repulsion, induction/polarization, and London dispersion interactions) and how they vary as a function of the helical parameters Rise, Twist, and Slide. Computed average values of Rise and Twist agree well with experimental data for B-form DNA from the Nucleic Acids Database, even though the model computations omitted the backbone atoms (suggesting that the backbone in B-form DNA is compatible with having the bases adopt their ideal stacking geometries). London dispersion forces are the most important attractive component in base stacking, followed by electrostatic interactions. At values of Rise typical of those in DNA (3.36 Å), the electrostatic contribution is nearly always attractive, providing further evidence for the importance of charge-penetration effects in π-π interactions (a term neglected in classical force fields). Comparison of the computed stacking energies with those from model complexes made of the "parent" nucleobases purine and 2-pyrimidone indicates that chemical substituents in DNA and RNA account for 20-40% of the base-stacking energy. A lack of correspondence between the SAPT results and experiment for Slide in RNA base-pair steps suggests that the backbone plays a larger role in determining stacking geometries in RNA than in B-form DNA. In comparisons of base-pair steps with thymine versus uracil, the thymine methyl group tends to enhance the strength of the stacking interaction through a combination of dispersion and electrosatic interactions.

Exact tensor hypercontraction: a universal technique for the resolution of matrix elements of local finite-range N-body potentials in many-body quantum problems.

Configuration-space matrix elements of N-body potentials arise naturally and ubiquitously in the Ritz-Galerkin solution of many-body quantum problems. For the common specialization of local, finite-range potentials, we develop the exact tensor hypercontraction method, which provides a quantized renormalization of the coordinate-space form of the N-body potential, allowing for a highly separable tensor factorization of the configuration-space matrix elements. This representation allows for substantial computational savings in chemical, atomic, and nuclear physics simulations, particularly with respect to difficult "exchangelike" contractions.