When Gold Is Not Enough: Platinum Standard of Quantum Chemistry with N7 Cost.

In this paper, we extend the rank-reduced coupled-cluster formalism to the calculation of non-iterative energy corrections due to quadruple excitations. There are two major components of the proposed formalism. The first is an approximate compression of the quadruple excitation amplitudes using the Tucker format. The second is a modified functional used for the evaluation of the corrections which gives exactly the same results for the exact amplitudes, but is less susceptible to errors resulting from the aforementioned compression. We show, both theoretically and numerically, that the computational cost of the proposed method scales as the seventh power of the system size. Using reference results for a set of small molecules, the method is calibrated to deliver relative accuracy of a few percent in energy corrections. To illustrate the potential of the theory, we calculate the isomerization energy of ortho/meta benzyne (C6H4) and the barrier height for the Cope rearrangement in bullvalene (C10H10). The method retains a near-black-box nature of the conventional coupled-cluster formalism and depends on only one additional parameter that controls the accuracy.

Explicitly Correlated Electronic Structure Calculations with Transcorrelated Matrix Product Operators.

In this work, we present the first implementation of the transcorrelated electronic Hamiltonian in an optimization procedure for matrix product states by the density matrix renormalization group (DMRG) algorithm. In the transcorrelation ansatz, the electronic Hamiltonian is similarity-transformed with a Jastrow factor to describe the cusp in the wave function at electron-electron coalescence. As a result, the wave function is easier to approximate accurately with the conventional expansion in terms of one-particle basis functions and Slater determinants. The transcorrelated Hamiltonian in first quantization comprises up to three-body interactions, which we deal with in the standard way by applying robust density fitting to two- and three-body integrals entering the second-quantized representation of this Hamiltonian. The lack of hermiticity of the transcorrelated Hamiltonian is taken care of along the lines of the first work on transcorrelated DMRG [ J. Chem. Phys. 2020, 153, 164115] by encoding it as a matrix product operator and optimizing the corresponding ground state wave function with imaginary-time time-dependent DMRG. We demonstrate our quantum chemical transcorrelated DMRG approach at the example of several atoms and first-row diatomic molecules. We show that transcorrelation improves the convergence rate to the complete basis set limit in comparison to conventional DMRG. Moreover, we study extensions of our approach that aim at reducing the cost of handling the matrix product operator representation of the transcorrelated Hamiltonian.

Quintic-scaling rank-reduced coupled cluster theory with single and double excitations.

We consider the rank-reduced coupled-cluster theory with single and double (RR-CCSD) excitations introduced recently [Parrish et al., J. Chem. Phys. 150, 164118 (2019)]. The main feature of this method is the decomposed form of doubly excited amplitudes, which are expanded in the basis of largest magnitude eigenvectors of MP2 or MP3 amplitudes. This approach enables a substantial compression of amplitudes with only minor loss of accuracy. However, the formal scaling of the computational costs with the system size (N) is unaffected in comparison with the conventional CCSD theory (∝N6) due to the presence of some terms quadratic in amplitudes, which do not naturally factorize to a simpler form even within the rank-reduced framework. We show how to solve this problem, exploiting the fact that their effective rank increases only linearly with the system size. We provide a systematic way to approximate the problematic terms using the singular value decomposition and reduce the scaling of the RR-CCSD iterations down to the level of N5. This is combined with an iterative method of finding dominant eigenpairs of MP2 or MP3 amplitudes, which eliminates the necessity to perform the complete diagonalization, making the cost of this step proportional to the fifth power of the system size, as well. Next, we consider the evaluation of perturbative corrections to CCSD energies resulting from triply excited configurations. The triply excited amplitudes present in the CCSD(T) method are decomposed to the Tucker-3 format using the higher-order orthogonal iteration procedure. This enables us to compute the energy correction due to triple excitations non-iteratively with N6 cost. The accuracy of the resulting rank-reduced CCSD(T) method is studied for both total and relative correlation energies of a diverse set of molecules. Accuracy levels better than 99.9% can be achieved with a substantial reduction of the computational costs. Concerning the computational timings, the break-even point between the rank-reduced and conventional CCSD implementations occurs for systems with about 30-40 active electrons.

