Progress in Nuclear Magnetic Resonance Spectroscopy

(1)

The quantum-chemical calculation of NMR indirect spin–spin coupling constants

Trygve Helgaker

^a,^*

, Michał Jaszun´ski

^b

, Magdalena Pecul

^c

aCentre for Theoretical and Computational Chemistry, Department of Chemistry, University of Oslo, P.O. Box 1033 Blindern, N-0315 Oslo, Norway

bInstitute of Organic Chemistry, Polish Academy of Sciences, Kasprzaka 44, 01224 Warsaw, Poland

cDepartment of Chemistry, University of Warsaw, Pasteura 1, 02-093 Warsaw, Poland

a r t i c l e i n f o

Article history:

Received 7 January 2008 Accepted 19 February 2008 Available online 29 February 2008

Keywords:

Spin-spin coupling constants Electronic-structure theory Quantum chemistry Ab initio calculations Density-functional theory

Contents

1. Introduction . . . 250

2. Theory . . . 251

2.1. Ramsey’s expression: the mechanism of nuclear spin–spin coupling . . . 251

2.2. The wave-function approach to indirect spin–spin coupling constants. . . 252

2.2.1. CI and MCSCF theories . . . 252

2.2.2. Coupled-cluster theory. . . 253

2.2.3. SOPPA theory . . . 253

2.2.4. FCI theory . . . 253

2.2.5. Comparison of results obtained using different wave-function methods . . . 254

2.3. The DFT approach to indirect spin–spin coupling constants . . . 255

2.3.1. Approximate exchange–correlation functionals . . . 255

2.3.2. Computer implementations of DFT for indirect spin–spin coupling constants. . . 255

2.3.3. The performance of DFT for indirect spin–spin coupling constants . . . 256

2.4. One-electron basis sets for indirect spin–spin coupling constants . . . 258

2.5. Relativistic effects . . . 259

2.5.1. Four-component methods . . . 259

2.5.2. Two-component methods . . . 259

2.5.3. Perturbation theory . . . 259

2.5.4. Basis-set considerations . . . 260

2.6. Rovibrational and temperature effects. . . 260

2.7. Solvent effects . . . 260

2.8. Available models and programs . . . 261

3. Interpretations . . . 262

3.1. The Dirac vector model and the Karplus relation . . . 262

3.2. Three-dimensional pictures . . . 262

3.3. Orbital contributions. . . 263

doi:10.1016/j.pnmrs.2008.02.002

* Corresponding author.

E-mail address:[email protected](T. Helgaker).

Contents lists available atScienceDirect

Progress in Nuclear Magnetic Resonance Spectroscopy

j o u r n a l h o m e p a g e :w w w . e l s e v i e r . c o m / l o c a t e / p n m r s

(2)

4. Applications . . . 263

4.1. Configurational and conformational analysis . . . 263

4.1.1. Conformational studies . . . 263

4.1.2. Modeling of the effects of molecular motion . . . 264

4.2. Investigation of intermolecular interactions . . . 264

4.2.1. Through-space indirect nuclear spin–spin coupling constants . . . 264

4.3. Studies of inorganic and organometallic compounds . . . 265

5. Conclusions. . . 265

Acknowledgements . . . 266

References . . . 266

1. Introduction

We present here a review of currently available methods for the calculation of indirect nuclear spin–spin coupling constants, discussing some of the questions and difficulties connected with their calculation. The main purpose of our review is to serve as an intro- ductory guide to the literature for nonspecialists beginning their work in the field of computational NMR spectroscopy rather than to provide an overview of all existing literature. Advice and sugges- tions concerning the choice of appropriate tools for a given prob- lem are given, including the selection of wave function or density-functional treatment. We discuss also the selection of an atomic basis set, the treatment of relativistic effects, vibrational corrections and solvent effects. Methods for the interpretation and visualization of the calculated spin–spin coupling constants are also briefly examined. Finally, we illustrate some of the described theoretical approaches by discussing the results obtained by means of them, limiting ourselves to those we consider representative of the field of computational NMR spectroscopy or relevant to the analysis of the outlined theoretical methods.

Many reviews of the calculation of indirect nuclear spin–spin coupling constants have appeared over the years. In particular, theoretical aspects of such calculations are reviewed annually in the Specialist Periodical Reports on NMR, lately by Fukui (see Ref.

[1], for the most recent review). In the Annual Reports on NMR Spectroscopy, several reviews by Contreras et al. have appeared over the years[2–4], while the progress from May 2003 to July 2006 has more recently been reviewed by Krivdin and Contreras in the same journal, citing more than 450 articles, a large proportion of which were published during these few years[5].

Autschbach and Ziegler [6–9] have discussed the calculation of spin–spin coupling constants in systems containing heavy atoms, for which relativistic effects are important.

In 2003, a special issue of the International Journal of Molecular Sciences edited by Sauer was devoted to the calculation of nuclear spin–spin coupling constants[10]. In the following year, a collection of articles discussing theory and applications was published, edited by Kaupp et al.[11]. Other reviews, devoted to specific aspects of theoretical studies of spin–spin coupling constants[12], to calculations of NMR parameters[13–16], or to specific applica- Nomenclature

General acronyms

CASSCF complete-active-space SCF [model]

cc-pCVXZ correlation-consistent core–valenceX-tuple-zeta [basis set]

cc-pVXZ correlation-consistent valenceX-tuple-zeta [basis set]

CCSD coupled-cluster singles-and-doubles [model]

CCSDT coupled-cluster singles–doubles–triples [model]

CCSDTQ coupled-cluster singles–doubles–triples–quadruples [model]

CCSD(T) CCSD [model] with a noniterative perturbative triples correction

CC3 CCSD [model] with an iterative perturbative triples correction

CI configuration-interaction [model]

COSMO conductor-like screening model DFT density-functional theory DHF Dirac–Hartree–Fock [model]

DSO diamagnetic spin–orbit

EOM-CCSD equation-of-motion CCSD [model]

FC Fermi-contact

FCI full configuration-interaction [theory]

FPT finite-perturbation theory HF Hartree–Fock [model]

IEF-PCM integral equation formalism polarizable continuum model

IORAmm infinite-order regular approximation with modified metric

