Domain-Based Local Pair Natural Orbital Version of Mukherjee's State

Subscriber access provided by READING UNIV

Article

Domain-based Local Pair Natural Orbital version of Mukherjee’s state specific coupled cluster method Jiri Brabec, Jakub Lang, Masaaki Saitow, Jiri Pittner, Frank Neese, and Ondrej Demel J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.7b01184 • Publication Date (Web): 18 Jan 2018 Downloaded from http://pubs.acs.org on January 20, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Domain-based Local Pair Natural Orbital version of Mukherjee’s state specific coupled cluster method Jiri Brabeca , Jakub Langa , Masaaki Saitowb , Jiˇr´ı Pittnera , Frank Neeseb , and Ondˇrej Demela∗∗ a

J. Heyrovský Institute of Physical Chemistry, v.v.i., Academy of Sciences of the Czech Republic, Dolejˇskova 3, 18223 Prague 8, Czech Republic

b

Max Planck Institute for Chemical Energy Conversion, M¨ ulheim an der Ruhr, Germany E-mail: [email protected] Abstract This article reports development of a local variant of Mukherjee’s state specific multireference coupled cluster method based on the domain-based pair natural orbital approach (DLPNO-MkCC). The current implementation is restricted to connected single and double excitations and model space with up to biexcited references. The performance of the DLPNO-MkCCSD was tested on calculations of tetramethyleneethane. The results show that above 99.9% of correlation energy was recovered, with respect to the conventional MkCC method. To demonstrate the applicability of the method to large systems, singlet-triplet gaps of triangulene and bis(1-(2,6-diisopropylphenyl)3,3,5,5-tetramethylpyrrolidine-2-ylidene)beryllium complex. For the last system (105 atoms), we were able to perform a calculation in cc-pVTZ with 2158 basis functions on a single CPU in less than 9 days.

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

The coupled cluster method 1 is one of the most successful approaches for the description of dynamical correlation. Among its advantages are compactness of the exponential parameterization of wave function, size-extensivity, invariance to orbital rotation within orbital subspaces, and a systematic hierarchy of approximations converging to the full configuration interaction limit.

However, the traditional single reference CC approach has also several disadvantages. One of them is a relatively poor performance at low truncation of cluster operator for systems where non-dynamical correlation is important. However, the multireference generalization of coupled cluster methods is not unique and many types of approaches have been developed over years. These methods can be divided into several classes: Hilbert space methods, 2–7 that are based on CAS-like reference and use a distinct cluster operator for each reference; Fock space and internally contracted methods 8–17 that also use a CAS-like reference but the cluster operator is common for all references; and methods that keep the framework of single reference methods. 18–37

In Hilbert space approaches, the wave operator is assumed in the form of JeziorskiMonkhorst ansatz. 2 For state universal methods, all the H eff eigenvalues are physical, whereas for state specific methods 5–7 only one has a physical meaning. Among the state-specific Hilbert space approaches, the Mukherjee’s variant (MkCC) 5,38 is most commonly used, due to its rigorous size-extensivity. Other advantages include insensivity to intruder problem, which plagues the state universal approaches, and a simplicity of structure of the working equations, which enables reusing much of single-reference code in its implementation. This method was developed in many groups, including Mukherjee, Schaefer, Gauss, Kállay, and Pittner. 39–43 Nowadays codes with iterative and noniterative triple excitations, 41,44 higher excitations, 42 analytic gradient, 45 F12 variant, 46 and uncoupled approximation, 47,48 as well 2


Page 2 of 42

Page 3 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


as massively parallel codes 49 are available.

Another drawback of the canonical CC method is the high computational cost due to the unfavorable scaling (O(n6 ) for CCSD, O(n7 ) for CCSD(T)). This represents a severe limitation for applicability to large molecules. One possibility how to overcome the high cost of calculation is massive parallelization. CC codes using distributed memory tensor networks – eg. Super Instruction Assembly Language, 50,51 Tensor Contraction Engine, 49,52–54 TiledArray, 55 and Cyclops Tensor Framework 56 – have been developed. However, this does not remove the steep scaling and the methods are still applicable only to medium sized molecules, even on the most powerful machines available today. Furthermore, only a fraction of computational chemistry groups have an unlimited access to clusters with tens or hundreds of thousands of cores, which are needed for these calculations. A well established way to reduce the scaling are the local approaches introduced by Pulay 57,58 which are based on the short range (R−6 ) character of electron correlation. In the basis of localized orbitals, the matrix of HN operator is sparse and non-zero elements are only encountered when all the indices correspond to spatially closed orbitals. For occupied orbitals, standard Foster-Boys and Pipek Mezey procedures 59,60 are commonly used. For virtual space, there are several possibilities. In the original local coupled cluster method of Werner and Sch¨ utz, 61–63 the projected atomic orbitals (PAOs) are used. These are a priori local and can be easily used to form domains corresponding to local occupied orbitals. Pairs of occupied orbitals are classified according to the real space distance of the respective centers. The strong pairs are treated at the full level of theory (CC), weak and distant pairs only approximately (MP2), while very distant pairs are neglected altogether. In this way, linear scaling methods were developed. 64,65 The remaining disadvantage is the sensitivity of the results to the values of cut-off parameters.

Many other local approaches have been developed based on various variants of divide-and-

3



conquer strategy: fragment molecular orbital method, 66 the divide-and-conquer method, 67 divide-expand-consolidate method, 68,69 the incremental scheme, 70 and the cluster-in-molecule approach, 71 among many others. Here, the system is divided into fragments. Series of separate calculations are performed for individual fragments or their pairs, and therefore the methods are embarrassingly parallel. The drawback of these method are redundancies originating from overlapping fragments and potentially limiting high cost of individual fragment calculations for large fragments. An alternative approach is to use a virtual space spanned by pair natural orbitals (PNOs). 72,73 The PNOs were introduced in the 1960’s by Edmiston and Krause and were shown to provide a compact parameterization of virtual space, yielding a accelerated convergence of the CI expansion. The CISD and CEPA code of Meyer based on PNOs 74,75 provided state-of-art description by that time for a very modest computational cost. The PNO approach was also used for multireference MR CISD method by Taylor 76 and MR CEPA by Fink and Staemmler. 77 The full potential of PNOs could only be realized in combination with local correlation approaches, modern hardware and modern integral transformation technologies, especially the density fitting (DF) or resolution of the identity (RI) methodology. 78 In late 2000’s, the PNO technique was used in the LPNO CEPA and LPNO-CCSD. 79,80 Within this scheme, each pair of localized occupied orbitals is assigned a separate set of PNOs. These orbitals provide the most compact description of the virtual space of the corresponding pair an are as local (or nonlocal) as the physical situation requires. In this approach, only a limited number of cut-off parameters are used: T CutP airs for contributing electron pairs as the minimum MP2 pair energy for the pair to be kept; T CutP N O as the minimum occupation number of a PNO to be kept, and T CutM KN for expansion of the auxiliary basis. None of these cut-off parameters is based on real space distance. The advantages of this approach include compactness of the LPNO-based wave function, high accuracy measured as the percentage of correlation energy relative to calculation with canonical orbitals, smooth dependence of

4


Page 4 of 42

Page 5 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the remaining error on the cut-off parameters, and black box character of the method. In order to achieve near linear scaling, the above mentioned scheme was improved by Riplinger et. al., yielding the domain based local pair natural orbital coupled cluster method (DLPNO-CCSD). 81 The main differences compared to LPNO-CCSD are i) construction of the PNOs as linear combinations of PAOs, ii) setting up the connected singles in the PNO basis, iii) prescreening of the contributing pair energies using dipole approximation, and iv) efficient integral transformation based on the PAO basis. In this way, the bottlenecks of the LPNO-CCSD method are removed and the resulting method can be applied to systems with over 500 atoms and 10 000 basis of basis set functions, 81–83 where the cost of correlation treatment is comparable or even smaller than the cost of preceding Hartree Fock calculation. At the same time, the robustness of the PNO approaches is preserved. Nowadays, the PNO based approaches are investigated by several groups, including Werner and Hättig 83–91 and are widely applied for studies of real chemical interest. 92–101 In this article we report development of the domain based local pair natural orbital variant of Mukherjee’s state specific multireference method (MkCC). The MkCC was selected because of the simple structure of working equations closely resembling the single reference CCSD, as well as size-extensivity which is crucial for large molecules. The DLPNO-MkCC method has been implemented in the ORCA package. To test the performance of the DLPNO-MkCCSD method, benchmark calculations of isomers of naphtynes, rotational barrier of tetramethyleneethane, [Be(Me L)2 ] complex, where Me

L =1-(2,6-diisopropylphenyl)-3,3,5,5-tetramethylpyrrolidine-2-ylidene, and triangulene were

carried out. The TME and naphthyne systems were used also in the previous LPNOMkCCSD study 102 and are thus suitable to test the accuracy and efficiency of the DLPNOMkCCSD method.

5



2

Page 6 of 42

Theory

This section is divided into three parts. First, we review basics of MkCC theory. Next we describe how to obtain PNOs for subsequent MkCC calculation. Finally we use these PNOs to formulate the DLPNO-MkCC method.

2.1

Basics of MkCC method

The reference function is assumed to be a linear combination of reference configurations spanning a model space |ΨPω i =

M X

Cµω |Φµ i.

