Shedding Light on the Accuracy of Optimally Tuned Range-Separated

May 17, 2017 - Citation data is made available by participants in Crossref's Cited-by ... with multi-reference character from optimally tuned range-se...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/JPCA

Shedding Light on the Accuracy of Optimally Tuned RangeSeparated Approximations for Evaluating Oxidation Potentials Mojtaba Alipour* and Soheila Mohseni Department of Chemistry, College of Sciences, Shiraz University, Shiraz, Iran ABSTRACT: There is a surge in the literature on the development of exchange− correlation density functionals for a wide variety of physical and chemical properties. As a recent endeavor toward the systematic and nonempirical design of density functional approximations, optimally tuned range-separated hybrid (OT-RSH) models have been introduced. In this work, we propose novel OT-RSH density functionals for predicting the oxidation potentials of organic compounds from different categories. In this regard, detailed analysis of the role of nonempirical optimization of the range separation parameter and importance of short- and longrange exact-like exchange in OT-RSH calculations of the oxidation potential has also been done. It is shown that the newly developed OT-RSH approximations not only perform better than other standard long-range corrected functionals but also in many cases outperform other conventional hybrid functionals with a fixed amount of exactlike exchange. Plus, we find that the proposed functionals describe well the oxidation potentials of compounds for which the tuning of the range separation parameter was not performed. From a different perspective, accountability of the computed frontier orbital energies from the OT-RSH density functionals for estimation of oxidation potentials has also been evaluated. Our results reveal that the negative highest occupied molecular orbital energies of molecules and the negative lowest unoccupied molecular orbital energies of their cations correlate remarkably with the observed oxidation potentials. Admittedly, with more efforts along this line, modern OT-RSH functionals with broader applicability can be released for computational electrochemistry.

1. INTRODUCTION Nowadays, density functional theory (DFT)1−3 has become a workhorse of quantum chemistry to bypass solving the Schrödinger equation. From a theoretical viewpoint, DFTbased calculations provide important quantities among which, and of particular relevance to the present work, molecular orbital (MO) eigenvalues are of interest.4−8 Introducing occupation number ni for each of the electronic states ψi and defining an electron density ρ(r) by the relation ρ(r) = ∑i ni |ψi(r)|2 , Janak proved that variation of the total energy E with respect to the orbital occupation ni is equal to the state eigenvalue εi9 ∂E = εi(ni) ∂ni

orbital ni is equal to 0, and this orbital turned out to be the lowest unoccupied molecular orbital (LUMO) of the cation M+. Therefore, using the approximate treatments on the right side of eq 2, this equation can be simplified as follows7

(1)

∫0

1

εi(ni) dni

(2)

From another perspective, if ni is equal to 1, the ith orbital, namely, the space part of the ith spin−orbit, may be considered as the highest occupied molecular orbital (HOMO) of the neutral molecule M, while when removing an electron from the ith © XXXX American Chemical Society

(3)

E(M+) − E(M) ≈ −LUMO(M+)

(4)

In the gas phase, the quantity on the left side of these equations is the ionization potential (IP) of the molecule M. Accordingly, the IP of a molecule can be estimated by either the negative HOMO eigenvalue of molecule M or the negative LUMO eigenvalue of its cation. Consider now a solution as the medium of concern. On the one hand, the left side of the above equations can be considered as the required energy for removing an electron from the molecule in solution, which in turn is a measure of the oxidation potential. On the other hand, on the basis of eqs 3 and 4, we expect that the frontier orbital energies in solution can also correlate to the oxidation potentials. Put together, we anticipate that besides directly using the E(M+) − E(M) in solution the frontier orbital energies obtained from DFT calculations can also be employed as measures of experimental oxidation potentials. Hereby, analysis of this issue and proposing novel approximations without invoking any empirical fitting

The ith spin−orbit is occupied if ni = 1, and it is unoccupied if ni = 0, and consequently, this equation can relate the ground-state energy of an N-electron molecule M to the ground-state energy of its (N−1)-electron cation, M+. Considering E(M) as the ground-state energy of molecule M for ni = 1 and E(M+) as the ground-state energy of M+ for ni = 0, we have E(M+) − E(M) = −

E(M+) − E(M) ≈ −HOMO(M)

Received: April 23, 2017 Revised: May 7, 2017

A

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

1 1 1 = {1 − [α + β erf(μr12)]} + [α + β erf(μr12)] r12 r12 r12

within a class of DFT in which the frontier orbital energies come into play, namely, the optimally tuned range-separated hybrid (OT-RSH) functionals,10−13 for predicting oxidation potentials constitute the subject of this contribution. To the best of our knowledge, no report in this regard has been appeared in the literature as of yet. The present paper has been organized in the following manner. In the next section, the theoretical foundation of the OT-RSH density functionals under study is overviewed. Then, technical details of our calculations are explained. Thereafter, the results of OT-RSHs on oxidation potentials, their discussion, and comparisons with other functionals from various rungs are provided. Last, we finalize the paper by some concluding remarks.

(6)

where α and α + β define the exact exchange percentage at r12 = 0 and ∞, respectively, and the inequalities 0 ≤ α + β ≤ 1, 0 ≤ α ≤ 1, and 0 ≤ β ≤ 1 should also be satisfied. In this case, the LR and SR contributions are expressed respectively as occupied occupied

ExLR = αExHF − β

ExSR = −

occupied occupied





i

j

ϕϕ i j

j

erf(μr12) ϕϕ j i r12

1 2

⎧ ⎪

∑ ∫ ρσ 4/3kσ ⎨(1 − α) − ⎪

σ



8 βaσ 3