Near-Exact CCSDT Energetics from Rank-Reduced Formalism Supplemented by Non-iterative Corrections.

We introduce a non-iterative energy correction, added on top of the rank-reduced coupled-cluster method with single, double, and triple substitutions, that accounts for excitations excluded from the parent triple excitation subspace. The formula for the correction is derived by employing the coupled-cluster Lagrangian formalism, with an additional assumption that the parent excitation subspace is closed under the action of the Fock operator. Owing to the rank-reduced form of the triple excitation amplitudes tensor, the computational cost of evaluating the correction scales as N7, where N is the system size. The accuracy and computational efficiency of the proposed method is assessed for both total and relative correlation energies. We show that the non-iterative correction can fulfill two separate roles. If the accuracy level of a fraction of kJ/mol is sufficient for a given system, the correction significantly reduces the dimension of the parent triple excitation subspace needed in the iterative part of the calculations. Simultaneously, it enables reproducing the exact CCSDT results to an accuracy level below 0.1 kJ/mol, with a larger, yet still reasonable, dimension of the parent excitation subspace. This typically can be achieved at a computational cost only several times larger than required for the CCSD(T) method. The proposed method retains the black-box features of the single-reference coupled-cluster theory; the dimension of the parent excitation subspace remains the only additional parameter that has to be specified.

A systematic construction of Gaussian basis sets for the description of laser field ionization and high-harmonic generation.

A precise understanding of mechanisms governing the dynamics of electrons in atoms and molecules subjected to intense laser fields has a key importance for the description of attosecond processes such as the high-harmonic generation and ionization. From the theoretical point of view, this is still a challenging task, as new approaches to solve the time-dependent Schrödinger equation with both good accuracy and efficiency are still emerging. Until recently, the purely numerical methods of real-time propagation of the wavefunction using finite grids have been frequently and successfully used to capture the electron dynamics in small one- or two-electron systems. However, as the main focus of attoscience shifts toward many-electron systems, such techniques are no longer effective and need to be replaced by more approximate but computationally efficient ones. In this paper, we explore the increasingly popular method of expanding the wavefunction of the examined system into a linear combination of atomic orbitals and present a novel systematic scheme for constructing an optimal Gaussian basis set suitable for the description of excited and continuum atomic or molecular states. We analyze the performance of the proposed basis sets by carrying out a series of time-dependent configuration interaction calculations for the hydrogen atom in fields of intensity varying from 5 × 1013 W/cm2 to 5 × 1014 W/cm2. We also compare the results with the data obtained using Gaussian basis sets proposed previously by other authors.

Analysis of QED and non-adiabaticity effects on the rovibrational spectrum of H3 + using geometry-dependent effective nuclear masses.

The influence of QED effects (including one- and two-electron Lamb-shift, Araki-Sucher term, one-loop self-energy, and finite nuclear size correction) together with non-adiabatic effects on the rovibrational bound states of H3 + has been investigated. Non-adiabaticity is modeled by using geometry-dependent effective nuclear masses together with only one single potential energy surface. In conclusion, for rovibrational states below 20 000 cm-1, QED and relativistic effects do nearly compensate, and a potential energy surface based on Born-Oppenheimer energies and diagonal adiabatic corrections has nearly the same quality as the one including relativity with QED; the deviations between the two approaches for individual rovibrational states are mostly below 0.02 cm-1. The inclusion of non-adiabatic effects is important, and it reduces deviations from experiments mostly below 0.1 cm-1.

A straightforward a posteriori method for reduction of density-fitting error in coupled-cluster calculations.