IPPP-CLOPPA inner projection of the polarization-propagator contributions from localized orbitals within the polarization-propagator approach

MCSCF multiconfigurational self-consistent field [model]

MD molecular dynamics MO molecular orbital

MP2 Møller–Plesset second-order [perturbation theory]

MPE multipole-expansion [method]

NBO natural bond orbital NMR nuclear magnetic resonance PSO paramagnetic spin–orbit

RASSCF restricted-active-space SCF [model]

RHF restricted Hartree–Fock [model]

SCF self-consistent field [theory]

SD spin–dipole

SOPPA second-order polarization-propagator approximation SOPPA(CCSD) SOPPA with CCSD correlation amplitudes ZORA zero-order regular approximation

ZPV zero-point vibrations

Acronyms for exchange–correlation functionals BLYP Becke–Lee–Yang–Parr

B3LYP Becke–3-parameter–Lee–Yang–Parr B97-1, B97-2, B97-3 Becke

CAMB3LYP Coulomb-attenuated-method B3LYP GGA generalized gradient approximation KT1, KT2 Keal–Tozer

LDA local-density approximation PBE, PBE0 Perdew–Burke–Ernzerhof PW91 Perdew–Wang 91

(3)

tions[17]have appeared during the last few years; we shall refer to some of them in what follows.

In this journal, the calculation of spin–spin coupling constants has been discussed on several occasions, from the review by Kowa- lewski more than 30 years ago[18]to more recent articles[19,20].

Likewise, the tensor properties of spin–spin coupling, as derived experimentally or theoretically[21], and the effects of solvation and weak interactions[22], have been discussed in this journal.

In the present review, we focus on recent work related to the development and analysis of theoretical calculations, paying less attention to the many recent articles devoted to the study of special spin–spin coupling constants of interest in specific chemical compounds, noting that these have been adequately reviewed elsewhere.

2. Theory

We present here only the essential features of the theory underlying the calculation of indirect nuclear spin–spin coupling constants, referring to textbooks for a detailed description of the various methods of quantum chemistry – see, for example, Ref.

[23] for wave-function methods and Ref. [24] for density-functional theory (DFT). A more detailed review of the theory of NMR parameters is given in Ref.[13].

In the theoretical description of molecular electronic structure, the nuclei are usually treated as point charges. To describe the NMR spectrum, we must also take into account the fact that, with each nuclear spinIPin a molecule, there is an associated magnetic moment

MP¼cP!hIP; ð1Þ

where cPis the gyromagnetic ratio. These nuclear magnetic moments interact with one another as well as with the electrons and with externally applied fields. As the changes in the electronic structure due to the external magnetic field and to the nuclear magnetic moments are very small, the NMR parameters can be adequately analyzed using perturbation theory. We shall not discuss here the nuclear shielding constants, which represent the interaction of these moments with an external magnetic field in the pres- ence of electrons. We neglect also the direct interaction between nuclear moments, described by a classical dipole mechanism, since this direct coupling is anisotropic and vanishes in isotropic media (gases and liquids). We are thus interested only in the small indirect contribution to the nuclear spin–spin coupling, which arises from hyperfine interactions with the intervening electrons and does not vanish in isotropic media.

2.1. Ramsey’s expression: the mechanism of nuclear spin–spin coupling

Consider a molecular system, where the nuclear magnetic moments are given byMP. Associated with these moments, there is a vector potentialAnuc(r) and a magnetic inductionBnuc(r), given by (in atomic units)

AnucðrÞ ¼a²X

P

MP$rP

r³_P ; ð2Þ

BnucðrÞ ¼$$AnucðrÞ: ð3Þ

In these equations, the fine-structure constanta=c^%1(where the speed of lightc&137 in atomic units) andrPis the position relative to nucleus P. The surrounding electrons will interact with these magnetic moments, both since the electrons are charged particles in motion relative to the nuclei, and since the electrons themselves have magnetic moments. However, these interactions are tiny rela-

tive to the electrostatic interactions between the electrons and nuclei and are therefore adequately described by perturbation theory, because the magnetic interactions will modify the electronic energy only slightly. For closed-shell molecules, there is no first-order change in the electronic energy. To second order, the change is described in terms of the reduced indirect spin–spin coupling tensors KPQ

EðMÞ ¼Eð0Þ þX

P>Q

M^T_PKPQMQþOðM³Þ; ð4Þ

whereMdenotes the collection of all magnetic momentsMPin the molecule. Therefore, theKPQare simply the second derivatives of E(M) at zero-magnetic momentsM=0

KPQ¼ d²EðMÞ dMPdMQ

!!

!M¼0

: ð5Þ

Thus, the coupling constants may be calculated by standard tech- niques for second-order molecular properties, such as those used for calculating, for example, quadratic force constants (which are the second derivatives with respect to nuclear coordinates rather than magnetic dipole moments). The spectroscopically observed indirect spin–spin coupling tensor is proportional to the reduced tensor

J_PQ¼hcP

2p cQ

2pKPQ: ð6Þ

The corresponding scalar spin–spin coupling constant J_PQ is one- third of the trace of this tensor – that is, the average of the diagonal elements ofJPQ. All higher-order terms in Eq.(4)are very small and may be safely neglected.

The basic equations underlying the nonrelativistic approach to the calculation of spin–spin coupling constants were first derived by Ramsey[25]. In nonrelativistic theory, there are several distinct contributions to the spin–spin coupling constants, arising from the magnetic hyperfine coupling of the spin of the nuclei to the orbital motion of the electrons and to the spin of the electrons.

First, the orbital hyperfine or spin–orbit (SO) coupling represents the interaction of the nuclei with charged particles (i.e., the electrons) moving in the vector potential Anuc(r) generated by the nuclei. There are two such orbital hyperfine operators – namely, the diamagnetic SO (DSO) operator and the paramagnetic SO (PSO) operator

h^DSO_PQ ¼a⁴X

i

r^T_iPriQI%riPr^T_iQ

r³_iPr³_iQ ; ð7Þ

h^PSO_P ¼a²X

i

riP$p_i

r³_iP ; ð8Þ

where piis the linear-momentum operator of electron i, Iis the 3$3 unit matrix, and the summations are over all electrons. Next, the spin hyperfine interactions are determined by the nuclear magnetic fieldBnuc(r). This field, which interacts with spin of the elec- tronss_i, gives rise to two distinct first-order triplet operators h^FCP ¼8pa²

