The Calculation of NMR Parameters by Density-Functional Theory

Chapter 23

The Calculation of N M R Parameters by Density-Functional Theory An Approach Based on Gauge Including Atomic Orbitals Georg Schreckenbach, Ross M . Dickson, Yosadara Ruiz-Morales, and Tom Ziegler

Downloaded by UNIV ILLINOIS URBANA on March 17, 2013 | http://pubs.acs.org Publication Date: May 5, 1996 | doi: 10.1021/bk-1996-0629.ch023

1

Department of Chemistry, University of Calgary, 2500 University Drive Northwest, Calgary, Alberta T2N 1N4, Canada

Schemes for calculating nuclear magnetic resonance (NMR) shielding tensors and spin-spin coupling constants have been implemented in the Amsterdam Density Functional program system (ADF). The shielding tensors are calculated by the Gauge Including Atomic Orbitals (GIAO). This method and the calculation of the coupling constants are tested for a number of smaller molecules and it is shown that calculations of couplings to transition-metal nuclei and shielding tensors in metal complexes are feasible. 1. Introduction. Nuclear magnetic resonance (NMR) is used extensively (1 ) as a practical tool in chemical research. Many of its applications can be carried out based on a simple effective Hamiltonian in which the observed shifts and spin-spin coupling constants are used as parameters without any further interpretation. However, an understanding of how electronic and geometrical effects influence these parameters has not been established in detail except for a few classes of compounds (la-c), although such an understanding might enhance the amount of useful information obtainedfromNMR experiments. Computational methods based on molecular orbital theory can in principle provide the required insight (la-c), and the comparison between calculated and observed NMR spectra might further help in the identification of new species. With this in mind, several first principle methods capable of calculating NMR parameters have appeared over the last decade (la-c). Density functional theory (DFT) (2) forms the basis for some of the approaches used in computational studies of the shielding tensor σ. Recent advances in DFT have made it possible to use this approach for shielding calculations. Malkin et al. have published a series of pioneering papers on the calculation of NMR properties, including shielding (3a-g) and spin-spin coupling (3h). To calculate the shielding, they combine modern DFT with the "individual gauge for localized orbitals" (IGLO) method (4). We have recently presented a method in which the NMR shielding tensor is calculated by combining the "Gauge Including Atomic Orbital" (GIAO) (5) approach 1

Corresponding author 0097-6156/96/0629-0328$15.00/0 © 1996 American Chemical Society In Chemical Applications of Density-Functional Theory; Laird, B., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1996.

23. SCHRECKENBACH ET AL.

Calculation

of NMR

Parameters by DFT

329

with DFT. Our implementation makes full use of the modem features of DFT in terms of accurate exchange-correlation (XC) energy functionals and large basis sets


(6).

Of particular interest - and still a challenge -- are applications in multi-nuclear NMR (lc, d), e.g., transition metal NMR. We present here results for carbonyl complexes M(CO)6 (M= Cr, Mo, W). In the present paper we address as well the calculation of nuclear spin-spin coupling constants. Spin-spin coupling constants are very difficult to compute. In the non-relativistic formulation of Ramsey (7) mere are four terms, each of which has different requirements with respect to correlation and basis set. Malkin et al. (3a,g) provided the first practical implementation of these calculations in a density functional program, using a basis of Gaussian-type basis functions. In this paper we illustrate an implementation using a basis of Slater-type rather than Gaussian functions, and examine the feasibility of calculating spin-spin coupling constants in transition-metal systems. 2. The Calculation of Shielding Tensors.

Shielding tensors and the G I A O - D F T method. The details of the GIAO-DFT method have already been described previously (6a,b). However, we will have to stress a few points about NMR in general and the GIAO formalism in particular to facilitate the discussion in the next sections. In NMR one considers the interaction energy between a nuclear magnetic moment fi in an electronic system and an external homogeneous magnetic field N

B. 0

The presence of B will induce an internal magnetic field B^ -Β^Λ-Β 0

ni

ρ

in the