There has been tremendous interest in the literature to apply the DFT approximations based on RS and CAM approaches for different properties.18−27 In more detail, some remarkable points are of concern. Setting α = 0 and consequently β = 1 in eq 6, eq 5 is reproduced. However, the presence of the two extra parameters α and β in eq 6 provides more flexibility for analyzing the effects of SR and LR exchange on properties under study. In earlier investigations, it has also been shown that the SR HF exchange (α ≠ 0) can lead to improved results.28−33 On the other hand, it should be noted that in some cases the resulting functionals based on eq 6 do not incorporate a full range separation in their composition, while there are also other density functionals including a full 100% contribution of asymptotic HF exchange (α + β = 1). As a matter of fact, α and β dictate the limiting behavior of the HF exchange, which in turn tends to α/r for r → 0 and to (α + β)/r for r → ∞. Accordingly, the correct LR asymptotic behavior of the exchange−correlation, 1/r, is achieved simply by enforcing α + β = 1.28,34,35 2.2. Optimally Tuning Scheme in Range-Separated Hybrid Functionals. In this approach, parameter tuning is performed without invoking any empirical fitting process. For each range separation parameter μ, the IP and electron affinity (EA) can be calculated as the differences between ground-state energies of systems with N and N ± 1 electrons

(5)

erf(μr12) ϕϕ j i r12

IP μ(N ) = Eμ(N − 1) − Eμ(N )

(7)

EAμ(N ) = Eμ(N ) − Eμ(N + 1)

and ExSR

i

ϕϕ i j

⎡ ⎤⎫ ⎛ 1 ⎞ ⎪ ⎢π 1/2 erf⎜ dr ⎟ + 2aσ (bσ − cσ )⎥⎬ ⎪ ⎢⎣ ⎝ 2aσ ⎠ ⎦⎥⎭

where the erf term stands for the standard error function and μ is the range separation parameter. The first term describes the SR effect and is modeled through DFT exchange, whereas the second term accounts for the LR contribution, leading to longrange-corrected (LC) approximations. The second term is of particular importance because it enforces a rigorously correct 100% contribution of asymptotic HF exchange, which in turn is essential for accurately describing many properties. Accordingly, the total exchange term is decomposed into the SR LR LR and SR contributions, Ex = ELR x + Ex , where the terms Ex and ESR are expressed respectively as x ExLR = −



where EHF x is the HF exchange and

2. THEORETICAL FRAMEWORK As noted earlier, the class of DFT under consideration is the OTRSH functionals. The related theoretical backgrounds that will be employed for interpretation of our results come in the following. Further details can be found in the original references cited below. 2.1. Range-Separated Hybrid Functionals. The RS approach in DFT indicates the correction of exchange functionals for long-range (LR) electron−electron exchange interactions, which are insufficiently incorporated in conventional exchange functionals. This approach mixes short-range (SR) density functional exchange with LR Hartree−Fock (HF) exchange by partitioning the electron repulsion operator into SR and LR terms as follows14−16 1 1 1 = {1 − [erf(μr12)]} + [erf(μr12)] r12 r12 r12



= −∑ σ

∫ ρσ

4/3

⎧ ⎪ 8 ⎡ 1/2 ⎛ 1 ⎞ − kσ ⎨ 1 aσ ⎢π erf⎜ ⎟ ⎪ 3 ⎢⎣ ⎝ 2aσ ⎠ ⎩

= IP μ(N + 1)

(8)

For the exact DFT, the IP of a neutral molecule is identical to the minus HOMO energy of the neutral molecule, εHOMO(N), and the EA of a neutral molecule is identical to the minus HOMO energy of the anion, εHOMO(N + 1). Moreover, by definition, EA(N) = IP(N + 1). However, due to the approximate nature of DFT functionals, these identities are not held and there are differences between IP/EA and the corresponding frontier orbital energies. The remaining deviations from the exact state for each case can be considered separately as the following functions

⎤⎫ ⎪ + 2aσ (bσ − cσ )⎥⎬ dr ⎥⎦⎪ ⎭

Here, ϕ is the MO, kσ is the enhancement factor, and ρσ represents the density of the electron with spin σ at point r. The explicit relations of aσ, bσ, and cσ can also be found in refs 14−16. Aiming to generalization of eq 5, Yanai et al. used the two extra parameters, α and β, and proposed their approach called the Coulomb Attenuating Method (CAM)17

μ J0 (μ) = εHOMO (N ) + IP μ(N )

B

(9) DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 1. Geometrical structures and their names for the molecules in the first benchmarked set. μ J1(μ) = εHOMO (N + 1) + IP μ(N + 1)

in sign to the IP and EA, respectively. However, the fundamental gap as the difference between IP and EA differs notably from the orbital energy difference between HOMO and LUMO energies. Therefore, to have a proper description of the fundamental gap,

(10)

Note that, as expected from the exact functional, setting the defined functions in eqs 9 and 10 equal to zero implies that the HOMO energies of the neutral and anion are equal and opposite C

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 2. continued

D

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 2. Tuning curves of the range separation parameter for (a) 2,5-dimethyl-1,4-dimethoxybenzene, (b) 4-fluoro-1,2-dimethoxybenzene, (c) 2bromo-1,4-dimethoxybenzene, and (d) 4-t-butylpyridine-N-oxide using various density functional approximations (BLYP, PBE, and TPSS). For each molecule, the figure on the left side corresponds to the case of α = 0.2 and β = 0.8 and the figure on the right side corresponds to the case of α = 0.0 and β = 1.0.

where the three terms Egas, δGsolvation, and δGthermal denote, respectively, the total electronic energy in the gas phase, the solvation contribution to the Gibbs free energy, and the thermal contribution to the Gibbs free energy. The sum of the first two terms, Egas and δGsolvation, can be considered as the total electronic energy in solution, Esolution. Moreover, it has earlier been shown that the thermal contributions to the oxidation potentials are small, ranging around 0.04 V in comparison to 3−5 V, and therefore, neglecting thermal contributions in the calculations of oxidation potentials does not lead to significant deviations.37 Therefore, the oxidation potential of a molecule M can be calculated by the equation below