3 X

i

dðriPÞsi; ð9Þ

h^SD_P ¼a²X

i

3riPr^T_iP%r²_iPI

r⁵_iP si; ð10Þ

where the Fermi contact (FC) operatorh^FCP represents interaction at the position of the nucleus and the spin–dipole (SD) operatorh^SDP

interaction at a distance.

Following Ramsey’s nonrelativistic theory[25], the final expression for the indirect nuclear spin–spin coupling tensorKPQcan be written in the sum-over-states formulation as

(4)

KPQ¼D0!!!h^DSO_PQ !!!0E

þ2X

n_S6¼0

0!!!h^PSO_P !!!n_S

D E

n_S!"h^PSO_Q #T

!!

!

!!

!!0

$ %

En_S%E0

þ2X

nT

0!!!h^FC_P þh^SD_P !!!n_T

D E

n_T"h^FC_Q#T

þ"h^SD_Q#T

!!

!!0

$ %

En_T%E0 : ð11Þ

The first two terms in Eq.(11)represent the orbital contribution to the coupling tensor (with a summation over all singlet excited states labeled byn_S), whereas the third term represents the spin contribution (with a summation over all triplet states labeled by n_T). The mixed FC–SD and SD–FC terms contribute only to the anisotropic part ofKPQand can be omitted in the calculation of the isotropic couplingKPQ¼¹₃trKPQ, leading to

KPQ¼K^DSO_PQ þK^PSO_PQ þK^FC_PQþK^SD_PQ ð12Þ The first contributionK^DSO_PQ is easily evaluated, being the ground- state expectation value of the one-electron second-order DSO operator. The evaluation of the remaining three contributionsK^PSO_PQ,K^FC_PQ andK^SD_PQ is much more difficult. Direct use of the sum-over-state expressions in Eq.(11)requires a knowledge of all singlet and triplet excited states and cannot be performed in practice. Instead, these contributions to the spin–spin coupling constants are obtained from linear response theory[26,27], by solving three singlet response equations and seven triplet response equations for each nucleus. For nuclear shielding constants, by contrast, only three equations are solved for any molecule, independent of its size.

Indirect scalar nuclear spin–spin coupling constants are often dominated by the FC contribution in Eq.(12); sometimes even to the extent that the remaining three contributions (SD, DSO and PSO) may be omitted. However, there are many exceptions to this behavior. For example, the PSO or SD contributions dominate in ClF3[28]; in many other molecules, the PSO and SD contributions to fluorine–fluorine coupling constants are at least as important as the FC contribution[29]. Similarly, the SD contribution domi- nates¹J_CPin HCP[30]. It should also be noted that long-range coupling constants (i.e., coupling constants between nuclei separated by several bonds) are usually dominated not by the FC contribution but by the DSO and PSO contributions, which at large separations typically have opposite signs, the DSO contribution being negative and the PSO contribution positive[31].

Clearly, all four Ramsey contributions to the spin–spin coupling constants in Eq.(12)should be considered in any attempt at quan- titative accuracy – the PSO, DSO and SD contributions may often be small but can rarely be neglected. For a qualitative description, a practical compromise may sometimes be to calculate the non-FC contributions – in particular, the expensive SD contribution – at a lower level of theory than the FC contribution (e.g., in a smaller basis set). Finally, we note that the decomposition of nuclear spin–spin coupling constants into DSO, PSO, FC and SD contributions is only valid in the nonrelativistic limit, which may not al- ways be adequate for spin–spin coupling constants.

Restricting ourselves to the scalar coupling constants (i.e., to the traces of the coupling tensors), we shall not analyze here the full set of components of the spin–spin coupling tensors, referring instead to the review by Vaara et al.[21]for a discussion of theoretical and experimental issues related to the anisotropy and antisymmetry of such tensors. Let us only note that, for a diatomic molecule, we can often determine both the spin–spin coupling constant and the anisotropy from their relation to the parameters of molecular beam electric and magnetic resonance spectra [21,32].

In the following, we shall discuss the reduced coupling con- stantsK_PQof Eq.(5)rather than the isotope-dependent coupling constantsJ_PQof Eq.(6)whenever we compare the accuracy of var-

ious computational approaches – in such cases, the gyromagnetic ratios enteringJ_PQwould only obscure the analysis when different nuclei are considered. On the other hand, whenever a comparison is made with experiment, we shall discuss the results in terms of the observable coupling constants J_PQ. To simplify the notation, we shall often omit the isotope number. Unless otherwise stated, we discuss the coupling constants for the spin-1/2 isotopes such as¹H,¹³C and¹⁵N; if no such isotope exists, we consider instead the magnetically active isotopes of the largest natural abundance – for example,¹¹B and¹⁷O. For a fixed molecular geometry, the coupling constants for other isotopes may then be derived from the gyromagnetic ratios.

2.2. The wave-function approach to indirect spin–spin coupling constants

There are two complementary approaches to the calculation of indirect nuclear spin–spin coupling constants from electronic- structure theory – we may calculate the spin–spin coupling constants either from an (approximate) electronic wave function or from an (approximate) electronic density, using DFT. As we shall see in the following, both approaches have their advantages and disadvantages. In the present section, we discuss calculations based on the construction of an approximate wave function; in the subsequent section, we discuss calculations based on the electron density. For an overview of the different methods for the calculation of spin–spin coupling constants, seeFig. 1, where we have arranged these methods (including the relativistic methods discussed in Section2.5) schematically in relation to one another.

2.2.1. CI and MCSCF theories

The basic wave-function model ofab initioquantum chemistry – that is, the restricted Hartree–Fock (RHF) model – is unfortunately ill suited to the calculation of spin–spin coupling constants, giving a poor description of triplet perturbations. In the RHF model, triplet excitation energies are typically strongly underestimated and may sometimes even become negative (triplet instabilities).

Because of such (near) instabilities, the RHF model often overesti- mates the FC and SD contributions to the coupling constants – see Eq.(11), which shows that these contributions depend inversely on the triplet excitation energies. Thus, it was realized early on that electron correlation must be included to obtain reliable spin–spin coupling constants. Already in 1974, Kowalewski et al.[33]used the configuration-interaction (CI) model to obtain an accurate value of the coupling constant in HD.

