Fermi Contact Spin Density Calculations of Aromatic Radicals

Sep 17, 2010 - The Gaussian-weighted RC operator with small r0 has been used in DFT spin density calculations.13 First-row atoms, first- row diatomic ...
0 downloads 0 Views 1MB Size
Download PDF
20648

J. Phys. Chem. C 2010, 114, 20648–20658

Fermi Contact Spin Density Calculations of Aromatic Radicals† Liyuan Liang and Vitaly A. Rassolov* Department of Chemistry and Biochemistry, UniVersity of South Carolina, South Carolina 29208 ReceiVed: June 15, 2010; ReVised Manuscript ReceiVed: August 26, 2010

The Gaussian-weighted operator derived by Rassolov and Chipman (RC operator) is applied to the spin density of large organic π-radicals using B3LYP and PBE0 density functionals and compared to the performance of the delta-function operator. The influence of the basis sets is examined. It is shown that both the deltafunction operator and the RC operator provide comparable spin densities with the basis sets ext-6-311++G** and ext-6-311++G(2d,2p). An optimum range parameter of 0.25 or 0.35 a0 is appropriate for the RC operator with the extended basis sets. For standard basis sets 6-311++G** and 6-311++G(2d,2p), the delta-function operator may yield considerable errors for heavy atoms due to the deficiency of the basis sets near the nuclei. To a large extent, however, this deficiency can be cured with the RC operator. With the same range parameter of 0.25 or 0.35 a0 previously recommended, the errors are significantly reduced with respect to the deltafunction method. If a relatively larger range parameter of 0.45 a0 is used, a further improvement can be obtained and the spin densities at the B3LYP/6-311++G** or B3LYP/6-311++G(2d,2p) level of methods are similar to the ones at the B3LYP/ext-6-311++G** or B3LYP/ext-6-311++G(2d,2p) level of methods. Generally, the UB3LYP/6-311++G** and the UB3LYP/ext-6-311++G** methods are able to provide satisfactory spin densities. For some anions, the B3LYP/6-311++G(2d,2p) and the B3LYP/ext-6311++G(2d,2p) methods are needed to achieve better agreement with the experiment. It is demonstrated that the B3LYP functional combined with the standard Pople-type basis sets 6-311++G** and 6-311++G(2d,2p), which are designed for conventional electronic structure calculations, serves as an efficient and reliable method for the calculation of spin densities of large organic radicals when the Rassolov-Chipman operator is utilized, while the delta-function operator may lead to far less reliable results. An optimum range parameter of 0.45 a0 is recommended. I. Introduction 1

Electron spin resonance (ESR) spectroscopy is of central importance for the study of radicals. Radicals exhibit Zeeman splitting in the presence of an external magnetic field. The interaction between the spin of the unpaired electron and nearby nuclear spins endows most ESR spectra with a multilined spectral pattern. The phenomenon, known as hyperfine structure, contains information about geometric and electronic structure of radicals. Specifically, it provides a sensitive probe into the mapping of the distribution of the unpaired electron spin. Two mechanisms lead to hyperfine interaction. The dipolar, or anisotropic hyperfine interaction, arises from the interaction between an electron and a nuclear dipole some distance away. Most often the ESR spectra of radicals, especially of the larger ones, are recorded in solution. The random rotational molecular motion in solution is so rapid that this interaction averages in time to zero. In this case the only remaining hyperfine interaction term is Fermi contact, or isotropic hyperfine interaction, which is proportional to the electron spin density at the nucleus. The calculation of anisotropic hyperfine coupling constant poses no special difficulties. However, accurate calculation of Fermi contact hyperfine coupling constant remains a challenging problem in computational chemistry.2 The reason for the difficulty is mainly due to the deltafunction singularity present in the isotropic hyperfine interaction operator. While it is formally easy to compute, the spin densities at the nuclei are highly sensitive to the quality of the molecular †

Part of the “Mark A. Ratner Festschrift”.

wave function at the nuclear position. As a consequence, the spin densities at the nuclei, therefore the isotropic hyperfine coupling constants (HFCCs), are strongly influenced by the electron correlation, the one-electron basis set, and the radical geometry. This is especially true for π-radicals, where spinrestricted open shell uncorrelated wave function yields no hyperfine interaction. Reliable calculations of isotropic HFCC require the use of the most refined post Hartree-Fock ab initio approaches, such as multiconfiguration self-consistent field (MCSCF), multireference configuration interaction (MRCI), or coupled cluster (CC) techniques. Moreover, large and specifically tailored basis sets are needed to remedy the deficiency of the Cartesian Gaussian basis functions (CGTOs), which are commonly used in molecular quantum mechanical calculations but do not satisfy the desired Kato’s cusp condition. Such calculations are highly demanding in computational resources, and only applicable to very small radicals. It is desirable to devise a simpler and more efficient method to predict isotropic HFCC for large systems. Density functional theory (DFT) has emerged as an effective alternative to the traditional ab initio approaches for the electronic structure calculation. The main merit of DFT lies in its ability to incorporate electron correlation through the exchange-correlation functionals while retaining a favorable scaling with the system size in comparison with correlated ab initio methods. As a consequence, larger molecular systems are accessible. In particular, Becke’s 3-parameter exchange functional3 in conjunction with the nonlocal Lee-Yang-Parr electron correlation functional4 (B3LYP) has been shown to provide accurate descriptions of isotropic and anisotropic

10.1021/jp105511a  2010 American Chemical Society Published on Web 09/17/2010

Spin Density Calculations of Aromatic Radicals

J. Phys. Chem. C, Vol. 114, No. 48, 2010 20649