the functions defined in eqs 9 and 10 for IP and EA should be minimized simultaneously, leading to minimization of a general function J(μ) with the form below J(μ) = J0 (μ) + J1(μ)

(11)

Although this function can be considered for minimization of the range separation parameter, it has been shown that, depending on the particular form of the dependencies of εHOMO and the total energy of the range separation parameter, the values of μ are very close to minima of functions defined in eqs 9 and 10.13 To avoid such a biased behavior, minimization of an objective function, J2(μ), with the following form has been proposed and used in our calculations

U (M) =

μ J 2 (μ) = [εHOMO (N ) + IP μ(N )]2 μ + [εHOMO (N + 1) + IP μ(N + 1)]2

[G(M+) − G(M)] − 1.46 e

(13)

The Gibbs free energy in solution can approximately be expressed as follows G = Egas + δGsolvation + δGthermal

(15)

As systems under study to compute their oxidation potentials, the sets of organic compounds from different categories including aromatic molecules, pyridine-N-oxide-like molecules, and 10-methylphenothiazine(MPT)-like molecules for which experimental data were available7,37 have been considered. Shown in Figure 1 are the geometrical structures of the corresponding molecules along with their names. To prevent problems arising from different choices of geometries on functional benchmarking, the optimized geometries by the widely used B3LYP38,39 functional and the standard Pople’s basis set 6-311++G(d,p) were used in all calculations. Next, all of the subsequent calculations on range separation parameter tuning, oxidation potentials, and other related quantities have been carried out using various functionals with the same basis set, as was advocated as a good enough choice in earlier efforts.7 For computations in solution, we employed the Polarizable Continuum Model (PCM)40 as a solvation paradigm whose reasonable results for describing solvation effects on redox potentials in electrochemical reactions have been confirmed.36,41 To be in accordance as much as possible with the electrolyte environment in experiment, the two key parameters, the dielectric constant and solvent radius, have been adjusted to describe the used solvent. Following previous analyses in this context, we considered an aqueous solution with the solvent

(12)

3. COMPUTATIONAL DETAILS The oxidation potential U (in V) of a neutral molecule M corresponds to the negative change in Gibbs free energies G (in eV) between the molecule M and its cation M+ divided by the electron charge e, [G(M+) − G(M)]/e. In order to obtain the values of oxidation potentials relative to the Li/Li+ reference electrode, one subtracts from this value the negative change in Gibbs free energies between a lithium atom in the bulk metal phase and a lithium ion in electrolyte solution divided by the electron charge e, 1.46 V.36 Accordingly, we arrive at the following relation for calculating the oxidation potential of molecule M, U(M) U (M) =

[Esolution(M+) − Esolution(M)] − 1.46 e

(14) E

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

obtained using the tuning processes of the range separation parameter for the other organic compounds considered. Having determined the optimal values for parameter μ, we are now in a position to use them for subsequent calculations of oxidation potentials. Table 1 lists the corresponding statistical

radius equal to 5.0 Å and dielectric constant equal to 60.0 in our calculations.37 However, the solvent effects on the computed oxidation potentials have also been taken into account by testing other solvents such as tetrahydrofuran and cyclohexane. It should also be taken into account that the use of solvation models such as the PCM imposes, in turn, other approximations on the results of OT-RSH functionals, and Janak theorem (frontier orbital energies). Tuning processes of the range separation parameter μ based on minimization of eq 12 have been performed using two widely used generalized gradient approximations (BLYP39,42 and PBE43) and one meta-generalized gradient approximation (TPSS44) as underlying density functional approximations. In the parameter tuning process, we adopted two options as follows: (i) without any SR exchange (α = 0) and (ii) with some amounts of exchange over the entire range (α ≠ 0). In the latter case, we mention that, although the value of α can be chosen from firstprinciples arguments,34,35 the value of 0.2 has been found to be satisfactory in previous studies,28−33,45−47 and consequently, this value is used herein. Last, it should be noted that, although the tuning processes include different exchange contributions at SR (α = 0 and α = 0.2), both of them recover the full 100% exchange at asymptotic distance (α + β = 1). All of the OT-RSH density functional implementations and standard DFT computations were carried out using the Gaussian09 suite of programs.48 As metrics for comparisons between experimental and calculated data, we used usual statistical measures as the mean signed deviation (MSD), mean absolute deviation (MAD), maximum absolute deviation (MaxAD), and root-mean-square deviation (RMSD). With Δi as the deviation of computed data from the corresponding reference values and n as given data, these measures are 1 n 1 n respectively defined as MSD = n ∑i Δi , MAD = n ∑i |Δi |, MaxAD = max|Δi|, and RMSD =

1/2

( 1n ∑in=1 Δi2)

Table 1. Statistical Measures of the Performance of the Proposed OT-RSH Functionals for Oxidation Potentials of the Molecules in the First Benchmarked Seta OT-RSH functionals

MSD

MAD

α = 0, β = 1, and μ = 0.25 Bohr−1 LC-BLYP 0.47 0.47 LC-PBE 0.62 0.62 LC-TPSS 0.58 0.58 α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 LC-BLYP 0.44 0.44 LC-PBE 0.60 0.60 LC-TPSS 0.56 0.56 a

MaxAD

RMSD

0.56 0.75 0.70

0.48 0.63 0.59

0.54 0.71 0.69

0.45 0.60 0.57

All values are in V.