(1)

µ=1

ˆ ω as Furthermore, let us define the wave operator Ω ˆ ω |ΨPω i. |Ψω i = Ω

(2)

ˆ eff = Pˆ H ˆΩ ˆ ω Pˆ , H

(3)

and the effective Hamiltonian H eff as

where Ψω stands for the exact wave function, P for projection operator on the model space, and index ω to the studied state.

ˆ eff |ΨPω i = Eω |ΨPω i. H

(4)

The matrix representation within the model space reads as X

eff ω Hµν Cν = Eω Cµω .

ν

6


(5)

Page 7 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The wave operator is taken in the form of Jeziorski-Monkhorst ansatz ˆω = Ω

M X

ˆ

eT (µ) |Φµ ihΦµ |,

(6)

µ=1

where T (µ) is a cluster operator assigned to the µ-th reference configuration. Assuming a complete model space (CMS) is employed, the amplitudes corresponding to internal excitations are zero due to intermediate normalization. The remaining cluster amplitudes are obtained by solving the cluster equations hΦ(µ)ϑ |e−T (µ) HeT (µ) |Φ(µ)iCµω +

X

eff ω Hµν Cν hΦ(µ)ϑ |e−T (µ) eT (ν) |Φ(µ)i = 0,

(7)

ν6=µ

where index ϑ corresponds to an excitation outside of the model space. Here, the first term is identical to single reference CC residual calculated with a set of cluster amplitudes corresponding to the µ-th reference configuration and only the second term, referred as coupling term, connects cluster amplitudes of different references. The matrix elements of the effective Hamiltonian are given as eff Tˆ(ν) ˆ eff |Φ(ν)i = hΦ(µ)|H|Φ(µ)iδ ˆ ˆ Hµν = hΦ(µ)|H |Φ(ν)i. µν + hΦ(µ)|HN (ν)e

2.2

(8)

Prescreening and PNO generation for DLPNO-MkCCSD

In this section we describe how to generate a common set of PNOs for all reference configurations used in subsequent DLPNO-MkCCSD calculation. In this paragraph indices i, j, . . . stand for inactive (occupied) orbitals, p, q, . . . for active, a, b, . . . for external orbitals, i′ , j ′ , . . . for orbitals occupied in a given reference, and a′ , b′ , . . . for orbitals unoccupied in a given reference. The DLPNO-MkCC calculation starts with a CAS SCF calculation. In the next step, the inactive orbitals are localized using standard Foster-Boys procedure. The orbital domains of projected atomic orbitals (PAOs) are set up in the same way as in the

7



Page 8 of 42

single reference DLPNO-CCSD. 81 The PAOs are canonicalized with respect to CAS SCF Fock matrix. In the next step, the prescreening of pair energies is performed. At this step, orbital pairs involving active orbitals, diagonal pairs, as well as pairs with differential overlap greater than T CutDO cut-off parameter (default 1.0 × 10−2 ) are always included. For the remaining pairs, the pair energy is estimated at the strongly contracted NEVPT2 level 103–105 using the dipole approximation. The contributions Vij (i.e. class including ij → pq excitations, the active indices are omitted in the name of the class), Vija (i.e. class including ij → ap excitations), and Vijab (i.e. class including ij → ab excitations) are all included into the pair energy. Details can be found in reference, 86 e.g. in equation (24). Pairs with estimated pair-energy greater than or equal to T CutP re (default value 1.0 × 10−6 Hartree) are kept for further calculation, while the remaining pairs (i.e. off-diagonal ij pairs with low pair energy and differential overlap) are neglected and only their pair energies are collected. As was established earlier, the dipole approximation works sufficiently well for those pairs with small pair energy contributions. Here, the number of active orbitals is assumed to be independent of the system size, the number of ip and pq pairs is thus negligible compared to ij pairs, and therefore the prescreening for ip and pq pairs is neither required nor performed. The next step is the evaluation of perturbation energies of the “surviving” pairs and subsequent construction of PNOs. We have devised two schemes, one based on MP2 and one based on NEVPT2. In the first scheme, both the pair energies and the PNOs are constructed at the MP2 level. For each pair of occupied/active spatial orbitals surviving the prescreening step we calculate first order amplitudes ′ ′

′ ′ tia′ jb′ (µ)

Kia′ jb′ = ǫa′ (µ) + ǫb′ (µ) − ǫi′ (µ) − ǫj ′ (µ)

(9)

where indices a′ ,b′ correspond to PAOs and are restricted to the domain of the i′ j ′ pair,

8


Page 9 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


′ ′

Kia′ jb′ is the exchange integral and ǫa′ (µ) is the diagonal element of Fock matrix of the µ-th reference. Then we use these amplitudes to calculate the unrestricted MP2 energy. The UMP2 energy is summed over spin and averaged over M ′ references where i′ j ′ are both occupied. If the average pair energy is greater than or equal to T CutP airs cut-off parameter (default value 1.0×10−4 ), the pair is kept, otherwise it is discarded. No pq pairs are neglected in this step, since the number of these pairs is negligible and all these pairs are assumed to be important. For the surviving pairs, the RHF MP2 densities for the αβ spin case are constructed for each reference Di′ j ′ (µ) = T˜i†′ j ′ (µ)Ti′ j ′ (µ) + T˜i′ j ′ (µ)Ti†′ j ′ (µ),

(10)

where T˜i′ j ′ (µ) =

1 (4Ti′ j ′ (µ) − 2Tj ′ i′ (µ)) 1 + δ i′ j ′

(11)

The densities are then averaged over M ′ references M 1 X = ′ Di′ j ′ (µ) M µ ′

Di ′ j ′

(12)

where both i′ and j ′ are occupied. ′ ′

Subsequent diagonalization of Di j matrix yields the PNOs

Di′ j ′ dxi′ j ′ = nxi′ j ′ dxi′ j ′

(13)

that are used for all spin cases of the pair. Only those PNOs with occupation number nxi′ j ′ greater or equal to T CutP N O (default value 3.33 × 10−7 ) are kept. The energy correction for neglected PNOs is calculated for each reference separately as the difference between unrestricted MP2 energy in the PAO basis and the truncated PNO basis. This correction is added to pair energies of pairs neglected in previous steps, for each reference. The PNOs for connected singles are constructed as PNOs of diagonal i = j or p = q pairs with a reduced 9



value of cut-off parameter T CutP N OSingles = 0.03 × T CutP N O analogously to Eq. 13. In the alternative PNO construction scheme, the semicanonical partially contracted Nelectron valence perturbation theory (PC-NEVPT2) energies are calculated for the surviving pairs. For the ij pairs, the energy contributions from Vijab , Vija , and Vij classes are added together. Detailed expression are presented in reference. 86 Similarly, the Viab , Via , and Vi contributions are added for ip pairs, as are Vab and Va added for pq terms. Among ij and ip pairs, only those with PC-NEVPT2 pair energies greater or equal to T CutP airs are included for the subsequent coupled cluster step. Again, no pq pairs are neglected. For the surviving pairs, the PC-NEVPT2 pair densities are constructed. For ij pairs, this includes contributions from classes Vijab and Vija , for ip pairs contributions from classes Viab and Via , and for pq pairs Vab and Va contributions. The densities are diagonalized, in order to obtain the PNOs. Only those PNOs with occupation number greater or equal to T CutP N O are kept. The energy contribution of the remaining PNOs is calculated as the difference between the pair energy calculated in PAO basis and truncated PNO basis. It was expected that the latter scheme will be preferable, due to more satisfactory treatment of the multireference nature of the studied systems. However, as detailed below, this expectation was not confirmed by benchmark calculations

2.3

DLPNO-MkCCSD working equations

In the next step, we make use of the fact that we have a common set of PNOs for all the references. This means that we can construct a common set of molecular integrals for all references, in order to save computational effort for the rather costly integral transformation as well as to cut disk storage requirement. For the same reason, we want to have the same PNOs for all the spin cases (alpha-alpha, alpha-beta, and beta-beta). However, the active orbitals are occupied for some references and unoccupied for others. It is therefore convenient to include the active orbitals to both occupied and virtual orbital window. The virtual orbital space is spanned by PNOs and these are in principle different for each pair 10


Page 10 of 42

Page 11 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of occupied orbitals. To keep the scheme simple, we keep the active orbitals unchanged and include them to the set of PNOs of all pairs for occupied orbitals. Therefore, we also adapt the notation in the subsequent text: indices i, j, . . . stand for both internal and active orbitals; a, b, . . . for active and external in MO basis; and a ¯ij , ¯bij , . . . for both active and external in the PNO basis of ij pair. The cluster operators are truncated at the connected singles and doubles level and are constructed in the PNO basis set as described above.

T (µ) = T1 (µ) + T2 (µ) XX T1 (µ) = t(µ)ai¯i a†a¯i ai T2 (µ) =

XX

a ¯i a ¯ij ¯bij † t(µ)ij aa¯ij a¯†bij aj ai

(14) (15)

i

(16)

ij a ¯ij ¯bij

The DLPNO-MkCCSD cluster equation (7) becomes hΦ(µ)ai¯ii |e−T (µ) HeT (µ) |Φ(µ)iCµα +

X

eff α Hµν Cν hΦ(µ)ai¯i i |e−T (µ) eT (ν) |Φ(µ)i = 0

(17)

ν6=µ