We present a simple method for a posteriori removal of a significant fraction of the density-fitting error from the calculated total coupled-cluster energies. The method treats the difference between the exact and density-fitted integrals as a perturbation, and simplified response-like equations allow us to calculate improved amplitudes and the corresponding energy correction. The proposed method is tested at the coupled-cluster singles and doubles level of theory for a diverse set of moderately-sized molecules. On average, error reductions by a factor of approximately 10 and 20 are observed in double-zeta and triple-zeta basis sets, respectively. Similar reductions are observed in calculations of interaction energies of several model complexes. The computational cost of the procedure is small in comparison with the preceding coupled-cluster iterations. The applicability of this method is not limited to the density-fitting approximation; in principle, it can be used in conjunction with an arbitrary decomposition scheme of the electron repulsion integrals.

Implementation of the Coupled-Cluster Method with Single, Double, and Triple Excitations using Tensor Decompositions.

We report a complete implementation of the coupled-cluster method with single, double, and triple excitations (CCSDT) in which tensor decompositions are used to reduce scaling and overall computational costs. For the decomposition of the electron repulsion integrals the standard density fitting (or Cholesky decomposition) format is used. The coupled-cluster single and double amplitudes are treated conventionally, and for the triple amplitudes tensor we employ the Tucker-3 compression formula, tijkabc ≈ tXYZ UaiX UbjY UckZ. The auxiliary quantities UaiX come from singular value decomposition (SVD) of an approximate triple amplitudes tensor based on perturbation theory. The efficiency of the proposed method relies on an observation that the dimension of the "compressed" tensor tXYZ sufficient to deliver a constant relative accuracy of the correlation energy grows only linearly with the size of the system, N. This fact, combined with proper factorization of the coupled-cluster equations, leads to practically N6 scaling of the computational costs of the proposed method, as illustrated numerically for linear alkanes with increasing chain length. This constitutes a considerable improvement over the N8 scaling of the conventional (uncompressed) CCSDT theory. The accuracy of the proposed method is verified by benchmark calculations of total and relative energies for several small molecular systems and comparison with the exact CCSDT method. The accuracy levels of 1 kJ/mol are easily achievable with reasonable SVD subspace size, and even more demanding levels of accuracy can be reached with a considerable reduction of the computational costs. Extensions of the proposed method to include higher excitations are briefly discussed, along with possible strategies of reducing other residual errors.

Complete Basis Set Extrapolation of Electronic Correlation Energies Using the Riemann Zeta Function.

In this article, we demonstrate the effectiveness of the method of complete basis set (CBS) extrapolation of correlation energies based on the application of the Riemann zeta function. Instead of fitting the results obtained with a systematic sequence of one-electron bases with a certain functional form, an analytic resummation of the missing contributions coming from higher angular momenta, l, is performed. The assumption that these contributions vanish asymptotically as an inverse power of l leads to an expression for the CBS limit given in terms of the zeta function. This result is turned into an extrapolation method that is very easy to use and requires no "empirical" parameters to be optimized. The performance of the method is assessed by comparing the results with very accurate reference data obtained with explicitly correlated theories and with results obtained with standard extrapolation schemes. On average, the errors of the zeta-function extrapolation are several times smaller compared with the conventional schemes employing the same sequence of bases. A recipe for the estimation of the residual extrapolation error is also proposed.

Symphony on strong field approximation.