electronic system so that the total interaction energy is given by E = -fi (B +B N

d

+ B)

p

0

= - / 2 ( l + d)-B N

(1)

0

= -ji -B s d

0

Βρ = - σ

Here σ = σ

ρ

+¥

ρ

Β

(2)

0

are referred to as the shielding tensors; one third of the trace of σ

is the shielding constant. The vector B represents the diamagnetic component of the d

induced field. It is in most cases opposite to B with B -B 0 for any 0

p

o

p

orientation of B . In this case the paramagnetic shielding tensor σ must have a negative trace and negative symmetrical diagonal components according to eq. 2. ρ

0

The diamagnetic shielding ¥ depends only on the unperturbed electron density, p°, while the paramagnetic shielding, σ , respect to die external magnetic field. ρ

contains the density up to first order with

Then, the st-component of the diamagnetic tensor ¥ formalism as (6)

is given in our GIAO

In Chemical Applications of Density-Functional Theory; Laird, B., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1996.

330

CHEMICAL APPLICATIONS OF DENSITY-FUNCTIONAL THEORY

of, = α Σ J ^ p ^ ^ r ^ J V - M a i - ' i v / ^ ^ v 2

(3) λ,ν Here the occupied molecular orbital (MO) Ψ,· has been expanded in terms of atomic functions χ

and the coefficients c^-. The zero superscript indicates that Ψ,·

ν

is calculated with zero magnetic field strength, B = 0, as an eigenfunction to the Downloaded by UNIV ILLINOIS URBANA on March 17, 2013 | http://pubs.acs.org Publication Date: May 5, 1996 | doi: 10.1021/bk-1996-0629.ch023

Q

Kohn-Sham operator h(o) with the eigenvalue ef.

The vector r = r -R v

v

denotes

the position of the electron relative to the nucleus at which the atomic orbital χ

ν

centered, and r

N

is

is the position of the electron relative to the nuclear magnetic

moment under consideration (relative to the NMR nucleus N). Further α is the dimensionless fine structure constant given as 1/137. The diamagnetic shielding tensor â of eq. 3 is gauge invariant and an expectation value of Hermitian operators. It depends only on the unperturbed occupied orbitals for which nj * 0; this has been pointed out before. é

The st-component of the paramagnetic tensor σ formalism as

ρ

is given according to the GIAO

occ Σ 4 Α ( ζ λ,ν

Γ

λ

unocc

-^χν

(4)

-it • «

.

AU) ε

_J0) {U) s

(°)- (°)

with i ^ a

(5),

ε

Χν

(6),

and xV

Χν)+{Χλ

| ( ^ v - ^ ) ] A(0)b x

(7),

The index J runs over orbitals occupied in the field free ground state and the index a runs over the corresponding unoccupied orbitals. The paramagnetic shielding tensor σ of eq. 4 is also gauge invariant by itself and an expectation value of Hermitian operators. The leading contribution to the paramagnetic shielding is the last term in ρ


SCHRECKENBACHETAL.

23.

Calculation of NMR

Parameters by DFT

331

eq. 4. It represents the first order magnetic coupling between an occupied molecular orbital, z, and a virtual orbital, a. This coupling is facilitated by way of the first order coefficient u^ \ which is inversely proportional to the difference of the eigenvalues, eq. 5. It is worth noting that the shielding is completely formulated in terms of the (occupied and virtual) MO's, eqs. 3-7. This is an advantage of the GIAO formalism as it allows for the detailed analysis of the different orbital contributions to the shielding tensors.


s

Computational details. The above established formalism has been implemented into the Amsterdam Density Functional Package (ADF) (8). A l l properties are evaluated using the given numerical integration scheme of ADF. We will use non-local X C energy functionals (9), unless otherwise stated. We employ uncontracted Slater Type Orbitals as basis functions. Our basis sets are generally of triple ζ quality in the valence region, and of double ζ quality for core MO's. They are augmented by two sets of (p or d) polarization functions per atomic center. The experimental geometries are the basis for all calculations. G I A O - D F T calculation of shielding constants for simple molecules. Calculated