hyperfine coupling constants.5,6 These computational advantages have made DFT the most widely used method for calculating hyperfine coupling constants and even the only option when large radicals are tackled. However, the calculation of spin densities with DFT, like its ab initio counterparts, strongly depends on the quality of the basis sets. One approach to overcome the basis set dependence is to substitute the delta-function operator with another global operator which samples the wave function at all points in space. An operator of this type was first proposed by Hiller, Sucher, and Feinberg (HSF).7 HSF operator and the delta-function operator give identical results for the spin density if the exact wave function is applied. It was later proven that the HSF operator and delta-function have the same expectation values in the complete basis set limit for UHF8 and any MCSCF9 wave functions. With other approximate wave functions they generally give different results, and an HSF operator may be preferable to the delta-function operator. However, the HSF operator displays incorrect long-range asymptotic behavior of the density with most approximate wave functions and therefore limits the accuracy of spin density calculations. An additional serious problem for implementation of the HSF operator in molecular systems is that the evaluation of many difficult two-electron integrals is needed. This calculation may require more computer resources than the calculation of the wave function itself. Rassolov and Chipman10,11 have developed a new class of alternative operators (RC operators) that are free of many of the drawbacks of both the delta-function and the HSF formulations. Its derivation is based on the hypervirial theorem.12 It ˆ , then states that if Ψ is an eigenfunction of the Hamiltonian H ˆ , W], where W stands the expectation value of the commutator [H for any time-independent operator, is equal to zero, i.e.,

ˆ ,W]|Ψ〉 ) 0 〈Ψ|[H

(1)

Setting W ) ΣiF(ri)(∂/∂ri)2Sˆzi, where summation runs over all electrons in the system and F(r) is an arbitrary radial weight function yields the expression

1 q (0) ) 〈Ψ| 2πF(0) F

[(

(

Lˆi2

∑ F(ri) - r 3 i

)

i

+

)

∂U ∂V + + ∂ri ∂ri

]

2 ∂F(ri) ∂2 1 ∂ 1 1 ∂ F(ri) ∂ ˆ + + 2S |Ψ〉 (2) ∂ri ri ∂ri 2 ∂r 2 ∂ri zi ∂ri2 2ri2 i

Here Lˆ is the total angular momentum operator for the ith electron while U and V designate the one- and two-electron parts of the potential energy operator, respectively. The expectation value of this operator is equal to that of a delta-function for any fully optimized (in a complete basis set) MCSCF wave function.9 The term dependent on V involves complicated twoelectron integrals, which can be made negligible by setting the radial weight function F(r) to be of sufficiently short-range. In the present work we use a Gaussian weight function F(r) ) 2 e-(r/r0) centered at the nucleus of interest with a sufficiently small parameter r0. This parameter should be rescaled by the value of a nuclear charge Z to achieve uniform accuracy at different nuclei. The RC operator leads to the results equivalent to the delta-function in the limit of r0 f 0. At the opposite extreme, when r0 f ∞, the new operator becomes identical with the HSF operator.

The Gaussian-weighted RC operator with small r0 has been used in DFT spin density calculations.13 First-row atoms, firstrow diatomic hydrides, and several small polyatomic radicals were considered. With first-row atoms it was found that the best DFT results are obtained with a value of 0.35 a0 for Zr0, which is a little larger than the 0.25 a0 value recommend by Rassolov and Chipman for ab initio studies. The study concluded that with moderate-sized and carefully selected Gaussian basis sets, the RC operator combined with DFT has clear advantage for the calculation of nuclear spin densities in molecules. The aim of this paper is to further explore the performance of the RC operator combined with the DFT approach in predicting the isotropic hyperfine coupling constants of mediumsized organic radicals. II. Computational Methods The spin-unrestricted formalism was employed for all radicals. All geometries were optimized with the 6-311++G** basis set and B3LYP hybrid functional. Subsequently, computations of single-point spin densities were performed using four basis sets of contracted Gaussian functions: (i) 6-311++G**, (ii) ext-6311++G**, (iii) 6-311++G(2d,2p), (iv) ext-6-311++G(2d,2p). The two extended basis sets, ext-6-311++G** and ext-6311++G(2d,2p), are constructed from the corresponding standard basis sets by adding a tight s-function to the core s-orbitals and uncontracting two core s-functions with the largest and the second largest exponents. The exponent of the tight s-function ζ is determined by the ratio relation ζ:ζ1 ) ζ1:ζ2, where ζ1 and ζ2 are the largest and the second largest exponents of the core s-orbitals. Such a procedure was used by Chipman in deltafunction-based studies of spin densities.14 The exponents of the tight s-function are listed as the following: H ) 225.10, C ) 30 531.417, N ) 41 734.500, O ) 56 861.400. The pure spherical harmonic Gaussians are employed in the extended basis sets. The choice of the basis sets in this study is based on the following considerations. It is known that basis sets of triple-ζ quality plus multiple polarization functions and diffuse functions are required for accurate spin density calculation. Furthermore, correct core-valence balance is needed.14 Pople-style 6311++G** and 6-311++G(2d,2p) basis sets, which are designed for conventional electronic structure calculations, satisfy these basic requirements and provide a good starting point for basis sets with the additional flexibility in the tight region. The effect of addition of tight s-functions and uncontraction of innermost s-function has proved to be crucial in the prediction of spin densities within the conventional delta-function formalism. In this paper we explore the need for tight region basis set flexibility with the new operators. A remarkable feature of Pople’s basis sets 6-311++G** and 6-311++G(2d,2p) is that they are smaller in magnitude than Dunning’s augmented correlation-consistent polarized valence triple-ζ basis set aug-cc-pVTZ. Although Dunning’s basis sets have the advantage to extrapolate to the basis set limit, this procedure is hardly applicable to a large organic radical. In this study we will compare the performance of the two sets of standard basis sets for HFCCs calculations. All of the calculations were performed using modified version of Q-Chem.15 The geometrical structures of the radicals are shown in Figures 1 and 2. III. Results and Discussion A. Benzyl. Experimental HFCCs for 13C (1), 13C (7) and all of the 1H nuclei are available.16 This radical is neutral with very