This paper has been prepared by the Symphony collaboration (University of Warsaw, Uniwersytet Jagiellonski, DESY/CNR and ICFO) on the occasion of the 25th anniversary of the 'simple man's models' which underlie most of the phenomena that occur when intense ultrashort laser pulses interact with matter. The phenomena in question include high-harmonic generation (HHG), above-threshold ionization (ATI), and non-sequential multielectron ionization (NSMI). 'Simple man's models' provide both an intuitive basis for understanding the numerical solutions of the time-dependent Schrödinger equation and the motivation for the powerful analytic approximations generally known as the strong field approximation (SFA). In this paper we first review the SFA in the form developed by us in the last 25 years. In this approach the SFA is a method to solve the TDSE, in which the non-perturbative interactions are described by including continuum-continuum interactions in a systematic perturbation-like theory. In this review we focus on recent applications of the SFA to HHG, ATI and NSMI from multi-electron atoms and from multi-atom molecules. The main novel part of the presented theory concerns generalizations of the SFA to: (i) time-dependent treatment of two-electron atoms, allowing for studies of an interplay between electron impact ionization and resonant excitation with subsequent ionization; (ii) time-dependent treatment in the single active electron approximation of 'large' molecules and targets which are themselves undergoing dynamics during the HHG or ATI processes. In particular, we formulate the general expressions for the case of arbitrary molecules, combining input from quantum chemistry and quantum dynamics. We formulate also theory of time-dependent separable molecular potentials to model analytically the dynamics of realistic electronic wave packets for molecules in strong laser fields. We dedicate this work to the memory of Bertrand Carré, who passed away in March 2018 at the age of 60.

Ab initio Potential Energy Curve for the Ground State of Beryllium Dimer.

This work concerns ab initio calculations of the complete potential energy curve and spectroscopic constants for the ground state X1Σ g+ of the beryllium dimer, Be2. High accuracy and reliability of the results is one of the primary goals of the paper. To this end, we apply large basis sets of Slater-type orbitals combined with high-level electronic structure methods including triple and quadruple excitations. The effects of the relativity are also fully accounted for in the theoretical description. For the first time the leading-order quantum electrodynamics effects are fully incorporated for a many-electron molecule. Influence of the finite nuclear mass corrections (post-Born-Oppenheimer effects) turns out to be completely negligible for this system. The predicted well-depth ( De = 934.6 ± 2.5 cm-1) and the dissociation energy ( D0 = 807.7 cm-1) are in a very good agreement with the most recent experimental data. We confirm the existence of the weakly bound twelfth vibrational level [Patkowski et al. Science 2009, 326, 1382] that it lies just below the onset of the continuum.

Efficient singular-value decomposition of the coupled-cluster triple excitation amplitudes.

We demonstrate a novel technique to obtain singular-value decomposition (SVD) of the coupled-cluster triple excitations amplitudes, t ijk abc . The presented method is based on the Golub-Kahan bidiagonalization strategy and does not require t ijk abc to be stored. The computational cost of the method is comparable to several coupled cluster singles and doubles (CCSD) iterations. Moreover, the number of singular vectors to be found can be predetermined by the user and only those singular vectors which correspond to the largest singular values are obtained at convergence. We show how the subspace of the most important singular vectors obtained from an approximate triple amplitudes tensor can be used to solve equations of the CC3 method. The new method is tested for a set of small and medium-sized molecular systems in basis sets ranging in quality from double- to quintuple-zeta. It is found that to reach the chemical accuracy (≈1 kJ/mol) in the total CC3 energies as little as 5 - 15% of SVD vectors are required. This corresponds to the compression of the t ijk abc amplitudes by a factor of about 0.0001 - 0.005. Significant savings are obtained also in calculation of interaction energies or rotational barriers, as well as in bond-breaking processes. © 2019 Wiley Periodicals, Inc.

Exploring Chemical Space with Alchemical Derivatives: BN-Simultaneous Substitution Patterns in C60.

