Testing Exact Upper Bounds to Exact Exchange - Journal of Chemical

Department of Chemistry and Biochemistry, Texas Christian University, Fort Worth, Texas ... For a more comprehensive list of citations to this article...
0 downloads 0 Views 720KB Size
Subscriber access provided by UB + Fachbibliothek Chemie | (FU-Bibliothekssystem)

Article

Testing exact upper bounds to exact exchange Emil Proynov, and Benjamin G. Janesko J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.7b00095 • Publication Date (Web): 31 Mar 2017 Downloaded from http://pubs.acs.org on April 1, 2017

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 34

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

Testing exact upper bounds to exact exchange Emil Proynov∗ and Benjamin G. Janesko∗ Department of Chemistry and Biochemistry, Texas Christian University, Fort Worth, TX 76110, USA E-mail: [email protected]; [email protected]

Abstract The exact exchange energy and its energy density are useful but computationally expensive ingredients in density functional approximations for Kohn-Sham density functional theory. We present detailed tests of some exact nonempirical upper bounds to exact exchange. These "Rung 3.5" upper bounds contract the Kohn-Sham one-particle density matrix with model density matrices used to construct semilocal model exchange holes, and invoke the Cauchy-Schwarz inequality. The contraction automatically eliminates the computationally expensive long-range component of the exact exchange hole. Numerical tests show that the exchange upper bounds underestimate total exchange energies while predict other properties with accuracy approaching standard hybrid approximations.

1 1.1

Introduction Density functional theory

Kohn-Sham density functional theory (DFT) has become one of the most widely used approximations in computational chemistry and solid-state physics. 1,2 Simple "semilocal" density 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

Page 2 of 34

functional approximations (DFAs) prove to capture some of the nondynamical left-right correlation critical to chemical bonding, by replacing the delocalized exact exchange hole with fairly localized model DFT exchange hole. 3,4 The nonsystematic nature of these approximations prompted a search for improvement via an admixture with exactly computed exchange. DFAs incorporating exact exchange range from the global hybrids standard in computational chemistry, 5–7 to the screened hybrids that are increasingly popular for semiconductors, 8,9 to long-range-corrected hybrids accurate for anions and spectroscopy, 10,11 to the latest more sophisticated local hybrid and hyper-GGA (Generalized Gradient Approximation) approaches, 12,13 some of which are capable of explicitly simulating left-right correlation. 14–16

1.2

Exact exchange

One stumbling block to continued DFA development is the computational expense of exactly computed exchange. Modern hyper-GGAs rely on the conventional-gauge 17 σ-spin exchange energy density at point ⃗r, 1 eXσ (⃗r) = − 2 ′



γσ (⃗r, ⃗r ) = Nσ

∫ ds1

∫ 3

d ⃗x2 ..

∫ d3⃗r

r, ⃗r ′ |γσ (⃗



)|2 , |⃗r − ⃗r ′ |

(1)

d3⃗xNσ σ(s1 )Φ(⃗rs1 , ⃗x2 ..⃗xNσ )Φ∗ (⃗r ′ s1 , ⃗x2 ..⃗xNσ )

(2)

Here γσ (⃗r, ⃗r ′ ) is the σ-spin one-particle density matrix, Φ is the ground-state wavefunction of the Kohn-Sham system of N noninteracting Fermions, ⃗xi denotes a combined space ⃗r and spin s coordinate of electron i, and σ(s1 ) projects onto the desired σ spin. Like any energy density, the exchange energy density is nonunique. 17 The delocalization of γσ (⃗r, ⃗r ′ ) over |⃗r − ⃗r ′ | makes Eq.1 computationally expensive in systems such as bulk metals. 8 Semilocal DFAs instead use a model density matrix localized about the reference point ⃗r:

r) emod Xσ (⃗

1 = − 2

∫ d3⃗r

r), |∇ρσ (⃗r)|, ⃗r ′ |γmodσ (ρσ (⃗ |⃗r − ⃗r ′ |

ACS Paragon Plus Environment

− ⃗r ′ )|2

(3)

Page 3 of 34

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

Here γmodσ is a function of the electron density and its gradient norm at reference point ⃗r. Localization of γmodσ about the reference point ⃗r helps make Eq.3 computationally tractable, and accounts for DFAs’ aforementioned crude but effective simulation of nondynamical correlation. (Briefly, the localized exchange hole is incorrect in stretched H+ 2 where the σ bond orbital and γσ delocalize over both centers, but simulates the localizing effects of left-right correlation in stretched singlet H2 . 4,18 ) Hybrid DFAs are based on various combinations of Eq.1 and Eq.3. Here we focus particularly on local-hybrid (hyper-GGA) forms, 19–21 where the gain of substituting for the exact exchange energy density may be significant:

ehyb r) = a(⃗r) exσ (⃗r) + (1 − aσ (⃗r)) emod r) , xσ (⃗ xσ (⃗

(4)

The bottleneck of using the exact-exchange energy density is in the cost and in the numerical problems that are met in extended systems. For the past several years, we have worked to develop "Rung 3.5" DFAs based on a density matrix product ("Π") ansatz:

eΠ−mod (⃗r) Xσ

1 = − 2

∫ d3⃗r

r, ⃗r ′ γσ (⃗



)γmodσ (ρσ (⃗r), ∇ρσ (⃗r), ⃗r − ⃗r ′ ) |⃗r − ⃗r ′ |

(5)

We have presented nonempirical derivations of these DFAs 22 and a practical auxiliary basis set implementation up to analytic second derivatives. 23 Minimally empirical versions provide accuracy between semilocal and hybrid DFAs for a broad range of properties. 23

1.3

An exact upper bound to exact exchange

The present work extends a useful property of Eq.5 first pointed out in ref. 24 The CauchySchwarz inequality ensures that Eqs.3,5, provide a negative-semidefinite upper bound to the conventional-gauge exact exchange energy density |eΠ−mod (⃗r)|2 Xσ ≡ eupb−mod (⃗r) . eXσ (⃗rACS ) ≤Paragon Xσ modPlus Environment eXσ (⃗r)

(6)

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

Page 4 of 34

It is important to stress that in these new exchange energy densities no exact exchange is actually involved and no calculation of four-center exchange integrals is needed. Being a strict upper bound to the (finite) exact exchange assures regularity and finiteness of eupb−mod (⃗r) Xσ for any choice of the semi-local model density matrix. ref. 24 used this upper bound in nonempirical DFAs based on the "local" Lieb-Oxford bound. 25–27 Here we explore its use more generally in constructing modern DFAs. An exact upper bound to the low-density uniform electron gas exchange-correlation energy was very recently reported in ref. 28 In this work we do not consider how different gauge terms added to eupb−mod would affect Xσ the results, a point that we will address in the future in view of the recent evidence for its importance in local hybrids design. 21

1.4

Tested DFAs

We employ two choices of model density matrix in Eq.6, first presented in ref. 23 Using a spherical-symmetric Gaussian approximation to the Local Density Approximation (LDA) exchange hole leads to ( ) γσLDA (ρ(⃗r), ⃗s) = ρσ (⃗r) exp −bρ2/3 r)|⃗s|2 σ (⃗

(7)

Here ⃗s is the electron-electron distance vector ⃗s = ⃗r2 − ⃗r1 , and b = 1.483869 ensures Eq.5 is exact in the uniform electron gas. Using this LDA form of γ as zeroth order in a GGA type expansion, the following nonspherical form of γ was also suggested: 23 ) ( 1 ⃗ σ (⃗r)) exp −f (Sσ ) ρ2/3 γσGGA (ρσ (⃗r), ∇ρσ (⃗r), ⃗s) = γσLDA (ρσ (⃗r), ⃗s) + (⃗s · ∇ρ r) s2 σ (⃗ 2

(8)

Here f (Sσ ) is a function of the reduced density gradient −4/3

Sσ (⃗r) = 0.5|∇ρσ (⃗r)|ρσ

(⃗r)(6π 2 )−1/3

used to adjust the model so as to recover the Perdew-Burke-Ernzerhof generalized gradient 29 Plus Environment Paragon approximation (PBE GGA) forACS exchange. Nonspherical forms of the one-electron density

Page 5 of 34

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

matrix and the corresponding exchange hole, Eq.3, arguably describe more realistically the electron interactions in nonuniform systems. It is particularly so in the context of Rung 3.5 density matrix product used in Eq.5. The result of the integration of the Rung 3.5 energy density depends sensibly on the sphericity or nonspericity of each density matrix involved in the product. We compare the Rung 3.5 DFAs constructed from Eqs.7-8 based on Eq.5, with the corresponding upper bound DFAs constructed form Eq.6. In the latter case the parameter b in the expression for γσLDA in Eqs.7,8, is rescaled to b = 2.24470. This is required for the upper-bound exchange, Eq.6, to about recover the uniform-gas limit. The asymptotic decay of the exchange energy density in Eqs.5,6 at long range was analysed in ref. 30 It was shown that the Rung 3.5 exchange energy density decays faster than either the Hartre-Fock (HF) or LDA one. Because of this, the upper-bound exchange energy density, Eq.6, underestimate substantially the exact exchange in density tails. Correcting these functionals at long range is next on our agenda. We consider local and global hybrid forms in which the exact-exchange energy density Π−mod 31,32 in Eq.4 is replaced by either the 3.5 Rung form exσ or by the corresponding upper

bound eupb−mod , Eq.6: xσ

eΠ−hyb (⃗r) = aσ (⃗r) eΠ−mod + (1 − aσ (⃗r)) emod r) , xσ xσ xσ (⃗

(9)

eupb−hyb (⃗r) = aσ (⃗r) eupb−mod + (1 − aσ (⃗r)) emod r) , xσ xσ xσ (⃗

(10)

We also test the 3.5 Rung exchange and the corresponding exact upper-bound exchange on stand alone (aσ (⃗r) = 1) to examine similarities and differences with the Hartree-Fock/exact exchange. Various choices of the local mixing function aσ (⃗r) in Eqs.9,10 can be constructed using either already published forms, 21,31–35 or trying new ones. Local mixing functions that were shown to work well with the exact exchange energy density may not necessarily do so with the Rung 3.5 model energy densities used in Eqs.9,10. We explore here a very simple local 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

Page 6 of 34

mixing function of the form:

aσ (⃗r) = α exp(−β Sσ2 ) ,

(11)

where Sσ is the reduced density gradient and α and β are empirical parameters that are optimized separately for each given 3.5 Rung local hybrid combination. A simpler version with α = 1 was tested in ref. 23 and it was found quite suitable for this kind of local hybrids. This mixing function may look somewhat unusual but it turned out quite useful in the present context. It scales down the eupb contributions from regions with large reduced x density gradient S, where the Rung 3.5 and upper bound exchange energies significantly underestimate accurate values.

2

Computational Methods

Efficient implementation of the Rung 3.5 DFAs was reported in ref. 23 using RI fitting with a finite set of auxiliary Gaussian functions {χ˜λ (⃗s)} to separate the ⃗r and ⃗s variables in γmod (spin-indices are omitted below for simplicity but implied):

γmod (ρ(⃗r), ∇ρ(⃗r), ⃗s) =



−1 mod χ˜λ (⃗s)Sλδ cδ (ρ(⃗r), ∇ρ(⃗r)) ,

(12)

λ,δ −1 Sλδ =



1 d3 ⃗s χ˜λ (⃗s) χ˜δ (⃗s) . s

(13)

For the spherically symmetric LDA-based model density matrix the fitting coefficients in Eq.12 with s-type Gaussian auxiliary functions were obtained in analytic form in ref. 23

(ρ(⃗r)) = ρ(⃗r) cLDA λ

2π αλ + b ρ2/3 (⃗r)

(14)

where b = 1.483869 (b = 2.24470 for the upper-bound as explained above), αλ is the orbital exponent of uncontracted s-type auxiliary function: ACS Paragon Plus Environment

Page 7 of 34

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

χ˜λ,s (⃗s) = exp(−αλ s2 ) The nonspherical component of γGGA , Eq.8 is adjusted to recover the enhancement factor of the PBE exchange functional, and is expanded in a set of p-type Gaussian functions: χ˜λ,px (⃗s) = sx exp(−αλ s2 ). The corresponding fitting coefficients have the form 23 (see Eq.8): ( BE cPλ,px

(ρ(⃗r), ∇ρ(⃗r)) =

∂ρ(⃗r) ∂x

)

π 3 (αλ + f (S) ρ2/3 (⃗r))

2

(15)

Using auxiliary-basis fitting allows for analytical separation of the ⃗r and ⃗s dependence of the model density matrix γmod . This significantly boosts the efficiency of the implementation that otherwise would be impeded by the ρ(⃗r) dependence of the exponents in Eqs.7,8. 23 The convergence of the results with respect to the fitting was thoroughly checked. In all cases an even-tempered set of about 13 to 17 functions was sufficient.

3

Results and discussion

Our main goal is to provide an alternative to the exact-exchange energy density. To this end we focus on how well the Rung 3.5upper bound, Eq.6, can be a substitute for exact exchange in different situations and for different properties. We try several functional combinations labeled as follows: ’ΠLDAup’ denotes the upper-bound exchange originating from the 3.5 Rung LDA-based model ’ΠLDA’, Eqs.7,6, 23 and combined with Perdew-Burke-Ernzerhof (PBE) correlation; ’ΠPBEup’ is the upper-bound exchange related to the Rung 3.5 exchange ΠPBE 23 using the GGA model density matrix, Eq.8 and combined with PBE correlation; ’ΠGG2’ and ’ΠGG2Fx’ are 3.5 Rung Π exchange functionals developed in ref. 36 that are not explicitly described here and used just for the sake of comparison; ’ΠPBE1’ is a global hybrid combination between Rung 3.5 exchange ΠPBE and PBE exchange, combined with PBE correlation; ’ΠPBEup1’ is the analogous global hybrid combination between the ΠPBEup upper-bound 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

exchange and the PBE exchange, combined with PBE correlation. ’ΠPBE1w’ is a local hybrid combination, Eqs.9,11, between Rung 3.5 exchange ΠPBE and PBE exchange, combined with PBE correlation; ’ΠPBEup1w’ is the analogous local hybrid combination, Eq.10, involving the ΠPBEup upperbound exchange and PBE exchange combined with PBE correlation; ’PBE0’ is the popular nonempirical global hybrid based on exact exchange and the PBE exchange and correlation.

3.1

Exchange energy of atoms

Table 1 reports various DFAs’ predictions for the total exchange energies of first and second row atoms. Our Hartree-Fock /cc-pVQZ values are quite close to the Optimized Effective Potential (OEP) estimates of ref. 37 The accuracy of the upper bound depends sensibly on the type of model density matrix used in the corresponding similarity metric 5. The mean unsigned error (MUE) with respect to the HF estimate decreases in the order:

ΠLDAup > ΠPBEup > ΠLDA > ΠGG2 > ΠPBE > ΠGG2Fx

The Rung 3.5 functional ΠGG2Fx has the smallest MUE here. It has no upper-bound counterpart because it involves GGA enhancement factor. 36 The ΠPBE Rung 3.5 exchange is the next closest to the HF exchange for these atoms. The upper-bound versions give somewhat larger MUE here because these functionals are less accurate in the region of large reduced density gradient which is especially pronounced in atoms.

3.2

Atomization energies and reaction barriers

Next we consider the small, representative AE6 set of molecular atomization energies and BH6 set of reaction barriers. 38 These test sets lead to parameter values that are close to 23,38 Environment those optimized on larger test ACS setsParagon like G2.Plus Table 2 reports the mean unsigned error

Page 8 of 34

Page 9 of 34

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

MUE from various DFAs. The results with most of the Rung 3.5 functionals are obtained in two ways (i) with post-PBE0 single iteration using the PBE0 SCF density, and (ii) with SCF solution, which gives an insight on the influence of the input electron density and the SCF potential. ΠPBE (100% Rung 3.5 exchange) and ΠPBEup (100% upper bound exchange) outperform Hartree-Fock exchange for atomization energies and underestimate reaction barriers, indicating that both are more "hybrid-like" than Hartree-Fock theory itself. The one-parameter global-hybrid forms ΠPBE1 and ΠPBEup1 give performance approaching the PBE0 global hybrid. Figure 1 illustrates the performance as a function of α, while Table 2 gives the optimal results. The local hybrid schemes ΠPBE1w/PBE0 and ΠPBEup1w/PBE0 based on Eq.9-10 are the most accurate Rung 3.5 functionals in these tests. Local hybrids were introduced with the intent to improve over the global hybrids, which in practice has rarely been the case, with the important exception of the schemes developed recently in refs. 21,35 Here we show that using Rung 3.5 functionals in local-hybrid combinations with PBE exchange-correlation leads to improvement over the respective Rung 3.5 global hybrids as well.

3.3

Enthalpies of formation of transition-metal compounds: the TM193 set

After optimizing the Rung 3.5 functionals on the AE6 and BH6 sets, we assess further their performance on the much larger test set of 193 transition-metal (TM) compounds considered in ref. 39 This is a stringent test for any theoretical method due to possible multireference effects and instability of the ground state in many of these compounds. DFT methods are known to have problems in describing transition-metal-ligand bonding. 40 The calculated MUE values for the enthalpy of formation on the TM193 set are presented in Table 3. The one-parameter Rung 3.5 functionals ΠPBE, ΠPBE1 and ΠPBEup1 all outperform the PBE GGA, an important improvement achieved without exact exchange. The conventional ACS Paragon Plus Environment hybrid functionals PBE0 and B3LYP outperform the Rung 3.5 functionals at the expense

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

of exact exchange, which is potentially computationally expensive for large transition-metal compounds and solids. The pure upper-bound exchange ΠPBEup on stand-alone combination with PBE correlation gives somewhat larger MUE in these tests. Using just a 0.25 fraction of the upper-bound exchange as in ΠPBEup1 brings down the mean absolute error considerably.

3.4

Highest occupied orbital energies and HOMO-LUMO gap

Besides atomization energies and enthalpy of formation, useful characteristics of a given approximate functional are the energy of the highest occupied Kohn-Sham orbital(HOMO) and the gap to the lowest unoccupied orbital (LUMO). Recent studies reiterate the physics behind Kohn-Sham orbital energies, 41,42 confirming that the exact HOMO energy should equal the minus first vertical ionization potential (Ip). 43 DFAs without exact exchange tend to underestimate Ip, while hybrids tend to improve it. Table 4 and Table S1 of our Support Information (SI) contain the mean unsigned error for the HOMO energies of a sample set of atoms and molecules, including a few TM atoms and compounds. DFAs with 100% Rung 3.5 exchange ΠPBE and ΠPBEup give HOMO energies that are closer to Hartree-Fock estimates, and to experimental values, than the PBE GGA. Moreover, Rung 3.5 DFAs give HOMO-LUMO gaps approaching the PBE0 global hybrid (MUE3 in Table 4), and showing notable improvement over PBE. (The PBE0 HOMO-LUMO gap is a good approximation for the experimental optical gap, 41 while the Hartree-Fock gap better approximates the experimental fundamental gap. 41,42,44 )

3.5

Fractional charge and spin errors

The dissociation of H+ 2 and singlet H2 illustrates DFAs’ trade-off of "fractional charge" and "fractional spin" errors. 45 Hartree-Fock theory is exact for H+ 2 dissociation, while semilocal DFAs improve H2 dissociation due to their crude simulation of left-right correlation. 3,4 Table ACS Paragon 5 illustrates various DFAs’ fractional chargePlus andEnvironment spin errors at dissociation, while Figure

Page 10 of 34

Page 11 of 34

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

2 shows predicted H+ 2 and H2 dissociation curves. Table 5 shows that the Rung 3.5 upper bounds are relatively "Hartree-Fock-like" at dissociation, with small fractional charge errors and larger fractional spin errors. However, Figure 2 shows that the Rung 3.5 upper bounds inherit several aspects of GGAs, including the spurious maximum in H+ 2 dissociation curve and the relatively deep minimum in H2 dissociation. This confirms that the Rung 3.5 upper bound incorporates some characteristic properties of semilocal exchange. Table 5 also confirms the relation of fractional charge and fractional spin errors seen in ref.: 46 DFAs with small fractional charge error tend to have large fractional spin error, and vice versa. Fractional charge and fractional spin errors are connected to the derivative discontinuity of the energy upon variation of electron number Nel . Ref. 46 presented a useful illustration: minimal-basis calculations on eight widely separated H atoms in a cube with N = 2, 4, . . . 16 total electrons. Figure 3 presents our predictions for this system. The full configuration interaction (FCI) solution shows piecewise linear E(N ) with a clear derivative discontinuity at N = 8. In contrast, Hartree-Fock theory and all tested DFAs show significant nonlinearities and no derivative discontinuity. However, a notable result is that the LSDA (denoted as HFS), Xα, ΠPBEup and ΠLDA correctly recover an energy minimum at N = 8, while the rest of approximate functionals (and HF) have their minimum at N = 10.

3.6

Magnetic coupling constants of molecular magnets

Magnetic exchange couplings are an important area of research in the design of molecular magnets and for assessing new methods. Calculations are carried out almost exclusively with DFT due to the size and complexity of the magnetic systems. In most cases their low-spin states have strong multireference character similar to open-shell singlet diradicals. Accurate estimates of the singlet-triplet energy gap and the localized magnetic moments in such systems are challenging. Low-spin sates are usually calculated by enforcing spinsymmetry breaking in order toACS obtain localized magnetic moments on certain ’magnetic’ 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

Page 12 of 34

atoms. The main characteristic here is the magnetic exchange coupling parameter JAB entering the Heisenberg spin-Hamiltonian. 47–49 Following ref. 48 we calculate JAB using the spin projected formula ELS − EHS = 4⟨SA ⟩⟨SB ⟩JAB

(16)

where LS,HS denote low-spin,high-spin state respectively, ⟨SA ⟩ is localized atomic spin on magnetic center A given in atomic units by the calculated atomic spin-population divided by two: 1 ⟨SA ⟩ = qS (A) . 2

(17)

We have calculated seven molecular magnets from the set considered in ref. 48 (Structures (1) to (7) in Table 1 of ref. 48 ). Figure 4 shows the structure and the Hartree-Fock spin density distribution of a representative di-copper species (Table 6 Structure (1)). The geometry of these molecular magnets are from the support information of ref. 48 To obtain localized atomic magnetic moments often an orbital localization is performed. 50 Alternatively, one can still use the canonical molecular orbitals but add to the spin population of a given magnetic center the population on its neighbours if any. 48 We follow the second approach using the ’magnetic neighbors’ for each structure following the support information of ref. 48 The composite atomic spin populations of the magnetic molecules obtained with different methods are given in our SI, Table S3. The way of obtaining broken-symmetry solutions for low-spin states is not unique. We follow here a procedure somewhat different from refs. 48,49 We use the option "Guess=mix" in Gaussian 51 which enforces spin-symmetry breaking, together with a stability analysis searching for the lowest energy broken-symmetry solution. Depending on the procedure used to obtain localized magnetic moments the results for JAB may somewhat vary. The magnetic couplings calculated with broken-symmetry unrestricted HF (UHF) are here relatively closer to the experimental data compared to the DFT estimates (except for Structure 5, where a proper broken-symmetry UHF solution was not found), tending to ACS Paragon Plus Environment

Page 13 of 34

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

underestimate the JAB absolute values for the first four structures (6). The GGA PBE tends to considerably overestimate the JAB absolute values for all structures and gives the magnetic coupling for Structure 6 with opposite sign. Mixing HF exchange as in B3LYP improves the results considerably perhaps due to a partial compensation of opposite trends. It is encouraging that the one-parameter Rung 3.5 hybrid functionals also improve the accuracy of the magnetic coupling constants compared to PBE and and are not too far behind B3LYP. Rung 3.5 functionals yield larger localized magnetic moments than PBE but smaller than those with HF (see SI, Table S3). We tried also DFT calculations of the magnetic couplings in a post-HF fashion. As it is seen from Table 6, there is a systematic difference in the calculated JAB values this way, compared to the SCF DFT results. Post-HF DFT calculations use the same HF HOMOLUMO gap and magnetic moments, which are generally different from the self-consistent DFT values. Finding the lowest-energy broken solution self-consistently with stability analysis is essential in some of these tests, like for Structure (6) where only the SCF energy lowering of the most stable state gives the correct sign of JAB as with B3LYP and UHF. In such cases the post-HF approach is less informative.

3.7

The geometry of 1-propynyl radical

We conclude by considering the upper bound functionals’ performance for an interesting case of density-driven error, 52 the geometry of 1-propynyl radical. Hartree-Fock theory, post-Hartree-Fock ab initio approximations, 53,54 and recent experiments 55 predict a C3v symmetric σ radical, while B3LYP DFT predicts a Cs geometry with significant C-C-C bond bending. 56,57 A systematic study by Oyeyemi and coworkers highlighted the strong DFA dependence of the predicted geometry. 58 We recently showed that even GGAs can predict linear σ radicals, bent σ − π radicals, and nearly linear π radicals with modest JahnTeller distortions, simply by changing the DFA form. 59 Here we extend these findings by presenting the Rung 3.5 upper ACS bound predictions. Figure 5 shows the geometries obtained 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

with ΠLDAup and ΠPBEup. ΠLDAup gives a linear σ radical similar to ΠLDA, while ΠPBEup gives an almost linear π radical structure similar to ΠPBE. 59 This again confirms that trends from the underlying Rung 3.5 DFAs often carry over into the upper bound functionals.

4

Conclusion

The Rung 3.5 functionals offer a new avenue in DFT development giving a good promise that the functional accuracy can be improved considerably without involving exact exchange. This has a practical and also methodological importance since DFT methods that involve exact exchange somewhat deviate from the original idea and goal of the DFT formulation. The avenue of 3.5 upper-bound exchange functionals has the advantage of having the approximations well controlled and gives the assurance that the model would not blow up out of bounds. The self-consistent implementation of the upper-bound exchange models presented here allows assessing these functionals on a broader scale for various properties where selfconsistency is essential (geometry optimization, orbital energies, TDDFT, NMR chemical shifts etc.) We find that the 3.5 run upper bound exchange functionals can be of more general purpose in practical applications, especially when used in global and local hybrid-like combinations with the PBE exchange-correlation. Such combinations closely approach the accuracy of conventional hybrids without using exact exchange. A problem that remains to resolve is the incorrect long range behavior of 3.5 Rung energy densities. We are currently working on implementing long-range correction to some of the 3.5 rung exchange forms.

References ACS Paragon Plus Environment (1) Burke, K. Perspective on density functional theory. J. Chem. Phys. 2012, 136, 150901.

Page 14 of 34

Page 15 of 34

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

(2) Yu, H. S.; Li, S. L.; Truhlar, D. G. Perspective: Kohn-Sham density functional theory descending staircase. J. Chem. Phys. 2016, 145, 130901. (3) Becke, A. D. Simulation of delocalized exchange by local density functionals. J. Chem. Phys. 2000, 112, 4020–4026. (4) Cremer, D. Density functional theory: coverage of dynamic and non-dynamic electron correlation effects. Mol. Phys. 2001, 99, 1899–1940. (5) Becke, A. D. A New Mixing of Hartee-Fock and Local Density-Functional Theories. J. Chem. Phys. 1993, 98, 1372–1377. (6) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648–5652. (7) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab Initio Calculation of Vibrational Absorption and Circular Dichroism Spectra Using Density Functional Force Fields. J. Phys. Chem. 1994, 98, 11623–11627. (8) Janesko, B. G.; Henderson, T. M.; Scuseria, G. E. Perspective: Screened hybrid density functionals for solid-state chemistry and physics. Phys. Chem. Chem. Phys. 2009, 11, 443–454. (9) Paier, J. Hybrid Density Functionals Applied to Complex Solid Catalysts: Successes, Limitations, and Prospects. Catalysis Letters 2016, 146, 861–885. (10) Henderson, T. M.; Janesko, B. G.; Scuseria, G. E. Range Separation and Local Hybridization in Density Functional Theory. J. Phys. Chem. A 2008, 112, 12530–12542. (11) Henderson, T. M.; Janesko, B. G.; Scuseria, G. E. Generalized Gradient Approximation Model Exchange Holes for Range-Separated Hybrids. J. Chem. Phys. 2008, 128, 194105. 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

(12) Schmidt, T.; Kraisler, E.; Makmal, A.; Kronik, L.; Kummel, S. A self-interactionfree local hybrid functional: accurate binding energies vis-a-vis accurate ionization potentials from Kohn-Sham eigenvalues. J. Chem. Phys. 2014, 140, 18A510. (13) de Silva, P.; Corminboeuf, C. Local hybrid functionals with orbital-free mixing functions and balanced elimination of self-interaction error. J. Chem. Phys. 2015, 142, 074112. (14) Becke, A. D. A real-space model of nondynamical correlation. J. Chem. Phys. 2003, 119, 2972–2977. (15) Proynov, E.; Shao, Y.; Kong, J. Efficient self-consistent DFT calculation of nondynamic correlation based on the B05 method. Chem. Phys. Lett. 2010, 493, 381–385. (16) Kong, J.; Proynov, E. Density functional model for nondynamic and strong correlation. J. Chem. Theory. Comput. 2016, 12, 133–143. (17) Tao, J.; Staroverov, V. N.; Scuseria, G. E.; Perdew, J. P. Exact-exchange energy density in the gauge of semilocal density functional approximation. Phys. Rev. A 2008, 77, 012509. (18) Baerends, E. J. Exact Exchange-Correlation Treatment of Dissociated H2 in Density Functional Theory. Phys. Rev. Lett. 2001, 87, 133004. (19) Jaramillo, J.; Scuseria, G. E.; Ernzerhof, M. Local hybrid functionals. J. Chem. Phys. 2003, 118, 1068–1073. (20) Perdew, J. P.; Staroverov, V. N.; Tao, J.; Scuseria, G. E. Density functional with full exact exchange, balanced nonlocality of correlation, and constraint satisfaction. Phys. Rev. A 2008, 78, 052513. (21) Theilacker, K.; Arbuznikov, A. V.; Kaupp, M. Gauge effects in local hybrid functionals evaluated for weak interactions and the GMTKN30 test set. Mol. Phys. 2016, 114, 1118–1127.

ACS Paragon Plus Environment

Page 16 of 34

Page 17 of 34

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

(22) Aguero, A.; Janesko, B. G. Nonspherical model density matrices for Rung 3.5 density functionals. J. Chem. Phys. 2012, 136, 024111. (23) Janesko, B. G.; Scalmani, G.; Frisch, M. J. Practical auxiliary basis implementation of Rung 3.5 functionals. J. Chem. Phys. 2014, 141, 034103. (24) Janesko, B. G. Nonempirical Rung 3.5 density functionals from the Lieb-Oxford bound. J. Chem. Phys. 2012, 137, 224110. (25) Lieb, E. H.; Oxford, S. Improved lower bound on the indirect Coulomb energy. Int. J. Quantum Chem. 1981, 19, 427–439. (26) Odashima, M. M.; Capelle, K. Nonempirical hyper-generalized-gradient functionals constructed from the Lieb-Oxford bound. Phys. Rev. A 2009, 79, 062515. (27) Haunschild, R.; Odashima, M. M.; Scuseria, G. E.; Perdew, J. P.; Capelle, K. Hypergeneralized-gradient functionals constructed from the Lieb-Oxford bound: Implementation via local hybrids and thermochemical assessment. J. Chem. Phys. 2012, 136, 184102. (28) M. Seidl and S. Vuckovic abd P. Gori-Giorgi, Challenging the LiebŰOxford bound in a systematic way. Mol. Phys. 2016, 114, 1076–1085. (29) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865–3868. (30) Janesko, B. G. A simple nonlocal model for exchange. J. Chem. Phys. 2009, 131, 234111. (31) Janesko, B. G.; Scuseria, G. E. Parameterized local hybrid functionals from densitymatrix similarity metrics. J. Chem. Phys. 2008, 128, 084111.

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

(32) Janesko, B. G.; Krukau, A. V.; Scuseria, G. E. Self-consistent generalized Kohn-Sham local hybrid functionals of screened exchange: Combining local and range-separated hybridization. J. Chem. Phys. 2008, 129, 124110. (33) Arbuznikov, A. V.; Bahmann, H.; Kaupp, M. Local Hybrid Functionals with an Explicit Dependence on Spin Polarization. J. Phys. Chem. A 2009, 113, 11898–11906. (34) Arbuznikov, A. V.; Kaupp, M. Advances in local hybrid exchange-correlation functionals: from thermochemistry to magnetic-resonance parameters and hyperpolarizabilities. Int. J. Quant. Chem. 2011, 111, 2625–2638. (35) Arbuznikov, A. V.; Kaupp, M. Towards improved local hybrid functionals by calibration of exchange-energy densities. J. Chem. Phys. 2014, 141, 204101. (36) Janesko, B. G. Rung 3.5 density functionals. J. Chem. Phys. 2010, 133, 104103. (37) Grabo, T.; Gross, E. Density-functional theory using an optimized exchange-correlation potential. Chem. Phys. Lett 1995, 240, 141–150. (38) Lynch, B. J.; Truhlar, D. G. Small Representative Benchmarks for Thermochemical Calculations. J. Phys. Chem. A 2003, 107, 8996, 108, 1460(E) (2004). (39) DeYonker, N. J.; Wiliams, T. G.; Imel, A. E.; Cundari, T. R.; Wilson, A. K. Accurate thermochemistry for transition metal complexes from first-principles calculations . J. Chem. Phys. 2009, 131, 024106. (40) Kepp, K. P. Consistent descriptions of metal–ligand bonds and spin–crossover in inorganic chemistry. Coord. Chem. Reviews 2013, 257, 196–209. (41) Baerends, E. J.; Gritsenko, O. V.; van Meer, R. The Kohn–Sham gap, the fundamental gap and the optical gap: the physical meaning of occupied and virtual KohnŰSham orbital energies. Phys. Chem. Chem. Phys. 2013, 15, 16408–16425. ACS Paragon Plus Environment

Page 18 of 34

Page 19 of 34

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

(42) Baerends, E. J.; Gritsenko, O. V.; van Meer, R. On the errors of local density (LDA) and generalized gradient (GGA) approximations to the Kohn-Sham potential and orbital energies. J. Chem. Phys. 2016, 144, 204114. (43) Almbladh, C.-O.; Pedroza, A. C. Density-functional exchange-correlation potentials and orbital eigenvalues for light atoms. Phys. Rev. A 1984, 29, 2322–2330. (44) Brosco, V.; Lorenzana, Y. J. Exact echange-correlation potential of an ionic Hubbard model with a free surface. Nature: Sci. Reps. (Open) 2013, 3, 2172. (45) Cohen, A. J.; Mori-Sanchez, P.; Yang, W. Insights into Current Limitations of Density Functional Theory. Science 2008, 321, 792–794. (46) Mori-Sánchez, P.; Cohen, A. J. Perspective: The derivative discontinuity of the exchangeŰcorrelation functional. Phys. Chem. Chem. Phys. 2014, 16, 14378–14387. (47) Noodleman, L.; D., P.; Baerends, E. J. Symmetry breaking and ionization from symmetry equivalent inner shells and lone pairs in Xa theory. Chem. Phys. 1982, 64, 159–166. (48) Peralta, J. E.; Melo, J. I. Magnetic exchange couplings with range-separated hybrid density functionals. J. Chem. Theory Comput. 2010, 6, 1894–1899. (49) Philips, J. J.; Peralta, J. E.; Janesko, B. J. Magnetic exchange couplings evaluated with Rung 3.5 density functionals. J. Chem. Phys. 2011, 134, 21401. (50) Ruiz, E.; Cano, J.; Alvarez, S.; Alemany, P. Broken symmetry approach to calculaiton of exchange coupling constants for homobinuclear and heteronuclear transition metal complexes. J. Comput. Chem. 1999, 20, 1391–1400. (51) Gaussian Development Version, Revision I.02+, M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, J. Bloino, ACS Paragon Plus Environment B. G. Janesko, A. F. Izmaylov, A. Marenich, F. Lipparini, G. Zheng, J. L. Sonnenberg,

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

W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, T. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, P. V. Parandekar, N. J. Mayhall, A. D. Daniels, O. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski, and D. J. Fox, Gaussian, Inc., Wallingford CT, 2014. (52) Kim, M.-C.; Sim, E.; Burke, K. Understanding and Reducing Errors in Density Functional Calculations. Phys. Rev. Lett. 2013, 111, 073003. (53) Mebel, A. M.; Jackson, W. M.; Chang, A. H. H.; Lin, S. H. Photodissociation Dynamics of Propyne and Allene: A View from ab Initio Calculations of the C3 Hn (n = 1 − 4) Species and the Isomerization Mechanism for C3 H2 . J. Am. Chem. Soc. 1998, 120, 5751–5763. (54) Wheeler, S. E.; Robertson, K. A.; Allen, W. D.; Schaefer III, H. F.; Bomble, Y. J.; Stanton, J. F. Thermochemistry of Key Soot Formation Intermediates: C3 H3 Isomers. J. Phys. Chem. A 2007, 111, 3819–3830. (55) Zhou, J.; Garand, E.; Eisfeld, W.; Neumark, D. M. Slow electron velocity-map imaging spectroscopy of the 1-propynyl radical. J. Chem. Phys. 2007, 127, 034304–. (56) Vereecken, L.; Pierloot, K.; Peeters, J. B3LYP-DFT characterization of the potential energy surface of the CH(X 2 Π) + C2 H2 Reaction. J. Phys. Chem. 1997, 108, 1068– 1080. ACS Paragon Plus Environment (57) Sreeruttun, R. K.; Ramasami, P.; Wannere, C. S.; Paul, A.; v. R. Schleyer, P.; Schae-

Page 20 of 34

Page 21 of 34

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

fer III, H. F. Effects of Fluorine on the Structures and Energetics of the Propynyl and Propargyl Radicals and Their Anions. J. Org. Chem. 2005, 70, 8676–8686. (58) Oyeyemi, V. B.; Pavone, M.; Carter, E. A. Accurate Bond Energies of Hydrocarbons from Complete Basis Set Extrapolated Multi-Reference Singles and Doubles Configuration Interaction. Chem. Phys. Phys. Chem. 2011, 12, 3354–3364. (59) Janesko, B. G.; Proynov, E. Accurate Alkynyl Radical Structures from Density Functional Calculations Without Hartree-Fock Exchange. J. Chem. Phys. 2017, 146, 054109.

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

Page 22 of 34

Table 1: Atomic exchange energy values (Hartree) with different functionals,self-consistent cc-pVQZ calculations.

Atom

-Ex (HF)

-ΠLDA

-ΠLDAup

-ΠPBEup

-ΠPBE

-ΠGG2-A9

-ΠGG2Fx

H He Li Be B C N O F Ne Na Mg Al Si P S Cl Ar

0.312 1.025 1.781 2.676 3.769 5.076 6.607 8.218 10.045 12.110 14.017 15.994 18.091 20.304 22.642 25.034 27.544 30.185

0.253 0.832 1.450 2.190 3.100 4.230 5.605 6.989 8.608 10.503 12.205 13.985 15.853 17.843 19.974 22.114 24.390 26.819

0.212 0.697 1.213 1.830 2.582 3.530 4.711 5.853 7.213 8.841 10.290 11.805 13.390 15.085 16.913 18.725 20.663 22.749

0.265 0.865 1.508 2.274 3.174 4.255 5.556 6.851 8.342 10.082 11.706 13.406 15.174 17.045 19.046 21.046 23.163 25.424

0.303 0.996 1.733 2.612 3.664 4.926 6.423 7.962 9.719 11.735 13.605 15.560 17.603 19.762 22.059 24.381 26.830 29.427

0.264 0.881 1.543 2.330 3.268 4.410 5.793 7.170 8.768 10.637 12.379 14.158 16.019 18.000 20.115 22.229 24.474 26.874

0.319 1.047 1.810 2.714 3.771 5.020 6.507 8.012 9.726 11.705 13.546 15.476 17.480 19.600 21.842 24.099 26.477 29.006

1.583

3.285

2.014

0.341

1.088

0.151

MUE (DFT-HF)

ACS Paragon Plus Environment

Page 23 of 34

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

Table 2: Mean unsigned error (MUE, kcal/mol) in AE6 atomization energies and BH6 reation barriers, 38 6-311++G(2d,2p) basis.

Method

AE6

BH6

αa

βb

HFPBE ΠPBEup ΠPBE

29.77 5.32 18.7 11.34 14.42 8.03

1.0 1.0 1.0

0.0 0.0 0.0

PBE ΠPBEup1 ΠPBE1

12.9 7.78 6.96

9.3 9.08 8.26

0.0 0.25 0.45

0.0 0.0 0.0

6.10 4.84 4.81

7.91 6.37 4.61

0.12 0.65 0.25

0.45 2.25 0.0

ΠPBEup1w ΠPBE1w c PBE0

c

a

b

Values of the α mixing parameter in Eq.11 Values of the β exponential factor in the local mixing function, Eq.11 c Post-PBE0 calculations

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

Table 3: Mean error (ME), mean unsigned error (MUE) and RMSD of enthalpy of formation (kcal/mol) on TM193 thermochemistry test 39 with different DFAs, cc-PVQZ basis

a

Method

MAE

ME

RMSD

PBE ΠPBE a ΠPBEup a ΠPBE1 a ΠPBEup1 a PBE0 B3LYP B3LYP a

29.0 26.0 40.2 22.6 25.9 13.6 12.5 13.2

28.4 5.5 17.9 19.0 24.4 -4.5 -7.6 -11.2

41.5 35.3 58.2 31.2 35.3 21.9 18.9 20.1

Without performing stability analysis for the molecules

ACS Paragon Plus Environment

Page 24 of 34

Page 25 of 34

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

Table 4: Mean unsigned errors of calculated HOMO energy and HOMO-LUMO gap (Hartree,6-311++G(2d,2p) basis). The sample test set is listed in Table S1 of SI.

system MUE1 MUE2 MUE3

a b c

HFa

ΠPBEup

ΠPBE

ΠPBEup1

0.055

0.100 0.137 0.035

0.090 0.137 0.037

0.134 0.170 0.053

ΠPBE1

PBE

-IP(expt)

0.119 0.096 0.142 0.158 0.133 0.178 0.044 0.00 0.072

a b

PBE0

Mean unsigned error for HOMO with respect to the experimental Ip Mean unsigned error for HOMO with respect to the HF HOMO values c MUE for HOMO-LUMO gap with respect to the PBE0 gap

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

Table 5: Fractional charge error E(H) − E(H2+ ) and fractional spin error E(H2 ) − 2E(H) (mH) of H+ 2 and singlet H2 at dissociation (100 Angstrom). Self-consistent calculations, aug-cc-PVQZ basis

Method PBE BLYP TPSS HF PBE0 M06L B3LYP ΠLDA ΠPBEup ΠPBE ΠPBE1 ΠPBEup1

Charge

Spin

105 109 101 0 78 99 85 50 27 39 75 83

82 71 83 282 130 101 108 160 202 194 133 117

ACS Paragon Plus Environment

Page 26 of 34

Page 27 of 34

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

Table 6: Magnetic exchange couplings of Peralta’s set of magnetic molecules (cm−1 ) with different DFAs using the basis set of ref. 48 Post-HF results are in parentheses a

System (1)CuII -CuII

B3LYP

PBE

ΠPBE1

ΠPBEup1

UHF

-100.8 (-13.5)

-591.7 (-14.0)

-193.0 (14.7)

-244.9 (-15.0)

-6.8

Expt -30.9

(2)CuII -CuII

-108.1 (-13.1)

-541.4 (-15.3)

-219.4 (-14.8)

-251.9 (-15.0)

-4.9

-37.4

(3)MnII -CuII

-38.3 (-5.4)

-173.0 (-6.1)

-85.3 (5.8)

-95.7 (-6.0)

-2.3

-15.7

(4)VIV -VIV

-100.2 (33.2)

-279.1 (+44.7)

-176.0 (+47.1)

-202.9 (+47.9)

-1.5

-107.0

(5)MnIII -MnIV

-210.0

-596.7

-215.9

-245.4

-

b

-110.0

(6)CuII -CuII

161.2 (-439.1)

-197.4 (-557.8)

-22.5 (-489.7)

-36.1 (-500.8)

+134.1

+84.0

(7)CuII -CuII

+202.3 +326.5 (+185.7) (+202.1) 70.4 347.4

+244.0 (+198.4) 126.0

+265.8 (+204.5) 152.7

+71.3

+57.0

MUE

From reference. 48 We were not able to find proper broken-symmetry UHF singlet solution for Structure 5. a

b

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

PiPBE−AE6 PiPBEUP−AE6 PiPBE−BH6 PiPBEUP−BH6

17.95 16.95 15.95 14.95 MUE (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 28 of 34

13.95 12.95 11.95 10.95 9.95 8.95 7.95 6.95

0

0.1

0.2

0.3

0.4

0.5

α

0.6

0.7

0.8

0.9

1

Figure 1: Mean unsigned error on the AE6-BH6 set as a function of α. SCF results, 6311++G(2d,2p) basis.

ACS Paragon Plus Environment

Page 29 of 34

PBE PiLDAup PiLDA PiPBEup PiPBE GG2UP HF

0.22

−0.015 De (a.u)

0.12

De (a.u)

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

−0.045

−0.075

FCI HF PBE PiPBEPBE PiPBE1PBE PiPBEupPBE PiPBEup1PBE PBE0

0.02

−0.08

−0.105 1.5

2.5

3.5

4.5

5.5

6.5

7.5

8.5

9.5 10.5 11.5

0.6

1.6

2.6

3.6

4.6

5.6

6.6

7.6

R(H−H) (Å)

R(H−H) (Å)

Figure 2: Dissociation curves of doublet H+ 2 (top) and singlet symmetric H2 (bottom), 6311++G(3df,3pd) basis.

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

0

HF FCI PiLDA PiPBEUP HFS PBE0 Xalpha GG2-Fx PiPBE

-0.5 -1 -1.5 E (a.u)

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 30 of 34

-2 -2.5 -3 -3.5 0

2

4

6

8

10

12

14

16

Nel

Figure 3: Energy of H8 cube as a function of the number of electrons. Full CI results are from ref. 46

ACS Paragon Plus Environment

Page 31 of 34

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

Figure 4: Structure and broken-symmetry UHF singlet spin density of magnetic molecule Structure (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

Figure 5: ΠLDAupPBE (top) and ΠPBEuPBE (bottom) structure and spin density of 1propyne radical.

ACS Paragon Plus Environment

Page 32 of 34

Page 33 of 34

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

ΠLDAupPBE gives the correct structure and spin density of 1-propyne radical. For Table of Contents Use Only

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

Graphical TOC Entry

ACS Paragon Plus Environment

Page 34 of 34