measures on the performance of OT-RSH functionals for the calculations of oxidation potentials of all of the considered systems. Moreover, the linear correlations between the experimental oxidation potentials and computed values using these OT-RSH functionals are shown in Figure 3. According to our numerical data and correlation plots, the following ordering from lower to higher deviations is found for both cases α = 0.0, β = 1.0, and μ = 0.25 Bohr−1 and α = 0.2, β = 0.8, and μ = 0.20 Bohr−1: LC-BLYP > LC-TPSS > LC-PBE. Altogether, the new OT-RSH functionals based on BLYP density functional approximation with smaller statistical measures and larger regression coefficients with respect to others are shown to be superior for overall performance and will be employed in subsequent comparisons. Concerning the role of continuum solvent on the computed oxidation potentials, we have also performed the calculations in tetrahydrofuran and cyclohexane solvents. The corresponding results on the performance of the OT-RSHs are collected in Table 2. We find that the deviations of the computed oxidation potentials with respect to experiment in these solvents are larger than that provided using the solvent with modified parameters. Moreover, similar classifications on the performance of OT-RSH functionals are obtained by employing all of the solvents. It should be noted however that there are discrepancies between the experimental and computed data, which in turn can be attributed to various factors. For instance, among other theoretical problems worthy of investigation (methods, basis sets, and so on), the different conditions of experimental measurements and theoretical predictions, the employed solvation models, and so forth can be mentioned. Up to now, comparisons have been made among the developed OT-RSH functionals. However, from a practical standpoint, it is also important to further assess the applicability of these models with respect to other approximations from various rungs of Jacob’s ladder.49 To this end, we have included different density functionals with both interelectronic distancedependent HF exchange and a fixed amount of HF exchange that is constant in space. The considered functionals are the standard LC approximations LC-BLYP, LC-PBE, and LC-TPSS with α = 0.0, β = 1.0, and μ = 0.47 Bohr−1, the RS functionals CAMB3LYP (α = 0.19, β = 0.46, and μ = 0.33 Bohr−1),17 ωB97 (α = 0 and μ = 0.4 Bohr−1) and ωB97X (α ≈ 0.16 and μ = 0.3 Bohr−1),20 and the hybrid functionals M05-2X,50 M06-2X,51 M06-HF,52

.

4. RESULTS AND DISCUSSION We have performed the tuning process of the range separation parameter μ based on minimization of eq 12 for all of the molecules under study. As illustrative examples, shown in Figure 2 are the corresponding plots of J2 as a function of μ for four molecules from different categories, 2,5-dimethyl-1,4-dimethoxybenzene, 4-fluoro-1,2-dimethoxybenzene, 2-bromo-1,4-dimethoxybenzene, and 4-t-butylpyridine-N-oxide, using BLYP, PBE, and TPSS as density functional approximations with the two options α = 0.0, β = 1.0 and α = 0.2, β = 0.8. Perusing the curves we find that for all of the systems the optimally tuned μ values nearly appeared at 0.25 Bohr−1 for α = 0.0, β = 1.0, while in the case of α = 0.2, β = 0.8 the minima of curves are located roughly at 0.20 Bohr−1. Furthermore, it can be seen from the figures that the same optimally tuned values for parameter μ are obtained for all of the tested approximations irrespective the type of exchange and correlation functionals used. This may indicate more dependency of the range separation parameter on the nature of the systems. However, the obtained optimal parameters are expected to be nearly the same for similar compounds in each set.29−33 On the other hand, we can also observe that the optimally tuned values of range separation parameters decreased when switching from the case α = 0.0, β = 1.0 to α = 0.2, β = 0.8, highlighting the role of the contribution of HF exchange (about 20%) over the entire range. Almost similar findings were F

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table 2. Comparisons of the Performance of the Proposed OT-RSH Functionals for Oxidation Potentials in Different Solventsa OT-RSH functionals

MSD

MAD

MaxAD

RMSD

Tetrahydrofuran α = 0, β = 1, and μ = 0.25 Bohr−1 LC-BLYP 0.67 0.67 0.78 LC-PBE 0.81 0.81 0.94 LC-TPSS 0.78 0.78 0.91 α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 LC-BLYP 0.64 0.64 0.76 LC-PBE 0.79 0.79 0.93 LC-TPSS 0.76 0.76 0.90 Cyclohexane α = 0, β = 1, and μ = 0.25 Bohr−1 LC-BLYP 1.29 1.29 1.45 LC-PBE 1.42 1.42 1.60 LC-TPSS 1.39 1.39 1.58 α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 LC-BLYP 1.26 1.26 1.43 LC-PBE 1.40 1.40 1.59 LC-TPSS 1.37 1.37 1.57 a

0.67 0.82 0.94 0.65 0.79 0.76

1.07 1.43 1.39 1.26 1.40 1.37

All values are in V.

Table 3. Comparisons of the Performance of the Best OTRSHs in This Work with Respect to Other DFT Approximations for the Oxidation Potentials of the Molecules in the First Benchmarked Setsa

Figure 3. Least-square fittings of the correlations between experimental oxidation potentials and computed data using OT-RSH functionals with various density functional approximations (BLYP, PBE, and TPSS) for the cases (a) α = 0.0 and β = 1.0 and (b) α = 0.2 and β = 0.8.

functionals

MSD

MAD

MaxAD

RMSD

LC-BLYPb LC-BLYPc LC-BLYPd LC-PBEd LC-TPSSd CAM-B3LYPe ωB97 ωB97X M05-2X M06-2X M06-HF BMK B3LYP B3P86

0.44 0.47 0.68 0.84 0.81 0.51 0.54 0.54 0.69 0.68 1.03 0.50 0.38 1.00

0.44 0.47 0.68 0.84 0.81 0.51 0.54 0.54 0.69 0.68 1.03 0.50 0.38 1.00

0.54 0.56 0.87 1.04 1.01 0.61 0.75 0.71 0.82 0.79 1.27 0.60 0.45 1.07

0.45 0.48 0.70 0.85 0.82 0.52 0.56 0.55 0.70 0.68 1.05 0.51 0.38 1.00

All values are in V. bα = 0.2, β = 0.8, and μ = 0.20 Bohr−1. cα = 0.0, β = 1.0, and μ = 0.25 Bohr−1. dα = 0.0, β = 1.0, and μ = 0.47 Bohr−1. eα = 0.19, β = 0.46, and μ = 0.33 Bohr−1.