RHF

CORELA

MCSCFTION

CC

SOPPA DFT

2-comp. DFT FCI

4-comp. HF ZORA-DFT Douglas-Kroll DFT IORAmm DFT CCSD

CC3

RASSCF CASSCF

LDA BLYP B3LYP other functionals...

numerical cost

Fig. 1.Electronic-structure methods for the calculation of indirect spin–spin coupling constants, arranged according to their treatment of electron correlation and relativity. The RHF method, which is not recommended for the calculation of spin–

spin coupling constants, has been crossed out.

(5)

The triplet instabilities of RHF theory arise from the use of a single Slater determinant in situations where instead several degenerate or nearly degenerate determinants should be included. In CI theory, the wave function is expressed as a linear combination of Slater determinants constructed from (occupied and virtual) RHF orbitals, but without a reoptimization of the orbitals. In the more flexible multiconfigurational self-consistent field (MCSCF) model, the wave function is constructed in the same manner but the orbitals are variationally optimized simul- taneously with the expansion coefficients of the determinants.

This simultaneous optimization makes the MCSCF model well suited to treat the static correlation that arises from the near degeneracy of several configurations.

In MCSCF theory, two models are typically used to compute spin–spin coupling constants: the complete-active-space SCF (CASSCF) model, in which the determinants of the wave function are generated by allowing all possible virtual electronic excitations within a set of ‘‘active” orbitals; and the restricted-active-space (RASSCF) model, where additional low-order excitations are allowed into an extended set of orbitals. Only static electron correlation (arising from near degeneracies) can be treated within the CASSCF model, whereas also dynamical correlation can be (at least partly) recovered by the RASSCF model. An implementation of MCSCF theory for the calculation of spin–spin coupling constants has been presented by Vahtras et al.[27]. For small molecules such as C2H2, highly accurate results may be obtained by applying MCSCF theory[34]; for larger systems, it becomes increasingly difficult to choose an adequate RASSCF active space[35]. For a discussion of MCSCF theory and its application to the calculation of indirect spin–spin coupling constants, see Ref.[13].

2.2.2. Coupled-cluster theory

Nowadays, coupled-cluster theory is the most popular method for recovering correlation energy by means of the traditional tech- nique of quantum chemistry – namely, by expanding the wave function in a linear combination of Slater determinants. In coupled-cluster theory, the coefficient of each determinant is not sep- arately optimized as in CI theory. Instead, the weight of each determinant is obtained by adding together the probabilities of all virtual excitation mechanisms that lead to this determinant from the Hartree–Fock determinant. To illustrate, consider the basic coupled-cluster singles–doubles (CCSD) model, where all possible combinations of single and double virtual excitations are included in the description. In the CCSD wave function, a given doubly excited determinant may be formed either as the result of a single ‘‘connected” double excitation, or as the result of two simultaneous ‘‘disconnected” single excitations. The total weight of each doubly excited determinant in the wave function is then the sum of the probabilities of both these excitation mechanisms.

At the coupled-cluster singles–doubles–triples (CCSDT) level of theory, we also include connected triple excitations, thereby improving our description by increasing the number of excitation mechanisms leading to a given determinant. Likewise, we include, at the coupled-cluster singles–doubles–triples–quadruples (CCSDTQ) level of theory, connected quadruple excitations in our description, and so on.

The coupled-cluster method is a convenient black-box method, where we may systematically approach the exact solution to the Schrödinger equation by including higher and higher connected excitations in our description, rapidly converging towards the exact solution. Indeed, for many purposes, the CCSD model yields satisfactory results, although triple excitations are often needed for agreement with experiment. However, since the CCSDT model is prohibitively expensive, the amplitudes of the connected triple excitations are usually estimated by perturbation theory, yielding the less expensive CCSD(T) (CCSD with a noniterative perturbative

treatment of the triples) method and the CC3 (CCSD with an iterative perturbative treatment of the triples) method.

For the accurate calculation of indirect spin–spin coupling constants, the CCSD model is often used. Its main drawbacks are high computational demands (it scales asN⁶, whereNis the number of basis functions) combined with a need for large basis sets to re- cover electron correlation (unlike in DFT). It should be also kept in mind that the CCSD description, being based on the RHF description, may also be affected by triplet instabilities or near instabilities. These problems are avoided by performing calculations in an ‘‘unrelaxed” manner, where the RHF orbitals are not allowed to relax in response to triplet perturbations. In this manner, the CCSD model is prevented from inheriting the unphysical RHF response to such perturbations. The effects of orbital relaxation are instead described by means of the coupled-cluster singles amplitudes[36].

There are two implementations of the CCSD model for spin–

spin calculations: the equation-of-motion CCSD (EOM-CCSD) implementation by Perera et al.[37]and the analytic second-derivative implementation by Stanton and Gauss[38]. The application of the latter derivative-based coupled-cluster approach to the evaluation of spin–spin coupling constants has been described in Ref.

[36]. A recent new implementation, which includes a parallel code for coupled-cluster analytic derivatives, should enable coupled- cluster calculations on larger systems[39].

The first implementations of coupled-cluster theory for the calculation of spin–spin coupling constants beyond the CCSD model are the CCSDT, CCSD(T) and CC3 implementations of Auer and Gauss [36]. These authors conclude that the CCSD(T) model – in most other respects the gold standard of quantum chemistry – is unsuitable for the calculation of spin–spin coupling constants, requiring a relaxed treatment of the orbitals[36]. If connected triples excitations are needed, the more expensive CC3 model is more appropriate, being suited to an unrelaxed orbital treatment. In- deed, more accurate results have been obtained with the CC3 method than with the CCSD(T) method. Very recently, the analytical-derivative approach has been extended to the calculation of indirect spin–spin coupling constants at the orbital-relaxed CCSDT and CCSDTQ levels of theory[40].

2.2.3. SOPPA theory