for the connected singles and a ¯ ¯b

hΦ(µ)ijij ij |e−T (µ) HeT (µ) |Φ(µ)iCµα +

X

a ¯ ¯b

eff α Hµν Cν hΦ(µ)ijij ij |e−T (µ) eT (ν) |Φ(µ)i = 0

(18)

ν6=µ

for the connected doubles. The first terms (direct terms) are calculated using slightly modified single reference open-shell DLPNO-CCSD code 83 in ORCA. The changes are needed to properly handle different occupation of active orbitals for the references. Since in the SR DLPNO-CCSD code all alpha-spin singly occupied molecular orbitals (SOMO’s) are occupied and all beta-spin SOMO’s unoccupied, it is efficient to generate amplitudes and intermediates within these specific orbital windows, which is different from the orbital layout used here (see above). Further changes are needed to zero out the amplitudes and residuals

11



Page 12 of 42

corresponding to excitations from unoccupied / to occupied spin orbitals for the respective reference. This technique was used previously in the LPNO-MkCCSD method. 102 The second terms in equations (17) and (18), the coupling term, can be evaluated as h(Φ(µ))ai¯i |e−T (µ) eT (ν) |Φ(µ)i = ∆tai¯i (ν/µ, µ)

(19)

for the connected singles and ¯b 1 a ¯ ¯b a ¯ ¯b a ¯ h(Φ(µ))ijij ij |e−T (µ) eT (ν) |Φ(µ)i = ∆tijij ij (ν/µ, µ) + P (ij)P (ab)∆ti ij (ν/µ, µ)∆tjij (ν/µ, ν) 2 (20) a... a... a... for the connected doubles. Here, ∆ta... i... (ν/µ, µ) = ti... (ν/µ) − ti... (µ) and ti... (ν/µ) are equal

to ta... i... (ν) if all i . . . are occupied and all a . . . unoccupied for both the ν-th and µ-th reference and zero otherwise. In contrast with the previous LPNO-MkCCSD variant, the singles couplings are formulated in the singles PNO basis. Since the PNOs are identical for all the references, this causes no complications. The doubles couplings in equation 20 are equivalent to their counterparts in LPNO-MkCCSD. The diagonal elements of the effective Hamiltonian are calculated as the CC energies corresponding to the respective references plus the correction for neglected pairs and truncated PNOs, denoted as EREST. If we restrict ourselves to active spaces with mutually eff up to biexcited references, the off-diagonal Hµν elements can be collected as the residuals

corresponding to the internal excitations.

3

Results and discussion

In this section, we present results of benchmark calculations on tetramethyleneethane (TME), isomers of naphtynes, [Be(Me L)2 ] complex, where

Me

L =1-(2,6-diisopropylphenyl)-3,3,5,5-

tetramethylpyrrolidine-2-ylidene, and triangulene. In all calculations we utilized CMS(2,2) 12


Page 13 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


orbitals obtained from CASSCF calculation.

3.1

Accuracy of the DLPNO-MkCCSD approach with respect to canonical MkCCSD

In order to compare the performance of our new method with the canonical version of MkCCSD, we performed benchmark calculations on two selected systems. The size of those systems cannot be too small, since the number of the pairs and the PNOs would not be appropriate for benchmarking purposes. On the other hand it cannot be too large, otherwise we could not perform canonical MkCCSD calculations. We have thus chosen tetramethyleneethane (TME) and naphthyne isomers, which have strong multireference character requiring CMS(2,2) model space and where the DLPNO-MkCCSD results can be compared with recently published LPNO-MkCCSD results (See Ref. 102 ). 3.1.1

Rotational barrier of TME

First, we calculated the rotational barrier of TME using the geometries from Ref. 44 For this system, a model space consisting of two references HOMO2 LUMO0 and HOMO0 LUMO2 was employed, the remaining two references do not contribute because of symmetry consideration. The standard correlation consistent cc-pVTZ basis set 106 in spherical representation was used throughout the calculations. The aim of these calculations was to answer several basic questions: (a) what percentage of correlation energy is recovered, how does that depend on cut-off parameters, (b) what is the accuracy for relative energies, and (c) how does it all depend on the method used for the PNO construction. As a first step, we calculated the rotational barriers as well as percentages of correlation energy recovered with the default cut-off settings with both NEVPT2 and MP2 PNOs. The results are shown in Fig.1(a) and Fig.1(b). In both cases, approximately 99.9 % of the canonical correlation energy was recovered. However, for relative energies the two curves show significant differences especially in the region of 30◦ -50◦ where the multireference char13



acter is more complicated (See 107 ). Clearly, DLPNO-MkCCSD based on MP2 PNOs gives a better description, which can be also seen from the values of nonparallelity error (NPE) – 0.17 kcal/mol for MP2 PNOs vs. 0.66 kcal/mol for NEVPT2 PNOs. 3.0

MkCCSD NEVPT2 PNOs PQ MP2 PNOs PQ+IP MP2 PNOs MP2 PNOs

2.0

100.00

Correlation energy [%]

2.5

Relative energy [kcal/mol]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 42

1.5 1.0 0.5

99.95

99.90

99.85

0.0 -0.5 0.0

NEVPT2 PNOs PQ MP2 PNOs PQ+IP MP2 PNOs MP2 PNOs

10.0

20.0

30.0

40.0

50.0

60.0

70.0

80.0

90.0

0.0

Angle

10.0

20.0

30.0

40.0

50.0

60.0

70.0

80.0

90.0

Angle

(a)

(b)

Figure 1: (a) Rotational barrier of TME calculated by DLPNO-MkCCSD method with different source of PNOs, (b) Percentage of correlation energy recovered for TME by DLPNOMkCCSD with different source of PNOs.

To understand the origin of these differences, we checked individual pair energies at the MP2 and NEVPT2 level and the corresponding number of PNOs. As a reference geometry, we took the angle 40◦ . The most interesting difference concerned the energy of the pq pairs, where the NEVPT2 pair energies were almost two orders of magnitude smaller than their MP2 counterparts (0.2 mHartree vs. 14.8 mHartree). Correspondingly, the number of PNOs obtained for these pairs was at NEVPT2 level smaller by a factor of four. To illustrate these findings, we have plotted diagrams showing ip and pq pair energies for NEVPT2, MP2, for MkCCSD at first iteration, and for converged MkCCSD. We also plot NEVPT2 pair energies summed up over the active index, in order to obtain an adequate quantity comparable with MP2 and MkCCSD results (MP2 or MkCCSD pair energies are reference-specific. For example of TME and CMS(2,2) containing only two closed shell references, only one of the ip pairs contributes to each reference for a given inactive orbital i, while for NEVPT2 all the pair energies are non-zero and these are not reference-specific.). The

14


Page 15 of 42

results are presented in Fig.2. In all cases MP2 ip pair energies are consistently closer their

2.5e-02

NEVPT2 NEVPT2 (sum) MP2 MkCCSD (1. iter) MkCCSD (converged)

2.0e-02

Pair energy [au]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


IP Pairs

PQ Pairs

1.5e-02

1.0e-02

5.0e-03

0.0e+00

Pair

Figure 2: The pair energies of ip (on the left) and pq (on the right) pairs for TME (α = 40◦ ) estimated at NEVPT2, MP2, and MkCCSD level, sorted by their magnitude. NEVPT2 NEVPT2(sum) = (sum) stands for NEVPT2 pair energies summed up over active label: Ei P NEVPT2 E . Default thresholds were used. p ip MkCCSD counterparts, slightly overestimating MkCCSD. NEVPT2 energies are smaller than MP2 energies in all cases, NEVPT2(sum) energies are smaller in all cases except for four largest contributions, where they are much larger than the ones of other methods. For pq pairs, the difference between NEVPT2 and MP2 or MkCCSD is dramatically larger. Overall, MP2 approach gives a much better description of pairs involving active orbitals, which is crucial for a proper construction of ip and pq PNOs. The same discrepancy in pq pair energies was observed for other systems we studied - naphtynes as well as for the beryllium complex. Furthermore, significant differences were observed in the energy of the neglected pairs and PNOs (EREST) for the ip pairs. For NEVPT2 PNOs, the EREST contribution was ten times larger than for MP2 PNOs, which indicates a worse description in the PNO basis for given pairs. The same was not observed for the pq pairs, where the EREST is approximately 0.03 mHartee. However, that was only due to the very small NEVPT2 pair energy. When the MP2 energy is calculated for the pq pairs using NEVPT2 PNOs, EREST rises by two orders of magnitude which lead to serious overestimation of the correlation energy. Again, this shows that, surprisingly, the NEVPT2 PNOs are not well suited for Hilbert space MRCC 15



methods. To verify this further, we tried to selectively tighten the T CutP N O cut-off parameter by a factor 0.01 for pq pairs and 0.1 for ip pairs. This has lead to a decrease in NPE from 0.66 kcal/mol to approximately 0.33 kcal/mol. Furthermore, we performed a series of calculations with mixed PNOs – a) MP2 PNOs for pq pairs, NEVPT2 for ip and ij pairs; and b) MP2 PNOs for pq and ip pairs, NEVPT2 for ij pairs. These curves are shown in Fig.1(a). The corresponding NPE’s of 0.30 kcal/mol or 0.35 kcal/mol, respectively, lie between the NPE for NEVPT2 and MP2 PNOs. It should be also pointed out that part of the relative success of calculations with MP2 PNOs may be connected with the relatively small active space employed. In the CMS(2,2) space, only one reference has nonzero contributions for a given pq pair since at least one of the orbitals p, q is unoccupied for the remaining references. Therefore, the PNOs for these pairs are optimized for the only contributing reference. When only two references of the CMS(2,2) space are used, the same is true also for the ip pairs. Based on the results we obtained in these benchmark calculations, we employed MP2 PNOs in all subsequent calculations. 3.0