shielding constants for a representative set of small molecules are collected in Table I.

T A B L E I. Calculated and Experimental Shielding Constants for a Number of Small Molecules

^ ""*"*" ^^^^^^^^^ Tsôm>pî^

M

MM

a

M

Molecule

Atom

DFT-GIAOa"

DFT-IGLO^

Experiment

CH4

C H

191.2 31.4

187.7 31.2

195.1 30.6

C F H

111.4 462.3 27.2

101.4 450.7 26.7

116.8 471.6 26.6

H2O

0 H

331.5 31.2

324.3 31.1

344.0 30.1

N2

Ν

-72.9

-78.9

-61.6

C2H2

C Η

110.4 30.4

108.9 30.0

117.2 29.3

Benzene

C

50.0

48.8C

57.2

H2CO

C 0 H

-15.7 -418.8 20.7

-26.6 -455.6 20.8

-8.4 -312.1 18.3

F

-282.7

-250.6

-232.8

CH3F

F2

Reference 6a. Uncoupled DFT-IGLO: Ref. 3b. We cite the results for the same L D A / N L functional as in our DFT-GIAO method. Ref. 3e.

a

b

c

We compare our results with those obtained by the "uncoupled" DFT-IGLO of Malkin and co-workers (3a-c) as well as with experiment. The agreement with In Chemical Applications of Density-Functional Theory; Laird, B., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1996.

332


experiment is generally satisfactory. The GIAO shieldings are of about the same quality as those obtained by the IGLO method, in most cases even slightly better, Table I. A direct comparison is however not possible since this would require to use exacdy the same basis sets. The agreement with experiment is not as good for non-hydrogen shifts in molecules like H2CO or F2 (this applies to both the GIAO and the uncoupled IGLO methods). These molecules are difficult cases for DFT in general.


Calculated G I A O - D F T shielding anisotropics for simple molecules.

The

shielding anisotropy, i.e., the individual tensor components of the shielding tensor, should be even more sensitive to the quality of the quantum-chemical method used than the (averaged) shielding constant. Here, we define the shielding anisotropy Ac as the difference between the parallel and the orthogonal principle components Δσ = σι - σ

(8)

±

In Table I I , we compare shielding constants and anisotropics for another representative set of small molecules with experimental results. We can see that the quality of the averaged shielding constant and the tensor components is comparable for any given molecule. Thus, we get in general good agreement between theory and experiment. However, we note large deviations between the calculated and experimental shielding anisotropics for those molecules (notably CO out of the given list) where the calculated shielding constant isn't reliable either. Note that the magnitude of the anisotropy can exceed the shielding constant considerably, Table Π. Table IL Calculated and Experimental Shielding Constants and Anisotropics

Moîëcûîë^tôîî^

-

Isotropic Shielding DFT-GIAO* Experiment

Shielding Anisotropy DFT-GIAO Experiment 3

H2

H

26.46b

26.26±1.5

1.64b

2.0C

HF

F

412.5

410

104.2

108

NH3

Ν

262.0

264.5

-48.1

-40

CO2

C

56.1

58.8

345.9

335

HCN

C N

91.5 8.4

82.1 -20.4

286.1 502.9

284.6±20 563.818

CO

C Ο

-9.3 -68.4

1.0 -42.3

424.1 718.7

406 676.1

Reference 6a. h LDA. Calculated with the Coupled Hartree-Fock method c

XC Functionals. In Table III we look at the influence of the XC functional on the calculated shielding constants. By comparing results of first and second generation DFT (LDA and LDA/NL, respectively) with experiment, we note a remarkable improvement for the latter method for all non-hydrogen nuclei (up to 57 ppm change for H2CO). Malkin et al. had observed a similar strong influence for their DFTIGLO method (3a,b), Table ΠΙ. We have also included into Table m the results of the "coupled DFT-IGLO" approach of the same authors (3a,c,e). The idea of this approach is to model the current dependency of the XC potential (which is neglected in the "uncoupled" methods) by introducing a first order change into it. The authors In Chemical Applications of Density-Functional Theory; Laird, B., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1996.