Although coupled-cluster theory provides the most satisfactory framework for including dynamical correlation in calculations of indirect spin–spin coupling constants, it is often too expensive for large systems. A useful, less expensive alternative is the second-order polarization-propagator approximation (SOPPA), where the relevant terms in the calculation of spin–spin coupling constants are obtained from perturbation theory. This method was first applied to the study of spin–spin coupling constants by Geertsen and Oddershede[41]. In the more recent SOPPA(CCSD) model[42], the Møller–Plesset (MP2) correlation amplitudes are replaced by the corresponding CCSD amplitudes. (This approach thus requires a CCSD calculation but only for the unperturbed system.) Today, the SOPPA and SOPPA(CCSD) models are both successfully used to compute coupling constants in relatively large molecules. It should be noted that the SOPPA model differs from the MP2 model, which has also been used (less successfully) for the calculation of spin–spin coupling constants – see, for example, Refs.[43,44].

2.2.4. FCI theory

The method that accounts for electron correlation in the most complete fashion (within a given basis set) is full configuration- interaction (FCI) theory, which is equivalent to coupled-cluster theory with all possible virtual excitations included. For a two- electron system such as the hydrogen molecule, it is equivalent to the CCSD model; for a four-electron system, it is equivalent to

(6)

the CCSDTQ model, and so on. Unfortunately, the number of FCI variables increases very rapidly with system size – namely, as

M N

& '

whereNandMare the numbers of electrons and basis functions, respectively. Systems containing more than ten active electrons (larger than, for instance, the nitrogen molecule with the core electrons frozen) therefore cannot be studied with FCI theory.

For spin–spin coupling constants, FCI theory has been applied to evaluate the coupling constant in the four-electron helium dimer [45]and to benchmark coupled-cluster results in the six-electron BH molecule[46].

2.2.5. Comparison of results obtained using different wave-function methods

InTable 1, we list the reduced indirect spin–spin coupling constants calculated using different wave-function models for a number of small molecules. All coupling constants have been calculated at the equilibrium geometry and are compared with experimental equilibrium values, obtained from the observed coupling constants by subtracting the vibrational contributions listed in the last column.

These vibrational corrections have been calculated theoretically at the DFT level[47]. The errors listed are the mean absolute errors jDjrelative to the empirical values (in the units of K and in percent).

The spin–spin coupling constants calculated using the RHF model are clearly unreliable. Whereas the coupling constants for HF, H2O

and NH3are too high but otherwise in reasonable agreement with experiment, those for CO and N2have wrong signs. For C2H4, the calculated constants are in error by an order of magnitude or more. As discussed in Section2.2.1, these errors arise from a poor description of triplet excitations, which are often strongly underestimated in RHF theory, giving rise to too large FC and SD contributions to the spin–spin coupling constants – see Eq.(11).

The triplet-instability problems are removed at the CASSCF level of theory, which includes a treatment of static (but not dynamical) correlation into the description. The spin–spin coupling constants calculated using simple full-valence CASSCF wave functions inTable 1are mostly in qualitative agreement with experiment, although the CO and N2 coupling constants are strongly underestimated. Since the CASSCF treatment ignores the effects of dynamical correlation we conclude that these are sometimes significant for indirect nuclear spin–spin coupling constants (although we note that static and dynamic correlation cannot be unambiguously separated).

In RASSCF theory, the description is further improved by incor- porating a simple treatment of dynamical correlation, leading to a better performance in most cases – in particular, for CO, whose calculated spin–spin coupling constant is now in good agreement with experiment. On the other hand, although the RASSCF treatment improves the spin–spin coupling constant of N2considerably, this constant is still in error by more than 50%.

Table 1

Reduced indirect nuclear spin–spin coupling constants calculated by different wave-function methods (10¹⁹kg m^%2s^%2Å^%2)

RHF CASSCF RASSCF SOPPA CCSD CC3 Exp. Vib.

HF ¹K_HF 59.2^a 48.0^b 48.1^f 46.8^j 46.1^m 46.1^m 47.6^q %3.4^w

CO ¹K_CO 13.4^a %28.1^c %39.3^c %45.4^j %38.3^m %37.3^m %38.3^r %1.7^w

N2 1K_NN 175.0^a %5.7^c %9.1^c %23.9^j %20.4^m %20.4^m %19.3^s %1.1^w

H2O ¹K_OH 63.7^a 51.5^d 47.1^g 49.5^j 48.4^m 48.2^m 52.8^t %3.3^w

2K_HH %1.9^a %0.8^d %0.6^g %0.7^j %0.6^m %0.6^m %0.7^t 0.1^w

NH3 1KNH 61.4â 48.7ê 50.2^h 51.0^k 48.1ⁿ 50.8û %0.3^w

2KHH %1.9â %0.8ê %0.9^h %0.9^k %1.0ⁿ %0.9û 0.1^w

C2H4 1K_CC 1672.0â 99.6ê 90.5ⁱ 92.5^l 92.3ô 87.8^v 1.2^w

1K_CH 249.7^a 51.5^e 50.2ⁱ 52.0^l 50.7^p 50.0^v 1.7^w

2KCH %189.3^a %1.9^e %0.5ⁱ %1.0^l %1.0^p %0.4^v %0.4^w

2KHH %28.7^a %0.2^e 0.1ⁱ 0.1^l 0.0^p 0.2^v 0.0^w

3K_cis 30.0^a 1.0^e 1.0ⁱ 1.0^l 1.0^p 0.9^v 0.1^w

3K_trans 33.3^a 1.5^e 1.5ⁱ 1.5^l 1.5^p 1.4^v 0.2^w

jDj Abs. 180.3 3.3 1.6 1.8 1.2 1.6

% 5709 60 14 24 23 6

The listed experimental equilibrium values have been obtained from the observed coupling constants by subtracting the calculated zero-point vibrational contributions shown in the last column.

aUnpublished results for this review.

b Ref.[48].

c Ref.[49].

d Ref.[50].

eRef.[13].

f Ref.[51].

gRef.[52].

h Ref.[53].

i Ref.[54].

j Ref.[55].

kRef.[56].

l Ref.[57].

mRef.[36].

n Ref.[37].

o Ref.[58].

p Ref.[59].

q Refs.[60,61].

rRef.[62].

s Ref.[63].

tRef.[64].

u Ref.[65].

vRefs.[48,50].

w Ref.[47].

(7)