0.60

2.0

Non-parallelity error [kcal/mol]

MkCCSD 1.0e-7 3.3e-7 1.0e-6

2.5

Relative energy [kcal/mol]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 42

1.5 1.0 0.5 0.0 -0.5 0.0

NPE, TCutPairs = 1.0e-4 NPE, TCutPairs = 5.0e-5

0.50

0.40

0.30

0.20

0.10

0.00 10.0

20.0

30.0

40.0

50.0

60.0

70.0

80.0

90.0

1.00e-05

Angle

1.00e-06

1.00e-07

TCutPNO

(b)

(a)

Figure 3: (a) Rotational barrier of TME with respect to T CutP N O, (b) corresponding NPE (in kcal/mol) as a function of T CutP N O.

In Fig.3(a), we compare the rotational barrier with respect to T CutP N O threshold, 16


Page 17 of 42

where all the PNOs were constructed from the MP2 amplitudes. These results show very similar curves for all the values of T CutP N O. For T CutP N O smaller than 1.0 × 10−6 , the relative DLPNO-MkCCSD energy w.r.t α = 90◦ is below the canonical MkCCSD energy at each point. The NPE associated with the T CutP N O parameter is shown in Fig.3(b). For default threshold, the NPE is only 0.17 kcal/mol, for lower values of T CutP N O it decreases to 0.12 kcal/mol. However, similar improvement can be obtained by tightening T CutP airs to 5.0 × 10−5 . This can be of practical importance since the computational cost for the two sets of calculations is very similar. With T CutP airs = 5.0 × 10−5 and T CutP N O = 1.0 × 10−7 the NPE drops below 0.1 kcal/mol. 100.50

100.8 α=0 α=40 α=90

α=0 α=40 α=90

100.7


100.40


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


100.30 100.20 100.10 100.00 99.90

100.6 100.5 100.4 100.3 100.2 100.1

99.80

100.0 1e-05

1e-06

1e-07

0.001

TCutPNO

0.0001

1e-05

TCutPairs

(a)

(b)

Figure 4: Percentage of canonical correlation energy recovered for three angles (α = 0◦ , α = 40◦ , and α = 90◦ ) with respect to (a) T CutP N O (b) T CutP airs threshold.

Furthermore, we have studied how the percentage of correlation energy changes with the cut-offs used in DLPNO-MkCCSD. The results are shown in Fig.4(a) for T CutP N O, Fig.4(b) for T CutP airs, and Fig.5 for T CutDO. In all cases, the dependency on cut-off parameter is smooth. The error is very small or negligible with T CutP airs smaller than 1 × 10−4 Hartree (default value), because the number of cut pairs is very small. As T CutP airs increases, the percentage of the recovered correlation energy (as well as the error) steeply increases to more than 100%, which is caused by overestimation of the pair energy of the neglected pairs at the perturbation theory level. Next, let’s investigate the T CutP N O curve. For very 17



100.0


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

α=0 α=40 α=90

99.8

99.6

99.4

99.2 0.1

0.01

0.001

TCutDO

Figure 5: Percentage of correlation energy recovered for three TME molecules for α = 0◦ , α = 40◦ and α = 90◦ , with respect to T CutDO threshold.

large values T CutP N O, only few PNOs per pair survive and the MP2 used for the neglected PNOs significantly overestimates of correlation energy. The curve has then a minimum near T CutP N O = 1.0 × 10−6 and then the percentage of recovered correlation grows ideally to 100%. The obtained curve is highly satisfactory since it is smooth, over 99.9% of correlation energy is recovered along all the curve, and for very small T CutP N O it is approaching 100% of correlation energy. We also calculated the energy with respect to T CutDO (Fig.5) threshold. The energy is almost constant with T CutDO up to 0.01, which is also the default value. With looser values of T CutDO, the amount of correlation energy T CutP N O quickly decreases as important atoms are cut out of the domains. 3.1.2

Naphthyne isomers

As the second benchmark, we have investigated four isomers – (1,7), (1,8), (2,6), and 2,7naphthyne, which exhibit multireference character. All calculations employed the cc-pVTZ basis set. For 1,7-naphthyne, also the two open-shell references of the CMS(2,2) space were included, while only the two closed-shell references are needed for the remaining three isomers due to the symmetry considerations. The geometries for single-point calculations were taken from Ref. 49 The percentage of correlation energies recovered by the DLPNO-MkCCSD method are 18


Page 18 of 42

Page 19 of 42

101.00 Naphtyne-2,7 Naphtyne-2,6 Naphtyne-1,8 Naphtyne-1,7


100.80 100.60 100.40 100.20 100.00 99.80 1.00e-04

1.00e-05

1.00e-06

1.00e-07

TCutPNO

Figure 6: Average percentage of correlation energy recovered for naphthyne isomers by DLPNO-MkCCSD as a function of T CutP N O

3.0

3.0 naphtyne-2,7 naphtyne-1,8 naphtyne-1,7

2.5

∆∆E wrt. naphtyne-2,6 [kcal/mol]

∆∆E wrt. naphtyne-2,6 [kcal/mol]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


2.0 1.5 1.0 0.5 0.0

naphtyne-2,7 naphtyne-1,8 naphtyne-1,7

2.5 2.0 1.5 1.0 0.5 0.0

1.00e-05

1.00e-06

1.00e-07

1.00e-05

TCutPNO

1.00e-06

1.00e-07

TCutPNO

(b)

(a)

Figure 7: Relative ∆∆E computed as ∆Ex − ∆E2,6 , where ∆Ex = ExDLPNO−MkCCSD − ExMkCCSD calculated with a) T CutP airs = 1.0 × 10−4 , b) T CutP airs = 5.0 × 10−5 .

19



shown in Fig.6. With the default T CutP N O threshold, approximately 99.90% to 99.93% of the correlation energy is recovered for all four isomers. This represents a significant improvement over LPNO-MkCCSD, which recovered less than 99.8% (Ref. 102 ). The minimum is between 1.0 × 10−6 and 1.0 × 10−5 ; however, DLPNO-MkCCSD still recovers more than 99.8% of the correlation energy. The smallest value 99.82% is for 1,8-naphthyne at T CutP N O = 1.0 × 10−5 , while 2,6-naphtyne has the minimum of 99.87% at T CutP N O = 1.0 × 10−6 . The minimum for LPNO-MkCCSD was approximately 0.2% lower. 102 Next, the isomerization energies with respect to 2,6-naphtyne were calculated. Canonical MkCCSD yields values of 0.57 kcal/mol for 2,7-naphthyne, 2.18 kcal/mol for 1,6-naphthyne, and 5.03 kcal/mol for 1,7-naphthyne. For DLPNO-MkCCSD, the errors of isomerization energies with respect to canonical MkCCSD are listed in Fig.7(a) and Fig.7(b). At the default value of T CutP N O, the deviations for 2,7-naphtyne is 0.06 kcal/mol and for 1, 8and 1, 7-naphtyne approximately 0.26 kcal/mol. With TCutPairs equal to 5 × 10−5 and default T CutP N O, i. e. without significant increase of computational cost the error is reduced to approximately 0.17 kcal/mol for both 1,7- and 1,8-naphtyne. Importantly, when T CutP N O is reduced to 1.0 × 10−7 and TCutPairs to 5 × 10−5 , the error is reduced to only 0.06 kcal/mol or less for all three isomers. For 1,7-naphtyne, this represents an improvement over LPNO-MkCCSD result 102 by approximately 0.1 kcal/mol from 0.36 kcal/mol to 0.26 kcal/mol. For the other two isomers the effect is smaller, the for 1,8-naphtyne the error is increased by 0.05 kcal/mol and for 2,7naphtyne reduced by 0.02 kcal/mol. However, the tightening of T CutP N O at the LPNOMkCCSD level had only a small effect. An explanation for the improvement for tight values of cut-offs is likely due to the inclusion of the terms in cluster equations that are neglected at the LPNO level. These consist mainly of 3-external dressings of 4-external doubles contributions and several terms in dressing of 2-external doubles contribution. It can be expected that these terms are larger in multireference case, where the T1 amplitudes are in general larger due to nonzero Fock el-

20


Page 20 of 42

Page 21 of 42

ements. Fig.8(a) and Fig.8(b) show the change of the percentage of recovered correlation 100.50 1,8-naphthyne

100.00


100.40


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


100.30

100.20

100.10

TCutDO TCutMKN

99.90

99.80

99.70

99.60 100.00 1.00e-03

1.00e-04

1.00e-05

1.00e-06

1.00e-07

1e-01

TCutPairs

1e-02

1e-03

1e-04

1e-05

1e-06

Threshold

(a)

(b)

Figure 8: Average percentage of correlation energy recovered for 1,8-naphthyne isomer by DLPNO-MkCCSD with respect to (a) T CutP airs, (b) T CutDO and T CutM KN thresholds.