20650

J. Phys. Chem. C, Vol. 114, No. 48, 2010

Liang and Rassolov

Figure 1. Sketch of benzyl, 1-pyrolyl, and fulvalene anions.

Figure 2. Sketch of p-benzosemiquinone anion.

weak polarity. The calculated dipole moment at the UB3LYP/ 6-311++G** level of theory is 0.1546 D. Previous B3LYP calculations have shown that vibrational averaging and solvent effects play only a negligible role in the values of HFCCs.17 Therefore only gas phase calculations are carried out to compare with the experiment. The spin density calculations are shown in Table 1. The effects of basis sets on HFCCs is investigated separately to assess the role of the (1) polarization functions by comparing

ext-6-311++G** with the ext-6-311++G(2d,2p) calculations and (2) of the tight s-region, by comparing 6-311++G(...) with the ext-6-311++G(...). First, let us examine the performance of the delta-function operator. For hydrogen atoms, the 6-311++G** basis set yields accurate HFCCs. Further extension or expansion of this basis set brings no significant changes in HFCCs. Thus for this radical, 6-311++G** may be close to the complete basis set limit. The calculated spin densities of 13C show stronger dependence on the basis sets used. UB3LYP/6-311++G** and UB3LYP/ext6-311++G** produce spin densities of -0.03514 and -0.03832, respectively, for C(1), which are close to the experimental value of -0.03603. UB3LYP/6-311++G** is slightly better due to, we suggest, fortuitous error cancellation. The difference in spin density of C(7) shows larger variation with UB3LYP/ext-6311++G**, yielding a result much closer to experiment. The relative errors for UB3LYP/6-311++G** and UB3LYP/ext6-311++G** are 26% and 6.6%, respectively. It is clear that the flexibility of the basis set in the vicinity of the nuclear position is necessary. The expansion of polarization function plays a minor role on the HFCCs of 13C. The performance of the RC operator depends on an appropriate choice of the range parameter r0. The dependence of HFCCs on r0 is shown in Figures 3, 4, and 5. We observe the following: (i) The standard basis sets show stronger r0 dependence. In the range of Zr0 ∈ [0.1, 0.8] the HFCCs of heavy atoms are approximately linear with an appreciable slope. The hydrogen atoms show a larger deviation from a linear relation. However, the HFCC still varies monotonously with r0. (ii) The extended basis set results are much less r0-dependent than their corresponding standard basis sets, in agreement with the original rationale for the RC operators.10,11 The plateau region of HFCC has Zr0 in the range [0.10, 0.50]. For heavy atoms, when Zr0 is in the range [0.20, 0.40], a relatively flat parabola-shaped curve is observed. The HFCCs are not sensitive to the variation of r0, and the optimum value of r0 is located in this approximately flat region. For hydrogen, this flat region is extended to a larger upper limit. At larger r0

TABLE 1: Spin Densities in a Benzyl Radical Obtained with B3LYP Functional and Different Basis Sets C(1)

C(7)

H(2,6)

H(3,5)

H(4)

H(7)

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

aug-cc-pVTZ

expta

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

-0.03514 -0.03777 -0.03821 -0.03870 0.04488 0.05197 0.05422 0.05658 -0.003315 -0.003592 -0.003564 -0.003530 0.001392 0.001507 0.001495 0.001479 -0.003832 -0.004147 -0.004112 -0.004071 -0.01006 -0.01086 -0.01075 -0.01063

-0.03832 -0.03872 -0.03870 -0.03903 0.05692 0.05636 0.05639 0.05802 -0.003403 -0.003479 -0.003495 -0.003497 0.001426 0.001457 0.001464 0.001464 -0.003916 -0.004005 -0.004025 -0.004027 -0.01019 -0.01041 -0.01047 -0.01048

-0.03531 -0.03801 -0.03847 -0.03898 0.04624 0.05342 0.05567 0.05803 -0.003377 -0.003620 -0.003572 -0.003519 0.001431 0.001531 0.001508 0.001483 -0.003941 -0.004222 -0.004165 -0.004101 -0.01020 -0.01089 -0.01073 -0.01055

-0.03862 -0.03900 -0.03896 -0.03931 0.05829 0.05774 0.05776 0.05938 -0.003392 -0.003466 -0.003477 -0.003469 0.001428 0.001458 0.001463 0.001459 -0.003948 -0.004035 -0.004049 -0.004040 -0.01015 -0.01036 -0.01040 -0.01037

-0.02963 -0.03172 -0.03214 -0.03259 0.02943 0.03618 0.03847 0.04080 -0.003188 -0.003463 -0.003442 -0.003420 0.001355 0.001472 0.001462 0.001451 -0.003721 -0.004044 -0.004019 -0.003992 -0.00957 -0.01039 -0.01033 -0.01026

-0.03603

Reference 16.

0.06096

-0.00323

0.00112

-0.00387

-0.0102

Spin Density Calculations of Aromatic Radicals

J. Phys. Chem. C, Vol. 114, No. 48, 2010 20651

Figure 3. Spin density on a C(1) nucleus of benzyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

Figure 4. Spin density on a C(7) nucleus of benzyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