In SOPPA theory, all constants are finally in a fairly good agreement with experiment, although the N2coupling is still overestimated by about 20%. As expected, coupled-cluster theory performs well, also for the difficult N2molecule. In the few cases where CC3 theory has been applied, the difference between the CCSD and CC3 results is small, indicating that the CCSD treatment may be sufficient for most purposes.

It is worth noting that the range of molecules to which wave- function methods can be applied may be soon extended with the emergence of local-correlation methods, which have already been implemented for some molecular properties[66–68]. To our best knowledge, no attempts have yet been made to apply these methods to the calculation of NMR properties, but it seems that local- correlation methods may be useful in this field.

2.3. The DFT approach to indirect spin–spin coupling constants The standard starting point for a wave-function treatment of molecular properties is the Hartree–Fock model, which for many properties yields results in qualitative agreement with experimental measurements. However, as discussed in Section2.2, the RHF model is not a good model for spin–spin coupling constants, because of its poor description of triplet perturbations. For spin–spin coupling constants, we must therefore either resort to MCSCF theory or use a method for dynamical correlation such as the CCSD method with unrelaxed orbitals. In short, there are no inexpensive wave-function methods available for spin–spin coupling constants.

For this reason, Kohn–Sham DFT, whose cost is similar to that of Hartree–Fock theory, plays an important role in the theoretical treatment of spin–spin coupling constants. Unlike correlated wave-function methods, it can be applied to large systems, at a fraction of their cost. Indeed, it was only with the introduction of Kohn–Sham theory that the quantum-mechanical evaluation of indirect nuclear spin–spin coupling constants became a routine procedure, regularly used in conjunction with experimental work.

The premise of Kohn–Sham theory is that the exact electron density can be generated from a single-determinant reference wave function, which plays the role of the exact wave function for a fictitious system of noninteracting electrons, whose effective Hamiltonian is adjusted in such a manner that its electron density becomes equal to the density of the interacting system. The advantage of this approach is that we may – in principle, at least – gen- erate the exact density of any electronic system at a cost broadly similar to that of Hartree–Fock theory. The validity of this approach is ensured by the Hohenberg–Kohn theorems, according to which the ground-state electron density determines the one-electron potential of the system and hence all its properties (including the ground-state energy). Moreover, from a knowledge of the depen- dence of the energy on the electron density, we may obtain the exact electron density and energy by application of the variation principle. In Kohn–Sham DFT, we apply this principle in a slightly modified form, calculating the bulk of the kinetic energy from an orbital expression in the same manner as in Hartree–Fock theory, while expressing the exchange and correlation energies (including corrections to the kinetic energy) directly in terms of the electron density, by means of a universal exchange–correlation functional.

2.3.1. Approximate exchange–correlation functionals

Unfortunately, the exact form of the universal exchange–correlation functional (from which we may calculate the exchange–correlation energy directly from the electron density for any molecule) is unknown. Over the years, many approximate exchange–correlation functionals have therefore been developed and tested. In the absence of a single, universal functional, some of these are better suited than others to the calculation of indirect nuclear spin–spin coupling constants. We shall here describe the

main classes of exchange–correlation functionals, discussing briefly the more common functionals and those functionals whose results are discussed in the present review.

In the simplest approximation, the exchange–correlation energy is calculated by integrating a local function of the electron density R

FðqðrÞÞdr over all three-dimensional space. This local- density approximation (LDA) has been very successful in solid- state physics but less so in chemistry, being less accurate than traditional wave-function theory. The LDA exchange–correlation functional is usually constructed by combining the Dirac–Slater exchange functional with the Vosko–Wilk–Nusair (VWN) correlation functional[69], a parameterization based on accurate simulations of the uniform electron gas[70].

In the generalized gradient approximation (GGA), the exchange–correlation energy is obtained by integrating over a function that depends not only on the electron density at each point in space but also on the density gradientRFðqðrÞ; $qðrÞÞdr. With the emergence of GGA and the development of useful gradient-cor- rected exchange–correlation functionals in the late 1980s, Kohn–

Sham theory became competitive with wave-function theory, pro- viding an accuracy often comparable with that of CCSD theory. Two commonly used GGA exchange–correlation functionals are the Becke–Lee–Yang–Parr (BLYP) functional [71,72] (fitted semi- empirically to noble-gas data) and the Perdew–Burke–Ernzerhof (PBE) functional [73] (nonempirical, based on the properties of the slowly varying electron gas). We shall later consider also the nonempirical Perdew–Wang 91 (PW91) functional [77] and the Keal–Tozer KT1 and semi-empirical KT2 functionals[75].

Often, even better results are obtained by combining GGA functionals with some proportion of exact (Hartree–Fock) exchange. In this manner, we obtain the popular semi-empirical Becke–3- parameter–Lee–Yang–Parr (B3LYP) hybrid functional [72,76]

(based on the BLYP functional, with 20% exact exchange, and fitted to thermochemical data) and the nonempirical PBE0 functional [74](based on the PBE functional, with 25% exact exchange). We shall also consider two hybrid modifications of the Becke 97 (B97) functional [78]: the 10-parameter semi-empirical B97-1 functional (fitted to thermochemical data) [79] and the B97-2 functional (fitted also to accurate exchange–correlation potentials) [80]. In Coulomb-attenuated hybrid functionals, the amount of exact exchange included in the functional depends on the interelec- tronic distance. One such functional is the Coulomb-attenuated- method B3LYP (CAMB3LYP) functional[81,82].

2.3.2. Computer implementations of DFT for indirect spin–spin coupling constants

The first successful (and still popular) implementation of Kohn–

Sham theory for the calculation of indirect nuclear spin–spin coupling constants was that of Malkin et al.[83–85], at the LDA and GGA levels of theory. In their code, the FC contribution to the spin–spin coupling constant was obtained by finite perturbation theory (FPT) and the SD term was omitted, although the important SD–FC contribution to the spin–spin coupling anisotropy was obtained by FPT. The PSO term was approximated by the sum-over- states approach.

Another early Kohn–Sham implementation of spin–spin coupling constants was the LDA implementation by Dickson and Zie- gler[86], also with the neglect of the SD term. Their code, which uses Slater rather than Gaussian atomic orbitals, has subsequently been made fully analytical at the GGA level of theory (with the SD term included) and extended to account for relativity in the zero- order regular approximation (ZORA)[87,88].