energy for 1,8-naphthyne as a function of T CutP airs, T CutDO, and T CutM KN thresholds. One can clearly see, that for T CutP airs the energy remains constant up to 5.0 × 10−5 . For larger thresholds, the energy increases monotonically as the number of neglected pairs grows significantly, and their pair energy is overestimated at the perturbative level. At the default value of 1.0 × 10−4 , the error is only 0.01%. T CutDO and T CutM KN have an opposite effect, because for larger values of cut-offs for the orbital domains (T CutDO) and fitting domains (T CutM KN ) start to lose essential components, which leads to substantial loss of the correlation energy contributions. The Fig.8(b) shows that for up to default values T CutDO = 1.0 × 10−2 and T CutDO = 1.0 × 10−3 the percentage is almost constant and for larger thresholds the correlation energy significantly decreases. 3.1.3

Potential curve for F2

In the process of generation of PNOs, the pair density matrices are averaged over contributing references. In principle, this could potentially lead to significant errors when the references have considerably different weights. In order to investigate this, we have calculated two

21



curves obtained by DLPNO-MkCCSD method in cc-pVQZ basis, one where all contributing references have the same weight and one where weights from preceding CAS SCF calculations are used. Here, the CAS SCF weights of the references change drastically along the curve. The difference between the two curves is shown in Fig.9. 100

80

∆E [µHartree]

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 42

60

40

20

0 1.5

2.0

2.5

3.0

3.5

4.0

4.5

Distance [Angstrom]

Figure 9: Difference between DLPNO-MkCCSD energy of F2 calculated in cc-pVQZ basis with a) all contributing references having the same weight in the construction of PNOs and b) reference weights obtained from preceding CAS SCF calculation.

As is clear from the graph, the difference between the two curves are very small. The largest difference of approximately 60 µHartree is observed at 1.6˚ A. Near the equilibrium (1.42 ˚ A), the differences are approximately 20 - 25 µHartree and in the dissociation limit, they disappear, since the weights are the same in both sets of calculations. Therefore, we conclude that the errors originating from the weighting method are significantly smaller than other errors encountered in DLPNO approaches.

22


Page 23 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(a)

(b)

Figure 10: The molecules (a) [Be(Me L)2 ], and (b) triangulene.

3.2

Application to larger systems

To illustrate the effectiveness of DLPNO-MkCCSD approach we have performed calculations on [Be(Me L)2 ] complex, where Me L =1-(2,6-diisopropylphenyl)-3,3,5,5-tetramethylpyrrolidine2-ylidene, and on triangulene diradical (See Fig.10). The [Be(Me L)2 ] complex was recently prepared as one of the first compounds containing Be(0) as a metal centre, where obtained results indicate that the stability is given by a strong three-centre two-electron bond C-Be-C (Ref. 108 ). The triangular Non-Kekule structures as triangulene have a high-spin state as ground state, in this particular case the triplet state. The first excited singlet state has an open-shell character, the published S-T energy gap is 0.6-0.8 eV at MR-CISD level (Ref. 109 ) or 0.16 eV at the spin-polarized DFT level (Ref. 110 ), respectively. The closed-shell singlet state has been predicted to be about 0.35 eV above the ground triplet state at spin-polarized DFT level (Ref. 110 ). In our calculations we focused at closed-shell singlet states. For a proper description of multireference character of singlet states we employed CMS(2,2) and utilize ccpVTZ basis set. The geometry of [Be(Me L)2 ] was taken from Ref., 108 for triangulene we used the geometry from Ref. 109 The triangulene has 828 basis functions, while [Be(Me L)2 ] complex were computed using 2158 basis functions. As the auxiliary basis set for RI projection we utilize cc-pVQZ/C giving 3564 or 9086 basis functions, respectively.

23



Page 24 of 42

Table 1: The total number of pairs (Pairs), number of surviving pairs, total number of amplitudes which would appear in canonical calculation, number of optimized PNO amplitudes in DLPNO-MkCCSD, average number of PNO per pair, the energy of the lowest singlet state, corresponding correlation energy, weak-pairs correction, and singlet-triplet gap energy. The S-T gap is defined as E singlet − E triplet . The triplet state energy was calculated using the single-reference openshell DLPNO-CCSD code (See Ref. 83 for method and implementation). System Number of atoms Number of basis functions Number of auxiliary basis functions Total number of RHF/UHF pairs Number of RHF/UHF pairs kept Total number of PNOs Number of canonical amplitudes Number of PNO amplitudes Average number of PNOs per (RHF) pair Total energy [au] Correlation energy [au] Right eigenvector of H eff Weak-pair energy correction [au] S-T gap [eV]

24

Triangulene 34 828 3564 2601/5151 1505/2959 59604 1807549275 4506028 (0.25%) 39 -843.767587 -3.493519 (-0.712,-0.703) -0.024 0.78


[Be(Me L)2 ] 105 2158 9086 13924/27730 3229/6322 99512 28474384260 5958516 (0.02%) 30 -1682.370676 -7.608054 (-0.937,0.350) -0.080 -0.42

Page 25 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The results obtained for default thresholds are shown in Tab.1. The number of RHF pairs for triangulene that were included in the coupled cluster treatment was 1505, i.e. 58% of the original 2601 pairs. Among the neglected pairs, 46 were excluded in the prescreening step. For [Be(Me L)2 ] the percentage of pairs kept was reduced to 23%. Also, the number of pairs dropped in prescreening has risen to 1536, i.e. 11% of the original pairs. In terms of the number of cluster amplitudes, for triangulene and DLPNO-MkCCSD method, 0.25% of the canonical amplitudes is needed. For larger [Be(Me L)2 ] the number is only 0.02%. The dropped pairs energy contribution represented approximately 70% of the weak pairs energy while the remaining 30% corresponds to neglected PNOs. The right eigenvectors of effective Hamiltonian are (-0.712,-0.703) for triangulene and (-0.937,0.350) for [Be(Me L)2 ], respectively, which confirms strong multireference character of both systems.

The singlet-triplet (S-T) gap obtained for triangulene is 0.78 eV, which is larger than recently published 0.35 eV in Ref. 110 This is in line with the situation for the lower open-shell singlet state, where the MR CISD value was also considerably higher than the spin-polarized DFT one (see above). Different situation was observed for [Be(Me L)2 ], where the singlet state is approximately 0.42 eV below the triplet state.

In order to quantify the influence of T CutP airs threshold (as discussed in previous section) on relative energy, we performed also [Be(Me L)2 ] calculations with T CutP airs = 5× 10− 5. The total DLPNO-MkCCSD energy increased for about 6.4 mHartree to -1682.364315 Hartree, whereas the single-reference triplet energy converged to -1682.349709 Hartree, which is for 5.6 mHartree higher than for the default cut-off. This results in a small change in S-T gap from -0.42 eV to -0.40 eV. On the other hand, the total number of UHF pairs treated at strongly-correlated level increased by 19% from 6322 to 7499, which negatively affects timings. We also performed a timing test on a single node (Tab.2), using 1 CPU core and 128 GB

25



Table 2: The computational timings for [Be(Me L)2 ]. The calculation was carried out with default thresholds, using 1 CPU core and 128 GB memory.

Fock Matrix Formation Organizing maps RI-PNO integral transformation Initial Guess State Vector Update Sigma-vector construction h0|H|Di h0|H|Si hD|H|Di(0-ext) hD|H|Di(2-ext Fock) hD|H|Di(2-ext) hD|H|Di(4-ext) hD|H|Di(4-ext-corr) CCSD doubles correction hS|H|Si hS|H|Di(1-ext) hD|H|Si(1-ext) hS|H|Di(3-ext) CCSD singles correction Fock-dressing Singles amplitudes (ik|jl)-dressing (ij|ab),(ia|jb)-dressing Pair energies

26

Time in sec 64540.3 6.9 195480.1 77744.7 73.1 421561.8 0.9 0.8 55771.4 296.2 145998.3 6913.1 64900.7 574.1 2163.4 4022.2 138.0 661.5 0.4 142704.3 587.3 15950.7 2680.3 20.9


Time in % 8.5 0.0 25.7 10.2 0.0 55.5 0.0 of sigma 0.0 of sigma 13.2 of sigma 0.1 of sigma 34.6 of sigma 1.6 of sigma 15.4 of sigma 0.1 of sigma 0.5 of sigma 1.0 of sigma 0.0 of sigma 0.2 of sigma 0.0 of sigma 33.9 of sigma 0.1 of sigma 3.8 of sigma 0.6 of sigma 0.0 of sigma

Page 26 of 42

Page 27 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of memory. As expected, the most expensive part of the calculation is an evaluation of sigma vectors followed by RI-PNO integral transformation. Prescreening and PNO construction took 10% of the total execution time, and 8.5% time was spent on Fock matrix formation. The total computation time was less than 9 days. As far as sigma construction, the most time was spent on linear two-external doubles contributions hD|H|Di(2-ext), dressing of Fock matrix, and three-external dressings of four-external doubles contributions hD|H|Di(4-extcorr), which is consistent with the single-reference DLPNO-CCSD code.

3.3

Conclusions