values the HFCCs are also linear, suggesting that the contributions of the two-electron integrals begin to increase. (iii) When r0 increases, the difference (gap) between the calculated HFCCs with the standard basis sets and their corresponding extended basis sets goes down and reaches a minimum. Afterward the HFCC-r0 relations for the extended basis sets are also linear with approximately same slope. An optimum range parameter of 0.25 or 0.35 is obtained for extended basis functions ext-6311++G** and ext-6-311++G(2d,2p). The two range parameters produce almost identical spin densities. While for standard basis sets 6-311++G** and 6-311++G(2d,2p), a larger opti-

mum Zr0, which locates in the range of [0.45, 0.50], is expected. We recommend 0.45 as optimum value. Comparing the performance of the delta-function and the RC operators, the extended basis sets ext-6-311++G** and ext-6311++G(2d,2p) yield very close and relatively accurate spin densities. With the standard basis sets 6-311++G** and 6-311++G(2d,2p), the delta-function operator method may produce large errors. In this radical, the most challenging nucleus is C(7). As indicated previously, the relative error for UB3LYP/ 6-311++G** is 26%. The RC-operator method, however, produces a result with a relative error of 7.2%.

20652

J. Phys. Chem. C, Vol. 114, No. 48, 2010

Liang and Rassolov

Figure 5. Spin density on a H(4) nucleus of benzyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

TABLE 2: Spin Densities in a Pyrolyl Radical Obtained with B3LYP Functional and Different Basis Sets N

H(3,4)

H(2,5)

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

aug-cc-pVTZ

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

-0.02135 -0.02411 -0.02493 -0.02580 -0.001960 -0.002123 -0.002107 -0.002090 -0.007988 -0.008640 -0.008568 -0.008485

-0.02561 -0.02557 -0.02564 -0.02627 -0.002019 -0.002066 -0.002074 -0.002074 -0.008164 -0.008347 -0.008387 -0.008394

-0.02189 -0.02469 -0.02550 -0.02637 -0.001983 -0.002133 -0.002111 -0.002088 -0.008102 -0.008691 -0.008582 -0.008464

-0.02613 -0.02610 -0.02617 -0.02680 -0.002018 -0.002065 -0.002071 -0.002068 -0.008161 -0.008343 -0.008369 -0.008354

-0.01412 -0.01652 -0.01732 -0.01815 -0.001882 -0.002046 -0.002036 -0.002026 -0.007660 -0.008335 -0.008293 -0.008250

expta -0.0252

-0.00223

-0.008314

Reference 18.

The aug-cc-pVTZ basis set yields hydrogen HFCCs comparable with 6-311++G** in accuracy. However, it gives much less satisfactory HFCCs for C(1) and reveals its deficiency in comparison with 6-311++G** for HFCCs calculations. B. 1-Pyrolyl. The radical formed by reaction of OH with pyrrole has been studied in aqueous solution. Experimental data for 1H and 14N are available.18 This radical is neutral. The calculated dipole moment at the UB3LYP/6-311++G** level of theory is 2.1847 D. Because of the strong dipole-dipole interaction between solute and solvent, the calculation based on SS(V)P model19–21 has been carried out to evaluate the solvent effects on HFCCs. The calculation, however, shows that the HFCC shifts due to solvent effects are negligible. The results are shown in Table 2. Starting the analysis with the delta-function operator, the dependence of the calculated spin densities on the basis sets used is shown in Table 2. UB3LYP/6-311++G** yields the spin density of -0.0214 for 14N, while UB3LYP/ext-6311++G** leads to a spin density of -0.0256, which is very close to the experimental value of -0.0252. The relative errors for the two levels of method are 15% and 1.6%, respectively.

The flexibility of the basis set is important. The expansion of polarization function space has little influence on 14N. Very close HFCCs are obtained with the expansion of basis sets from 6-311++G**/ext-6-311++G** to 6-311++G(2d,2p)/ext-6311++G(2d,2p). The extension of basis sets from 6-311++G** to ext-6311++G** yields slightly improved hydrogen HFCCs, especially for H(2,5). Further expansion the polarization functions of the basis sets from ext-6-311++G** to ext-6-311++G(2d,2p) produces almost identical HFCCs, suggesting that addition of multiple polarization functions is not necessary. The evaluation of the performance of the RC operator follows the same trend as in the benzyl radical. The dependence of HFCCs on r0 is shown in Figures 6, 7, and 8. We observe very similar trends in HFCC dependence on r0, and similar conclusions can be obtained for 14N. The same optimum range parameters, which are 0.25 or 0.35 for extended basis sets and 0.45 for the standard basis sets, are applicable to 14N. The comparison of the performance between the delta-function operator method and the RC-operator method also lead to the

Spin Density Calculations of Aromatic Radicals

J. Phys. Chem. C, Vol. 114, No. 48, 2010 20653

Figure 6. Spin density on nitrogen nucleus of pyrolyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

Figure 7. Spin density on H(2) and H(5) nuclei of pyrolyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

conclusion that the two methods yield very close and accurate spin densities with the extended basis sets ext-6-311++G** and ext-6-311++G(2d,2p). With the standard basis sets 6-311++G(2d,2p) and 6-311++G(2d,2p), the RC-operator method is systematically better than the delta-function operator method, especially for the heavy atoms. The 6-311++G** basis set yields fairly accurate spin densities even when the delta-function method is applied. The aug-cc-pVTZ basis set gives HFCCs with larger errors, espe-

cially for 14N. This observation confirms again the advantage of the 6-311++G** basis set over the aug-cc-pVTZ basis set for spin density calculations. C. Fulvalene Anion. The fulvalene radical anion has been prepared in tetrahydrofuran (THF) by the oxidation of cyclopentadienyllithium or dilithium fulvalenediide with oxygen, and by the reduction of fulvalene with sodium, or by electrolysis.22 Experimental data for all nuclei are available. THF is a moderately polar aprotic solvent with a dielectric constant of

20654

J. Phys. Chem. C, Vol. 114, No. 48, 2010