The first fully analytical Kohn–Sham implementations of indirect spin–spin coupling constants, including all four Ramsey terms, were those presented in 2000 by Sychrovsky´ et al.[89]and by Hel- gaker et al. [90]. Both implementations included hybrid DFT, in

(8)

addition to LDA and GGA. In 2004, Watson et al.[91]presented an implementation of spin–spin coupling constants that uses Slater orbitals, at the hybrid level of theory and with all four Ramsey terms included. Very recently, DFT-based perturbation theory was applied to compute spin–spin coupling tensors in extended systems, subject to periodic boundary conditions[92].

For reviews of Kohn–Sham theory for the calculation of indirect nuclear spin–spin coupling constants, see the one by Malkin et al.

[93]and the more recent ones by Alkorta and Elguero[94]and by Helgaker and Pecul[95].

2.3.3. The performance of DFT for indirect spin–spin coupling constants

With some notable exceptions, the accuracy of DFT for indirect nuclear spin–spin coupling constants is surprisingly good – seeTa- ble 2, where we have listed the coupling constants calculated using Kohn–Sham theory for the same molecules as inTable 1, which contains the corresponding wave-function results. Kohn–Sham theory represents a huge improvement on Hartree–Fock theory, being much less affected by triplet-instability problems.

First, at the simple LDA level of theory, the agreement with experiment is rather poor – in particular, for the difficult N2mole- cule, for which LDA gives an incorrect sign. Moreover, the HF coupling constant is underestimated by about 25%, while the CO constant is overestimated by about 70%. Comparing withTable 1, we note that, even though the LDA method represents a vast improvement on RHF theory, its performance is still poorer than that of the correlated methods. Clearly, to compete with these methods, we must use exchange–correlation functionals better than the LDA functional.

Next, proceeding to GGA theory, where the exchange–correlation treatment includes a gradient correction, we observe a significant improvement on LDA theory, with errors reduced by about a factor of two when the semi-empirical BLYP functional and the nonempirical PBE functional are used. However, even though the sign of the N2

coupling constant is now correctly predicted, the HF, CO and N2re-

sults are still in poor agreement with experimental measurements.

By contrast, the spin–spin coupling constants of C2H4are accurately predicted by the BLYP functional and, in particular, by the PBE functional. The excellent performance of the PBE functional for¹J_CHcou- pling constants has been noted by Maximoff et al.[98].

Finally, with the introduction of exact exchange in the B3LYP, B97-2 and B97-3 exchange–correlation functionals, errors relative to experiment are reduced further – indeed, the performance of these hybrid methods is better than that of the CASSCF method inTable 1. Recalling that DFT is applicable to a much broader range of molecules than CASSCF theory, this observation illustrates the usefulness of DFT for the prediction of spin–spin coupling constants[90]. Thus, in Ref.[99], it was found that DFT performs better than MCSCF theory but less well than CCSD theory for small rigid hydrocarbons. For further comparisons with wave-function methods, see Refs.[100–102].

In general, the performance of DFT for indirect nuclear spin–

spin coupling constants improves from LDA to GGA functionals and then from GGA to hybrid functionals. However, the difference in performance between the different classes of functionals is not clear-cut. Indeed, Maximoff et al.[98]found no such correlation in their study of¹J_CHcoupling constants, for which a GGA functional (PBE) performs best. In a follow-up study, Keal et al.[96]

compared the PBE, B3LYP, B97-2 and B97-3 results for ¹J_CH and for other types of spin–spin coupling constants, concluding that the best overall performance is that of the hybrid B97-3 functional.

The performance of the CAMB3LYP functional (which greatly improves Rydberg and charge-transfer excitation energies as well as polarizabilities) was in Ref.[82]found to be similar to that of the B3LYP functional for spin–spin coupling constants.

For typical organic molecules, good results are usually obtained with the standard DFT exchange–correlation functionals [96–

98,103]. InFig. 2, we have plotted the experimental spin–spin coupling constants against the B3LYP results for a set of small organic molecules [103], illustrating the usefulness of this popular exchange–correlation functional for the prediction of spin–spin Table 2

Reduced indirect nuclear spin–spin coupling constants calculated using different exchange–correlation functionals (10¹⁹kg m^%2s^%2Å^%2)

LDA BLYP B3LYP PBE B97-2 B97-3 Exp. Vib.

HF ¹K_HF 35.0^a 34.5^a 38.9^c 32.6^c 39.3^c 40.5^c 47.6^d %3.4^j

CO ¹K_CO %65.4^b %55.7^b %47.4^c %62.0^c %44.9^c %43.4^c %38.3^e %1.7^j N2 1K_NN 32.9^b %46.6^b %20.4^c %43.2^c %22.7^c %12.5^c %19.3^f %1.1^j

H2O ¹K_OH 40.3^b 44.6^b 47.2^c 41.2^c 45.7^c 46.3^c 52.8^g %3.3^j

2K_HH %0.3^b %0.9^b %0.7^c %0.5^c %0.7^c %0.6^c %0.7^g 0.1^j

NH3 1KNH 41.0^b 49.6^b 52.3^c 47.0^c 48.7^c 50.1^c 50.8^h %0.3^j

2KHH %0.4^b %0.7^b %0.9^c %0.7^c %0.8^c %0.8^c %0.9^h 0.1^j

C2H4 1K_CC 66.6^b 90.3^b 96.2^c 83.4^c 94.5^c 92.9^c 87.8ⁱ 1.2^j

1K_CH 42.5^b 55.3^b 55.0^c 50.0^c 50.6^c 51.4^c 50.0ⁱ 1.7^j

2KCH 0.4^b 0.0^b %0.5^c %0.2^c %0.7^c %0.3^c %0.4ⁱ %0.4^j

2KHH 0.4^b 0.4^b 0.3^c 0.3^c 0.2^c 0.3^c 0.2ⁱ 0.0^j

3K_cis 0.8^b 1.1^b 1.1^c 1.0^c 0.9^c 1.0^c 0.9ⁱ 0.1^j

3K_trans 1.2^b 1.7^b 1.7^c 1.6^c 1.5^c 1.5^c 1.4ⁱ 0.2^j

jDj Abs. 11.2 5.9 3.1 6.4 2.7 2.6