a

BMK,53 B3LYP,38,39 and B3P86.38,54 Regarding the reasons for the selection of these approximations for comparison, we note that from the viewpoint of HF exchange amount at SR the considered functionals cover both categories of methods with α = 0 and α ≠ 0. From a different perspective, it should be mentioned that, unlike the standard LC functionals, there are other functionals that do not incorporate a full range separation. For instance, the CAM-B3LYP functional has only 65% HF exchange at LR. Some hybrid functionals have also been included in our calculations because of their widespread use and good performances for a wide range of properties. Finally, despite their different parametrizations, some of the specific hybrid functionals such as B3LYP and M06-HF have been considered because they represent different extremes of hybrids where the contribution of HF exchange ranges from a low value for the former (20%) to a very high value for the latter (100%). Table 3 compares the performance of the best OT-RSH functionals in this work with respect to other approximations. First, it is found

that the values of signed deviations using all of the tested approximations are positive, leading to the more reinforced positive MSDs in Table 3 and consequently a tendency of overestimating the oxidation potentials. We can also see from these data that the LC-BLYP functional with the ingredients α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 and α = 0.0, β = 1.0, and μ = 0.25 Bohr−1(MADs of 0.44 and 0.47 V, respectively) outperforms other functionals with both fixed and interelectronic distancedependent HF exchange contributions. With more details, comparing the results of the new OT-RSH functionals with other RS functionals and standard LC approximations, it can be concluded that not only is the lower value of parameter μ an important factor for oxidation potential calculations but also imposing the correct 100% asymptotic HF exchange with an G

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 4. Geometrical structures and their names for the compounds in the second benchmarked set.

appropriate proportion of range separation parameter can improve the results. However, the results of the developed OT-RSH functionals are also comparable to those provided by other approximations (such as BMK, CAM-B3LYP, etc.) and even in some cases can be outperformed by some other parametrized functionals for which the fitting processes against various data sets have been carried out in their construction (such as B3LYP). Finally, the parametrized M06-HF functional with full 100% HF exchange over the entire range has poor performance among others. As further evaluation of the proposed OT-RSH functionals, we have also benchmarked the performance of these functionals for predicting the oxidation potentials of molecules for which the tuning processes of the range separation parameter have not been performed. In fact, it would be interesting to see whether or not there exist similar good performances for these approximations with new parameters to estimate the oxidation potentials for some other compounds. To this end, we have considered several molecules whose geometrical structures and names are presented in Figure 4. The statistical measures of the performance of the OT-RSH functionals in comparison to the standard LC functionals for predicting the oxidation potentials of these molecules are gathered in Table 4. Upon examination of the table, we observe that not only the OT-RSHs outperform the

Table 4. Statistical Measures of the Performance of the Proposed OT-RSH Functionals and Other Standard LC Approximations for Estimation of Oxidation Potentials of the Molecules in the Second Benchmarked Seta functionals

MSD

MAD

MaxAD

RMSD

0.51 0.63 0.61

0.40 0.54 0.52

0.60 0.68 0.70

0.40 0.55 0.53

0.96 1.06 1.07

0.69 0.85 0.82

−1

α = 0, β = 1, and μ = 0.25 Bohr LC-BLYP 0.39 0.39 LC-PBE 0.54 0.54 LC-TPSS 0.51 0.51 α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 LC-BLYP 0.39 0.39 LC-PBE 0.54 0.54 LC-TPSS 0.53 0.53 α = 0, β = 1, and μ = 0.47 Bohr−1 LC-BLYP 0.67 0.67 LC-PBE 0.84 0.84 LC-TPSS 0.81 0.81 a

All values are in V.

standard LC functionals but also the two newly LC-BLYP functionals with MADs of roughly 0.39 V are the top performers. Moreover, there is still the same ordering as that in previous cases for the OT-RSH functionals. Finally, from a more general H

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 5. Correlation plots of the negative highest occupied molecular orbital of molecules in solution, −HOMO(M), and the negative lowest unoccupied molecular orbital of cations in solution, −LUMO(M+), with the oxidation potentials computed using OT-RSH functionals with various density functional approximations (BLYP, PBE, and TPSS) for the cases (a) α = 0.0 and β = 1.0 and (b) α = 0.2 and β = 0.8.

energies. Also, from the regression coefficients, in all cases, we see that the best correlations appeared for the OT-RSH functionals based on BLYP. On the other hand, considering the correlations between computed and experimental oxidation potentials, it is expected that there should be an approximate quantitative description of the oxidation potentials based on frontier orbital energies. As a piece of supporting evidence to these arguments, Figure 6 exhibits the correlation plots between experimental oxidation potentials and the OT-RSH functional energies of −HOMO(M) and −LUMO(M+) in solution. It can be seen from the figures that the experimental data and the frontier orbital energies in solution correlate with each other, confirming that the −HOMO(M) and −LUMO(M+) energies can be good approximations to the oxidation potentials. In addition, comparing the correlations related to the −HOMO(M) energies with the −LUMO(M+) counterparts, we find that much stronger correlations are obtained when using the energies of −HOMO(M). However, it is noted here that if there are no self-interaction errors in the used approach and the MOs do not relax when an electron is removed the eigenvalue εi will be a constant for ni between 1 and 0 (HOMO(M) equals LUMO(M+)), and