Liang and Rassolov

Figure 8. Spin density on H(3) and H(4) nuclei of pyrolyl radical as a function of Zr0, computed with the RC operator with the B3LYP functional.

TABLE 3: Spin Densities in Fulvalene in Vacuum Obtained with B3LYP Functional and Different Basis Sets C(1,1′)

C(2,2′,5,5′)

C(3,3′,4,4′)

H(2,2′,5,5′)

H(3,3′,4,4′)

rms % error

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

aug-cc-pVTZ

expta

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

0.002519 0.003950 0.004545 0.005159 -0.004136 -0.003906 -0.003712 -0.003518 0.001868 0.002802 0.003187 0.003584 -0.0009952 -0.001082 -0.001075 -0.001067 -0.002011 -0.002173 -0.002154 -0.002131 0.4246 0.3048 0.2540 0.2040

0.005394 0.005065 0.005083 0.005507 -0.003359 -0.003561 -0.003550 -0.003416 0.003730 0.003524 0.003537 0.003812 -0.001034 -0.001058 -0.001062 -0.001061 -0.002055 -0.002101 -0.002110 -0.002109 0.1873 0.2117 0.2102 0.1774

0.003685 0.005220 0.005830 0.006455 -0.003811 -0.003562 -0.003368 -0.003174 0.002412 0.003406 0.003802 0.004208 -0.001004 -0.001080 -0.001068 -0.001056 -0.002078 -0.002224 -0.002193 -0.002158 0.3357 0.2120 0.1649 0.1251

0.006719 0.006389 0.006395 0.006822 -0.003011 -0.003209 -0.003199 -0.003065 0.004371 0.004160 0.004167 0.004445 -0.001021 -0.001044 -0.001047 -0.001045 -0.002083 -0.002129 -0.002134 -0.002127 0.1190 0.1284 0.1282 0.1093

-0.0003142 0.001118 0.001701 0.002297 -0.004663 -0.004306 -0.004101 -0.003897 -0.0004880 0.000530 0.000941 0.001354 -0.0009376 -0.001022 -0.001018 -0.001013 -0.001966 -0.002135 -0.002121 -0.002106

0.00723

-0.00349

0.00536

-0.000972

-0.00232

Reference 22.

7.6. The SS(V)P model19–21 is used to evaluate the solvent effects on HFCCs. The results are listed in Tables 3 and 4. Although the dielectric constant of the solvent is small, the SS(V)P model calculations show that the solvent effects on spin densities are not negligible. It is encouraging to notice that the SS(V)P corrections lead to systematically better spin densities as demonstrated by the root-mean-square deviations. It is also clear that the solvent effects do not change the relative performance of delta-function and RC operators. The deltafunction operator method and the RC-operator method give close

spin densities, if extended basis sets are used. For standard basis sets, the delta-function operator method gives spin densities with large errors for all of the three carbon atoms, while the RCoperator method gives spin densities close to extended basis sets results with much smaller errors. In this radical, the expansion of the polarization function from 6-311++G**/ext-6-311++G** to 6-311++G(2d,2p)/ext-6311++G(2d,2p) is important. The latter produces spin densities significantly closer to experimental values. This fact reflects the diffuse nature of the wave function for anions.

Spin Density Calculations of Aromatic Radicals

J. Phys. Chem. C, Vol. 114, No. 48, 2010 20655

TABLE 4: Spin Densities in Fulvalene in THF Obtained with U3LYP Functional and Different Basis Sets C(1,1′)

C(2,2′,5,5′)

C(3,3′,4,4′)

H(2,2′,5,5′)

H(3,3′,4,4′)

rms % error

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

0.002710 0.004160 0.004761 0.005380 -0.004502 -0.004324 -0.004146 -0.003969 0.002068 0.003041 0.003439 0.003850 -0.0009464 -0.001029 -0.001023 -0.001015 -0.002094 -0.002264 -0.002245 -0.002221 0.4153 0.2928 0.2384 0.1830

0.005613 0.005285 0.005303 0.005731 -0.003811 -0.004011 -0.003999 -0.003876 0.004001 0.003789 0.003802 0.004087 -0.0009842 -0.001007 -0.001011 -0.001010 -0.002142 -0.002190 -0.002199 -0.002199 0.1606 0.1924 0.1903 0.1522

0.003892 0.005449 0.006064 0.006695 -0.004183 -0.003987 -0.003809 -0.003632 0.002640 0.003676 0.004085 0.004506 -0.0009554 -0.001028 -0.001017 -0.001006 -0.002157 -0.002311 -0.002279 -0.002244 0.3210 0.1913 0.1367 0.0834

0.006957 0.006628 0.006635 0.007066 -0.003470 -0.003665 -0.003655 -0.003533 0.004674 0.004458 0.004464 0.004752 -0.0009729 -0.0009954 -0.0009980 -0.0009962 -0.002167 -0.002214 -0.002220 -0.002214 0.0666 0.0899 0.0889 0.0570

aug-cc-pVTZ

expta 0.00723

-0.00349

0.00536

-0.000972

-0.00232

Reference 22.

TABLE 5: Spin Densities p-Benzosemiquinone Anion Radical in Vacuum, Computed with B3LYP Functional with Different Basis Sets O(1,4)

C(1,4)

C(2,3,5,6)

H(2,3,5,6)

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

0.02805 0.03542 0.03836 0.04144 -0.009390 -0.009246 -0.008978 -0.008716 -0.001780 -0.001236 -0.0009516 -0.0006621 -0.001386 -0.001502 -0.001490 -0.001476

0.04053 0.03960 0.04022 0.04263 -0.008441 -0.008793 -0.008774 -0.008594 -0.0004927 -0.0007083 -0.0006970 -0.0004972 -0.001424 -0.001457 -0.001463 -0.001463