With the idea of using alchemical derivatives to explore in an efficient, computer- and cost-effective way Chemical Space was launched several years ago. In the context of Conceptual DFT response functions, these energies vs nuclear charge derivatives permit the estimatation of the energy of transmutants of a given starting or reference molecule showing different nuclear compositions. After an explorative study on small and planar molecules ( Balawender et al. J. Chem. Theory Comput. 2013 , 9 , 5327 ) by the present authors of this paper, the present study fully exploits the computational advantages of the alchemical derivatives in larger three-dimensional systems. Starting from a single reference calculation on C60, the complete BN substitution pattern, from single substituted C58BN via the belt (C20(BN)20 and the ball C12(BN)24 structures to the fully substituted (BN)30, is explored. Successive and simultaneous substitution strategies are followed and compared, indicating that both techniques yield identical results up to 13 substitutions but that for higher substitutions the simultaneous approach needs to be taken. Due to the cost-efficiency of the algorithm this path can indeed be followed as opposed to earlier work in the literature where for each step a full SCF calculation was at stake leading to prohibitively large computational demands for adopting the simultaneous approach. Previously formulated rules governing the substitution pattern by Kar and co-workers are scrutinized in this context and reformulated giving chemical insight in the gradual substitution process and the relative energies of the isomers. In its present form the method offers an interesting venue to study BN substitution patterns in higher fullerenes and graphene and in general paves the way for more efficient exploration of the Chemical Space.

Transition moments between excited electronic states from the Hermitian formulation of the coupled cluster quadratic response function.

We introduce a new method for the computation of the transition moments between the excited electronic states based on the expectation value formalism of the coupled cluster theory [B. Jeziorski and R. Moszynski, Int. J. Quantum Chem. 48, 161 (1993)]. The working expressions of the new method solely employ the coupled cluster operator T and an auxiliary operator S that is expressed as a finite commutator expansion in terms of T and T†. In the approximation adopted in the present paper, the cluster expansion is limited to single, double, and linear triple excitations. The computed dipole transition probabilities for the singlet-singlet and triplet-triplet transitions in alkali earth atoms agree well with the available theoretical and experimental data. In contrast to the existing coupled cluster response theory, the matrix elements obtained by using our approach satisfy the Hermitian symmetry even if the excitations in the cluster operator are truncated, but the operator S is exact. The Hermitian symmetry is slightly broken if the commutator series for the operator S are truncated. As a part of the numerical evidence for the new method, we report calculations of the transition moments between the excited triplet states which have not yet been reported in the literature within the coupled cluster theory. Slater-type basis sets constructed according to the correlation-consistency principle are used in our calculations.

Calculation of the molecular integrals with the range-separated correlation factor.

Explicitly correlated quantum chemical calculations require calculations of five types of two-electron integrals beyond the standard electron repulsion integrals. We present a novel scheme, which utilises general ideas of the McMurchie-Davidson technique, to compute these integrals when the so-called "range-separated" correlation factor is used. This correlation factor combines the well-known short range behaviour resulting from the electronic cusp condition, with the exact long-range asymptotics derived for the helium atom [Lesiuk, Jeziorski, and Moszynski, J. Chem. Phys. 139, 134102 (2013)]. Almost all steps of the presented procedure are formulated recursively, so that an efficient implementation and control of the precision are possible. Additionally, the present formulation is very flexible and general, and it allows for use of an arbitrary correlation factor in the electronic structure calculations with minor or no changes.

Reexamination of the calculation of two-center, two-electron integrals over Slater-type orbitals. I. Coulomb and hybrid integrals.

In this paper, which constitutes the first part of the series, we consider calculation of two-center Coulomb and hybrid integrals over Slater-type orbitals. General formulas for these integrals are derived with no restrictions on the values of the quantum numbers and nonlinear parameters. Direct integration over the coordinates of one of the electrons leaves us with the set of overlaplike integrals which are evaluated by using two distinct methods. The first one is based on the transformation to the ellipsoidal coordinates system and the second utilizes a recursive scheme for consecutive increase of the angular momenta in the integrand. In both methods simple one-dimensional numerical integrations are used in order to avoid severe digital erosion connected with the straightforward use of the alternative analytical formulas. It is discussed that the numerical integration does not introduce a large computational overhead since the integrands are well-behaved functions, calculated recursively with decent speed. Special attention is paid to the numerical stability of the algorithms. Applicability of the resulting scheme over a large range of the nonlinear parameters is tested on examples of the most difficult integrals appearing in the actual calculations including, at most, 7i-type functions (l=6).

Reexamination of the calculation of two-center, two-electron integrals over Slater-type orbitals. II. Neumann expansion of the exchange integrals.

In this paper we consider the calculation of two-center exchange integrals over Slater-type orbitals (STOs). We apply the Neumann expansion of the Coulomb interaction potential and consider calculation of all basic quantities which appear in the resulting expression. Analytical closed-form equations for all auxiliary quantities have already been known but they suffer from large digital erosion when some of the parameters are large or small. We derive two differential equations which are obeyed by the most difficult basic integrals. Taking them as a starting point, useful series expansions for small parameter values or asymptotic expansions for large parameter values are systematically derived. The resulting expansions replace the corresponding analytical expressions when the latter introduce significant cancellations. Additionally, we reconsider numerical integration of some necessary quantities and present a new way to calculate the integrand with a controlled precision. All proposed methods are combined to lead to a general, stable algorithm. We perform extensive numerical tests of the introduced expressions to verify their validity and usefulness. Advances reported here provide methodology to compute two-electron exchange integrals over STOs for a broad range of the nonlinear parameters and large angular momenta.

On the large interelectronic distance behavior of the correlation factor for explicitly correlated wave functions.

In currently most popular explicitly correlated electronic structure theories, the dependence of the wave function on the interelectronic distance rij is built via the correlation factor f(r(ij)). While the short-distance behavior of this factor is well understood, little is known about the form of f(r(ij)) at large r(ij). In this work, we investigate the optimal form of f(r12) on the example of the helium atom and helium-like ions and several well-motivated models of the wave function. Using the Rayleigh-Ritz variational principle, we derive a differential equation for f(r12) and solve it using numerical propagation or analytic asymptotic expansion techniques. We found that for every model under consideration, f(r12) behaves at large r(ij) as r12(ρ)e(Br12) and obtained simple analytic expressions for the system dependent values of ρ and B. For the ground state of the helium-like ions, the value of B is positive, so that f(r12) diverges as r12 tends to infinity. The numerical propagation confirms this result. When the Hartree-Fock orbitals, multiplied by the correlation factor, are expanded in terms of Slater functions r(n)e(-ßr), n = 0,...,N, the numerical propagation reveals a minimum in f(r12) with depth increasing with N. For the lowest triplet state, B is negative. Employing our analytical findings, we propose a new "range-separated" form of the correlation factor with the short- and long-range r12 regimes approximated by appropriate asymptotic formulas connected by a switching function. Exemplary calculations show that this new form of f(r12) performs somewhat better than the correlation factors used thus far in the standard R12 or F12 theories.

Molecular electrostatic potential at the atomic sites in the effective core potential approximation.

Considering calculations of the molecular electrostatic potential at the atomic sites (MEP@AS) in the presence of effective core potentials (ECP), we found that the consequent use of the definition of MEP@AS based on the energy derivative with respect to nuclear charge leads to a formula that differs by one term from the result of simple application of Coulomb's law. We have developed a general method to analytically treat derivatives of ECP with respect to nuclear charge. Benchmarking calculations performed on a set of simple molecules show that our formula leads to a systematic decrease in the error connected with the introduction of ECP when compared to all-electron results. Because of a straightforward implementation and relatively low costs of the developed procedure we suggest to use it by default.

Exploring Chemical Space with the Alchemical Derivatives.

In this paper, we verify the usefulness of the alchemical derivatives in the prediction of chemical properties. We concentrate on the stability of the transmutation products, where the term "transmutation" means the change of the nuclear charge at an atomic site at constant number of electrons. As illustrative transmutations showing the potential of the method in exploring chemical space, we present some examples of increasing complexity starting with the deprotonation, continuing with the transmutation of the nitrogen molecule, and ending with the substitution of isoelectronic B-N units for C-C units and N units for C-H units in carbocyclic systems. The basis set influence on the qualitative and quantitative accuracies of the alchemical predictions was investigated. The alchemical deprotonation energy (from the second order Taylor expansion) correlates well with the vertical deprotonation energy and can be used as a preliminary indicator for the experimental deprotonation energy. The results of calculations for the BN derivatives of benzene and pyrene show that this method has great potential for efficient and accurate scanning of chemical space.