viewpoint, our test calculations reveal that the optimally tuned values of the range separation parameter for the benchmarked compounds and other organic molecules similar to those considered here are smaller than those provided by other LC functionals for different properties.14,19,55−58 However, the obtained parameters are in reasonable agreement with the fixed range separation value recommended for the original LC-ωPBE0 functional, 0.2 Bohr−1.28 Putting all of the results together, we can recommend the LC-BLYP with α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 and/or α = 0.0, β = 1.0, and μ = 0.25 Bohr−1 for predicting the oxidation potentials. We now deal with the relationship between MO eigenvalues obtained from OT-RSHs and oxidation potentials. Looking at eqs 3 and 4, we infer that the frontier orbitals energies should be related to the computed oxidation potentials. To numerically verify this point, shown in Figure 5 are the correlations of computed oxidation potentials with −HOMO(M) and −LUMO(M+) energies in solution. As can be seen from the figures, it becomes unambiguously apparent that there exist remarkably strong linear relationships between the computed oxidation energies from the difference of total electronic energies of molecules and cations in solution and their frontier orbital I

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 6. Correlation plots of the experimental oxidation potentials with the negative highest occupied molecular orbital of molecules in solution, −HOMO(M), and the negative lowest unoccupied molecular orbital of cations in solution, −LUMO(M+), computed using OT-RSH functionals with various density functional approximations (BLYP, PBE, and TPSS) for the cases (a) α = 0.0 and β = 1.0 and (b) α = 0.2 and β = 0.8.

To wrap up, we are acutely aware that such functional developments and DFT benchmarking studies are not the whole story. Hopefully, the proposed OT-RSH models in this study can provide a balanced description for calculation of oxidation potentials and other related parameters in the research field of electrochemistry, at least for compounds with sizes and structures close to those of the considered systems herein.

consequently, the left side of eq 2 is exactly equal to either −HOMO(M) or −LUMO(M+). Finally, observing such correlations prompts us to check another proposal in our calculations, namely, using the MO eigenvalues in the gas phase. As a matter of fact, we found from the numerical data that the HOMO energies of molecules in the gas phase correlate with their counterparts in solution. This is because the effect of the surrounding solvent on the molecular geometry and MOs of a neutral molecule can be small. Accordingly, besides total electronic energies and MO eigenvalues in solution, the OT-RSH functional energies of −HOMO(M) in the gas phase may also be considered as another approximate route for estimation of the oxidation potentials, as shown from the correlations in Figure 7. However, although the correlations are not as strong as those found using the frontier orbital energies in solution, such relationships are encouraging and deserve additional analysis, wherein a computed quantity in the gas phase can be used to predict a molecular property in solution, the oxidation potential.

5. FINAL COMMENTS During this work, the applicability and accountability of optimally tuned range-separated hybrid functionals for predicting the oxidation potentials of organic compounds from different categories have been evaluated. More specifically, we dissected the role of nonempirical optimization of the range separation parameter and the importance of other factors such as SR and LR exact-like exchanges in calculations of the oxidation potential using the proposed optimally tuned functionals. It was shown that among the optimally tuned approximations, the top-ranked LC-BLYP with the ingredients α = 0, β = 1, and μ = 0.25 Bohr−1 and α = 0.2, β = 0.8, and μ = 0.20 Bohr−1 not only perform better J

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A



Article

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +98 71 36137160. Fax: +98 71 36460788. ORCID

Mojtaba Alipour: 0000-0003-3037-0232 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The authors are grateful to Shiraz University for providing the computational facilities for this project.



REFERENCES

(1) Hohenberg, P.; Kohn, W. Inhomogeneus Electron Gas. Phys. Rev. 1964, 136, B864−B871. (2) Kohn, W.; Sham, L. Self-consistent Equations Including Exchange and Correlation Effects. Phys. Rev. 1965, 140, A1133−A1138. (3) Tsuneda, T. Density Functional Theory in Quantum Chemistry; Springer: New York, 2014. (4) Albright, T. A.; Burdett, J. K.; Whangbo, M.-H. Orbital Interactions in Chemistry; John Wiley & Sons: Hoboken, NJ, 2013. (5) Stowasser, R.; Hoffmann, R. What Do the Kohn-Sham Orbitals and Eigenvalues Mean? J. Am. Chem. Soc. 1999, 121, 3414−3420. (6) Zhang, G.; Musgrave, C. B. Comparison of DFT Methods for Molecular Orbital Eigenvalue Calculations. J. Phys. Chem. A 2007, 111, 1554−1561. (7) Chen, J. H.; He, L. M.; Wang, R. L. Connection of DFT Molecular Orbital Eigenvalues with the Observable Oxidation Potentials/ Oxidation Energies. J. Phys. Chem. A 2013, 117, 5132−5139. (8) Liu, S.; Schauer, C. K. Origin of Molecular Conformational Stability: Perspectives from Molecular Orbital Interactions and Density Functional Reactivity Theory. J. Chem. Phys. 2015, 142, 054107−1− 054107−8. (9) Janak, J. F. Proof that ∂E/∂ni = ε in Density-functional Theory. Phys. Rev. B: Condens. Matter Mater. Phys. 1978, 18, 7165−7168. (10) Stein, T.; Kronik, L.; Baer, R. Reliable Prediction of Charge Transfer Excitations in Molecular Complexes Using Time-Dependent Density Functional Theory. J. Am. Chem. Soc. 2009, 131, 2818−2820. (11) Stein, T.; Kronik, L.; Baer, R. Prediction of Charge-Transfer Excitations in Coumarin-Based Dyes Using A Range-Separated Functional Tuned From First Principles. J. Chem. Phys. 2009, 131, 244119/1−244119/5. (12) Kronik, L.; Stein, T.; Refaely-Abramson, S.; Baer, R. Excitation Gaps of Finite-Sized Systems from Optimally Tuned Range-Separated Hybrid Functionals. J. Chem. Theory Comput. 2012, 8, 1515−1531. (13) Karolewski, A.; Kronik, L.; Kümmel, S. Using Optimally Tuned Range Separated Hybrid Functionals in Ground-State Calculations: Consequences and Caveats. J. Chem. Phys. 2013, 138, 204115/1− 204115/11. (14) Iikura, H.; Tsuneda, T.; Yanai, T.; Hirao, K. A Long-Range Correction Scheme For Generalized-Gradient-Approximation Exchange Functionals. J. Chem. Phys. 2001, 115, 3540−3544. (15) Tawada, Y.; Tsuneda, T.; Yanagisawa, S.; Yanai, T.; Hirao, K. A Long-Range-Corrected Time-Dependent Density Functional Theory. J. Chem. Phys. 2004, 120, 8425−8433. (16) Kamiya, M.; Sekino, H.; Tsuneda, T.; Hirao, K. Nonlinear Optical Property Calculations by the Long-Range-Corrected Coupled-Perturbed Kohn−Sham Method. J. Chem. Phys. 2005, 122, 234111/1− 234111/10. (17) Yanai, T.; Tew, D. P.; Handy, N. C. A New Hybrid Exchange− Correlation Functional Using the Coulomb-Attenuating Method (CAM-B3LYP). Chem. Phys. Lett. 2004, 393, 51−57. (18) Vydrov, O. V.; Scuseria, G. E. Assessment of a Long Range Corrected Hybrid Functional. J. Chem. Phys. 2006, 125, 234109/1− 234109/9.