0.02936 0.03679 0.03971 0.04279 -0.009069 -0.008915 -0.008648 -0.008386 -0.001430 -0.0008505 -0.0005596 -0.0002651 -0.001417 -0.001519 -0.001499 -0.001475

0.04186 0.04094 0.04156 0.04395 -0.008111 -0.008454 -0.008434 -0.008252 -0.0000841 -0.0003021 -0.0002948 -0.0000935 -0.001425 -0.001456 -0.001460 -0.001455

expt(DMSO)a

expt(H2O)a

-0.0437

0.0402

-0.00531

0.000598

-0.000299

-0.00177

-0.00152

-0.00148

References 32–34.

Although a range parameter of 0.45 yields spin densities with smaller root-mean-square deviation, an optimum range parameter of 0.25 or 0.35 is still appropriate for extended basis functions ext-6-311++G** and ext-6-311++G(2d,2p). Both range parameters of 0.25 or 0.35 located in the plateau region produce very close spin densities. D. p-Benzosemiquinone Anion. Semiquinone anion radicals play important roles in the electron transfer reactions of photosynthesis and respiration. Among these radicals, the p-benzosemiquinone anion has been extensively studied by ESR spectroscopy and detailed results have been documented for the 1 H, 13C, and 17O HFCCs in a variety of solvents. In addition, the influence of hydrogen bonding on the spin distribution can be observed. The quinone carbonyl oxygen atoms act as hydrogen bond acceptors. Hydrogen-bonding effects on the p-benzosemiquinone anion HFCCs have been investigated in several static DFT studies23–26 as well as the computationally demanding Car-Parrinello ab initio molecular dynamics methodology27 based on the delta-function operator method. Using

optimized structures, the hydrogen bonding interaction is described by a supermolecular complex with two water molecules coordinated at each carbonyl oxygen atom in such a manner as to maintain the same D2h symmetry as the nonhydrogen-bonded form (Figure 2).24,26 The calculations give satisfactory agreement with aqueous experimental data. In this paper, we adopt the same hydrogen bonding model and investigate the performance of RC-operator-based method. The results are listed in Tables 5 and 6. The delta-function approach shows strong basis set dependence both in the hydrogen-bonded and in non-hydrogen-bonded forms. For heavy atoms, the difference between the standard basis sets and the extended basis sets is significant. The extension of basis sets is required to produce more reliable HFCCs. The inspection of the dependence of HFCCs on r0 for the RC-operator method, both in the hydrogen-bonded and nonhydrogen-bonded forms, reveals that the plateau regions can be observed for all the nuclei if extended basis sets are used.

20656

J. Phys. Chem. C, Vol. 114, No. 48, 2010

Liang and Rassolov

TABLE 6: Spin Densities p-Benzosemiquinone Anion Radical Hydrogen Bonded to Four Water Molecules, Computed with B3LYP Functional with Different Basis Sets O(1,4)

C(1,4)

C(2,3,5,6)

H(2,3,5,6)

a

Zr0

6-311++G**

ext-6-311++G**

6-311++G(2d,2p)

ext-6-311++G(2d,2p)

0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45 0.0 0.25 0.35 0.45

0.02766 0.03439 0.03703 0.03982 -0.005067 -0.004331 -0.003889 -0.003442 -0.002678 -0.002236 -0.001978 -0.001718 -0.001325 -0.001437 -0.001427 -0.001414

0.03887 0.03809 0.03867 0.04084 -0.003160 -0.003541 -0.003522 -0.003214 -0.001541 -0.001758 -0.001747 -0.001568 -0.001366 -0.001398 -0.001403 -0.001403

0.02927 0.03611 0.03876 0.04156 -0.004730 -0.003985 -0.003545 -0.003101 -0.002304 -0.001826 -0.001563 -0.001298 -0.001356 -0.001455 -0.001437 -0.001416

0.04063 0.03987 0.04043 0.04261 -0.002811 -0.003183 -0.003165 -0.002859 -0.001113 -0.001331 -0.001324 -0.001143 -0.001369 -0.001400 -0.001404 -0.001399

expt(DMSO)a

expt(H2O)a

-0.0437

0.0402

-0.00531

0.000598

-0.000299

-0.00177

-0.00152

-0.00148

References 32–34.

Figure 9. Spin density on oxygen nuclei of p-benzosemiquinone anion radical, computed with the RC operator with the B3LYP functional, as a function of Zr0.

The dependence of 17O HFCCs on r0 is shown in Figures 9 and 10. The 17O has approximately the same plateau region as 13C. Same optimum range parameter 0.25-0.35 can be chosen for 17 O. The standard basis sets can produce close spin densities comparable with the standard basis sets results if the optimum range parameter 0.45 is chosen. Usually for extended basis sets, the delta-function operator and the RC-operator methods give close results. An anomaly can be observed for C(2,3,5,6) in the gas phase calculation. This may be caused by the numerical cancellation error associated with the delta-function method. The HFCCs of the two carbons are strongly influenced by the hydrogen bonding interaction. 1H and 17O HFCCs are less influenced. The combination of RC operator with standard basis sets provides a satisfactory description of the phenomenon. The expansion of the polarization function from 6-311++G**/ ext-6-311++G** to 6-311++G(2d,2p)/ext-6-311++G(2d,2p)