We have introduced a new DLPNO-MkCCSD approach, which is a Hilbert-space multireference method applicable to large systems (more than hundred atoms and thousands of basis functions). Even using only one CPU core, we can perform multireference CCSD calculation with 105 atoms and 2158 basis functions in less than 9 days, which is not feasible in canonical orbitals representation within any reasonable time. On two benchmark systems – TME and isomers of naphtynes – we have demonstrated the accuracy of the DLPNO-MkCCSD with respect to the canonical MkCCSD method. It is shown, that using default set of thresholds one could obtain more than 99.9% of the correlation energy. Also non-parallelity error tested on TME system is 0.17 kcal/mol and errors for isomerization energies of naphtynes are smaller than 0.3 kcal/mol. We also studied two larger systems - triangulene diradical and [Be(Me L)2 ] complex. Obtained singlet-triplet gap 0.78 eV for triangulene is larger than the recently published DFT energy, while for [Be(Me L)2 ] we found the singlet state as ground state with the S-T energy gap -0.42 eV. Compared to the previous LPNO-MkCCSD method, we have shown that DLPNO-MkCCSD is more accurate in terms of percentage of the correlation energy recovered and is applicable to larger systems due to almost linearly scaling computational cost. Furthermore, it is possible to develop a perturbative triples correction within DLPNO-MkCCSD scheme, that 27



has a computational cost comparable with the corresponding singles and doubles method, unlike in the LPNO-MkCCSD framework. This is particularly important, since it is known that the triples correction is highly beneficial for the accuracy of the MkCC method, and it will be the next step of our further development.

4

Acknowledgement

This work was supported by Czech Science Foundation (project 15-00058Y), DAAD (DAAD16-17), and Ministry of Education, Youth and Sports from the Large Infrastructures for Research, Experimental Development and Innovations (project IT4Innovations National Supercomputing Center - LM2015070).

References ˇ ıˇzek, J. On the Correlation Problem in Atomic and Molecular Systems. Calcula(1) C´ tion of Wavefunction Components in Ursell-Type Expansion Using Quantum-Field Theoretical Methods. The Journal of Chemical Physics 1966, 45, 4256–4266. (2) Jeziorski, B.; Monkhorst, H. J. Coupled-cluster method for multideterminantal reference states. Phys. Rev. A 1981, 24, 1668–1681. (3) Kucharski, S. A.; Bartlett, R. J. Hilbert space multireference coupled-cluster methods. I. The single and double excitation model. J. Chem. Phys. 1991, 95, 8227–8238. (4) Balková, A.; Kucharski, S. A.; Bartlett, R. J. The multi-reference Hilbert space coupled-cluster study of the Li2 molecule. Application in a complete model space. Chem. Phys. Lett. 1991, 182, 511–518. (5) Mahapatra, U. S.; Datta, B.; Mukherjee, D. A size-consistent state-specific multiref-

28


Page 28 of 42

Page 29 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


erence coupled cluster theory: Formal developments and molecular applications. J. Chem. Phys. 1999, 110, 6171–6188. (6) Máˇsik, J.; Hubaˇc, I. Multireference Brillouin-Wigner coupled-cluster theory. Singleroot approach. Adv. Quant. Chem. 1998, 31, 75–104. (7) Hanrath, M. An exponential multireference wave-function Ansatz. The Journal of Chemical Physics 2005, 123, 084102. (8) Banerjee, A.; Simons, J. The coupled-cluster method with a multiconfiguration reference state. International Journal of Quantum Chemistry 1981, 19, 207–216. (9) Hoffmann, M. R.; Simons, J. Analytical energy gradients for a unitary coupled-cluster theory. Chemical Physics Letters 1987, 142, 451 – 454. (10) Mukherjee, D. Normal ordering and a Wick-like reduction theorem for fermions with respect to a multi-determinantal reference state. Chemical Physics Letters 1997, 274, 561 – 566. (11) Evangelista, F. A.; Gauss, J. An orbital-invariant internally contracted multireference coupled cluster approach. The Journal of Chemical Physics 2011, 134, 114102. (12) Hanauer, M.; Köhn, A. Pilot applications of internally contracted multireference coupled cluster theory, and how to choose the cluster operator properly. The Journal of Chemical Physics 2011, 134, 204111. (13) Hanauer, M.; Köhn, A. Pilot applications of internally contracted multireference coupled cluster theory, and how to choose the cluster operator properly. The Journal of Chemical Physics 2011, 134, 204111. (14) Nooijen, M. State Selective Equation of Motion Coupled Cluster Theory: Some Preliminary Results. International Journal of Molecular Sciences 2002, 3, 656–675.

29



(15) Datta, D.; Nooijen, M. Multireference equation-of-motion coupled cluster theory. The Journal of Chemical Physics 2012, 137, 204107. (16) Demel, O.; Datta, D.; Nooijen, M. Additional global internal contraction in variations of multireference equation of motion coupled cluster theory. The Journal of Chemical Physics 2013, 138, 134108. (17) Nooijen, M.; Demel, O.; Datta, D.; Kong, L.; Shamasundar, K. R.; Lotrich, V.; Huntington, L. M.; Neese, F. Communication: Multireference equation of motion coupled cluster: A transform and diagonalize approach to electronic structure. The Journal of Chemical Physics 2014, 140, 081102. (18) Oliphant, N.; Adamowicz, L. Multireference coupled-cluster method using a singlereference formalism. The Journal of Chemical Physics 1991, 94, 1229–1235. (19) Piecuch, P.; Oliphant, N.; Adamowicz, L. A state-selective multireference coupledcluster theory employing the single-reference formalism. The Journal of Chemical Physics 1993, 99, 1875–1900. (20) Piecuch, P.; Adamowicz, L. State-selective multireference coupled-cluster theory employing the single-reference formalism: Implementation and application to the H8 model system. The Journal of Chemical Physics 1994, 100, 5792–5809. (21) Adamowicz, L.; Piecuch, P.; Ghose, K. B. The state-selective coupled cluster method for quasi-degenerate electronic states. Mol. Phys. 1998, 94, 225. (22) Piecuch, P. Active-space coupled-cluster methods. Mol. Phys. 2010, 108, 2987–3015. ˇ ıˇzek, J.; Takahashi, M. Approximate account of the connected quadruply (23) Paldus, J.; C´ excited clusters in the coupled-pair many-electron theory. Phys. Rev. A 1984, 30, 2193–2209.

30


Page 30 of 42

Page 31 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(24) Piecuch, P.; Tobola, R.; Paldus, J. Approximate account of connected quadruply excited clusters in single-reference coupled-cluster theory via cluster analysis of the projected unrestricted Hartree-Fock wave function. Phys. Rev. A 1996, 54, 1210–1241. (25) Li, X.; Paldus, J. Reduced multireference CCSD method: An effective approach to quasidegenerate states. J. Chem. Phys. 1997, 107, 6257. (26) Paldus, J.; Li, X. Externally corrected coupled-cluster approaches: Energy versus Amplitude corrected CCSD. Collect. Czech. Chem. Commun. 2003, 68, 554–586. (27) Kowalski, K.; Piecuch, P. The method of moments of coupled-cluster equations and the renormalized CCSD[T], CCSD(T), CCSD(TQ), and CCSDT(Q) approaches. J. Chem. Phys. 2000, 113, 18. (28) Kinoshita, T.; Hino, O.; Bartlett, R. J. Coupled-cluster method tailored by configuration interaction. J. Chem. Phys. 2005, 123, 074106. ¨ Pittner, J. Coupled Cluster (29) Veis, L.; Antal´ık, A.; Brabec, J.; Neese, F.; Legeza, O.; Method with Single and Double Excitations Tailored by Matrix Product State Wave Functions. The Journal of Physical Chemistry Letters 2016, 7, 4072–4078, PMID: 27682626. (30) Piecuch, P.; Kowalski, K.; Pimienta, I. S. O.; McGuire, M. J. Recent advances in electronic structure theory: Method of moments of coupled-cluster equations and renormalized coupled-cluster approaches. Int. Rev. Phys. Chem. 2002, 21, 527. (31) Piecuch, P.; Kowalski, K.; Pimienta, I. S. O. Method of Moments of Coupled-Cluster Equations: Externally Corrected Approaches Employing Configuration Interaction Wave Functions. Int. J. Mol. Sci. 2002, 3, 475–497. (32) Kowalski, K.; Piecuch, P. Extension of the method of moments of coupled-cluster equations to excited states: The triples and quadruples corrections to the equation31



of-motion coupled-cluster singles and doubles energies. J. Chem. Phys. 2002, 116, 7411. (33) Wloch, M.; Lodriguito, M. D.; Piecuch, P.; Gour, J. R. Two new classes of noniterative coupled-cluster methods derived from the method of moments of coupledcluster equations. Mol. Phys. 2006, 104, 2149, Erratum Mol. Phys. 104, 2991 (2006). (34) Lodriguito, M. D.; Kowalski, K.; Wloch, M.; Piecuch, P. Non-iterative coupled-cluster methods employing multi-reference perturbation theory wave functions. J. Molec. Struc. (Theochem) 2006, 771, 89–104. (35) Piecuch, P.; Kowalski, K.; Pimienta, I. S. O.; Fan, P.-D.; Lodriguito, M.; McGuire, M. J.; Kucharski, S. A.; Ku´s, T.; Musial, M. Method of moments of coupledcluster equations: a new formalism for designing accurate electronic structure methods for ground and excited states. Theor. Chem. Acc. 2004, 112, 349. (36) Kowalski, K.; Piecuch, P. New type of noniterative energy corrections for excited electronic states: extension of the method of moments of coupled-cluster equations to the equation-of-motion coupled-cluster formalism. J. Chem. Phys. 2001, 115, 2966– 2978. (37) Piecuch, P.; Kowalski, K. In Computational Chemistry: Reviews of Current Trends; Leszczy´ nski, J., Ed.; World Scientific, Singapore, 2011; Vol. 5; pp 1–104. (38) Mahapatra, U. S.; Datta, B.; Mukherjee, D. Development of a size-consistent statespecific multireference perturbation theory with relaxed model-space coefficients. Chem. Phys. Lett. 1999, 299, 42–50. (39) Evangelista, F. A.; Allen, W. D.; Schaefer III, H. F. High-order excitations in stateuniversal and state-specific multireference coupled-cluster theories: Model systems. J. Chem. Phys. 2006, 125, 154113.