23.

SCHRECKENBACHETAL.

Calculation of NMR

Parameters by DFT

333


do this in a somewhat ad hoc fashion. The influence of the new coupling is negligible in many cases. However, it leads to a significant reduction of the absolute value of the shielding and therefore to better agreement with experiment for some of the "difficult" molecules like F2 or H2CO, Table ΙΠ. The influence of the XC functional can be traced back to well-known changes in the differences of the orbital energies, eq. 5, between occupied and low lying virtual MO's (6a). The influence of these small changes is therefore less pronounced in singly bound systems with large HOMO-LUMO separations, Table ΠΙ. The term "uncoupled" DFT is somewhat misleading, since even this level of theory accounts for correlation effects, and only contributions of the induced current to the correlation are neglected. Uncoupled Hartree-Fock theory on the other hand excludes correlation effects completely. TaWj^n^^Th^Ii^ Molecule/ Isotropic Shieldings (ppm) Atom DFT-IGLO DFT-GIAOa Uncoupled * LDA LDA/NL LDA LDA/NLC 1

NH3

Exp. Coupled > LDA/NL c

d

Ν Η

267.2 31.2

262.0 31.6

259.3 30.8

253.4 31.2

253.3 31.2

264.5 32.43

HF

F H

415.1 29.4

412.5 30.0

412.7 29.1

409.0 29.7

409.6 29.6

410 28.7

N2

Ν

-83.2

-72.9

-86.7

-78.9

-69.3

-61.6

C2H2

C

102.9

110.4

102.5

108.9

109.6

117.2

-310.2

-293.7

-271.8

-250.6

-197.8

-232.8

-31.7 -475.8

-15.7 -418.8

-40.4 -504.7

-26.6 -455.6

-12.3 -362.6

-8.4 -312.1

F F2 H2CO C

O

Reference 6a. Ref. 3a. method. Ref. 3b.

a

b

c

We cite the results for the same LDA/NL functional as in our DFT-GIAO

d

Frozen core approximation. For heavier elements one might wonder to what degree it is possible to make use of the frozen core approximation in which orbitals of lower energy are taken from atomic calculations. The extension of the shielding calculations (eqs. 3 to 7) to include the frozen core approximation has been discussed in detail earlier (6b). Here, we address this question in connection with calculated 7?Se shieldings in Table IV. The core at the selenium atom contains the Is, 2s, and 2p shells while the core of the second period atoms carbon and fluorine contains the Is shell only. This choice was taken according to the discussion in reference 6b. The deviation in the calculated shielding between the frozen core results and the all electron calculations (numbers in brackets) is always smaller than 10 ppm, Table IV. We note also from Table IV that some of the deviation between the frozen core and all electron cases cancels when relative shifts are considered; the deviation does not exceed 5 ppm in this case. Let us now compare the calculated results to experimentally obtained values. The deviation between theory and experiment is considerable for all the ^^Se shieldings, Table IV. However, we get a much better agreement between theory and experiment when we consider relative shifts instead of absolute shieldings, Table IV. The In Chemical Applications of Density-Functional Theory; Laird, B., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1996.

334



experimental accepted standard for ??Se shifts is liquid Dimethyl Selenide, (CH3)2Se. We have therefore included this compound into our investigation. The agreement between theory and experiment is good for 7?Se shifts for the few compounds that have been considered here. The experimental uncertainty is certainly large, as are gas-to-liquid shifts and solvation effects. The former are as big as 119 ppm for H2Se. Solvent and counterion effects amount to a shift range of 25 ppm in the example of the cyclic tetra selenium dication, Table IV. Calculated shieldings refer of course to a single molecule at zero temperature whereas all the experimental data is obtained at finite temperatures and pressures; most of the experiments were carried out in solution or neat liquids. All these effects can yield considerable shifts, and make a direct comparison between theory and experiment difficult. On the theoretical side, it is likely that our basis sets are not yet completely saturated. Table IV. Calculated and Experimental S e Shieldings and Shifts for a few Molecules (Numbers in ppm) 7 7