makes a difference for heavy atoms, especially for the two carbons. The 6-311++G(2d,2p)/ext-6-311++G(2d,2p) basis set gives better results for O(1,4) and C(1,4), while the 6-311++G**/ ext-6-311++G** basis set produces better C(2,3,5,6) spin density. E. PBE0 Functional. Recent work of Barone et al. demonstrates28 that PBE0 hybrid functional yields comparable accuracy to the B3LYP model29 in hyperfine structure calculations. The PBE0 functional30 uses 25% of the exact exchange with Perdew, Burke, and Ernzerhof’s generalized gradient functional.31 At the suggestion of a referee we have investigated the performance of the RC operator with PBE0 functional. The selected data are shown in Tables 7, 8, 9, and 10. The basis sets with single set of polarization functions (**) are used for neutral radicals, and split-polarized basis sets (2d,2p) for anions; all calculations use spin-unrestricted wave functions.

Spin Density Calculations of Aromatic Radicals

J. Phys. Chem. C, Vol. 114, No. 48, 2010 20657

Figure 10. Spin density on oxygen nuclei of p-benzosemiquinone anion radical hydrogen bonded to four water molecules, computed with the RC operator with the B3LYP functional, as a function of Zr0.

TABLE 7: Spin Densities of a Benzyl Radical Obtained with B3LYP and PBE0 Functionals B3LYP/6-311++G** PBE0/6-311++G** B3LYP/6-311++G** PBE0/6-311++G** expta a

Zr0

C(1)

C(7)

H(2,6)

H(3,5)

H(4)

H(7)

0.0 0.0 0.45 0.45

-0.03514 -0.04139 -0.03870 -0.04568 -0.03603

0.04488 0.04813 0.05658 0.06044 0.06096

-0.003315 -0.003835 -0.003530 -0.004121 -0.00323

0.001392 0.001862 0.001479 0.001993 0.00112

-0.003832 -0.004372 -0.004071 -0.004681 -0.00387

-0.01006 -0.01139 -0.01063 -0.01214 -0.0102

Reference 16.

TABLE 8: Spin Densities of a Pyrolyl Radical Obtained with B3LYP and PBE0 Functionals B3LYP/6-311++G** PBE0/6-311++G** B3LYP/6-311++G** PBE0/6-311++G** expta a

Zr0

N

H(3,4)

H(2,5)

0.0 0.0 0.45 0.45

-0.02135 -0.02582 -0.02580 -0.03101 -0.0252

-0.001960 -0.002095 -0.002090 -0.002255 -0.00223

-0.007988 -0.008941 -0.008485 -0.009577 -0.008314

the recommended range parameter value Zr0 ) 0.45 drops to 0.0016 au for B3LYP, and 0.0031 au for PBE0. Use of smaller values of the range parameter, but the extended core basis sets as described in section II gives results nearly identical to these. Rms deviations for B3LYP are 0.0016 for Zr0 ) 0.25 and Zr0 ) 0.35, and for PBE0 they are 0.0028 and 0.0030 for Zr0 ) 0.25 and Zr0 ) 0.35, respectively.

Reference 18.

IV. Conclusions

The PBE0 functional is slightly superior to B3LYP with deltafunction based HFCCs. The rms deviation from experiment across all 18 data points using standard 6-311++G basis sets with proper polarization functions described above is 0.0050 au for B3LYP and 0.0043 fro PBE0. The agreement with experiment improves with the use of RC operators, albeit less dramatically for the PBE0 functional. The rms deviation for

We have examined the performance of the Gaussian-weighted RC operator combined with unrestricted B3LYP and PBE0 density functionals in the calculation of spin densities of larger aromatic radicals. Based on the standard Pople-type basis sets with triple-ζ quality, the extended basis sets are constructed. The effect of extension of standard Pople-type basis sets on the calculated spin densities is investigated. The HFCCs of hydrogen are less influenced by the extension of the standard

TABLE 9: Spin Densities of Fulvalene Anion in Thf Obtained with B3LYP and PBE0 Functionals Zr0 B3LYP/6-311++G(2d,2p) PBE0/6-311++G(2d,2p) B3LYP/6-311++G(2d,2p) PBE0/6-311++G(2d,2p) expta a

Reference 22.

0.0 0.0 0.45 0.45

13

C(1,1′)

0.003892 0.004125 0.006695 0.006935 0.00723

13

C(2,2′,5,5′)

-0.004183 -0.004920 -0.003632 -0.004472 -0.00349

13

C(3,3′,4,4′)

H(2,2′,5,5′)

H(3,3′,4,4′)

0.002640 0.002459 0.004506 0.004332 0.00536

-0.0009554 -0.0009573 -0.001006 -0.001020 -0.000972

-0.002157 -0.002396 -0.002244 -0.002512 -0.00232

20658

J. Phys. Chem. C, Vol. 114, No. 48, 2010

Liang and Rassolov

TABLE 10: Spin Densities of Benzosemiquinone Anion with Four Explicit Water Molecules Obtained with B3LYP and PBE0 Functionals Zr0 B3LYP/6-311++G(2d,2p) PBE0/6-311++G(2d,2p) B3LYP/6-311++G(2d,2p) PBE0/6-311++G(2d,2p) expt(H2O)a a

0.0 0.0 0.45 0.45

17

O(1,4)

0.02927 0.03170 0.04156 0.04402 0.0402

13

C(1,4)

-0.004730 -0.005675 -0.003101 -0.004147 0.000598

13

C(2,3,5,6)

H(2,3,5,6)

-0.002304 -0.002767 -0.001298 -0.001835 -0.00177

-0.001356 -0.001410 -0.001416 -0.001486 -0.00148

References 32–34.