32


Page 32 of 42

Page 33 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(40) Evangelista, F. A.; Allen, W. D.; Schaefer III, H. F. Coupling term derivation and general implementation of state-specific multireference coupled cluster theories. J. Chem. Phys. 2007, 127, 024102. (41) Evangelista, F. A.; Simmonett, A. C.; Allen, W. D.; Schaefer III, H. F.; Gauss, J. Triple excitations in state-specific multireference coupled cluster theory. J. Chem. Phys. 2008, 128, 124104. (42) Das, S.; Mukherjee, D.; Kállay, M. Full implementation and benchmark studies of Mukherjee’s state-specific multireference coupled cluster ansatz. J. Chem. Phys. 2010, 132, 074103. (43) Bhaskaran-Nair, K.; Demel, O.; Pittner, J. Multireference State-Specific Mukherjee’s coupled cluster method with non-iterative triexcitations. J. Chem. Phys. 2008, 129, 184105. ˇ (44) Bhaskaran-Nair, K.; Demel, O.; Smydke, J.; Pittner, J. Multireference state-specific Mukherjee’s coupled cluster method with noniterative triexcitations using uncoupled approximation. The Journal of Chemical Physics 2011, 134, 154106. (45) Prochnow, E.; Evangelista, F. A.; III, H. F. S.; Allen, W. D.; Gauss, J. Analytic gradients for the state-specific multireference coupled cluster singles and doubles model. The Journal of Chemical Physics 2009, 131, 064109. (46) Demel, O.; Kedzuch, S.; Svana, M.; Ten-no, S.; Pittner, J.; Noga, J. An explicitly correlated Mukherjee’s state specific coupled cluster method: development and pilot applications. Phys. Chem. Chem. Phys. 2012, 14, 4753–4762. (47) Das, S.; Bera, N.; Ghosh, S.; Mukherjee, D. An externally-corrected size-extensive single-root MRCC formalism: Its kinship with the rigorously size-extensive statespecific MRCC theory. Journal of Molecular Structure: THEOCHEM 2006, 771, 79

33



– 87, Modelling Structure and Reactivity: the 7th triennial conference of the World Association of Theoritical and Computational Chemists (WATOC 2005). (48) Demel, O.; Bhaskaran-Nair, K.; Pittner, J. Uncoupled multireference state-specific Mukherjee’s coupled cluster method with triexcitations. J. Chem. Phys. 2010, 133, 134106. (49) Brabec, J.; Bhaskaran-Nair, K.; Kowalski, K.; Pittner, J.; van Dam, H. J. J. Towards large-scale calculations with State-Specific Multireference Coupled Cluster methods: Studies on dodecane, naphthynes, and polycarbenes. Chem. Phys. Lett. 2012, 542, 128–133. (50) Lotrich, V.; Flocke, N.; Ponton, M.; Yau, A. D.; Perera, A.; Deumens, E.; Bartlett, R. J. Parallel implementation of electronic structure energy, gradient, and Hessian calculations. The Journal of Chemical Physics 2008, 128, 194104. (51) Deumens, E.; Lotrich, V. F.; Perera, A.; Ponton, M. J.; Sanders, B. A.; Bartlett, R. J. Software design of ACES III with the super instruction architecture. Wiley Interdisciplinary Reviews: Computational Molecular Science 2011, 1, 895–901. (52) Kowalski, K.; Krishnamoorthy, S.; Olson, R. M.; Tipparaju, V.; Aprà, E. Scalable Implementations of Accurate Excited-state Coupled Cluster Theories: Application of High-level Methods to Porphyrin-based Systems. Proceedings of 2011 International Conference for High Performance Computing, Networking, Storage and Analysis. New York, NY, USA, 2011; pp 72:1–72:10. (53) Aprà, E.; Klemm, M.; Kowalski, K. Efficient Implementation of Many-body Quantum Chemical Methods on the Intel Xeon Phi Coprocessor. Proceedings of the International Conference for High Performance Computing, Networking, Storage and Analysis. Piscataway, NJ, USA, 2014; pp 674–684.

34


Page 34 of 42

Page 35 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(54) Brabec, J.; Krishnamoorthy, S.; van Dam, H.; Kowalski, K.; Pittner, J. Massively parallel implementation of the Brillouin-Wigner MRCCSD method. Chem. Phys. Lett. 2011, 514, 347–351. (55) See https://github.com/ValeevGroup/tiledarray for a massively-parallel block-sparse tensor network written in C++, accessed May 13, 2017. (56) Solomonik, E.; Matthews, D.; Hammond, J. R.; Stanton, J. F.; Demmel, J. A massively parallel tensor contraction framework for coupled-cluster computations. Journal of Parallel and Distributed Computing 2014, 74, 3176 – 3190, Domain-Specific Languages and High-Level Frameworks for High-Performance Computing. (57) Pulay, P. Localizability of dynamic electron correlation. Chemical Physics Letters 1983, 100, 151 – 154. (58) Sæbø, S.; Pulay, P. Local configuration interaction: An efficient approach for larger molecules. Chemical Physics Letters 1985, 113, 13 – 18. (59) Boys, S. F. Construction of Some Molecular Orbitals to Be Approximately Invariant for Changes from One Molecule to Another. Rev. Mod. Phys. 1960, 32, 296–299. (60) Pipek, J.; Mezey, P. G. A fast intrinsic localization procedure applicable for ab initio and semiempirical linear combination of atomic orbital wave functions. The Journal of Chemical Physics 1989, 90, 4916–4926. (61) Hampel, C.; Werner, H.-J. Local treatment of electron correlation in coupled cluster theory. The Journal of Chemical Physics 1996, 104, 6286–6297. (62) Sch¨ utz, M.; Werner, H.-J. Low-order scaling local electron correlation methods. IV. Linear scaling local coupled-cluster (LCCSD). The Journal of Chemical Physics 2001, 114, 661–681.

35



(63) Schutz, M. A new, fast, semi-direct implementation of linear scaling local coupled cluster theory. Phys. Chem. Chem. Phys. 2002, 4, 3941–3947. (64) Werner, H.-J.; Sch¨ utz, M. An efficient local coupled cluster method for accurate thermochemistry of large systems. The Journal of Chemical Physics 2011, 135, 144116. (65) Yang, J.; Chan, G. K.-L.; Manby, F. R.; Sch¨ utz, M.; Werner, H.-J. The orbital-specificvirtual local coupled cluster singles and doubles method. The Journal of Chemical Physics 2012, 136, 144105. (66) Fedorov, D. G.; Kitaura, K. Coupled-cluster theory based upon the fragment molecular-orbital method. The Journal of Chemical Physics 2005, 123, 134103. (67) Kobayashi, M.; Nakai, H. Extension of linear-scaling divide-and-conquer-based correlation method to coupled cluster theory with singles and doubles excitations. The Journal of Chemical Physics 2008, 129, 044103. (68) Kristensen, K.; Ziólkowski, M.; Jans´ık, B.; Kjærgaard, T.; Jørgensen, P. A Locality Analysis of the Divide–Expand–Consolidate Coupled Cluster Amplitude Equations. Journal of Chemical Theory and Computation 2011, 7, 1677–1694, PMID: 26596432. (69) Høyvik, I.-M.; Kristensen, K.; Jansik, B.; Jørgensen, P. The divide-expand-consolidate family of coupled cluster methods: Numerical illustrations using second order MøllerPlesset perturbation theory. The Journal of Chemical Physics 2012, 136, 014105. (70) Stoll, H. The correlation energy of crystalline silicon. Chemical Physics Letters 1992, 191, 548 – 552. (71) Li, S.; Ma, J.; Jiang, Y. Linear scaling local correlation approach for solving the coupled cluster equations of large systems. Journal of Computational Chemistry 2002, 23, 237–244.

36


Page 36 of 42