Figure 7. Correlation plots of the experimental oxidation potentials with the gas-phase negative highest occupied molecular orbital of molecules, −HOMO(M) in the gas phase, computed using OT-RSH functionals with various density functional approximations (BLYP, PBE, and TPSS) for the cases (a) α = 0.0 and β = 1.0 and (b) α = 0.2 and β = 0.8.

than their PBE and TPSS counterparts and other standard LC functionals but also in many cases outperform other conventional hybrid functionals with a fixed amount of exact-like exchange. Moreover, we observed that these new optimally tuned functionals describe well the oxidation potentials of compounds for which tuning of the range separation parameter was not carried out. From another perspective, the utility of frontier orbital energies computed using the optimally tuned functionals for estimating the oxidation potentials has also been explored. It was found that the −HOMO(M) in solution (and also in the gas phase) and the −LUMO(M+) correlate remarkably with the experimental oxidation potentials. Thus, in addition to the use of total electronic energies of optimally tuned approximations, these quantities can also be considered as approximate candidates for prediction of oxidation potentials. K

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

(39) Lee, C. T.; Yang, W. T.; Parr, R. G. Development of the ColleSalvetti Correlation- Energy Formula into a Functional of the ElectronDensity. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (40) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum Mechanical Continuum Solvation Models. Chem. Rev. 2005, 105, 2999−3094. (41) Namazian, M.; Norouzi, P.; Ranjbar, R. Prediction of Electrode Potentials of Some Quinone Derivatives in Acetonitrile. J. Mol. Struct.: THEOCHEM 2003, 625, 235−241. (42) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (43) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (44) Tao, J.; Perdew, J. P.; Staroverov, V. N.; Scuseria, G. E. Climbing the Density Functional Ladder: Nonempirical Meta-Generalized Gradient Approximation Designed for Molecules and Solids. Phys. Rev. Lett. 2003, 91, 146401/1−146401/4. (45) Refaely-Abramson, S.; Sharifzadeh, S.; Jain, M.; Baer, R.; Neaton, J. B.; Kronik, L. Gap Renormalization of Molecular Crystals From Density-Functional Theory. Phys. Rev. B: Condens. Matter Mater. Phys. 2013, 88, 081204/1−081204/5. (46) Refaely-Abramson, S.; Jain, M.; Sharifzadeh, S.; Neaton, J. B.; Kronik, L. Solid-state optical absorption from optimally tuned timedependent range-separated hybrid density functional theory. Phys. Rev. B: Condens. Matter Mater. Phys. 2015, 92, 081204/1−081204/4. (47) Sarkar, S.; Kronik, L. Ionisation and (de-)Protonation Energies of Gas-Phase Amino Acids from an Optimally Tuned Range-Separated Hybrid Functional. Mol. Phys. 2016, 114, 1218−1224. (48) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; et al. Gaussian 09, revision B.01; Gaussian, Inc.: Wallingford, CT, 2010. (49) Perdew, J. P.; Schmidt, K. Jacob’s Ladder of Density Functional Approximations for the Exchange-Correlation Energy. AIP Conf. Proc. 2001, 577, 1−20. (50) Zhao, Y.; Schultz, N. E.; Truhlar, D. G. Design of Density Functionals by Combining the Method of Constraint Satisfaction with Parametrization for Thermochemistry, Thermochemical Kinetics, and Noncovalent Interactions. J. Chem. Theory Comput. 2006, 2, 364−382. (51) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-Class Functionals and 12 Other Functionals. Theor. Chem. Acc. 2008, 120, 215− 241. (52) Zhao, Y.; Truhlar, D. G. Density Functional for Spectroscopy: No Long-Range Self-Interaction Error, Good Performance for Rydberg and Charge-Transfer States, and Better Performance on Average than B3LYP for Ground States. J. Phys. Chem. A 2006, 110, 13126−13130. (53) Boese, A. D.; Martin, J. M. L. Development of Density Functionals for Thermochemical Kinetics. J. Chem. Phys. 2004, 121, 3405−3416. (54) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B: Condens. Matter Mater. Phys. 1986, 33, 8822−8824. (55) Garrett, K.; Sosa Vazquez, X. A.; Egri, S. B.; Wilmer, J.; Johnson, L. E.; Robinson, B. H.; Isborn, C. M. Optimum Exchange for Calculation of Excitation Energies and Hyperpolarizabilities of Organic Electro-Optic Chromophores. J. Chem. Theory Comput. 2014, 10, 3821−3831. (56) Johnson, L. E.; Dalton, L. R.; Robinson, B. H. Optimizing Calculations of Electronic Excitations and Relative Hyperpolarizabilities of Electrooptic Chromophores. Acc. Chem. Res. 2014, 47, 3258−3265. (57) Nénon, S.; Champagne, B.; Spassova, M. I. Assessing Long-Range Corrected Functionals with Physically-Adjusted Range-Separated Parameters for Calculating the Polarizability and the Second Hyperpolarizability of Polydiacetylene and Polybutatriene Chains. Phys. Chem. Chem. Phys. 2014, 16, 7083−7088. (58) Oviedo, M. B.; Ilawe, N. V.; Wong, B. M. Polarizabilities of πConjugated Chains Revisited: Improved Results from Broken-