basis set. However, the extended basis sets yield significantly improved HFCCs for heavy atoms. It is found that for the extended basis sets, which are of high flexibility near the nuclear region, the RC operator and the delta-function operator provide close spin densities. For neutral radicals, the UB3LYP/ext-6311++G** level of method gives reliable spin densities. The UB3LYP/ext-6-311++G(2d,2p) method generally does not provide significant overall improvement. For anions, the UB3LYP/ ext-6-311++G(2d,2p) method yields quite accurate spin densities while the UB3LYP/ext-6-311++G** method is less satisfactory. Different atoms display very similar HFCC-r0 dependence behavior with an approximately same plateau region. Therefore a common optimum range parameter of 0.25 or 0.35 can be chosen and recommended. The delta-function operator may yield large errors, especially for the heavy atoms, when used with standard 6-311++G (2d,2p) and 6-311++G(2d,2p) basis sets. The RC operator, on the other hand, still can yield accurate results comparable to extended basis set results. This is especially significant for the calculation of large organic radicals. The variation of HFCCs with respect to the range parameter can be regarded as a measure of the completeness of the basis sets for the calculation of HFCCs. The optimum choice of range parameter can help to recover the basis set error. An optimum range parameter of 0.45 is recommended with these bases. The UB3LYP/6-311++G** method can provide satisfactory spin densities for neutral radicals, while the UB3LYP/6-311++G(2d,2p) method is more accurate when applied to anion radicals. We have demonstrated that the standard Pople-type basis sets 6-311++G(2d,2p) and 6-311++G(2d,2p), which are designed for conventional electronic structure calculations, combined with the RC-operator formalism and the UB3LYP functional serves as an efficient and reliable method for the calculation of spin densities of large organic radicals. The Dunning’s basis set augcc-pVTZ leads to large errors for heavy atoms. The deltafunction operator is unreliable when used with these bases. The PBE0 functional gives slightly better overall accuracy than B3LYP with the delta-function operator and standard Popletype basis sets. The RC operator offers some improvement over the delta-function, but this improvement is not as dramatic as for B3LYP. The overall r0 dependence of the HFCC is qualitatively similar for two functionals. Solvation or hydrogen bonding may exert a strong influence on spin density distribution. The calculation to include such effects requires more demanding computational resources. We

have shown that RC-operator-based method is superior to the delta-function operator method. Acknowledgment. Acknowledgment is made to the Donors of the American Chemical Society Petroleum Research Fund for partial support of this research. References and Notes (1) Carrington, A. McLachlan, A. D. Introduction to Magnetic Resonance; Harper and Row: New York, 1967. (2) Chipman, D. M. In Quantum Mechanical Electronic Structure Calculations with Chemical Accuracy; Langhoff, S. R., Ed.; Kluwer: Dordrecht, The Netherlands, 1995. (3) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (4) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. B 1988, 37, 785. (5) Cohen, M. J.; Chong, D. P. Chem. Phys. Lett. 1995, 234, 405. (6) Adamo, C.; Barone, V.; Fortunelli, A. J. Chem. Phys. 1995, 102, 384. (7) Hiller, J.; Sucher, J.; Feinberg, G. Phys. ReV. A 1978, 18, 2399. (8) Larsson, S. Chem. Phys. Lett. 1981, 77, 176. (9) Rassolov, V. A.; Chipman, D. M. Theor. Chim. Acta 1995, 91, 1. (10) Rassolov, V. A.; Chipman, D. M. J. Chem. Phys. 1996, 105, 1470. (11) Rassolov, V. A.; Chipman, D. M. J. Chem. Phys. 1996, 105, 1479. (12) Hirschfelder, J. O. J. Chem. Phys. 1960, 33, 1462. (13) Wang, B.; Baker, J.; Pulay, P. Phys. Chem. Chem. Phys. 2000, 2, 2131. (14) Chipman, D. M. J. Chem. Phys. 1989, 91, 5455. (15) Shao, Y.; Molnar, L.; Jung, Y.; Kussmann, J.; Ochsenfeld, C.; Brown, S.; Gilbert, A.; Slipchenko, L.; Levchenko, S.; O’Neill, D.; et al. Phys. Chem. Chem. Phys. 2006, 8, 3172. (16) Fischer, H. Z. Naturforsch. 1965, 20, 488. (17) Adamo, C.; Subra, R.; Matteo, A. D.; Barone, V. J. Chem. Phys. 1998, 109, 10244. (18) Samunl, A.; Neta, P. J. Phys. Chem. 1973, 77, 1629. (19) Chipman, D. M. J. Chem. Phys. 2000, 112, 5558. (20) Chipman, D. M. Theor. Chim. Acta 2000, 107, 80. (21) Chipman, D. M.; Dupuis, M. Theor. Chim. Acta 2002, 107, 90. (22) Davies, A. G.; Giles, J. R. M.; Lusztyk, J. J. Chem. Soc., Perkin Trans. 1981, 2, 747. (23) O’Malley, P. J. J. Phys. Chem. A 1997, 101, 6334. (24) O’Malley, P. J. J. Phys. Chem. A 1997, 101, 9813. (25) Eriksson, L. A.; Himo, F.; Siegbahn, P. E. M.; Babcock, G. T. J. Phys. Chem. A 1997, 101, 9496. (26) Chipman, D. M. J. Phys. Chem. A 2000, 104, 11816. (27) Asher, J. R.; Kaupp, M. ChemPhysChem 2007, 8, 69. (28) Barone, V.; Cimino, P. Chem. Phys. Lett. 2008, 454, 139. (29) Barone, V.; Cimino, P.; Stendardo, E. J. Chem. Theor. Comput. 2008, 4, 751. (30) Adamo, C.; Barone, V. J. Chem. Phys. 1999, 110, 6158. (31) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (32) Stone, A. W.; Maki, A. H. J. Am. Chem. Soc. 1965, 87, 454. (33) Das, M. R.; Fraenkel, G. K. J. Chem. Phys. 1965, 42, 1350. (34) Gulick, W. M., Jr.; Geske, D. H. J. Am. Chem. Soc. 1966, 88, 4119.

JP105511A