Page 37 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(72) Edmiston, C.; Krauss, M. Configuration-Interaction Calculation of H3 and H2. The Journal of Chemical Physics 1965, 42, 1119–1120. (73) Ahlrichs, R.; Kutzelnigg, W. Ab initio calculations on small hydrides including electron correlation. Theoretica chimica acta 1968, 10, 377–387. (74) Meyer, W. Ionization energies of water from PNO-CI calculations. International Journal of Quantum Chemistry 1971, 5, 341–348. (75) Meyer, W. PNO-CI and CEPA studies of electron correlation effects. Theoretica chimica acta 1974, 35, 277–292. (76) Taylor, P. R. A rapidly convergent CI expansion based on several reference configurations, using optimized correlating orbitals. The Journal of Chemical Physics 1981, 74, 1256–1270. (77) Fink, R.; Staemmler, V. A multi-configuration reference CEPA method based on pair natural orbitals. Theoretica chimica acta 1993, 87, 129–145. (78) Vahtras, O.; Almlöf, J.; Feyereisen, M. Integral approximations for LCAO-SCF calculations. Chemical Physics Letters 1993, 213, 514 – 518. (79) Neese, F.; Wennmohs, F.; Hansen, A. Efficient and accurate local approximations to coupled-electron pair approaches: An attempt to revive the pair natural orbital method. The Journal of Chemical Physics 2009, 130, 114108. (80) Neese, F.; Hansen, A.; Liakos, D. G. Efficient and accurate approximations to the local coupled cluster singles doubles method using a truncated pair natural orbital basis. The Journal of Chemical Physics 2009, 131, 064103. (81) Riplinger, C.; Neese, F. An efficient and near linear scaling pair natural orbital based local coupled cluster method. The Journal of Chemical Physics 2013, 138, 034106.

37



(82) Riplinger, C.; Sandhoefer, B.; Hansen, A.; Neese, F. Natural triple excitations in local coupled cluster calculations with pair natural orbitals. The Journal of Chemical Physics 2013, 139, 134101. (83) Saitow, M.; Becker, U.; Riplinger, C.; Valeev, E. F.; Neese, F. A new near-linear scaling, efficient and accurate, open-shell domain-based local pair natural orbital coupled cluster singles and doubles theory. The Journal of Chemical Physics 2017, 146, 164105. (84) Pinski, P.; Riplinger, C.; Valeev, E. F.; Neese, F. Sparse maps—A systematic infrastructure for reduced-scaling electronic structure methods. I. An efficient and simple linear scaling local MP2 method that uses an intermediate basis of pair natural orbitals. The Journal of Chemical Physics 2015, 143, 034108. (85) Riplinger, C.; Pinski, P.; Becker, U.; Valeev, E. F.; Neese, F. Sparse maps—A systematic infrastructure for reduced-scaling electronic structure methods. II. Linear scaling domain based pair natural orbital coupled cluster theory. The Journal of Chemical Physics 2016, 144, 024109. (86) Guo, Y.; Sivalingam, K.; Valeev, E. F.; Neese, F. SparseMaps—A systematic infrastructure for reduced-scaling electronic structure methods. III. Linear-scaling multireference domain-based pair natural orbital N-electron valence perturbation theory. The Journal of Chemical Physics 2016, 144, 094111. (87) Werner, H.-J.; Knizia, G.; Krause, C.; Schwilk, M.; Dornbach, M. Scalable Electron Correlation Methods I.: PNO-LMP2 with Linear Scaling in the Molecular Size and Near-Inverse-Linear Scaling in the Number of Processors. Journal of Chemical Theory and Computation 2015, 11, 484–507, PMID: 26580908. (88) Menezes, F.; Kats, D.; Werner, H.-J. Local complete active space second-order pertur-

38


Page 38 of 42

Page 39 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


bation theory using pair natural orbitals (PNO-CASPT2). The Journal of Chemical Physics 2016, 145, 124115. (89) Helmich, B.; Hättig, C. A pair natural orbital implementation of the coupled cluster model CC2 for excitation energies. The Journal of Chemical Physics 2013, 139, 084114. (90) Schmitz, G.; Helmich, B.; Hättig, C. A scaling PNO–MP2 method using a hybrid OSV–PNO approach with an iterative direct generation of OSVs. Molecular Physics 2013, 111, 2463–2476. (91) Schwilk, M.; Ma, Q.; Köppl, C.; Werner, H.-J. Scalable Electron Correlation Methods. 3. Efficient and Accurate Parallel Local Coupled Cluster with Pair Natural Orbitals (PNO-LCCSD). Journal of Chemical Theory and Computation 2017, 13, 3650–3675, PMID: 28661673. (92) Antony, J.; Grimme, S.; Liakos, D. G.; Neese, F. Protein–Ligand Interaction Energies with Dispersion Corrected Density Functional Theory and High-Level Wave Function Based Methods. The Journal of Physical Chemistry A 2011, 115, 11210–11220, PMID: 21842894. (93) Anoop, A.; Thiel, W.; Neese, F. A Local Pair Natural Orbital Coupled Cluster Study of Rh Catalyzed Asymmetric Olefin Hydrogenation. Journal of Chemical Theory and Computation 2010, 6, 3137–3144, PMID: 26616776. (94) Liakos, D. G.; Neese, F. Interplay of Correlation and Relativistic Effects in Correlated Calculations on Transition-Metal Complexes: The (Cu2O2)2+ Core Revisited. Journal of Chemical Theory and Computation 2011, 7, 1511–1523, PMID: 26610142. (95) Zade, S. S.; Zamoshchik, N.; Reddy, A. R.; Fridman-Marueli, G.; Sheberla, D.; Bendikov, M. Products and Mechanism of Acene Dimerization. A Computational

39



Study. Journal of the American Chemical Society 2011, 133, 10803–10816, PMID: 21710966. (96) Kubas, A.; Brse, S.; Fink, K. Theoretical Approach Towards the Understanding of Asymmetric Additions of Dialkylzinc to Enals and Iminals Catalysed by [2.2]Paracyclophane-Based N,O-Ligands. Chemistry A European Journal 2012, 18, 8377–8385. (97) Ashtari, M.; Cann, N. Proline-based chiral stationary phases: A molecular dynamics study of the interfacial structure. Journal of Chromatography A 2011, 1218, 6331 – 6347. (98) Zhang, J.; Dolg, M. Dispersion Interaction Stabilizes Sterically Hindered Double Fullerenes. Chemistry A European Journal 2014, 20, 13909–13912. (99) Minenkov, Y.; Chermak, E.; Cavallo, L. Accuracy of DLPNO–CCSD(T) Method for Noncovalent Bond Dissociation Enthalpies from Coinage Metal Cation Complexes. Journal of Chemical Theory and Computation 2015, 11, 4664–4676, PMID: 26574257. (100) Sparta, M.; Neese, F. Chemical applications carried out by local pair natural orbital based coupled-cluster methods. Chem. Soc. Rev. 2014, 43, 5032–5041. (101) Liakos, D. G.; Sparta, M.; Kesharwani, M. K.; Martin, J. M. L.; Neese, F. Exploring the Accuracy Limits of Local Pair Natural Orbital Coupled-Cluster Theory. Journal of Chemical Theory and Computation 2015, 11, 1525–1539. (102) Demel, O.; Pittner, J.; Neese, F. A Local Pair Natural Orbital-Based Multireference Mukherjee’s Coupled Cluster Method. Journal of Chemical Theory and Computation 2015, 11, 3104–3114, PMID: 26575747. (103) Angeli, C.; Cimiraglia, R.; Evangelisti, S.; Leininger, T.; Malrieu, J.-P. Introduction

40


Page 40 of 42

Page 41 of 42 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of n-electron valence states for multireference perturbation theory. The Journal of Chemical Physics 2001, 114, 10252–10264. (104) Angeli, C.; Cimiraglia, R.; Malrieu, J.-P. n-electron valence state perturbation theory: A spinless formulation and an efficient implementation of the strongly contracted and of the partially contracted variants. The Journal of Chemical Physics 2002, 117, 9138–9153. (105) Angeli, C.; Pastore, M.; Cimiraglia, R. New perspectives in multireference perturbation theory: the n-electron valence state approach. Theoretical Chemistry Accounts 2007, 117, 743–754. (106) Dunning, T. H. The Cl + Bk extrapolation method. Application to hydrogen fluoride11Work performed under the auspices of the Office of Basic Energy Sciences of the U.S. Department of Energy. Chemical Physics 1979, 42, 249 – 258. (107) Pozun, Z. D.; Su, X.; Jordan, K. D. Establishing the Ground State of the Disjoint Diradical Tetramethyleneethane with Quantum Monte Carlo. Journal of the American Chemical Society 2013, 135, 13862–13869, PMID: 23947763. (108) Arrowsmith, M.; Braunschweig, H.; Celik, M. A.; Dellermann, T.; Dewhurst, R. D.; Ewing, W. C.; Hammond, K.; Kramer, T.; Krummenacher, I.; Mies, J.; Radacki, K.; Schuster, J. K. Neutral zero-valent s-block complexes with strong multiple bonding. Nat Chem 2016, 8, 890–894. (109) Das, A.; M¨ uller, T.; Plasser, F.; Lischka, H. Polyradical Character of Triangular NonKekulé Structures, Zethrenes, p-Quinodimethane-Linked Bisphenalenyl, and the Clar Goblet in Comparison: An Extended Multireference Study. The Journal of Physical Chemistry A 2016, 120, 1625–1636, PMID: 26859789. (110) Pavliˇcek, N.; Mistry, A.; Majzik, Z.; Moll, N.; Meyer, G.; Fox, D. J.; Gross, L. Synthesis and characterization of triangulene. Nat Nano 2017, 12, 308–311. 41



Figure 11: For Table of Contents Only

42


Page 42 of 42

Domain-Based Local Pair Natural Orbital Version of Mukherjee's State

Recommend Documents