% 72 48 14 33 14 14

The listed experimental equilibrium values have been obtained from the observed coupling constants by subtracting the calculated zero-point vibrational contributions shown in the last column.

aRef.[96].

b Unpublished results for this review.

c Ref.[97].

d Refs.[60,61].

eRef.[62].

f Ref.[63].

gRef.[64].

h Ref.[65].

i Refs.[48,50].

j Ref.[47].

(9)

coupling constants. In a similar study of more than 250 carbon–

carbon coupling constants across one, two and three bonds in substituted benzenes[104], an excellent agreement is observed between B3PW91/6-311++G(d,p) results and experiment. In this investigation, the molecular geometries were optimized using the same exchange–correlation functional and basis set as in the subsequent calculation of the spin–spin coupling constants, which ap- pears to be the best approach[102]. In another B3LYP study, more than 240 coupling constants in pyrazoles were analyzed [105].

Again, the results are satisfactory, with a linear regression of J_exp= (1.042 ± 0.003)Jcalcbetween the calculated and experimental values. Likewise, a total of 73 B3LYP coupling constants were com-

pared with available experimental data and EOM-CCSD results in Ref.[101].

Some care should be exercised when selecting the exchange–

correlation functional for the calculation of spin–spin coupling constants in highly unsaturated organic molecules, especially those containing triple bonds: the quality of calculated DFT coupling constants depends sensitively on the density of the electrons, deteriorating with increasing CAC bond order[15]. Theo-benzyne results inTable 3illustrate the performance of different exchange–

correlation functionals for the calculation of spin–spin coupling constants in a highly strained biradical system [102]. Because of near-equilibrium triplet instabilities that follow from its biradical character, the calculated spin–spin coupling constants of this molecule are unusually sensitive to details of the exchange–correlation functional. At the experimental equilibrium geometry, the exchange–correlation functionals give widely different results, nearly all far from experiment. By contrast, more uniform and better results are obtained when the spin–spin coupling constants are evaluated at the geometry optimized using the same exchange–

correlation functional as for the coupling constants. The reason for this difference in behavior is that the geometry optimization moves the system away from the triplet instabilities, thereby improving the results.

In general, DFT performs poorly for indirect spin–spin coupling constants involving lone-pair atoms – in particular, fluorine as illustrated for the HF molecule inTable 2. As a further illustration, we present inTable 4results for the T-shaped ClF3molecule. For this molecule, the Ramsey sum-over-states contributions FC, SD and PSO depend sensitively on the exchange–correlation functional. Moreover, their values vary strongly with the geometry and do not improve at the optimized geometries, suggesting that the poor performance is not related to triplet instabilities. The agreement of the DFT and wave-function results for¹J(Cl–Fax) is fortuitous, the individual Ramsey contributions differing signifi- cantly. The only available experimental values are the weighted average of the chlorine–fluorine coupling constants ¹J_ave(Cl–F) = 260 Hz and ²J(Feq–Fax) = 403 Hz. Further examples of the erratic behavior of DFT for spin–spin coupling constants to the fluorine atom are given in a recent review[15].

An obvious advantage of DFT is that it can be routinely applied to large molecules such as C60[109]and valinomycin[31], provid- ing results of the same quality as for small molecules – see Section 4for other examples and references.

50 100 150 200 250 J calc.

50 100 150 200 250

J exp.

J_pred=0.909J_calc−0.067 Hz

Fig. 2.DFT/B3LYP indirect spin–spin coupling constants plotted against empirical coupling constants for ethyne, ethene, allene, cyclopropene, cyclopropane, cyclobutene, pyrrole, furan, thiophene and benzene (Hz). The DFT constants have been calculated using the B3LYP functional in the HIV-su5 basis set[103]. The empirical coupling constants have been obtained from the experimental ones by subtracting theoretically calculated zero-point vibrational corrections[103]. The full line represents the least-squares fitJpred= 0.909Jcalc%0.067 Hz, which is compared with a straight (dashed) line of slope one.

Table 3

One-bond indirect nuclear spin–spin coupling constants ino-benzyne (Hz)^a

BLYP PW91 PBE KT1 KT2 B3LYP B97-1 B97-2 CCSD

Experimental geometry^b

1J(C1„C2) 236.4 664.9 225.9 213.0 434.1 577.3 448.5 566.4 210.2

1J(C2AC3) 75.4 121.9 69.4 69.8 96.7 108.9 109.3 110.1 83.4

1J(C3@C4) 59.7 68.8 54.1 55.8 63.9 68.8 74.3 69.4 63.3

1J(C4AC5) 77.9 90.0 71.8 73.8 83.9 90.1 93.2 89.2 77.2

1J(C3H3) 166.5 174.2 152.4 159.4 180.0 170.4 166.6 157.9 155.8

1J(C4H4) 147.5 146.3 134.6 141.0 155.4 144.6 143.9 134.7 137.4

BLYP PW91 PBE^c KT1 KT2 B3LYP B97-1 B97-2 Exp.^d

Optimized geometry

1J(C1„C2) 222.3 389.2 208.7 205.7 309.7 322.2 310.3 334.0 177.9 ± 0.7

1J(C2AC3) 80.6 103.7 74.7 74.0 92.4 96.8 104.4 99.0 75.7 ± 0.9

1J(C3@C4) 54.2 60.0 48.9 49.9 57.3 61.5 67.8 62.4 50.9 ± 0.8

1J(C4AC5) 76.0 84.1 71.3 72.9 81.6 82.6 88.7 83.6 71.0 ± 0.8

1J(C3H3) 167.4 172.2 153.7 162.9 177.1 166.6 165.2 154.4

1J(C4H4) 149.1 148.6 135.3 142.8 153.7 145.6 144.5 134.1

a Calculated in the cc-pCVTZ basis, see Ref.[102].

b From microwave measurements[106].

c Very similar results are obtained in the cc-pCVQZ basis – for instance,¹J(C1„C2) = 205.4 Hz,¹J(C2–C3) = 74.6 Hz,¹J(C3@C4) = 49.2 Hz and¹J(C4AC5) = 71.4 Hz; and in the uncontracted cc-pCVQZ basis.

d Experimental data of Ref.[107](o-benzyne incarcerated in a host molecule).