(19) Song, J.-W.; Hirosawa, T.; Tsuneda, T.; Hirao, K. Long-Range Corrected Density Functional Calculations of Chemical Reactions: Redetermination of Parameter. J. Chem. Phys. 2007, 126, 154105/1− 154105/7. (20) Chai, J.-D.; Head-Gordon, M. Systematic Optimization of LongRange Corrected Hybrid Density Functionals. J. Chem. Phys. 2008, 128, 084106/1−084106/15. (21) Chai, J.-D.; Head-Gordon, M. Long-Range Corrected Hybrid Density Functionals with Damped Atom-Atom Dispersion Corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (22) Peverati, R.; Truhlar, D. G. Improving the Accuracy of Hybrid Meta-GGA Density Functionals by Range Separation. J. Phys. Chem. Lett. 2011, 2, 2810−2817. (23) Lin, Y.-S.; Tsai, C.-W.; Li, G.-D.; Chai, J.-D. Long-Range Corrected Hybrid Meta-Generalized-Gradient Approximations with Dispersion Corrections. J. Chem. Phys. 2012, 136, 154109/1−154109/ 12. (24) Kar, R.; Song, J.-W.; Hirao, K. Long-Range Corrected Functionals Satisfy Koopmans’ Theorem: Calculation of Correlation and Relaxation Energies. J. Comput. Chem. 2013, 34, 958−964. (25) Alipour, M. How Does LCDFT Compare to SAC-CI for the Treatment of Valence and Rydberg Excited States of Organic Compounds? J. Phys. Chem. A 2014, 118, 1741−1747. (26) Tsuneda, T.; Hirao, K. Long-Range Correction for Density Functional Theory. WIREs Comput. Mol. Sci. 2014, 4, 375−390. (27) Song, J.-W.; Hirao, K. Long-Range Corrected Density Functional Theory with Accelerated Hartree-Fock Exchange Integration Using a Two-Gaussian Operator [LC-ωPBE(2Gau)]. J. Chem. Phys. 2015, 143, 144112/1−144112/8. (28) Rohrdanz, M. A.; Martins, K. M.; Herbert, J. M. A Long-RangeCorrected Density Functional That Performs Well for Both GroundState Properties and Time-Dependent Density Functional Theory Excitation Energies, Including Charge-Transfer Excited States. J. Chem. Phys. 2009, 130, 054112/1−054112/8. (29) Egger, D. A.; Weissman, S.; Refaely-Abramson, S.; Sharifzadeh, S.; Dauth, M.; Baer, R.; Kummel, S.; Neaton, J. B.; Zojer, E.; Kronik, L. Outer-valence Electron Spectra of Prototypical Aromatic Heterocycles from an Optimally Tuned Range-Separated Hybrid Functional. J. Chem. Theory Comput. 2014, 10, 1934−1952. (30) Raeber, A. E.; Wong, B. M. The Importance of Short- and LongRange Exchange on Various Excited State Properties of DNA Monomers, Stacked Complexes, and Watson-Crick Pairs. J. Chem. Theory Comput. 2015, 11, 2199−2209. (31) Bokareva, O. S.; Grell, G.; Bokarev, S. I.; Kühn, O. Tuning RangeSeparated Density Functional Theory for Photocatalytic Water Splitting Systems. J. Chem. Theory Comput. 2015, 11, 1700−1709. (32) Alipour, M.; Fallahzadeh, P. First Principles Optimally Tuned Range-Separated Density Functional Theory for Prediction of Phosphorus-Hydrogen Spin-Spin Coupling Constants. Phys. Chem. Chem. Phys. 2016, 18, 18431−18440. (33) Alipour, M.; Fallahzadeh, P. Nonempirically Tuning RangeSeparated Functionals for Dipole Polarizabilities of Nanostructures Containing Hydrogen Bonds. Theor. Chem. Acc. 2017, 136, 22/1−22/8. (34) Refaely-Abramson, S.; Sharifzadeh, S.; Govind, N.; Autschbach, J.; Neaton, J. B.; Baer, R.; Kronik, L. Quasiparticle Spectra from a Nonempirical Optimally Tuned Range-Separated Hybrid Density Functional. Phys. Rev. Lett. 2012, 109, 226405/1−226405/6. (35) Srebro, M.; Autschbach, J. Does a Molecule-Specific Density Functional Give an Accurate Electron Density? The Challenging Case of the CuCl Electric Field Gradient. J. Phys. Chem. Lett. 2012, 3, 576−581. (36) Vollmer, J. M.; Curtiss, L. A.; Vissers, D. R.; Amine, K. Reduction Mechanisms of Ethylene, Propylene, and Vinylethylene Carbonates A Quantum Chemical Study. J. Electrochem. Soc. 2004, 151, A178−A183. (37) Wang, R. L.; Buhrmester, C.; Dahn, J. R. Calculations of Oxidation Potentials of Redox Shuttle Additives for Li-Ion Cells. J. Electrochem. Soc. 2006, 153, A445−A449. (38) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. L

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A Symmetry Range-Separated DFT and New CCSD(T) Benchmarks. J. Chem. Theory Comput. 2016, 12, 3593−3602.

M

DOI: 10.1021/acs.jpca.7b03811 J. Phys. Chem. A XXXX, XXX, XXX−XXX