molecule

SeF6

absolute shielding exp. calculated 2Ô09 1,666 (1,673) 2,401 2,093 (2,096) 1,438 988 (992)

CSe2

1,738

lCH3)2Se, (C2V) H2Se

Se4

chemical shift experiment calculated*

3

15

TT^

-345 (g) -226 (1) 610(g) 631 (1) 331 (g) 299 (sol) 1,923-1,958 (sol)

1,441 (1,448) -170

2+

Ô (6)

-427 (-423) 678 (681) 225 (225) 1,836

c

Calculated shieldings from frozen core calculations and (in brackets) from all electron calculations. - gas phase; 1 - liquid; sol - solution. Result depending on solvent and counterion.

a

c

A special case is the Se4^+ ion. This is a highly correlated molecule, and traditional Hartree-Fock based methods are unable to predict the chemical shift for this ion. However, the DFT result compares well with experiment; DFT is indeed capable of handling such systems. Let us now come back to the absolute shieldings. We note that the calculated absolute shieldings seem to be uniformly too small by about 300 ppm, Table IV. The experimental absolute shielding scale is based on the absolute shielding of SeF6 that was found to be 1,438164 ppm. However, this value is based on a theoretically predicted (diamagnetic) shielding value of the free selenium atom. This theoretical value has been corrected explicitly for relativistic effects, in particular the relativistic contraction of the core density. Other relativistic effects are probably not yet important for ??Se chemical shifts. The magnitude of the necessary correction is estimated at 300 ppm. Therefore, we find that our calculated shielding values are uniformly too small by about 300 ppm, Table IV. This uniform error of 300 ppm cancels of course when (relative) chemical shifts are calculated. This point illustrates the importance of absolute shielding scales for the test of theoretical methods. Transition Metal Complexes. The Example of C r ( C O ) , M o ( C O ) , and W ( C O ) . 6

6

6

Table V displays calculated (6c) and experimental absolute C NMR shielding tensor components for the three hexacarbonyls M(CO)6 (M=Cr, Mo, and W). A similar compilation is given in Table V I for the O shielding. The calculated and 13

n

observed tensor components a

(s = x, y, z) compare well for both C and 0 , with 13

ss

differences within the experimental error limits. The a

ss

1 7

tensor components in


23.

Calculation

SCHRECKENBACHETAL.

of NMR

Parameters by DFT

335

Tables V and V I correspond to an orientation in which the C O ligand probed by NMR has the 1 3 C 1 7 0 bond vector along the z-axis.

Table V .

A Comparison between Experimental and Calculated

Absolute

C

for the C O Molecule and G r o u p 6

Metal Carbonyls System

ο^ΡΡ™

Anisotropy Δ σ

(«q*f

ppm

OyyPP

(exptf

CO


1

Chemical Shielding Tensor Components

(expt)*

Isotropic shielding (expt)'

-149.4 8

-149.4

273.6

423.0

•8.4

(-132.3)

(273.4)

(406(8)11.4)

(1.0 (sr))

-167.9 8 (-167.6115 ) e

-167.9 (-167.6±15)

4323

-23.9

(423130)

(-26.6115(8,1))

-164.0 8

-164.0

266.4

430.4

-203

(-157.6115)

(260.4115)

(417130)

(-17.6115)

W(CO)6

-161.5 8

-161.5

272.1

433.6

-16.9

ReL

-149.6

-149*

266.7

416.4

-10.9

(-138.6115)

(256.4115)

(395130)

(-6.6115(8.1))

e

h

h

(-138.6*15°)

-20.4

264.2 (255.4115)

(-157.6115 )

MoiCOte

Other woric* (Absolute Scale")

0

σι ppm

(-1323 ) d

CXCOfc

01

(expt)'

-19.4

-6.2

The Calculation of NMR Parameters by Density-Functional Theory

Recommend Documents