Accurate prediction of hyperfine coupling tensors for main group

Jan 30, 2019 - ... a unitary group based rigorously spin-adapted coupled-cluster theory ... for the previously reported unitary group based spin-adapt...
0 downloads 0 Views 1MB Size
Subscriber access provided by UB der LMU Muenchen

Quantum Electronic Structure

Accurate prediction of hyperfine coupling tensors for main group elements using a unitary group based rigorously spin-adapted coupled-cluster theory Dipayan Datta, and Jürgen Gauss J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b01048 • Publication Date (Web): 30 Jan 2019 Downloaded from http://pubs.acs.org on February 2, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Accurate prediction of hyperfine coupling tensors for main group elements using a unitary group based rigorously spin-adapted coupled-cluster theory Dipayan Datta∗,†,‡ and Jürgen Gauss∗,† †Institut für Physikalische Chemie, Johannes Gutenberg-Universität Mainz, Duesbergweg 10-14, 55128 Mainz, Germany ‡Current address: Max-Planck-Institut für Kohlenforschung, Kaiser-Wilhelm-Platz 1, 45470 Mülheim an der Ruhr, Germany E-mail: [email protected]; [email protected]

Abstract We present the development of a perturbative triples correction scheme for the previously reported unitary group based spin-adapted combinatoric open-shell coupledcluster (CC) singles and doubles (COS-CCSD) approach and report on the applications of the newly developed method, termed “COS-CCSD(T)”, to the calculation of hyperfine coupling (HFC) tensors for radicals consisting of hydrogen, second- and third-row elements. The COS-CCSD(T) method involves a single noniterative step with N 7 scaling of the computational cost for the calculation of triples corrections to the energy. The key feature of this development is the use of spatial semicanonical orbitals generated from standard restricted open-shell Hartree-Fock (ROHF) orbitals, which allows to

1

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

define the unperturbed Hamiltonian operator in terms of a diagonal spin-free Fock operator. The HFC tensors are computed as first-order property via implementation of an analytic derivative scheme. The required one-particle spin density matrix is computed by using one- and two-particle spin-free density matrices that are obtained from the analytic derivative implementation, in this way avoiding the use of any spin-dependent operator and maintaining spin adaptation of the CC wavefunction. Benchmark calculations of HFC tensors for a set of 21 radicals indicate reasonably good agreement of the COS-CCSD(T) results with experiment and a consistent improvement over the COSCCSD method. We demonstrate that the accuracies of the isotropic hyperfine coupling constants obtained in unrestricted HF (UHF) reference based spin-orbital CCSD(T) calculations deteriorate when spin contamination in the UHF wavefunction is large and the results may even become qualitatively incorrect when spin polarization is the driving mechanism. Within a similar noniterative perturbative treatment of triple excitations, the spin-adapted COS-CCSD(T) approach produces accurate results, thus ensuring cost-effectiveness together with reliability.

1

Introduction

The magnetic hyperfine interaction between the unpaired electrons in a molecule with nuclei having nonzero spins is described in terms of the hyperfine coupling (HFC) tensor, 1–3 which is an essential parameter for the elucidation of the electronic structure in open-shell molecules. In electron paramagnetic resonance (EPR) spectroscopy, the experimental spectra are interpreted in terms of the phenomenological spin Hamiltonian, 4,5 in which the HFC tensors appear as fitting parameters. An accurate theoretical prediction of the HFC tensor components thus is the key for obtaining detailed information about the electronic structure from the observed spectra. However, the complexity of the physical effects that lead to the hyperfine coupling interactions in molecules poses tremendous challenges to an accurate calculation of HFC tensors. In computational chemistry, the common approach is to employ

2

ACS Paragon Plus Environment

Page 2 of 61

Page 3 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

spin-unrestricted formalisms, e.g., the unrestricted Hartree-Fock (UHF) theory, which necessarily introduces spin contamination into the electronic wavefunction and thus leads to large errors. 6 In applications of density functional theory 7 (DFT) as well, spin contamination is known to introduce large inaccuracies 8,9 to the computed HFC tensors. Although coupledcluster (CC) theory is known to provide accurate results, 10–15 CC calculations for open-shell systems usually employ a spin-orbital ansatz, which introduces spin contamination. Whether spin contamination allows the standard spin-orbital CC approach to predict HFC tensors accurately in all cases, is an open question. The availability of efficient computational tools based on a rigorously spin-adapted CC theory is immensely helpful in this context, which we put forth in this paper. Spin contamination in correlated wavefunctions arising in spin-orbital based CC treatments can be reduced via spin projection or by imposing spin restrictions on the CC wavefunction. 16–19 However, achieving rigorous spin adaptation requires the use of spatial orbital parameters (cluster amplitudes). A spin-adapted CC theory for open-shell states using spatial orbital amplitudes is formulated by (a) using a spin-free unitary group adapted (UGA) linear combination of determinants, e.g., a configuration state function (CSF) as the reference state 20,21 and (b) defining excitations in the cluster operator Tˆ in terms of the spin-free unitary group generators. 22 The unitary group generators commute with the spin operators Sˆ2 and Sˆz , which ensures that the CC wavefunction is a spin eigenfunction. The pioneering work on the UGA-CC approach was due to Li and Paldus 23–27 as well as Jeziorski et al. 28,29 Further work in this field includes the development of the combinatoric open-shell CC (COS-CC) approach by Datta and Mukherjee 30–32 and the more recent developments of the spin-adapted state-specific 33 and state-universal 34 multireference CC (SS/SU-MRCC) approaches by Mukherjee and co-workers. We shall use the generic term “unitary group based CC methods” to refer to the UGA-CC, COS-CC, and SU/SS-MRCC approaches as a group of rigorously spin-adapted CC methods. The unitary group based CC methods offer the advantage of requiring lower computa-

3

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

tional costs than spin-orbital based CC theory for open-shell states by reducing the number of cluster amplitudes by about a factor of three at the singles and doubles (SD) level and about four when triple excitations are included. Furthermore, the use of a CSF instead of a single determinant as reference state enables the treatment of low-spin open-shell states, which are not accessible via standard spin-orbital based single-reference CC theory. However, a general implementation for open-shell states with arbitrary spin multiplicities cannot be achieved within the framework of unitary group based CC theory. The available implementations are able to handle doublet 24,25,30–35 (one unpaired electron), triplet, and low-spin singlet 24,25,33,34 (two unpaired electrons) states. Another drawback, which applies to the UGA-CC and the spin-adapted SU/SS-MRCC approaches, is the non-availability of analytic derivative implementations for property calculations, while for spin-orbital based CC theory such techniques 36–42 have long been used as standard and powerful computational tools. Applications of the UGA-CC approach to compute static electrical properties via finite-field technique and to compute equilibrium geometries and vibrational frequencies for diatomic molecules have been reported. 25,43,44 However, no application of the above two methods to the calculation of spin-dependent properties can so far be found, which in our opinion is the real testing ground for the performance of a spin-adapted CC theory. Having an analytic derivative implementation is indispensable for this purpose, as numerical differentiation is complicated in the case of spin-dependent perturbations. 45 Our recent work on unitary group based CC theory aimed at bridging the above gap by developing analytic derivative schemes for the COS-CC approach. We accomplished an efficient implementation of the COS-CCSD ansatz by developing programming tools for the automated derivation and implementation of general many-body methods. 35 This was followed by an implementation of analytic first derivatives of the COS-CCSD energy, 46 which was applied to the calculations of spin-independent properties, e.g., dipole moments and nuclear quadrupole coupling constants. These applications have so far been limited to dou-

4

ACS Paragon Plus Environment

Page 4 of 61

Page 5 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

blet systems. Before extending our implementation to open-shell states with higher spin multiplicities, e.g., triplets and low-spin singlets, we considered it worthwhile to assess the performance of a unitary group based CC theory, viz., the COS-CC approach, in the calculation of spin-dependent properties, since a spin-adapted treatment of electron correlation is expected to be important for accuracy in these calculations. The calculation of spin-dependent properties, e.g., the HFC tensors, requires the oneparticle spin density matrix. The unitary group based CC theory uses a purely spatial orbital based, i.e., spin-free formulation, while the one-particle spin density matrix is a ˆ pq = spin-dependent quantity. If we introduce the one-particle spin density operator ∆ ˆqβ , ... denote spin-orbital creation and annihilation operators) directly ˆqβ (ˆ a†pα , a ˆ†pβ a ˆ qα − a a ˆ†pα a into the spin-free formulation, the spin symmetry of the correlated wavefunction will be broken, i.e., it will no longer remain an eigenfunction of the Sˆ2 operator. The unitary group based CC theory offers the opportunity to compute energies and molecular properties in a rigorously spin-adapted manner. We aim at exploiting the full advantage of this feature and compute spin-dependent properties without severing the rigorous spin adaptation of the CC wavefunction. This choice makes the computation of the one-particle spin density matrix conceptually somewhat involved. This goal can nevertheless be accomplished and possible strategies have been discussed in the literature. 47–50 Among various procedures, the most pragmatic strategy is the one proposed by Luzanov et al. 51–54 (see Sec. 2.4). We implemented the latter procedure within our analytic derivative scheme 55 and demonstrated the applicability of the COS-CCSD method in the calculation of isotropic hyperfine coupling constants for 2 Π radicals. Our previous applications 46,55 pointed to the intriguing fact that at the SD truncation level, spin-adapted CC theory does not offer the adequate flexibility required for the accurate prediction of molecules properties. This occasionally leads to larger errors in COS-CCSD results than in spin-orbital based CCSD results. The source of this problem is as follows. At the single excitation level, the spin-adapted CC theory entails two types of amplitudes, viz.,

5

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the pure singles tia (i, j, .. denote doubly occupied spatial orbitals and a, b, .. spatial virtual orbitals) and the exchange-spectator type excitations tui au that correspond to a spin flip in the singly occupied orbital u. The full space spanned by the spin-adapted singly excited configurations can be reached from the reference state only if both amplitudes are present. 24 This is true also for the doubles. Within the CCSD approximation, the three-body exchangeij spectator counterpart tiuj abu of the pure doubles tab is usually missing, as they have been shown

to make insignificant contributions to the energy. 56 However, these “pseudo-doubles” might be important for accurate property calculations. This led us to the question as to whether the desired accuracy may be achieved by including triple excitations in the COS-CC ansatz, such that we arrive at a computational method that performs in a consistent manner in property calculations for a broad variety of molecules and therefore can be used as a reliable tool similar to the analytic derivative techniques for spin-orbital based CC theory. This sets the background for the present work, where we report on the development of a perturbative triples correction scheme on top of the COS-CCSD ansatz. We resort to a computationally efficient noniterative triples correction scheme with an N 7 scaling of the computational cost similar to the standard CCSD(T) approach. 57 Special care is exercised in the formulation, such that the orbital transformation properties of the underlying COSCCSD ansatz are preserved. This is crucial for achieving an efficient analytic derivative implementation. 40 Importantly, perturbative triples correction schemes were reported for the UGA-CC ansatz, 26,27 but their applications were limited only to energy calculations. Here, we present applications of the newly developed “COS-CCSD(T)” approach to compute HFC tensors for radicals consisting of main group elements and provide an assessment of its performance relative to COS-CCSD and spin-orbital based electron correlation methods. We then address the question concerning the importance of spin adaptation in CC theory.

6

ACS Paragon Plus Environment

Page 6 of 61

Page 7 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

2

Theory

2.1

Brief review of the COS-CC approach

We first briefly recall the essential features of the COS-CC approach. For details, we refer to previous papers. 30–32,35 The correlated wavefunction in this formalism is parametrized as

|ΨCC i = {{exp(Tˆ)}}|Ri,

(1)

where |Ri is a unitary group adapted spin-free open-shell reference state, e.g., a Gel’fandTsetlin state, 20,58 which is described solely in terms of spatial orbitals. For general open-shell states, |Ri is essentially a linear combination of Slater determinants, i.e, a CSF. For doublet states considered in this work, |Ri is a single determinant obtained from a restricted openshell Hartree-Fock (ROHF) calculation. The excitations in the cluster operator Tˆ are defined in terms of spatial orbital replacements. Defining such excitations requires the use of a purely spatial orbital based function as vacuum state for classifying the orbitals as occupied and virtual. The open-shell reference state |Ri being a spin-free function, can in principle itself serve as vacuum state. However, in the formulation of the COS-CC approach, 30–32,35 we chose to use a closed-shell determinant Φ0 as vacuum state, in which the orbital singly occupied in |Ri is doubly filled. While this choice is not mandatory, we show in this paper that this offers a simple route to the formulation of a noniterative triples correction scheme without destroying the invariance of the COS-CC ansatz under transformations with respect to doubly occupied (virtual) orbitals. The (spatial) orbitals are classified into holes (occupied, denoted by I, J, K, ...) and particles (virtual, denoted by a, b, c, ...) with respect to the vacuum state Φ0 . The hole orbitals I, J, K, ... are further classified into a doubly occupied set i, j, k, ... and the singly occupied orbital u. The indices p, q, r, ... are used for denoting general spatial orbitals. The excitations in Tˆ are defined in terms of the generators of the orbital unitary group 24

7

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 61

U (n) (n denotes the number of spatial orbitals)

E˜qp =

X

pq {ˆ a†pσ a ˆqσ }, E˜rs =

X

{ˆ a†pσ a ˆ†qρ a ˆ sρ a ˆrσ },

(2)

σ,ρ

σ=α,β

which are spin-summed normal ordered products of spin-orbital creation and annihilation operators. In our formulation, the normal ordering is taken with respect to the closedshell vacuum state Φ0 . This choice of the vacuum is indicated by the symbol {..} and the ˜ Within the singles and doubles corresponding unitary group generators are denoted by E. (SD) approximation, the excitations in Tˆ are defined as     (0) (1) ˆ ˆ ˆ ˆ ˆ ˆ ˆ ˆ T = T + T = T1 + T2 + T1e + T2e + T2x

(3)

with

P Tˆ1 = tIa E˜Ia , Tˆ2 = a,I

1 2!

P

P i ˜u ˜ ab ˆ tIJ tu Ei , ab EIJ , T1e = i

a,b, I,J

P P ij ˜ au ˆ ˜ au Tˆ2e = tau Eij , T2x = tui au Eui .

(4)

a,i

a,i,j

The Tˆ(0) operator is defined in terms of hole-particle excitations and when the singly occupied ˆ(1) orbital u is involved, Tˆ(0) contains excitations out of u only. Note that tuu ab = 0. The T operator involves excitations into the singly occupied orbital as well as “exchange spectator” type single excitations, i.e., Tˆ2x . The above choice of excitations within the SD approximation is equivalent to that of the UGA-CCSD(is) approach of Li and Paldus, 56 which includes only those excitations that appear at first order in perturbation theory. Importantly, the excited configurations {|L0 i} reached by the action of the excitations defined in eq 4 on |Ri are neither linearly independent nor orthogonal. The spin-free nature of these excitations only ensures that the configurations {|L0 i} are also spin-free CSFs, but they are not Gel’fand states. The UGA-CC approach, 24,25 on the other hand, uses unitary

8

ACS Paragon Plus Environment

Page 9 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

group adapted tensor operators to define excitations in the cluster operator, which ensure that the excited configurations are both linearly independent and orthogonal. In COS-CC approach 30–32,35 (as well as in spin-adapted SU/SS-MRCC methods 33,34 ), a nonorthogonal yet linearly independent excitation manifold {|Li} is employed. Such an excitation manifold can be obtained by taking simple linear combinations of the excited configurations {|L0 i} as 2 1 2 1 au 0 au 0a ab 0 ab 0 ba |Lui i = |L0 ui i, |Lai i = 23 |L0 ai i − 13 |L0 au ui i, |Lui i = 3 |L ui i − 3 |L i i, |Lij i = 3 |L ij i + 3 |L ij i, and 2 1 0 au 0 au |Lau ij i = 3 |L ij i + 3 |L ji i, without a rigorous unitary group adaptation.

A special feature of the above spin-free excitations is their non-commuting nature, e.g., [Tˆ(0) , Tˆ(1) ] 6= 0, [Tˆ(1) , Tˆ(1) ] 6= 0, which leads to Tˆ − Tˆ contractions. These contractions make the CCSD amplitude equations cumbersome when the standard exponential ansatz, exp(Tˆ), is used. These equations contain terms with up to the 8th power of cluster amplitudes, while standard spin-orbital based CC amplitude equations (naturally) terminate at the 4th power of cluster amplitudes. In UGA-CC theory, the amplitude equations are simplified via truncation at the quadratic power, 24 while in SU/SS-MRCC methods, 33,34 the terms with Tˆ − Tˆ contractions are altogether avoided by using a normal ordered exponential ansatz. 59 The Tˆ − Tˆ contractions mimic the powers of SD excitations, the presence of which is the source of the robustness of CC theory. For a core electron property such as the HFC tensors, the accurate prediction of which is tremendously challenging (see Sec. 3), a persistent success of these approximate schemes cannot be ensured. The primary idea behind the COS-CC approach was to develop a spin-adapted polynomial ansatz that captures the same physical effects as the standard exp(Tˆ) ansatz of CC theory, yet leads to amplitude equations as compact as in standard CC theory without the necessity of approximations. Such an ansatz was developed by modifying the standard exp(Tˆ) ansatz containing Tˆ − Tˆ contractions such that the following polynomial, 30,31 i.e.,

9

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 61

the COS-CC ansatz, which is denoted in eq 1 as {{exp(Tˆ)}}, is obtained 1 1 1 {{exp(Tˆ)}} = 1 + {Tˆ(0) } + {Tˆ(1) } + {Tˆ(0) }2 + ... + {Tˆ(0) Tˆ(1) } + {Tˆ(1) Tˆ(0) } + ... 2! 2! 2! X 1 X 1 (1) (1) (1) (1) (1) {TˆA TˆB } + {TˆA TˆB TˆC } + ..., + (5) f f AB ABC A,B A,B,C with an overline denoting a contraction. The COS-CC ansatz differs in two aspects from the standard exp(Tˆ) ansatz. First, it allows only Tˆ(1) − Tˆ(1) contractions, while Tˆ(1) − Tˆ(0) contractions are excluded. 31 By definition, the excitations in Tˆ(0) commute among themselves, i.e., [Tˆ(0) , Tˆ(0) ] = 0 and hence the powers of Tˆ(0) are automatically included in eq 5. The second and most important difference lies in the prefactors accompanying the contracted composites of Tˆ(1) operators. The prefactor for the nth power of the cluster operator in the exp(Tˆ) ansatz is 1/n! irrespective of whether Tˆ − Tˆ contractions exist or not. When developing the COS-CC ansatz, a careful analysis was made 30–32 to investigate whether the prefactor 1/n! is appropriate for capturing the same physical effects arising from n contracted Tˆ(1) operators as in the case with no Tˆ − Tˆ contractions. This analysis led to the conclusion that the prefactors for the powers of Tˆ(1) operators need to be replaced with the quantities fAB... , which were termed the “automorphic factors”. 30–32,35 These factors ensure that the powers of the Tˆ(1) operator (with contractions) are included in the COS-CC ansatz with appropriate weights and that the physical contents of such powers are not undermined. The interested reader is referred to Ref., 32 where a perturbative analysis for the emergence of these automorphic factors was presented. It was furthermore shown that the use of these automorphic factors is the key for arriving at a compact structure of the COS-CCSD amplitude equations via simplification, 30–32 which contain terms with at most quartic power of cluster amplitudes similar to the spin-orbital based CCSD amplitude equations. The COS-CCSD amplitude equations for the doublet case 30,35 read

ˆ exp(Tˆ)}}s − {exp(Tˆ(1) )}h H ˆ ef f )|Ri = 0, ∀i, hLi |({{H

10

ACS Paragon Plus Environment

(6)

Page 11 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

where the subscript “s” refers to strongly connected terms 30,31 that contain Tˆ(1) − Tˆ(1) contractions, but each t amplitude is contracted to the Hamiltonian operator through at least ˆ ef f denotes the effecone doubly occupied or virtual orbital index. In the above equation, H tive Hamiltonian operator defined via

ˆ ef f |Ri = 0, ∀i hLi |H

(7)

ˆ ef f |Ri = hR|{{H ˆ exp(Tˆ)}}s |Ri = ECOS−CCSD . hR|H

(8)

and

ˆ ef f is defined as the “closed” component of the strongly connected operator {{H ˆ exp(Tˆ)}}s Thus, H that is labelled solely by the singly occupied orbital index u, while all external indices (douˆ ef f with bly occupied and virtual) are fully contracted. Hence, the expectation value of H respect to the reference state |Ri leads to the COS-CCSD energy. In eq 6, the subscript “h” furthermore implies that each Tˆ(1) operator in {exp(Tˆ(1) )}h is individually contracted to the ˆ ef f operator through the u index. 30,31 H Let us recall that our formulation and implementation uses the Fock operator

F

(0)

no i X Xh pK pK hpq + (2vqK − vKq ) E˜qp = p,q

(9)

K=1

defined with respect to the closed-shell vacuum state Φ0 , with no denoting the number of pq (doubly) occupied orbitals in Φ0 , hpq denoting the one-electron integrals, and vrs = hpq|rsi

denoting the two-electron integrals in Dirac notation.

2.2

Perturbative treatment of triple excitations within the COS-CC ansatz

For developing a perturbative triples correction scheme for the COS-CC ansatz, we intend to satisfy the following criteria: (a) the method should be noniterative, (b) should be invariant 11

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

under rotations among doubly occupied (virtual) orbitals similar to COS-CCSD, and (c) the scaling of the computational cost of the triples corrections should be N 7 akin to the standard CCSD(T) approach. 40,57 ˆ 0 needs to be In order to satisfy criterion (a), the unperturbed Hamiltonian operator H defined in terms of a diagonal Fock matrix. However, the Fock matrix corresponding to the F (0) operator defined in eq 9 is nondiagonal in the ROHF orbital basis. In principle, a noniterative triples correction scheme can be developed by choosing only the diagonal part ˆ 0 . However, the resulting method will not satisfy criterion of this Fock matrix to define H (b), as the off-diagonal part of F (0) is ignored within a noniterative scheme. 60 The loss of orbital invariance would complicate the implementation of analytic derivatives. Moreover, the errors arising from the neglected off-diagonal Fock matrix elements are not guaranteed to be small. In the formulation of the spin-orbital based ROHF-CCSD(T) approach, 40 the problem due to a nondiagonal Fock matrix is resolved by rotating the standard ROHF orbitals to the semicanonical basis, 61 which diagonalizes the occupied-occupied (virtual-virtual) block of the Fock matrix. However, this transformation renders the spatial parts of the α and β spin orbitals nonidentical. These semicanonical orbitals are thus incompatible with our requirement that the reference state |Ri should be a spin-free function. Thus for a unitary group based CC theory that uses a ROHF reference state, it is not straightforward to find a diagonal Fock matrix. The noniterative triples correction schemes developed for the UGACC ansatz 26 bear the same problem and lack the orbital transformation properties of the underlying UGA-CCSD ansatz. In the COS-CC approach, our choice of the closed-shell vacuum state Φ0 in the formulation provides an indirect but straightforward solution to this problem, since we have the option to use the closed-shell Fock operator F (0) . While the ROHF orbitals cannot be directly rotated to a semicanonical basis without making the spatial parts of the α and β spin orbitals different, it is possible to first construct the nondiagonal Fock matrix for the

12

ACS Paragon Plus Environment

Page 12 of 61

Page 13 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

closed-shell vacuum state Φ0 in the standard ROHF orbital basis and then rotate the orbitals such that the occupied-occupied (virtual-virtual) block of the transformed closed-shell Fock matrix becomes diagonal. The transformed closed-shell Fock operator reads

F˜ (0) =

no X I

(0)

with ˜p

(0) ˜I E˜II +

nv X

˜a ˜(0) a Ea +

a

X

(0) f˜Ia E˜aI

(10)

I,a

denoting the pth spatial orbital energy in the transformed orbital basis and nv

denoting the number of the virtual orbitals. It is then straightforward to define the spin-free reference state |Ri for the target open-shell state in the transformed spatial orbital basis. Following the standard definition, 61 we name these transformed orbitals “semicanonical”. Importantly, the difference of our semicanonical orbitals with those used in the spin-orbital based ROHF-CCSD(T) approach 40 is implied. We emphasize that our formulation does not involve the quasi-restricted HF (QRHF) orbitals optimized for the closed-shell vacuum state, nor do we use them in any calculation. The underlying orbital basis is defined by standard ROHF orbitals and their transformation to the semicanonical basis is used as an intermediate step, such that the implemented triples correction scheme satisfies the criteria (a)-(c) above. We have numerically verified that the COS-CCSD energy and the first-order properties remain invariant under this orbital rotation. The spin-free triple excitation operators for the doublet case are defined as follows 1 X ijk ˜ abu ˆ 1 X IJK ˜ abc ˆ 1 X iuj ˜ abu tabc EIJK , T3e = tabu Eijk , T3x = t E . Tˆ3 = 3! I,J,K 2! i,j,k 2! i,j abu iuj a,b,c

a,b

(11)

a,b

The Tˆ3 operator represents hole-particle triple excitations from |Ri with at the most one of the indices I, J, K being equal to the singly occupied orbital index u and the Tˆ3e operator represents “semi-internal” triple excitations involving excitation into u. The amplitudes corresponding to I = J = K (or i = j = k) and a = b = c are excluded. The Tˆ3x operator in eq 11 indicates a three-body pseudo-double excitation, which is actually the

13

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 61

exchange spectator counterpart for Tˆ2 . Following the suggestion of Li and Paldus, 26,27 we have decided to resort to a perturbative treatment of the Tˆ3x operator. Likewise, we do not ˜ abc include the four-body exchange spectator type excitations corresponding to tijk abc Eijk in our perturbative triples correction scheme. It is important to note that the spin-free triple excitations defined in eq 11 do not generate an orthonormal excitation manifold when acting on |Ri and that the number of triples amplitudes is greater than those defined in Ref. 26 For example, there are six triples ampliijk ijk ijk ijk ijk tudes in our formulation, viz., tijk abc , tbac , tacb , tcba , tcab , tbca , in contrast to five corresponding

orthogonally spin-adapted amplitudes. 26 However, as pointed out by Matthews et al., 62 using these nonorthogonally spin-adapted amplitudes is rather advantageous for the derivation of equations and for achieving an efficient implementation. The equations can be derived by employing common tools of second quantization (e.g., the spin-summed form of Wick’s theorem 63 ). For the full CCSDT implementation, the use of a larger number of nonorthogonally spin-adapted triples amplitudes might lead to an unnecessary increase in the computational cost. However, this is not a major problem for the present noniterative triples correction scheme. For the formulation of a perturbative triples correction scheme, we employ the CC Lagrangian technique. 64,65 For this, we define an energy functional by augmenting the COS-CC ˆ ef f operator, with the amplitude residuals ri (eq 6) as conenergy, given in terms of the H straints

ˆ ef f |Ri + E = hR|H

X

λi ri ,

i

  ˆ exp(Tˆ)}}s − {exp(Tˆ(1) )}h H ˆ ef f |Ri ˆ ef f |Ri + hR|Λ ˆ {{H = hR|H

(12)

ˆ denotes a de-excitation operator that via a set of Lagrange multipliers {λi }. In eq 12, Λ is defined in terms of the {λi } multipliers as the amplitudes. While the COS-CC energy is nonvariational in the ti amplitudes, the energy functional E is made stationary with respect

14

ACS Paragon Plus Environment

Page 15 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

to both the ti and the λi parameters. This implies that the ti amplitudes are determined from the stationarity of E with respect to the variation of the λi parameters and vice versa. Since we are interested in corrections to the COS-CCSD energy due to triple excitations ˆ operators in eq 12 to at most three. only, we restrict the excitation ranks of both Tˆ and Λ ˆ operator is defined The excitations in Tˆ have been defined already in eqs 3 and 11. The Λ in terms of de-excitation components corresponding to the excitations in Tˆ as

ˆ = Λ ˆ (0) + Λ ˆ (1) = (Λ ˆ1 + Λ ˆ2 + Λ ˆ 3 ) + (Λ ˆ 1e + Λ ˆ 2e + Λ ˆ 3e + Λ ˆ 2x + Λ ˆ 3x ) Λ

(13)

with

ˆn = Λ

X 1 X ab.. ˜ IJ.. ˆ 1e = λIJ.. Eab.. , n = 1 − 3, Λ λui E˜ui , n! I,J,.. i

ˆ 2e = Λ

X

ˆ 2x = Λ

X

a,b,..

X ˜ ij , Λ ˆ 3e = 1 ˜ ijk λau E λabu ij au ijk Eabu , 2! i,j,a i,j,k a,b

˜ ui ˆ λau ui Eau , Λ3x

i,a

1 X abu ˜ iuj λ E . = 2! i,j iuj abu

(14)

a,b

ˆ operators we Splitting eq 12 into terms with various excitation ranks of the Tˆ and Λ obtain   (1) SD ˆ ˆ ˆ ˆ ˆ ˆ E = hR|Hef f |Ri + hR|ΛSD {{H exp(TSD )}}s − {exp(TSD )}h Hef f |Ri   ˆ SD {{H ˆ exp(TˆT )}}s |Ri + hR|Λ ˆ T {{H ˆ exp(Tˆ)}}s − {exp(Tˆ3e )}h H ˆ ef f |Ri +hR|Λ ˆ D {{H( ˆ Tˆ1 TˆT + Tˆ2 Tˆ3x + {Tˆ(1) Tˆ3x })}}s |Ri, +hR|Λ SD

(15)

ˆ ef f ), D, and T imply SD, double, and triple where the subscripts SD (superscript for H ˆ excitations in Tˆ, respectively, as well as the corresponding de-excitation components in Λ. In the fourth term, Tˆ = TˆSD + TˆT . As in the standard CCSD(T) approach, 40,57 we first iteratively solve the COS-CCSD 15

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 61

amplitude equations and keep the SD amplitudes fixed when triple excitations are included. Hence, the second term in eq 15 can be eliminated once the converged SD amplitudes are available and the first term can be substituted with eq 8 for the COS-CCSD energy. Equation 15 then simplifies to

ˆ SD {{H ˆ exp(TˆT )}}s |Ri E = ECOS−CCSD + hR|Λ   ˆ exp(Tˆ)}}s − {exp(Tˆ3e )}h H ˆ ef f |Ri ˆ T {{H +hR|Λ ˆ D {{H( ˆ Tˆ1 TˆT + Tˆ2 Tˆ3x + {Tˆ(1) Tˆ3x })}}s |Ri. +hR|Λ SD (16)

Now we partition the Hamiltonian into a zeroth-order component and a perturbation as

ˆ =H ˆ0 + W ˆ H

(17)

with ˆ0 = H

X

˜p ˆ ˜ (0) ˜(0) p Ep , W = FOV + V =

p

X I,a

1 X pq ˜ pq (0) v E f˜Ia E˜aI + 2! p,q, rs rs

(18)

r,s

ˆ operators in orders of the perturbation and expand the Tˆ and Λ

Tˆ = Tˆ[1] + Tˆ[2] + Tˆ[3] + ...,

(19)

ˆ = Λ ˆ [1] + Λ ˆ [2] + Λ ˆ [3] + ... Λ

(20)

With the above partitioning of the Hamiltonian, the first-order cluster amplitudes are given by [1] (0) [1] D1 Tˆ1 = F˜OV , Tˆ1e = 0, [1] [1] [1] D2 Tˆ2 = V, D2e Tˆ2e = V, D1 Tˆ2x = V.

16

ACS Paragon Plus Environment

(21)

Page 17 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

[1] with Dk denoting the orbital energy denominators. Note that Tˆ1e = 0 in the semicanonical

orbital basis, since the diagonalization of the occupied-occupied block of the closed-shell (0) (0) Fock matrix sets f˜iu = f˜ui = 0.

The lowest order of perturbation theory at which triple excitation amplitudes appear is two. These second-order triples amplitudes are obtained by setting the first derivatives of E ˆ T to zero. Clearly, with respect to the λi parameters that correspond to the amplitudes of Λ only the third term in eq 16 contributes to the corresponding triples amplitude residuals. Restricting all terms in these residuals to second order, we obtain [2]

ˆ 0 Tˆ3 } + {V Tˆ2 } = 0, {H [2] [2] ˆ 0 Tˆ3e ˆ 0 } = 0, {H } + {V Tˆ2 } + {V Tˆ2e } − {(Tˆ3e )h H [2] ˆ 0 Tˆ3x } + {V Tˆ2 } + {V Tˆ2e } + {V Tˆ2x } = 0, {H

(22)

where we have substituted the first-order double excitation amplitudes with their converged values obtained by solving the COS-CCSD amplitude equations. The subscript “s” has been dropped from the contracted terms for brevity. Furthermore, we have replaced the composites involving double curly braces with simple normal ordered products with contractions indicated by {...}. Since the above equations do not contain any term with Tˆ − Tˆ conˆ 0 operator in the tractions, consideration of the automorphic factors is not required. The H [2] (0) ˆ 0 } denotes the zeroth-order contribution to H ˆ ef f , which is simply f˜uu term {(Tˆ3e )h H = ˜0u .

Introducing the orbital energy denominators D3 = (˜0I + ˜0J + ˜0K − ˜0a − ˜0b − ˜0c ), D3e = (˜0i + ˜0j + ˜0k − ˜0a − ˜0b − ˜0u ), and D2 = (˜0i + ˜0j − ˜0a − ˜0b ), we rewrite eq 22 as [2]

= {V Tˆ2 },

(23)

[2]

= {V Tˆ2 } + {V Tˆ2e },

(24)

D3 Tˆ3

D3e Tˆ3e

[2] D2 Tˆ3x = {V Tˆ2 } + {V Tˆ2e } + {V Tˆ2x }.

17

ACS Paragon Plus Environment

(25)

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 61

Substitution of the above second-order triples amplitudes in eq 16 eliminates the third term through second order of perturbation theory. The correction to the COS-CCSD energy due to triple excitations then arise from the second and the fourth terms. Denoting the energy correction by δE = E − ECOS−CCSD and including only those terms in δE that are correct through fourth-order, we obtain from eq 16 [1]

(0)

[2]

ˆ (F˜ + V )T }|Ri. δE [4] = hR|{Λ SD OV T

(26)

In the above equation, we have dropped the explicit overline to indicate a contraction, though it is implied that the quantity within the braces represents a connected normalordered product. We follow the same notational simplification in the following equations as well. Importantly, all terms in eq 26 arise from the second term of eq 16. The fourth term contributes earliest at the fifth-order in perturbation theory. Thus, there are no terms involving Tˆ − Tˆ contractions in δE [4] . Since we aim at developing a triples correction scheme that is analogous to the CCSD(T) ˆ SD operator in eq 26 with Tˆ† . Substituting eqs 17 and approach, 40,57 we substitute the Λ SD 20 in the COS-CCSD lambda equations (see Ref. 46 ) and collecting the first-order terms, one obtains ˆ [1] ˜ (0) ˆ [1] D1 Λ 1 = FOV , Λ1e = 0, ˆ [1] ˆ [1] ˆ [1] D2 Λ 2 = V, D2e Λ2e = V, D1 Λ2x = V,

(27)

ˆ [1] ˆ[1]† ˆ [1] ˆ[1]† ˆ [1] ˆ[1]† ˆ [1] ˆ[1]† Λ 1 = T1 , Λ2 = T2 , Λ2e = T2e , Λ2x = T2x .

(28)

which leads via eq 21 to

Using eq 28 and substituting the first-order SD cluster amplitudes with the values ob-

18

ACS Paragon Plus Environment

Page 19 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

tained by solving the COS-CCSD amplitude equations, we obtain our final expression for the COS-CCSD(T) energy from eq 26 as [4]

[4]

[4]

ECOS−CCSD(T) = ECOS−CCSD + ET + EST + EDT

(29)

with [4]

ET

[2]† ˆ ˆ [2] ˆ[2]† ˆ ˆ[2] ˆ[2]† ˆ ˆ[2] = hR|{Tˆ3 H 0 T3 } + {T3e H0 T3e } + {T3x H0 T3x }|Ri,

[4] [2] [2] [2] EST = hR|{Tˆ1† V Tˆ3 } + {Tˆ1† V Tˆ3e } + {Tˆ1† V Tˆ3x }|Ri, [4] (0) [2] † ˜ (0) ˆ [2] † ˜ (0) ˆ [2] † ˜ (0) ˆ [2] EDT = hR|{Tˆ2† F˜OV Tˆ3 } + {Tˆ2e FOV T3e } + {Tˆ2e FOV T3x } + {Tˆ2x FOV T3x }|Ri.

(30)

[4]

In arriving at the expression for ET from eq 26, we have used the adjoints of the expressions in eq 22, e.g., [2]† [2]† ˆ ˆ 0 (Tˆ3e ˆ† ˆ† {H )h } − {Tˆ3e H 0 } = {(T2 + T2e )V }.

(31)

We use the acronym “COS-CCSD(T)” for the triples correction scheme described above because of its similar characteristics to the spin-orbital based ROHF-CCSD(T) approach. 40 This means in particular that (a) this method is noniterative, (b) the computational cost of the terms in eqs 23-25 scales as N 7 , (c) the triples contributions to the energy are correct through fourth-order perturbation theory, and (d) importantly, this method preserves the orbital invariance properties of the underlying COS-CCSD ansatz.

2.3

Analytic first derivatives of the COS-CCSD(T) energy

For developing an analytic derivative scheme for COS-CCSD(T), we employ the standard Lagrangian technique 66–71 and rewrite the energy functional E defined in Sec. 2.2 in terms

19

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 61

of the COS-CCSD(T) energy and the amplitude residuals as     X X E = ECOS−CCSD + λi ri + δE [4] + λi ri . i∈SD

(32)

i∈T

Differentiation of the above equation with respect to the cluster amplitudes yields the lambda equations for determining the corresponding Lagrange multipliers. It can be shown that at second order of perturbation theory, the three-body Lagrange multipliers corresponding to the triples amplitudes are given entirely in terms of the SD cluster amplitudes, F˜OV , and V . Therefore, no storage for these three-body Lagrange multipliers is required. The first two terms in eq 32 define the COS-CCSD Lagrangian and the differentiation of these terms with respect to the SD cluster amplitudes would yield the COS-CCSD lambda equations. 46 However, there will be additional contributions arising from the triples amplitudes. These contributions arise solely from the δE [4] term. The contributions from the last term in eq 32 vanish due to eqs 23-25. The triples contributions to the SD lambda equations are given in a compact form as [2]†

hR|{(Tˆ3

lai : ij lab :

[2]†

hR|{(2Tˆ3

[2]†

[2]†

+ Tˆ3e + Tˆ3x )V }|Lai i

[2]† (0) [2]† [2]† † † + Zˆ3† )V } + {Tˆ3 F˜OV } + {(2Tˆ3e + Zˆ3e )V } + {(2Tˆ3x + Zˆ3x )V }|Lab ij i,

ij lau :

[2]† [2]† [2]† (0) [2]† † † hR|{(2Tˆ3e + Zˆ3e )V } + {(Tˆ3e + Tˆ3x )F˜OV } + {(2Tˆ3x + Zˆ3x )V }|Lau ij i

ui lau :

[2]† [2]† (0) † hR|{(2Tˆ3x + Zˆ3x )V } + {Tˆ3x F˜OV }|Lau ui i.

(33)

ij In the above equation, lai , lab , ... denote lambda residuals that determine the λai , λab ij , ... pa-

rameters, respectively. The various intermediates are defined via D3 Zˆ3† = {Tˆ1† V } + {Tˆ2† F˜OV }, (0)

† † ˜ (0) D3e Zˆ3e = {Tˆ1† V } + {Tˆ2e FOV }, † † † ˜ (0) D2 Zˆ3x = {Tˆ1† V } + {(Tˆ2e + Tˆ2x )FOV }.

20

ACS Paragon Plus Environment

(34)

Page 21 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

The SD λi parameters are determined by iteratively solving the modified lambda equations. All terms in eqs 33 and 34 are independent of the SD lambda parameters and therefore need to be constructed only once before the iterative solution of the lambda equations. Once the converged cluster and lambda amplitudes are available, the first derivative of the energy with respect to an external perturbation parameter χ can be expressed in terms of one- and two-particle spin-free (charge) density matrices as (0) (0) (0) X ∂ f˜aa X X ∂ f˜ dE X ∂ f˜II ∂hpq|rsi = DI + Da + DIa Ia + . Dpq,rs dχ ∂χ ∂χ ∂χ ∂χ p,q, a I I,a

(35)

r,s

Since the semicanonical orbital based Fock matrix is diagonal in the occupied-occupied (virtual-virtual) block, we compute only the diagonal elements of the one-particle density matrix DI and Da as indicated in eq 35. The contributions to these density matrices arising from the SD cluster and lambda amplitudes were implemented in our previous work. 46 The triples contributions to the one-particle density matrices are given as †[2]

[2] †[2] †[2] [2] † ˜ i [2] † + Zˆ3† )E˜ii T3 } + {(Tˆ3e + Zˆ3e )Ei T3e } + {(Tˆ3x + Zˆ3x )E˜ii T3x }|Ri,

†[2]

† ˜u + Zˆ3† )E˜uu T3 }|Ri + hR|{(Tˆ3e + Zˆ3e )Eu T3e }|Ri,

DiT

= −hR|{(Tˆ3

DuT

= −hR|{(Tˆ3

DaT

= hR|{(Tˆ3

†[2]

†[2]

[2]

[2]

[2] †[2] †[2] [2] † † ˜ a [2] )E˜aa T3x }|Ri, )Ea T3e } + {(Tˆ3x + Zˆ3x + Zˆ3† )E˜aa T3 } + {(Tˆ3e + Zˆ3e

[2] † ˜ i [2] † ˜ i [2] T Dia = hR|{Tˆ2† E˜ai T3 } + {Tˆ2e Ea T3e } + {Tˆ2x Ea T3x }|Ri, [2] † ˜ u [2] T Ea T3x }|Ri. Dua = hR|{Tˆ2† E˜au T3 } + {Tˆ2e

(36)

Expressions for the two-particle density matrices have been similarly obtained in terms of the intermediates defined in eq 34. Since the standard ROHF orbitals are employed in the calculation, we finally need the relaxed one-particle density matrix in this orbital basis and not in terms of semicanonical orbitals. This in turn requires the calculation of the full one-particle density matrix, since the Fock matrix defined in eq 9 is nondiagonal in the standard ROHF orbital basis. While

21

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the diagonal elements of the one-particle density matrix are obtained from eq 36, the offdiagonal elements are computed in an indirect way by requiring that the orbitals remain ˜ ij (i > j) semicanonical under the external perturbation. 40,72 The off-diagonal elements D ˜ ab (a > b) are then computed as orbital response contributions. and D One could alternatively compute the off-diagonal elements directly via equations similar to eq 36 in the semicanonical basis. The final relaxed density matrix is invariant to the choice of orbitals. We have numerically verified this. This is where we benefit from maintaining the orbital transformation properties of the underlying COS-CCSD ansatz when formulating the perturbative triples correction scheme. Since the COS-CCSD(T) ansatz is also invariant under rotations among doubly occupied (virtual) orbitals, the diagonal and the off-diagonal elements of the one-particle density matrix can be computed in different ways. A direct computation of the triples contributions to the off-diagonal density matrix elements would require the full storage of the triples amplitudes, which is highly expensive and is otherwise not required in our noniterative scheme. Once the full one-particle density matrix is available, the final step is to account for the response of the molecular orbitals to the external perturbation via the Z-vector technique. 73 The Z-vector equations for ROHF orbitals were implemented in our previous work. 46 A slight modification is now necessary, 40 since the A matrices and the X intermediates involved in the Z-vector equations are available in the present case only in the semicanonical basis and we need the solution vector in the standard ROHF orbital basis. During the iterative solution of Z-vector equations, the approximate solution vector at the ith iteration is first constructed in semicanonical basis by contracting the A matrices with the Z-vector from the previous iteration and then transformed to the standard ROHF orbital basis for the next update. The proposed method has been implemented within an in-house spin-adapted CC program developed in this work. 35 This program works in conjunction with the CFOUR program package, 74 which provides the ROHF orbitals. The QRHF option in CFOUR is used for defining the closed-shell vacuum state Φ0 in the ROHF orbital basis. The orbitals are

22

ACS Paragon Plus Environment

Page 22 of 61

Page 23 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

subsequently rotated to the semicaonical basis in the case of COS-CCSD(T) calculations. This is followed by a transformation of the one- and two-electron integrals into the molecular orbital basis within CFOUR. The next steps, viz., the iterative solution of the COS-CC amplitude and lambda equations, construction of charge and spin density matrices, and the treatment of orbital-relaxation effects, are performed within the spin-adapted CC program (see Refs. 35,46 for details). The final output from this program are the orbital-relaxed oneparticle charge and spin density matrices, which are contracted with the requisite property integrals within CFOUR for the calculation of the desired properties. The triples contributions are computed within restricted loops over occupied indices. For 3 the tIJK abc amplitudes, for example, a vector of length nv is stored in the main memory for every

set of indices {I, J, K} with the restriction I ≥ J ≥ K. The triples amplitudes for a specific set {I, J, K} are first computed, which is followed by the evaluation of their contributions to the energy and within the same loop, to the SD lambda equations and the (charge) density matrices. All contributions are computed via vectorized matrix-matrix multiplications. The correctness of the implementation of the analytic derivative scheme has been confirmed by comparing the analytically computed orbital-relaxed COS-CCSD(T) dipole moments with those obtained in finite-field calculations.

2.4

Calculation of spin densities within the spin-free framework of the COS-CC ansatz

As discussed in Sec. 1, we aim at computing spin-dependent molecular properties without destroying the rigorous spin adaptation of the CC wavefunction, which arises from the explicitly spin-free formulation of the unitary group based CC theory. However, this makes the computation of the required one-particle spin density matrix nontrivial in comparison to the standard and direct approach in spin-orbital based CC theory. One possible strategy for computing spin densities within a spin-free quantum chemical method is to employ the “spin-dependent UGA” formalism. 47 This requires computing 23

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 61

matrix elements of the U (2n) group generators, viz., operators of type a ˆ†pσ a ˆqρ with unequal spin indices σ and ρ, or alternatively the matrix elements of a U (n) adjoint tensor operˆ which is given as a polynomial of degree two in the orbital U (n) generators. 47–49 ator ∆, Unfortunately, the procedure for computing such matrix elements 49,50 is complicated and applications of the spin-dependent UGA to modern quantum chemical methods are far from being standard. An alternative procedure for computing spin densities within the unitary group approach was proposed by Luzanov et al., 51–54 who provided the following compact expression for the one-particle spin density matrix for an N -electron state |Ψi with spin S and MS = +S n

∆pq =

i X N 1 h Γpr,rq . (2 − )Γpq − S+1 2 r

(37)

The one- and two-particle charge density matrices in eq 37 are defined as

pq Γpq = hΨ|Eˆqp |Ψi, Γpq,rs = hΨ|Eˆrs |Ψi

(38)

with the spin-free operators Eˆ defined in normal order with respect to the bare vacuum |0i. This relation was derived within the density matrix formulation rather than by using the spin-dependent UGA. This permits its application even within those spin-adapted quantum chemical methods, which only use the orbital U (n) generators and not the full machineries of the UGA, in particular the unitary group algebra for the evaluation of matrix elements. Furthermore, eq 37 can be alternatively derived from simple considerations of anticommutation relations between spin-orbital creation and annihilation operators and by using definitions of the spin operators Sˆ± and Sˆz in second quantization. 75 We recently reported the first implementation of eq 37 for computing the one-particle spin density matrix at the CC level using the unitary group based COS-CCSD method. 55 A subsequent implementation of this formula was reported for internally contracted complete active space second-order perturbation theory (CASPT2). 76

24

ACS Paragon Plus Environment

Page 25 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Since we use the closed-shell determinant Φ0 as vacuum state in the COS-CC formulation ˆ application of eq 37 within the and employ the unitary group generators E˜ rather than E, COS-CCSD(T) analytic derivative scheme requires a renormalization of the density matrices Dpq and Dpq,rs to obtain Γpq and Γpq,rs as described in Ref. 55 Note that we use eq 37 only to compute unrelaxed spin density matrices and add the orbital response contributions later on. For COS-CCSD(T) calculations, the spin density matrices are computed only after the full one-particle charge density matrices are available.

3

Results and discussions

The hyperfine coupling (HFC) tensor for the N th magnetic nucleus AN in a molecule is composed of two contributions 1–3

N A N = aN iso 1 + Adip ,

(39)

where the first term is the isotropic HFC tensor with aN iso denoting the isotropic hyperfine coupling constant (iHFCC) and the second term denotes the anisotropic or spin-dipolar HFC tensor, which is a traceless 3 × 3 tensor arising from dipole-dipole interaction between the spin magnetic moments of the electrons and the nuclei. In absence of spin-orbit coupling, the iHFCC is identical to the Fermi contact term aN F , which is proportional to the spin density at the nucleus and is given (in MHz) as

aN F =

X 4π α−β PN hSˆz i−1 Pµν hφµ |δ(rnN )|φν i 3 µ,ν

(40)

with PN = βe βN ge gN , where ge is the free-electron g-factor, gN is the nuclear g-factor, βe and βN denote the Bohr magneton and the nuclear magneton, respectively, rnN = |rnN | = |rn − RN | is the magnitude of the position vector of the electron n relative to the nucleus N , and δ(rnN ) is the Dirac delta function operator. In eq 40, {φµ } denote atomic basis 25

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 61

functions and Pα−β denotes the one-particle spin density matrix in the atomic orbital basis as obtained from an electronic structure calculation. The components of AN dip are given (in MHz) as X 1 −5 α−β 2 ˆ −1 Pµν hφµ |rnN (rnN δij − 3rnN ;i rnN ;j )|φν i, AN ij = PN hSz i 2 µ,ν

(41)

where rnN ;i and rnN ;j denote the Cartesian components of rnN . The HFC tensor is usually expressed in the principal axis system, in which it is diagonal. In our calculations, the HFC tensor is expressed in the principal axes of the moment of inertia tensor with the components of AN dip being Taa , Tbb , and Tcc (Taa + Tbb + Tcc = 0). The iHFCCs in a molecule indicate the distribution of unpaired electron spins among various atoms, while the spin-dipolar terms provide information about the bonding environments of the spin centres. An accurate theoretical prediction of HFC tensors is a challenging task with the Fermi contact term (iHFCC) being particularly difficult. An accurate description of the spin density near the nucleus is required for this purpose, which may arise from two major sources: (a) the direct contribution from an unpaired electron in a σ or s orbital, and (b) the indirect or spin polarization contribution. The latter usually arises when the unpaired electron is in a p or π orbital, although this may arise also for radicals in Σ states when the unpaired electron is localized on an atom different from the target nucleus (see Sec. 3.2). The accurate prediction of iHFCCs becomes highly challenging when spin polarization is the driving mechanism. The ROHF wavefunction does not account for this effect, while the UHF wavefunction describes spin polarization by introducing spin contamination and is known to largely overestimate the iHFCCs. 6 A proper treatment of electron correlation is thus mandatory for accuracy in these calculations. The early applications of wavefunction based post-HF methods employed the spin-adapted single- and multi-reference (MR) configuration interaction (CI) singles-doubles approach. 77–82 These calculations demonstrated that higher than double excitations are important for an accurate description of iHFCCs. In single-reference SD-CI calculations based on a ROHF reference, higher excitations are to be included explicitly. In MRCI calculations, on the other 26

ACS Paragon Plus Environment

Page 27 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

hand, this is achieved by employing sufficiently large active spaces. For simple doublet radicals, however, an MR treatment is unnecessary in most cases. The MR wavefunction rather serves to include triply and quadruply excited configurations by means of SD excitations from active space configurations other than the ROHF determinant. The direct contributions from higher excitations, e.g., triples, are not significant, but their indirect effect on the singles coefficients is important for the accuracy (see Ref. 83 for a numerical illustration). Single-reference CC theory achieves higher accuracy than CI in the computation of HFC tensors 10–15 already at the CCSD level. This is due to the use of an exponential ansatz, which incorporates higher excitations as products of SD excitations in contrast to the linear SD-CI ansatz. Higher excitations are nevertheless important for achieving a quantitative agreement with experiment. 13 Most CC applications employ the spin-orbital based CC ansatz, which necessarily introduces spin contamination into the CC wavefunction. Whether the undesirable spin contaminants in the CC wavefunction indeed make legitimate contributions to the computed spin densities, however, is not clear. A more practical question is, whether the spin-orbital based CC treatment provides the expected accuracy in all cases. The first application of a unitary group based spin-adapted CC approach, viz., the COSCCSD method, to the calculation of iHFCCs was recently reported by us. 55 The only other application of a spin-adapted CC approach, viz., the internally contracted MRCC (icMRCC) theory, 84 was presented more recently by Samanta and Köhn. 85 Importantly, their approach to spin adaptation is different from the one in COS-CC theory. While the COS-CC approach uses an explicitly spin-free formulation, the spin-adapted variant of the icMRCC approach uses spin projection operators to project out spin contaminants from the spin-orbital based wavefunction parameters. Interestingly, the results for the HFC tensors from ic-MRCC calculations pointed to exactly the same problem of spin-adapted CC theory in the description of spin polarization as our earlier work, 55 i.e., the inadequate flexibility of this ansatz at the SD truncation level. Samanta and Köhn demonstrated that either an enlargement of the active space or the inclusion of pseudo-triple excitations (denoted as icMRCCSDpsT) is

27

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

necessary in order to alleviate this problem. This observation together with insights gained from applications of the CI method provide a further incentive to investigate to what extent the inclusion of triple excitations in the COS-CC ansatz serves to redeem the accuracy issues.

3.1

Comparison with full CI results

In order to gain insights into the effect of triple excitations on the HFC tensor components computed using the COS-CC method, it is highly instructive to compare the COS-CC results with full CI (FCI) calculations. Here we consider the CH (2 Π) radical as test case, for which FCI results have been recently reported by Samanta and Köhn. 85 We use a C-H bond length of 1.12 Å and the modified EPR-III basis set 86 as described in Ref. 85 This is an interesting system, as spin densities at both C and H nuclei arise from spin polarization. In Table 1, we compare COS-CC results with FCI, icMRCC, and spin-orbital based CC results. The icMRCC results considered for this comparison were obtained with the (1e, 1o) active space. Hence, they correspond essentially to spin-adapted single-reference CC results. Let us note that the pseudo-triples in the icMRCCSDpsT ansatz are identical to the exchange-spectator type pseudo-doubles (i.e., the Tˆ3x operator) in our formulation, where we follow the terminology used by Li and Paldus in the context of UGA-CC theory. 24,56 In the following discussion, we shall use the term “pseudo-doubles” to refer to these excitations without distinguishing between the COS-CCSD(T) and icMRCCSDpsT methods. To assess the individual contributions from different types of triple excitations defined in Sec. 2.2, we report results obtained by separately including each type of triples amplitudes in the COS-CCSD(T) ansatz alongside the results obtained with all types of triple excitations. Table 1 indicates that the iHFCCs obtained at the COS-CCSD level substantially differ from the FCI values. The carbon iHFCC deviates from FCI by 15 MHz, which is similar to the 14 MHz deviation of the icMRCCSD result. The inclusion of only the Tˆ3e operator in the COS-CCSD(T) ansatz does not have a noticeable effect. The major improvement to the carbon iHFCC comes from the pseudo-doubles, which amounts to 5.7 MHz. This observation 28

ACS Paragon Plus Environment

Page 28 of 61

Page 29 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

is also in line with the icMRCC approach. The effect of the pseudo-doubles on the proton iHFCC is also large, though in this case the deviation of the COS-CCSD result from FCI is further increased by 1.4 MHz. The correction to the carbon iHFCC due to the hole-particle triple excitations (i.e., the Tˆ3 operator) is smaller than that due to the pseudo-doubles, while these two contributions are identical but opposite for the proton iHFCC. This observation is somewhat surprising. The CH radical is a very small system and so is the basis set used in this calculation. For larger molecules, the number of amplitudes corresponding to the Tˆ3 and Tˆ3e operators far exceeds the number of pseudo-doubles amplitudes. These genuine triples are then expected to make stronger contributions to the spin densities. The deviation of the carbon iHFCC from FCI reduces by about 9 MHz when all three types of triple excitations are included in the COS-CCSD(T) ansatz. The effect of triples on the proton iHFCC is smaller (∼1 MHz), which may be attributed to the fact that the effects of the hole-particle triples and the pseudo-doubles are opposite in this case. Importantly, the effects of different types of triples on the computed properties are not additive, as their contributions get coupled through the solution of the lambda equations. The errors in the Taa values relative to FCI are much smaller even at the COS-CCSD level and are further reduced upon inclusion of triple excitations. In contrast to iHFCCs, the accuracy of the spin-dipolar terms depends more strongly on the choice of basis set than on the treatment of electron correlation. 77,82 We note in passing that the EPR-III basis was developed mainly for DFT calculations 86 and it does not offer enough flexibility for the electron correlation treatment using wavefunction based methods. The results reported in Table 1 are thus far from being converged with respect to the basis set. In view of similar trends of the COS-CCSD and icMRCCSD approaches, we conclude that the inadequate flexibility of the spin-adapted CCSD wavefunction in describing spin polarization is a general feature. The results in addition indicate that the major source of this problem is the absence of the pseudo-doubles within the SD approximation. Although these excitations are not important for the energy, 56 their effect on the iHFCCs is significant

29

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(see also Ref. 85 ). The fact that the pseudo-doubles are the legitimate counterpart of the pure doubles (see Sec. 1) is a special feature of the spin-adapted CC formulation. In spin-orbital based CC theory, by contrast, the analogous excitations correspond to independent triple excitations. Hence, the spin-orbital based CC methods do not suffer from an incomplete treatment of double excitations at the CCSD level and in turn perform much better in this case (see Table 1). It is evident from Table 1 that the flexibility of the spin-adapted CC wavefunction can be effectively restored via the inclusion of triple excitations and a perturbative treatment of the triples is already able to resolve much of the accuracy problem. Let us note an operational difference between the COS-CC and icMRCC approaches. The pseudo-doubles in the icMRCCSDpsT ansatz are treated in a non-perturbative manner in contrast to their perturbative treatment within the COS-CCSD(T) ansatz. The smaller corrections to the iHFCCs introduced by the pseudo-doubles alone in the COS-CCSD(T) ansatz relative to those in the case of icMRCCSDpsT might be attributed to this difference. We also note that the effects of the pseudo-doubles on the carbon and proton iHFCCs are opposite in the case of COS-CCSD(T). However, it is not conclusive from the present results on a small system whether a non-perturbative treatment is really necessary. We think that a perturbative treatment should suffice as long as the genuine triples are also included in the COS-CCSD(T) ansatz. Moreover, the genuine triples are important for the description of spin polarization in complex situations (see Secs. 3.2, 3.3), in which the pseudo-doubles alone may not be sufficient.

3.2

Benchmark calculations

In order to assess the performance of the COS-CC approach in the calculation of HFC tensors, we consider a set of 21 radicals consisting of main group elements. Although these are small systems, they nevertheless represent a broad spectrum of situations, e.g., spin polarization vs direct σ or s orbital contribution as the origin of iHFCC, thus allowing a thorough assessment of the performance of the COS-CC method, in particular the effect of 30

ACS Paragon Plus Environment

Page 30 of 61

Page 31 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

triple excitations. The molecular geometries and the values for the element-specific constant PN employed in our calculations are reported in the Supporting Information. We aim at using sufficiently flexible basis sets for the calculation of the HFC tensors, such that the results are well converged, in this way allowing a meaningful comparison with experiment. Keeping this in mind, we employed the fully uncontracted aug-cc-pCVQZ basis set 87,88 (denoted as ACVQZ-unc) for B-F and Si-Cl. We used the uncontracted aug-cc-pV5Z basis set for hydrogen, and following previous suggestions, 13,89 added five even-tempered steep s functions, 55 whose exponents were obtained by multiplying the largest s exponent in the basis set with 5n (n=1-5) (denoted as APV5Z-unc5). Smaller basis sets have been used only for C3 N, F2 NO, SiF3 , and CH3 SO, for which we employed the corresponding uncontracted triple- and quadruple-zeta bases: ACVTZ-unc for C-F, Si, S, and the APVQZunc5 basis for H. All electrons were correlated in all calculations, as core correlation is important in the calculation of HFC tensors. The iHFCCs for the benchmark set are reported in Tables 2-4, 6 and the components of the spin-dipolar HFC tensors are reported in Table 7. The computed HFC tensors can be compared most reliably with gas phase experimental data. Hence, we compare our results whenever possible with data obtained from microwave spectroscopy in gas phase. However, in the case of BO, CO+ , FCO, NO2 , and SiF3 , gas phase values are available only for B, C, N, and F and in the case of HCO only for H. For other cases, we resort to a comparison with the available EPR data recorded in rare gas matrices (O in BO, 90 C in HCO and FCO, 91 O in CO+ , 92 N, F in F2 NO 93 ), in a methane matrix (O in HCO 94 ), or in a SF6 matrix (N, F in F2 NO, 95 Si in SiF3 96 ), in solid solution of SF6 (O in NO2 97 ), or in single crystals (CH3 SO 98 ). These results differ from the gas phase values and in some cases the differences are substantial, e.g., the fluorine iHFCC in FCO measured in the gas phase 99 and in solid carbon monoxide 100 differ by 80 MHz (see Table 4). Since our aim is to assess the performance of the rigorously spin-adapted COS-CC approach in the calculation of a spin-dependent property, viz., HFC tensors, it is interesting

31

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

to identify whether a molecule under consideration poses a difficult case, e.g., whether the corresponding UHF wavefunction has a moderate or large spin contamination. It is also interesting to investigate how the standard spin-orbital based CC method performs in such cases. For this purpose it is most reasonable to compare COS-CC results with UHF-CC results, which has been employed in most CC investigations of HFC tensors. 10–12,15,89 Although the CC treatment corrects for spin contamination in the UHF reference, it still remains an interesting question as to whether a non-iterative triples correction scheme as in the UHFCCSD(T) approach is adequate for the accuracy of the iHFCCs for difficult cases. For this reason, in Tables 2-4, 6 and 7 we compare our spin-adapted COS-CC results with UHF-CC results. In addition, we report results from the second-order Møller-Plesset perturbation theory (UHF-MP2) calculations. All CC and MP2 calculations incorporate orbital-relaxation effects due to the external perturbation.

3.2.1

Proton HFC tensors

Table 2 indicates that the proton iHFCCs computed at the COS-CCSD(T) level agree reasonably well with experiment, except for the large (∼20 MHz) deviations observed in the case of HCO, CH2 N and HSiS. Previous calculations have shown that corrections to the proton iHFCCs due to vibrational averaging are quite large for HCO and CH2 N, viz., 8 MHz for HCO 77 and 6 MHz for CH2 N. 13 Our calculations do not account for vibrational corrections. However, if the previously suggested values for these corrections are added to the COS-CCSD(T) results, the deviations from experiment become somewhat smaller. For NH2 , PH2 , and CH2 P, on the other hand, vibrational corrections are only of the order of 1-1.5 MHz. 13 Clearly, the agreement between COS-CCSD(T) and experiment is much better in these cases. The improvements to the COS-CCSD results due to triple excitations is obvious from Table 2 apart from the case of OH. For HCO, CH2 N, HSiS, CH2 P, and CH3 SO, the triples corrections are large, about 9-12 MHz, while for NH2 , PH2 , and CH3 , about 3 MHz corrections are observed. In the case of NH2 , PH2 , and CH2 P, the deviations of the

32

ACS Paragon Plus Environment

Page 32 of 61

Page 33 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

COS-CCSD(T) results from experiment are of the same order as the vibrational corrections. Therefore, the computed iHFCCs can be considered well converged with respect to the basis set (see also Table 5). The spin-adapted COS-CCSD(T) and the spin-orbital based UHF-CCSD(T) methods show a similar performance for NH2 and PH2 , while COS-CCSD(T) performs slightly better for HCO, HSiS, and CH3 SO. On the other hand, UHF-CCSD(T) performs better for OH, CH2 N and in particular for CH3 , for which the COS-CCSD(T) result is ∼5 MHz away from the experiment. Interestingly, spin contamination in the UHF reference state (indicated by the values of hSˆ2 iUHF in Tables 2-4 and 6) is not severe for these radicals. The exception is CH2 P, for which hSˆ2 iUHF > 1. Clearly, COS-CCSD(T) performs better than UHF-CCSD(T) in this case. Moreover, the UHF-CCSD(T) result shows a larger deviation from experiment than UHF-CCSD in this case. Table 7 indicates that reasonably good results for the spin-dipolar terms are already obtained at the COS-CCSD level and that the effect of triple excitations is much less pronounced than for the iHFCCs. This is because electron correlation plays a less important role for the accuracy of the spin-dipolar terms 77,82 and therefore, the COS-CCSD results are less affected by the inadequate flexibility of the spin-adapted CC wavefunction. The agreement of the COS-CCSD(T) results with experiment is also quite satisfactory, except for the Tbb − Tcc component in CH2 N, for which UHF-CCSD(T) seems to work slightly better. 3.2.2

Boron, carbon, and nitrogen HFC tensors

Table 3 indicates that the COS-CCSD(T) result for the boron iHFCC in BO, which is the only boron containing system in the present benchmark set, agrees very well with experiment and has only 0.5% deviation from the gas phase value of 1027.4 MHz, while the COS-CCSD result has an error of 27 MHz. This implies that in this case the improvement due to triple excitations is more significant than for the proton iHFCCs. The UHF-CCSD(T) result, on the other hand, is 20 MHz off from experiment and interestingly UHF-CCSD seems to perform

33

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

slightly better. This is another case in which spin contamination in the UHF wavefunction is large (hSˆ2 iUHF = 0.8) and COS-CCSD(T) performs better than UHF-CCSD(T). For the carbon iHFCCs in CH3 , HCO, and CH2 N, the COS-CCSD(T) results agree reasonably well with experiment. In the case of FCO, however, a deviation of 43 MHz from the EPR measurement in argon matrix 91 is observed. The corresponding UHF-CCSD(T) result shows an even larger deviation from experiment. Cochran et al. 91 gave a rough estimate of 110◦ for the FCO bond angle based on the observed 13 C HFC tensor, which is much smaller than our optimized value of 127.7◦ (see Supporting Information). However, our optimized value agrees very well with the experimental bond angle 111 of 127.3◦ . In this case, two factors might be responsible for the observed 43 MHz difference, viz., the matrix effects that are missing in our calculations and the geometric effects, which can also be severe for the iHFCCs. 77 While it is unknown as to which effect is more pronounced, our theoretical results provide a more reliable benchmark than the experimental value for gas phase measurements and for other theoretical methods. For CO+ , the carbon iHFCC values obtained from EPR measurements in neon matrix 92 and gas phase microwave spectroscopy 112 differ by 67 MHz. The gas phase result has a large uncertainty of 15 MHz, 112 which makes it in this case less reliable for comparison. The COS-CCSD(T) result is close to the one from the matrix EPR measurement, though far from the gas phase value. The UHF-CCSD(T) result, although far from both experimental values, nevertheless is closer to the gas phase result than COS-CCSD(T). Interestingly, the UHF-CCSD and COS-CCSD(T) results are close and the inclusion of triple excitations in the UHF-CC ansatz leads to a large deviation between the CCSD and CCSD(T) results (see also Refs. 12,15 ). While large deviations in COS-CCSD results upon inclusion of triple excitations are expected due to the insufficient flexibility of the spin-adapted CCSD wavefunction, a deviation of 88 MHz in the case of UHF-CC (only 36 MHz for COS-CC) is somewhat surprising. Large cluster amplitudes are observed neither in the COS-CC nor in the UHF-CC calculations, which implies that single-reference CC theory should perform well. Considering

34

ACS Paragon Plus Environment

Page 34 of 61

Page 35 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

the large spin contamination in the UHF reference state (hSˆ2 iUHF ≈ 1), the accuracy of the UHF-CCSD(T) result in this case is questionable. Let us now turn to the nitrogen iHFCCs. The best agreement between the COS-CCSD(T) results and experiment is observed for this nucleus except in the case of F2 NO, for which both COS-CCSD(T) and UHF-CCSD(T) results show a deviation of about 10 MHz. For this radical, both experimental results correspond to matrix EPR measurements and the missing matrix effects in our calculations might be responsible for the observed deviation. The triples corrections are roughly of the same order as for the proton iHFCCs except for F2 NO, for which the difference between COS-CCSD and COS-CCSD(T) values is about 25 MHz. The UHF-CCSD(T) approach performs slightly better than COS-CCSD(T) in the case of NH2 . The results are otherwise nearly identical, except for NO and NO2 . In these cases, the UHF-CCSD(T) results show larger deviations from experiment, even though spin contamination in the UHF wavefunction is not large. For the spin-dipolar terms, the effect of triple excitations is, as already observed for the protons, much less pronounced. The COS-CCSD(T) and UHF-CCSD(T) results have much smaller differences compared to the iHFCCs, although COS-CCSD(T) works slightly better in some cases, e.g., C in CO+ , N in NO, C in HCO, etc. The Tcc term for the carbon atom in CH3 as obtained in COS-CCSD(T) calculation differs by 34 MHz from the experimental value. The experimental data were recorded in sodium acetate trihydrate crystals. 113 The computational results thus lack the lattice effects. However, the COS-CCSD(T) and UHFCCSD(T) results compare well and agree well with the icMRCCSDpsT(1e, 1o) result of Samanta and Köhn. 85

3.2.3

Oxygen and fluorine HFC tensors

Table 4 indicates that the oxygen iHFCC in OH and NO computed at the COS-CCSD(T) level agree well with gas phase experimental results. For the other molecules, the reported experimental values correspond to matrix EPR measurements. Considering the missing

35

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

matrix effects in our calculations, the COS-CCSD(T) results for oxygen iHFCCs in BO and CO+ are reasonable. The largest deviation of the COS-CCSD(T) result from experiment is observed in the case of NO2 . While the missing matrix effects account for some of the error, the COS-CCSD(T) result is still somewhat large. The UHF-CCSD(T) result is better in this case as well as for HCO, while the deviations from experiment are larger than for COS-CCSD(T) in the case of NO and CO+ . For the fluorine iHFCCs, we note from Table 4 that the differences between the COSCCSD and COS-CCSD(T) results are larger than those for the other second-row elements and that COS-CCSD(T) overshoots the experimental values, in particular for SiF3 , NF2 and F2 NO. For SiF3 , the COS-CCSD result agrees better with the gas phase experimental data, while for NF2 , FCO the COS-CCSD results are too low. For the second-row atoms, contributions to the spin density arising from spin polarization of the 1s and 2s shells are nearly equal but have opposite signs. 79,80 A proper balance between these two contributions is the key to the accuracy of the computed iHFCC. Let us examine how COS-CC performs in the case of the fluorine atom. The results are reported in Table 5. Clearly, COS-CCSD(T) shows the best performance, although COS-CCSD predicts a by far too small value. This implies that while COS-CCSD suffers in balancing the above two spin polarization contributions, triple excitations are able to restore this balance effectively. Figure 1 shows the spin polarization contributions from the 1s and 2s shells to the net spin density at the F nucleus and indicates that the source of the too small COS-CCSD value for the iHFCC is the underestimation of the 2s shell contribution. Importantly, the 0.02 a.u. difference in spin density amounts to a difference of 84 MHz in the iHFCC, since the constant PN =502.23 MHz/au3 is large for

19

F. This explains the reason for large differences

between the COS-CCSD and COS-CCSD(T) results for fluorine iHFCCs as reported in Table 4. In summary, there is no specific issue with balancing the 1s and 2s spin polarization contributions to the iHFCC at the COS-CCSD(T) level in the case of F atom. Certain molecules, however, may require describing spin polarization in more complex

36

ACS Paragon Plus Environment

Page 36 of 61

Page 37 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

environments than in the atom. In SiF3 , for example, the unpaired electron is localized on the Si atom 121 as indicated by the large silicon iHFCC. The spin density at the F nucleus then arises from spin polarization due to the unpaired electron on Si. The description of spin polarization is further complicated due to the highly ionic nature of the Si-F bonds, 96 which stems from the large electronegativity difference between Si and F. An MR treatment seems to be more suitable in this case for the proper description of spin polarization across charge separated bonds and the single-reference framework of the COS-CC approach seems to fall short. It is gratifying though that COS-CCSD(T) performs much better for the silicon iHFCC (see below). The UHF-CC approach works surprisingly well for the fluorine iHFCC in SiF3 . The difference between UHF-CCSD and CCSD(T) results is much smaller. The spin-orbital formulation seems to favour the description of spin polarization across charge separated bonds even within the single-reference framework. The reason for the overestimation of fluorine iHFCCs in NF2 and F2 NO at the COSCCSD(T) level is less clear, as charge separation in the N-F bonds is expected to be smaller. As already mentioned, the values obtained for the fluorine iHFCC in FCO from microwave and EPR experiments differ by 80 MHz. Since the gas phase result 99 is more reliable for our comparison, the performance of COS-CCSD(T) is much better in this case than for SiF3 , NF2 and F2 NO. The UHF-CCSD(T) result, on the other hand, is far off from the gas phase value. The difference between the UHF-CCSD and UHF-CCSD(T) results is also larger than for SiF3 and NF2 . As with the other elements, the COS-CCSD(T) results show much smaller deviations from experiment for the

17

O spin-dipolar terms. A slightly better performance in comparison to

UHF-CCSD(T) is observed for CO+ as well as for the Tbb − Tcc term in NO. Similar to the case of the fluorine iHFCCs, the COS-CCSD(T) results deviate largely from experiment also for the spin-dipolar terms, though the result is quite satisfactory in the case of SiF3 .

37

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

3.2.4

HFC tensors for the third-row elements

The accurate prediction of iHFCCs for the third-row atoms is much more challenging than for the second-row atoms for both wavefunction based methods as well as DFT. 122 In this case, a subtle balance of spin polarization contributions from 1s and 2s shells as well as from the additional 3s shell is required, which is extremely difficult to achieve. Figure 1 shows the various spin polarization contributions in the case of Cl atom. The contributions from the 1s and 2s shells are of opposite signs similar to the second-row atoms, although their magnitudes are quite different. The 3s shell contribution is again positive. Figure 1 indicates that COS-CCSD accumulates errors in all three contributions relative to COSCCSD(T), although the 0.04 a.u. difference in the net spin density at the Cl nucleus (a 17.5 MHz difference in the iHFCC) arises mainly from the underestimation of the positive 2s and 3s shell contributions. For molecules containing third-row elements, the accurate prediction of iHFCC is somewhat simpler when the unpaired electron is in a σ orbital with substantial 3s character, while the case of π radicals is more challenging. 123 Apart from SiF3 and HSiS, in all other radicals considered in our benchmark study, the singly occupied molecular orbital (SOMO) is a π orbital. Table 6 indicates that COS-CCSD suffers severely in these cases in balancing the above spin polarization contributions. The P, S, and Cl iHFCCs are exceedingly small. However, triple excitations bring in tremendous improvements. It is evident from Table 6 that for all π radicals the COS-CCSD(T) results agree very well with experiment. These results demonstrate the effectiveness of the triples correction in amending the shortcoming due to the inadequate flexibility of the spin-adapted CCSD wavefunction. For the silicon iHFCC in SiF3 , the COS-CCSD result is already reasonable and the improvement due to triples correction is relatively small. The somewhat large deviation from experiment may be attributed to the missing matrix effect in our calculation. The result is nevertheless much more satisfactory than the previous DFT result. 124 For HSiS, the sulfur iHFCC obtained at the COS-CCSD level is again somewhat low. The silicon iHFCC 38

ACS Paragon Plus Environment

Page 38 of 61

Page 39 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

is, on the other hand, much larger in comparison to the UHF-CCSD value. Here, COSCCSD seems to struggle in describing the distribution of the unpaired electron spin among two third-row atoms. The triples correction improves the sulfur iHFCC, though its effect on the silicon iHFCC is small. Unfortunately, no experimental results are available for the Si and S iHFCCs in HSiS. For the spin-dipolar terms, the differences between the COS-CCSD and COS-CCSD(T) results are again much smaller. The UHF-CCSD(T) approach performs slightly better in certain cases, although for CH2 P the differences between the UHF-CCSD and CCSD(T) results are large. In general, the agreement between both COS-CCSD(T) and UHF-CCSD(T) results and experiment is less satisfactory than for the second-row elements. This might be attributed to deficiencies in the basis set. Further steep p (possibly also d) functions might be helpful for attaining a better agreement. The mean absolute errors (MAEs) of the iHFCCs for the various nuclei relative to experiment as obtained in spin-adapted COS-CC and the spin-orbital based UHF-MP2 and UHF-CC calculations are shown in Figure 2. For proton iHFCCs, all methods show a similar performance and the MAEs are rather small. The improvement of the COS-CC results upon inclusion of triple excitations is clear for all nuclei except

17

O and

19

F. Figure 2 also

indicates that, on an average, COS-CCSD(T) performs better than UHF-CCSD(T) in the case of boron, carbon, nitrogen, and in particular for the third-row elements. Let us also note from Tables 2-4, 6 that the errors in the UHF-MP2 results for the iHFCCs relative to experiment can be up to an order of magnitude larger than those for any CC method and the signs of the iHFCCs are also wrong in many cases.

3.3

HFC tensors for the 2 Σ state of C3 N molecule

The ground 2 Σ state of the linear C3 N molecule is an interesting test case for critically assessing the performance of the spin-adapted COS-CC approach vis-à-vis spin-orbital based CC theory because of the large spin contamination in the UHF wavefunction: hSˆ2 iUHF =1.623 39

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(see Table 8). Accurate experimental data for the HFC tensors at all four nuclei are available from gas phase Fourier transform microwave spectroscopy, 126–128 which makes a direct comparison with experiment feasible. The HFC tensors obtained in the COS-CC and UHF-CC calculations are reported in Table 8. The experimental iHFCCs indicate that the unpaired electron is localized mainly on the terminal C atom and to some extent on the neighbouring one. 126 The large positive iHFCCs at these two C nuclei indicate a large 2s character of the SOMO. The small iHFCC at the C nuclus next to N, which is almost an order of magnitude smaller than that at the preceding C nucleus, provides further evidence for the localization of the unpaired electron away from this C atom. The negative nitrogen iHFCC, which points to a negative spin density (the magnetic moment for

14

N is positive, see Supporting Information), results from

spin polarization effects. The unpaired electron accumulates α spin population around the two terminal C atoms, in turn accumulating β spin population at N. Table 8 indicates that the iHFCCs for all four nuclei calculated at both COS-CCSD and COS-CCSD(T) levels agree well with experiment except for the somewhat large value at the terminal C atom. The differences between the COS-CCSD and COS-CCSD(T) results are of the same order as observed for the carbon and nitrogen iHFCCs in other molecules (see Sec. 3.2). Let us specifically note the excellent agreement, both in magnitude and sign, of the COS-CCSD(T) result for nitrogen iHFCC with the experiment. For the UHF-CC results, one first notices that the large spin contamination in the UHF reference state is eliminated neither entirely at the CCSD nor at the CCSD(T) level. The reported UHF-CCSD value for hSˆ2 i is astonishingly large and the hSˆ2 i value deviates from the ideal value of 0.75 even at the UHF-CCSD(T) level. Nevertheless, UHF-CCSD manages to predict the iHFCCs reasonably well except for the value at the terminal C nucleus, which is overestimated to a somewhat larger extent than in the case of COS-CCSD. The inclusion of triple excitations in the UHF-CC ansatz, however, has an adverse effect. While the iHFCC at the terminal C nucleus is improved relative to UHF-CCSD, the iHFCC at the neighbouring

40

ACS Paragon Plus Environment

Page 40 of 61

Page 41 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

C nucleus becomes worse. The most severe impact is observed in the case of the nitrogen iHFCC, which not only changes sign but also becomes quite large in magnitude. This result is even qualitatively incorrect. The differences between the UHF-CCSD and UHF-CCSD(T) values are larger than those between COS-CCSD and COS-CCSD(T). This implies that triple excitations attempt to attenuate the residual spin contamination in the UHF-CCSD wavefunction. Clearly, the improvement in the quality of the wavefunction does not work in unison with the improvement of the iHFCCs. Although the relative signs of the iHFCCs cannot be experimentally determined, they nevertheless carry valuable information about the physical effects from which the spin density at a particular nucleus in a molecule originates. For example, a negative spin density is an unambiguous indication of spin polarization effects. Hence, a complete elucidation of the electronic structure from the observed spectra is not possible without the signs of the iHFCCs. As EPR experiments are unable to determine the signs, quantum-chemical calculations should provide accurate information about them. Interestingly, Fessenden has shown that it is indeed possible to determine the relative signs with fair certainty via a detailed analysis of the line positions in an EPR spectrum 129 when the iHFCCs are reasonably large. Since the spin-dipolar terms are less affected by the electron correlation treatment, COSCCSD(T) and UHF-CCSD(T) show a similar performance for the Taa term. While the role of triple excitations in the spin-orbital based CC ansatz is to eliminate errors due to spin contamination, the results for C3 N clearly demonstrate that when spin contamination in the UHF reference state is large, the effect of an approximate treatment of the triples may become problematic for spin-dependent properties. The rigorously spinadapted COS-CCSD(T) approach, by contrast, provides reliable results in these calculations.

3.4

Discussion on the importance of spin adaptation

The purpose of developing a rigorously spin-adapted CC theory is to avoid spin contamination in the electron correlation treatment and to gain as a consequence reliable accuracy in 41

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

the calculation of energies and molecular properties, the case of spin-dependent properties being particularly important. However, to arrive at a definite conclusion about the importance of spin adaptation is far from trivial. While extensive benchmarking and comparison of the results with experiment and other computational methods is the obvious way, special care should be exercised to avoid unwarranted conclusions. This is especially important for the iHFCCs, the accurate theoretical prediction of which poses different challenges based on the source of the spin density at a nucleus. For example, the accurate prediction of iHFCCs becomes far more challenging when spin polarization is the governing mechanism compared to the cases where the spin density at a nucleus arises from a direct σ or s orbital contribution. Hence, the consistent performance of a method for a wide range of elements and in diverse situations should be given precedence over the results for individual nuclei. Our results indicate that triple excitations are important in the COS-CC approach for an accurate description of spin polarization, while COS-CCSD often falls short. The crux of the problem of modelling spin polarization effects is the subtle balance of the contributions from various atomic shells, which usually have opposite signs, to the net spin density at a nucleus (see Figure 1). The situation becomes further complicated when the unpaired electron in a molecule is localized on an atom different from the target nucleus. It is not entirely surprising that a spin-adapted CC approach, which models spin polarization purely by means of electron correlation (i.e., configuration mixing), requires higher excitations, e.g., triples, in order to balance such intricate physical effects. The extensive literature on CI 77–83 also points to the importance of triple excitations for the accuracy of the computed iHFCCs. Spin polarization is related to the preferential alignment of electron spins in the vicinity of unpaired electron(s). The use of spin symmetry breaking excitations in spin-orbital based CC theory seems to facilitate the description of this effect owing to its very nature. This leads to a better performance of spin-orbital based CCSD in certain cases where COS-CC requires triple excitations. However, it is not guaranteed that spin-orbital based CC theory is always able to balance the underlying physical effects within the SD approximation, e.g.,

42

ACS Paragon Plus Environment

Page 42 of 61

Page 43 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

for third-row atoms the treatment of triple excitations becomes mandatory for attaining a discernible accuracy. 122 The same is expected also for transition metal elements. Thus, triple excitations are important in both cases. A matter of concern is then the computational cost that one has to pay for this purpose. Our benchmark calculations involving a wide range of main group elements demonstrate that for the spin-adapted COS-CC approach, (a) the triple excitations attain consistent improvements over the COS-CCSD results unless a single-reference description falls insufficient, (b) the results are qualitatively correct, i.e., the signs of iHFCCs are correct, and (c) the most widely used noniterative CCSD(T) approximation works well. With the UHF-CC approach, on the other hand, one occasionally encounters cases where the customary CCSD(T) scheme has an unreasonably adverse effect on the iHFCCs, e.g., the carbon iHFCC in CO+ (also the boron iHFCC in BO) and the nitrogen iHFCC in C3 N. In the former case, the iHFCC arises from a direct contribution from the σ SOMO, while in the latter case spin polarization is the mechanism. The UHF-CCSD(T) results are poor in both cases, although the effect is more severe for C3 N where the result is even qualitatively wrong. In all theses cases, the UHF reference is heavily spin contaminated. These results are suggestive of a more frequent occurrence of similar or more severe errors in UHF-CCSD(T) results for transition metal compounds, especially when spin polarization is the major source of the iHFCCs. The reasonable performance of UHF-CCSD in the case of C3 N and CO+ (also BO) implies that a full (iterative) treatment of triples within the UHF-CCSDT scheme might attain the desired accuracy. However, this makes the spin-orbital CC approach too expensive for routine applications in the calculation of HFC tensors. In summary, while triple excitations are important for an accurate prediction of iHFCCs, it is doubtful whether the computationally affordable CCSD(T) approximation is an infallible option for the UHF-CC method. The spin-adapted COS-CC approach, on the other hand, provides accurate results within the cost-effective, noniterative CCSD(T) approximation. From a computational standpoint, this is clearly a great advantage of rigorous spin

43

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

adaptation in CC theory. In view of the prevailing applications of the DFT and MP2 methods for the calculation of HFC tensors still today, which suffer due to spin contamination, 8,9 it is encouraging to note, as evident from Figure 2, that the spin-adapted COS-CCSD(T) approach performs significantly better than one of these methods, viz., UHF-MP2. Our results also argue in favour of a spin-adapted CC approach approach as the more reliable method for benchmarking DFT and MP2 results for HFC tensors than UHF-CC theory.

4

Conclusions

We have presented a perturbative triples correction scheme for the unitary group based spinadapted COS-CC ansatz, which has the following important properties: (a) computation of the triples correction to the energy requires a single noniterative step with N 7 scaling of the computational cost, and (b) the method conserves the orbital invariance properties of the underlying COS-CCSD ansatz. The first-order molecular properties at the COSCCSD(T) level are computed via implementation of an analytic derivative scheme. The one-particle spin density matrix required for the calculation of spin-dependent properties is obtained within our spin-free formulation by employing the highly useful spin density calculation procedure proposed by Luzanov and co-workers. 51–54 Although this procedure is somewhat involved in terms of implementation, as it requires the two-particle density matrix, this is straightforward for us due to the availability of a full-blown analytic derivative implementation. 46 As rigorous spin adaptation is our objective, we do not intend to introduce any spin-dependent operator. The reported benchmark calculations amply demonstrate the validity and merits of our strategy. We have presented an extensive benchmark study on the calculation of HFC tensors for molecules containing main group elements. The results clearly demonstrate that the shortcoming of the COS-CCSD method in describing spin polarization is effectively overcome

44

ACS Paragon Plus Environment

Page 44 of 61

Page 45 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

upon introducing perturbative triples corrections. Furthermore, a comparison of our results with benchmark FCI and icMRCCSD(1e, 1o) results 85 demonstrates that the insufficient flexibility of the COS-CCSD ansatz is not an idiosyncratic feature but is generally true for spin-adapted CC (also CI) approaches. The COS-CCSD(T) approach exhibits a consistently good performance in predicting the HFC tensors for a wide range of main group elements, the third-row elements in particular, for which the accurate calculation of iHFCCs poses high challenges. As with any other method, there are a few exceptions, viz., the fluorine iHFCCs in certain radicals, for which COS-CCSD(T) overshoots the experimental values. Importantly, however, COS-CCSD(T) predicts the HFC tensor components with correct signs and qualitatively wrong results were not encountered within the present benchmark study. In summary, COS-CCSD(T) provides reliable accuracy in these calculations as long as a single-reference description suffices to capture the essential spin polarization effect. The scope of a unitary group based CC theory, viz. the COS-CC method, has been extended in our work to the calculation of spin-dependent molecular properties via implementation of efficient analytic schemes similar to the ones available for spin-orbital based CC theory. Our results demonstrate that within the cost-effective, noniterative CCSD(T) approximation, which is the broadly acclaimed “gold standard” in quantum chemistry, the spin-adapted COS-CC approach offers higher accuracy and reliability in the calculation of HFC tensors than UHF-CC theory. Our work also raises the hope that, with the availability of an implementation for open-shell systems with arbitrary spin multiplicities, the unitary group based spin-adapted CC approach will excel in demanding applications such as the accurate prediction of magnetic properties.

Acknowledgement D.D. expresses his sincere gratitude to Prof. Frank Neese and the Max-Planck Gesellschaft for generous support during the final phase of this work. In Mainz, this work was supported

45

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

by funding from the Deutsche Forschungsgemeinschaft (grant no. GA 370/5-1) to J.G.

Supporting Information Available Tables containing optimized geometries of the radicals used in the benchmark calculations of hyperfine coupling tensors and the values of the constant PN = βe βN ge gN for various nuclei. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Atherton, N. M. Principles of Electron Spin Resonance; Prentice Hall: New York, 1993. (2) Weil, J. A.; Bolton, J. R.; Wertz, J. E. Electron Paramagnetic Resonance: Elementary Theory and Practical Applications; Willey & Sons: New York, 1994. (3) Weltner, W., Jr. Magnetic Atoms and Molecules; Van Nostrand: New York, 1983. (4) Abragam, A.; Pryce, M. H. L. Theory of the nuclear hyperfine structure of paramagnetic resonance spectra in crystals. Proc. R. Soc. Lond. A Math. Phys. Sci. 1951, 205, 135–153. (5) Barone, V.; Polimeno, A. Toward an integrated computational approach to CW-ESR spectra of free radicals. Phys. Chem. Chem. Phys. 2006, 8, 4609–4629. (6) Bagus, P. S.; Liu, B.; Schaefer, H. F., III. Study of the contact-term contribution to the hyperfine structure obtained from spin-unrestricted Hartree-Fock wave functions. Phys. Rev. A 1970, 2, 555–560. (7) Munzarová, M. In Calculation of NMR and EPR parameters, Theory and applications; Kaupp, M., Bühl, M., Malkin, V. G., Eds.; Wiley VCH: Weinheim, 2004; pp 463–482. (8) Munzarová, M.; Kaupp, M. A Critical validation of density functional and coupled-cluster approaches for the calculation of EPR hyperfine coupling constants in transition metal complexes. J. Phys. Chem. A 1999, 103, 9966–9983. (9) Schattenberg, C. J.; Maier, T. M.; Kaupp, M. Lessons from the Spin-Polarization/Spin-Contamination Dilemma of Transition-Metal Hyperfine Couplings for the Construction of Exchange-Correlation Functionals. J. Chem. Theory Comput. 2018, 14, 5653–5672. (10) Carmichael, I. Ab initio coupled-cluster calculations of isotropic hyperfine splitting in some diatomic hydrides. J. Phys. Chem. 1990, 94, 5734–5740. (11) Perera, S. A.; Watts, J. D.; Bartlett, R. J. A theoretical study of hyperfine coupling constants. J. Chem. Phys. 1994, 100, 1425–1434. Perera, S. A.; Salemi, L. M.; Bartlett, R. J. Hyperfine coupling constants of organic radicals. J. Chem. Phys. 1997, 106, 4061–4066. (12) Verma, P.; Perera, A.; Morales, J. A. Massively parallel implementations of coupled-cluster methods for electron spin resonance spectra. I. Isotropic hyperfine coupling tensors in large radicals. J. Chem. Phys. 2013, 139, 174103.

46

ACS Paragon Plus Environment

Page 46 of 61

Page 47 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(13) Puzzarini, C.; Barone, V. Toward spectroscopic accuracy for open-shell systems: Molecular structure and hyperfine coupling constants of H2 CN, H2 CP, NH2 , and PH2 as test cases. J. Chem. Phys. 2010, 133, 184301. (14) Puzzarini, C.; Barone, V. Toward spectroscopic accuracy for organic free radicals: Molecular structure, vibrational spectrum, and magnetic properties of F2 NO. J. Chem. Phys. 2008, 129, 084306. (15) Kossmann, S.; Neese, F. Correlated ab initio spin densities for larger molecules: Orbital-optimized spin-component-scaled MP2 Method. J. Phys. Chem. A 2010, 114, 11768–11781. (16) Knowles, P. J.; Hampel, C.; Werner, H.-J. Coupled cluster theory for high spin, open shell reference wave functions. J. Chem. Phys. 1993, 99, 5219–5227. (17) Neogrády, P.; Urban, M.; Hubač, I. Spin adapted restricted Hartree-Fock reference coupled cluster theory for open shell systems. J. Chem. Phys. 1994, 100, 3706–3716. (18) Szalay, P. G.; Gauss, J. Spin-restricted open-shell coupled-cluster theory. J. Chem. Phys. 1997, 107, 9028–9038. (19) Heckert, M.; Heun, O.; Gauss, J.; Szalay, P. G. Towards a spin-adapted coupled-cluster theory for high-spin open-shell states. J. Chem. Phys. 2006, 124, 124105. (20) Paldus, J. Group theoretical approach to the configuration interaction and perturbation theory calculations for atomic and molecular systems. J. Chem. Phys. 1974, 61, 5321–5330. Paldus, J. In Theoretical Chemistry: Advances and Perspectives, Eyring, H.; Henderson, D. J., Eds.; Academic Press: New-York, 1976; Vol. 2, pp 131-290. (21) Paldus, J. In Electrons in Finite and Infinite Structures; Phariseau, P., Scheire, L., Eds.; Plenum: New York, 1977; pp 411–429. (22) Matsen, F. A.; Pauncz, R. The Unitary Group in Quantum Chemistry; Elsevier: Amsterdam, 1986. (23) Paldus, J.; Li, X. In Symmetries in Science VI: From the Rotation Group to Quantum Algebras; Gruber, B., Ed.; Plenum: New York, 1993; pp 573–591. Li, X.; Paldus, J. Multiconfigurational spinadapted single-reference coupled cluster formalism. Int. J. Quantum Chem. 1993, 48, 269–285. (24) Li, X.; Paldus, J. Automation of the implementation of spin-adapted open-shell coupled-cluster theories relying on the unitary group formalism. J. Chem. Phys. 1994, 101, 8812–8826. (25) Li, X.; Paldus, J. In Recent Advances in Computational Chemistry; Bartlett, R. J., Ed.; World Scientific: Singapore, 1997; Vol. 3, pp 183-219. (26) Li, X.; Paldus, J. Unitary-group-based open-shell coupled-cluster method with corrections for connected triexcited clusters. I. Theory. Int. J. Quantum Chem. 1998, 70, 65–75. (27) Li, X.; Paldus, J. Unitary-group-based open-shell coupled-cluster method with corrections for connected triexcited clusters. II. Applications. Mol. Phys. 1998, 94, 41–54. (28) Jankowski, P.; Jeziorski, B. Unitary group based open-shell coupled cluster theory: Application to van der Waals interactions of high-spin systems. J. Chem. Phys. 1999, 111, 1857–1869. (29) Jeziorski, B.; Paldus, J.; Jankowski, P. Unitary group approach to spin-adapted open-shell coupled cluster theory. Int. J. Quantum Chem. 1995, 56, 129–155. (30) Datta, D.; Mukherjee, D. A compact spin-free combinatoric open-shell coupled cluster theory applied to single-reference doublets. Int. J. Quantum Chem. 2008, 108, 2211–2222.

47

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(31) Datta, D.; Mukherjee, D. An explicitly spin-free compact open-shell coupled cluster theory using a multireference combinatoric exponential ansatz: Formal development and pilot applications. J. Chem. Phys. 2009, 131, 044124. (32) Datta, D.; Mukherjee, D. The spin-free analogue of Mukherjee’s state-specific multireference coupled cluster theory. J. Chem. Phys. 2011, 134, 054122. (33) Maitra, R.; Sinha, D.; Mukherjee, D. Unitary group adapted state-specific multi-reference coupled cluster theory: Formulation and pilot numerical applications. J. Chem. Phys. 2012, 137, 024105. (34) Sen, S.; Shee, A.; Mukherjee, D. Formulation and implementation of a unitary group adapted state universal multi-reference coupled cluster (UGA-SUMRCC) theory: Excited and ionized state energies. J. Chem. Phys. 2012, 137, 074104. (35) Datta, D.; Gauss, J. A non-antisymmetric tensor contraction engine for the automated implementation of spin-adapted coupled cluster approaches. J. Chem. Theory Comput. 2013, 9, 2639–2653. (36) Scheiner, A. C.; Scuseria, G. E.; Rice, J. E.; Lee, T. J.; Schaefer, H. F., III. Analytic evaluation of energy gradients for the single and double excitation coupled cluster (CCSD) wave function: Theory and application. J. Chem. Phys. 1987, 87, 5361–5373. Scuseria, G. E. Analytic evaluation of energy gradients for the singles and doubles coupled cluster method including perturbative triple excitations: Theory and applications to FOOF and Cr2 . J. Chem. Phys. 1991, 94, 442–447. Gauss, J.; Stanton, J. F. Analytic gradients for the coupled-cluster singles, doubles, and triples (CCSDT) model. J. Chem. Phys. 2002, 116, 1773–1782. (37) Gauss, J.; Stanton, J. F.; Bartlett, R. J. Coupled-cluster open-shell analytic gradients: Implementation of the direct product decomposition approach in energy gradient calculations. J. Chem. Phys. 1991, 95, 2623-2638. (38) Watts, J. D.; Gauss, J.; Bartlett, R. J. Open-shell analytical energy gradients for triple excitation many-body, coupled-cluster methods: MBPT(4), CCSD+T(CCSD), CCSD(T),and QCISD(T). Chem. Phys. Lett. 1992, 200, 1–7. (39) Gauss, J.; Lauderdale, W. J.; Stanton, J. F.; Watts, J. D.; Bartlett, R. J. Analytic energy gradients for open-shell coupled-cluster singles and doubles (CCSD) calculations using restricted open-shell Hartree-Fock (ROHF) reference functions. Chem. Phys. Lett. 1991, 182, 207–215. (40) Watts, J. D.; Gauss, J.; Bartlett, R. J. Coupled-cluster methods with noniterative triple excitations for restricted open-shell Hartree-Fock and other general single determinant reference functions. Energies and analytical gradients. J. Chem. Phys. 1993, 98, 8718–8733. (41) Gauss, J.; Stanton, J. F.; Bartlett, R. J. Analytic evaluation of energy gradients at the coupled-cluster singles and doubles level using quasi-restricted Hartree-Fock open-shell reference functions. J. Chem. Phys. 1991, 95, 2639–2645. (42) Kállay, M.; Gauss, J.; Szalay, P. G. Analytic first derivatives for general coupled-cluster and configuration interaction models. J. Chem. Phys. 2003, 119, 2991–3004. (43) Paldus, J.; Li, X. Calculation of static molecular properties in the framework of the unitary group based coupled cluster approach. Can. J. Chem. 1996, 74, 918–930. (44) Li, X.; Paldus, J. A unitary group based open-shell coupled cluster study of vibrational frequencies in ground and excited states of first row diatomics. J. Chem. Phys. 1996, 104, 9555–9562. (45) Gauss, J. In Modern Methods and Algorithms of Quantum Chemistry; Grotendorst, J., Ed.; John von Neumann Institute for Computing: Jülich, Germany, 2000; Vol. 3, pp 541-592.

48

ACS Paragon Plus Environment

Page 48 of 61

Page 49 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(46) Datta, D.; Gauss, J. Analytic first derivatives for a spin-adapted open-shell coupled cluster theory: Evaluation of first-order electrical properties. J. Chem. Phys. 2014, 141, 104102. (47) Gould, M. D.; Paldus, J. Spin-dependent unitary group approach. I. General formalism. J. Chem. Phys. 1990, 92, 7394–7401. (48) Gould, M. D.; Paldus, J.; Chandler, G. S. Unitary group approach to reduced density matrices. J. Chem. Phys. 1990, 93, 4142–4153. (49) Gould, M. D.; Chandler, G. S. Unitary group approach to the many-electron problem. I. Matrix element evaluation and shift operators. Int. J. Quantum Chem. 1984, 25, 553–601. Gould, M. D.; Chandler, G. S. Unitary group approach to the many-electron problem. II. Adjoint tensor operators for U (n). Int. J. Quantum Chem. 1984, 25, 603–633. Gould, M. D.; Chandler, G. S. Unitary group approach to the many-electron problem. III. Matrix elements of spin-dependent Hamiltonians. Int. J. Quantum Chem. 1984, 25, 1089–1109. (50) Gould, M. D.; Battle, J. S. Spin-dependent unitary group approach. II. Derivation of matrix elements for spin-dependent operators. J. Chem. Phys. 1993, 99, 5961–5975. (51) Luzanov, A. V. Calculating spin density in the quantum-chemical unitary formalism. Theor. Exp. Chem. 1985, 21, 329–331. (52) Vaiman, G. E.; Luzanov, A. V.; Mestechkin, M. M. Separation of spin and the fock coordinate wave function method in the N-representability problem. Theor. Math. Phys. 1976, 28, 634–643. (53) Luzanov, A. V.; Whyman, G. E. Structure and spin-purity conditions for reduced density matrices of arbitrary order. Int. J. Quantum Chem. 1981, 20, 1179–1199. (54) Luzanov, A. V. Spin-free quantum chemistry: What one can gain from fock’s cyclic symmetry. Int. J. Quantum Chem. 2011, 111, 4042–4066. (55) Datta, D.; Gauss, J. Communication: Spin densities within a unitary group based spin-adapted openshell coupled-cluster theory: Analytic evaluation of isotropic hyperfine-coupling constants for the combinatoric open-shell coupled-cluster scheme. J. Chem. Phys. 2015, 143, 011101. (56) Li, X.; Paldus, J. Unitary group based state-selective coupled-cluster method: Comparison of the first order interacting space and the full single and double excitation space approximations. J. Chem. Phys. 1995, 102, 8897–8905. (57) Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. A fifth-order perturbation comparison of electron correlation theories. Chern. Phys. Lett. 1989, 157, 479–483. Bartlett, R. J.; Watts, J. D.; Kucharski, S. A.; Noga, J. Non-iterative fifth-order triple and quadruple excitation energy corrections in correlated methods. Chem. Phys. Lett. 1990, 165, 513–522. (58) Shavitt, I. Graph theoretical concepts for the unitary group approach to the many-electron correlation problem. Int. J. Quantum Chem. 1977, 12, 131–148. (59) Lindgren, I. A coupled-cluster approach to the many-body perturbation theory for open-shell systems. Int. J. Quantum Chem. 1978, 14, 33-58. (60) Scuseria, G. E. The open-shell restricted Hartree-Fock singles and doubles coupled-cluster method including triple excitations CCSD(T): application to C+ 3 . Chem. Phys. Lett. 1991, 176, 27–35. (61) Handy, N. C.; Pople, J. A.; Head-Gordon, M.; Raghavachari, K.; Trucks, G. W. Size-consistent Brueckner theory limited to double substitutions. Chem. Phys. Lett. 1989, 164, 185–192.

49

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(62) Matthews, D. A.; Gauss, J.; Stanton, J. F. Revisitation of nonorthogonal spin adaptation in coupled cluster theory. J. Chem. Theory Comput. 2013, 9, 2567–2572. (63) Wick, G. C. The Evaluation of the Collision Matrix. Phys. Rev. 1950, 80, 268–272. (64) Kucharski, S. A.; Bartlett, R. J. Noniterative energy corrections through fifth-order to the coupled cluster singles and doubles method. J. Chem. Phys. 1998, 108, 5243–5254. (65) Taube, A. G.; Bartlett, R. J. Improving upon CCSD(T): ΛCCSD(T). I. Potential energy surfaces. J. Chem. Phys. 2008, 128, 044110. (66) Helgaker, T. Simple derivation of the potential energy gradient for an arbitrary electronic wave function. Int. J. Quantum Chem. 1982, 21, 939–940. (67) Koch, H.; Jensen, H. J. Aa.; Jørgensen, P.; Helgaker, T.; Scuseria, G. E.; Schaefer, H. F., III. Coupled cluster energy derivatives. Analytic Hessian for the closed-shell coupled cluster singles and doubles wave function: Theory and applications. J. Chem. Phys. 1990, 92, 4924–4940. (68) Koch, H.; Jørgensen, P. Coupled cluster response functions. J. Chem. Phys. 1990, 93, 3333–3344. (69) Szalay, P. G. Analytic energy derivatives for coupled-cluster methods describing excited states: General formulas and comparison of computational costs. Int. J. Quantum Chem. 1995, 55, 151–163. (70) Arponen, J. Variational principles and linked-cluster exp S expansions for static and dynamic manybody problems. Ann. Phys. 1983, 151, 311-382. (71) Salter, E. A.; Trucks, G. W.; Bartlett, R. J. Analytic energy derivatives in many-body methods. I. First derivatives. J. Chem. Phys. 1989, 90, 1752–1766. (72) Lee, T. J.; Rendell, A. P. Analytic gradients for coupled-cluster energies that include noniterative connected triple excitations: Application to cis- and trans-HONO. J. Chem. Phys. 1991, 94, 6229–6236. (73) Handy, N. C.; Schaefer, H. F., III. On the evaluation of analytic energy derivatives for correlated wave functions. J. Chem. Phys. 1984, 81, 5031–5033. (74) Stanton, J. F.; Gauss, J.; Cheng, L.; Harding, M. E.; Matthews, D. A.; Szalay, P. G.; Auer, A. A.; Bartlett, R. J.; Benedikt, U.; Berger, C.; Bernholdt, D. E.; Bomble, Y. J.; Christiansen, O.; Engel, F.; Faber, R.; Heckert, M.; Heun, O.; Hilgenberg, M.; Huber, C.; Jagau, T.-C.; Jonsson, D.; Jusélius, J.; Kirsch, T.; Klein, K.; Lauderdale, W. J.; Lipparini, F.; Metzroth, T.; Mück, L. A.; O’Neill, D. P.; Price, D. R.; Prochnow, E.; Puzzarini, C.; Ruud, K.; Schiffmann, F.; Schwalbach, W.; Simmons, C.; Stopkowicz, S.; Tajti, A.; Vázquez, J.; Wang, F.; Watts, J. D. CFOUR: Coupled cluster techniques for computational chemistry, with the integral packages MOLECULE (Almlöf, J.; Taylor, P. R.), PROPS (Taylor, P. R.), ABACUS (Helgaker, T.; Jensen, H. J. Aa.; Jørgensen, P.; Olsen, J.), and ECP routines by Mitin, A. V.; van Wüllen, C. For the current version, see http://www.cfour.de. (75) Gidofalvi, G.; Shepard, R. The evaluation of spin-density matrices within the graphically contracted function method. Int. J. Quantum Chem. 2009, 109, 3552–3563. (76) Shiozaki, T.; Yanai, T. Hyperfine coupling constants from internally contracted multireference perturbation theory. J. Chem. Theor. Comput. 2016, 12 4347–4351. (77) Feller, D.; Davidson, E. R. Ab initio configuration interaction calculations of the hyperfine structure in small radicals. J. Chem. Phys. 1984, 80, 1006–1017; Feller, D.; Davidson, E. R. Difficulties in ab initio CI calculations of the hyperfine structure of small radicals. Theor. Chim. Acta 1985, 68, 57–67. (78) Feller, D.; Davidson, E. R. A multireference CI determination of the isotropic hyperfine constants for first row atoms B–F. J. Chem. Phys. 1988, 88, 7580–7587.

50

ACS Paragon Plus Environment

Page 50 of 61

Page 51 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(79) Engels, B.; Peyerimhoff, S. D.; Davidson, E. R. Calculation of hyperfine coupling constants. An ab initio MRD-CI study for nitrogen to analyse the effects of basis sets and CI parameters. Mol. Phys. 1987, 62, 109–127. (80) Engels, B.; Peyerimhoff, S. D. Study of the 1s and 2s shell contributions to the isotropic hyperfine coupling constant in nitrogen. J. Phys. At. Mol. Opt. Phys. 1988, 21, 3459–3471. (81) Engels, B.; Eriksson, L.; Lunell, S. Recent Developments in Configuration Interaction and Density Functional Theory Calculations of Radical Hyperfine Structure. Adv. Quantum Chem. 1996, 27, 297–369. (82) Funken, K.; Engels, B.; Peyerimhoff, S. D.; Grein, F. Study of the hyperfine coupling constants of the molecules NH2 , NHD and ND2 . Chem. Phys. Lett. 1990, 172, 180–186. (83) Engels, B. In Calculation of NMR and EPR parameters, Theory and applications; Kaupp, M., Bühl, M., Malkin, V. G., Eds.; Wiley VCH: Weinheim, 2004; pp 483–492. (84) Hanauer, M.; Köhn, A. Pilot applications of internally contracted multireference coupled cluster theory, and how to choose the cluster operator properly. J. Chem. Phys. 2011, 134, 204111. (85) Samanta, P. K.; Köhn, A. First-order properties from internally contracted multireference coupledcluster theory with particular focus on hyperfine coupling tensors. J. Chem. Phys. 2018, 149, 064101. (86) Barone, V. In Recent Advances in Density Functional Methods; Chong, D. P., Ed.; World Scientific; 1995; Part 1; pp. 287-334. (87) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. Electron affinities of the first-row atoms revisited. Systematic basis sets and wave functions. J. Chem. Phys. 1992, 96, 6796–6806. Woon, D. E.; Dunning, T. H., Jr. Gaussian basis sets for use in correlated molecular calculations. III. The atoms aluminum through argon. J. Chem. Phys. 1993, 98, 1358–1371. (88) Woon, D. E.; Dunning, T. H., Jr. Gaussian basis sets for use in correlated molecular calculations. V. Core-valence basis sets for boron through neon. J. Chem. Phys. 1995, 103, 4572–4585. Peterson, K. A.; Dunning, T. H., Jr. Accurate correlation consistent basis sets for molecular core–valence correlation effects: The second row atoms Al–Ar, and the first row atoms B–Ne revisited. J. Chem. Phys. 2002, 117, 10548–10560. (89) Fau, S.; Bartlett, R. J. Gaussian Basis Sets for Highly Accurate Calculations of Isotropic Hyperfine Coupling Constants at Hydrogen. J. Phys. Chem. A 2003, 107, 6648–6655. Al Derzi, A. R.; Fau, S.; Bartlett, R. J. Benchmark Study of Isotropic Hyperfine Coupling Constants for Hydrogen: Influence of Geometry, Correlation Method, and Basis Set. J. Phys. Chem. A 2003, 107, 6656–6667. (90) Knight, L. B., Jr.; Herlong, J. O.; Kirk, T. J.; Arrington, C. A. Laser vaporization generation of B14 NH, B15 NH, B14 ND, B16 O, and B17 O: Electron-spin-resonance investigation in neon matrices under ultracold trapping conditions. J. Chem. Phys. 1992, 96, 5604–5613. (91) Cochran, E. L.; Adrian, F. J.; Bowers, V. A. 13 C Hyperfine Splittings in the Electron Spin Resonance Spectra of HCO and FCO. J. Chem. Phys. 1966, 44, 4626–4629. (92) Knight, L. B., Jr.; Steadman, J. An ESR investigation of the C17 O+ cation radical isolated in a neon matrix at 4 K. J. Am. Chem. Soc. 1984, 106, 900–902. (93) Misochko, E. Ya.; Akimov, A. V.; Goldschleger, I. U.; Wight, C. A. Infrared and EPR Spectra of the Difluoronitroxide Radical. J. Am. Chem. Soc. 1998, 120, 11520-11521. (94) Adrian, F. J.; Kim, B. F.; Bohandy, J. Matrix isolation spectroscopy in methane. Isotropic ESR spectrum of HC17 O. J. Chem. Phys. 1985, 82, 1804–1809.

51

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 52 of 61

(95) Morton, J. R.; Preston, K. F. Statistical treatment of EPR data: The spectra of CF3 and F2 NO in SF6 . J. Chem. Phys. 1980, 73, 4914–4916. (96) Merritt, M. V.; Fessenden, R. W. ESR Spectra of the Fluorinated Silyl Radicals. J. Chem. Phys. 1972, 56, 2353–2357. (97) Morton, J. R.; Preston, K. F.; Strach, S. J. Paramagnetic oxides of nitrogen observed in a sulfur hexafluoride matrix. J. Phys. Chem. 1979, 83, 533–536. (98) Nishikida, K.; Williams, F. Angular dependence of proton hyperfine splittings in the electron spin resonance spectrum of the methylsulfinyl radical. J. Am. Chem. Soc. 1974, 96, 4781–4784. 0

(99) Habara, H.; Yamamoto, S. Microwave Spectrum of the FCO Radical in the 2 A Electronic Ground State. J. Mol. Spectrosc. 2001, 207, 238–242. (100) Adrian, F. J.; Cochran, E. L.; Bowers, V. A. Electron Spin Resonance Spectrum of FCO. J. Chem. Phys. 1965, 43, 462–464. (101) Leopold, K. R.; Evenson, K. M.; Comben, E. R.; Brown, J. M. The far-infrared Laser Magnetic Resonance spectrum of the 17 OH radical: Determination of nuclear hyperfine parameters. J. Mol. Spectrosc. 1987, 122, 440–454. (102) Tonooka, M.; Yamamoto, S.; Kobayashi, K.; Saito, S. The microwave spectrum of the NH2 radical: The hyperfine structure of the 2 B1 ground electronic state. J. Chem. Phys. 1997, 106, 2563–2568. (103) Chipman, D. Theoretical study of the properties of methyl radical. J. Chem. Phys. 1983, 78, 3112–3132. (104) Fessenden, R. W. Electron spin resonance spectra of some isotopically substituted hydrocarbon radicals. J. Phys. Chem. 1967, 71, 74–83. (105) Boland, B. J.; Brown, J. M.; Carrington, A. A microwave determination of some molecular parameters for HCO. Mol. Phys. 1977, 34, 453–464. (106) Adrian, F. J.; Cochran, E. L.; Bowers, V. A. ESR Spectrum and Structure of the Formyl Radical. J. Chem. Phys. 1962, 36, 1661–1672. (107) Margulès, L.; Herbst, E.; Ahrens, V.; Lewen, F.; Winnewisser, G.; Müller, H. S. P. The Phosphidogen Radical, PH2 : Terahertz Spectrum and Detectability in Space. J. Mol. Spectrosc. 2002, 211, 211–220. 0

(108) Brown, F. X.; Yamamoto, S.; Saito, S. The microwave spectrum of the HSiS radical in the 2 A ground electronic state. J. Mol. Struct. 1997, 413–414, 537–544. ˜ 2 B2 ground electronic (109) Yamamoto, S.; Saito, S. The microwave spectrum of the CH2 N radical in the X state. J. Chem. Phys. 1992, 96, 4157–4162. (110) Saito, S.; Yamamoto, S. The microwave spectrum of a new phosphorus-bearing radical CH2 P (2 B2 ). J. Chem. Phys. 1999, 111, 7916–7920. (111) Nagai, K.; Yamada, C.; Endo, Y.; Hirota, E. Infrared diode laser spectroscopy of FCO: The ν1 and ν2 bands. J. Mol. Spectrosc. 1981, 90, 249–272. (112) Piltch, N. D.; Szanto, P. G.; Anderson, T. G.; Gudeman, C. S.; Dixon, T. A.; Woods, R. C. The microwave spectrum of isotopically substituted CO+ ion. J. Chem. Phys. 1982, 76, 3385–3388. (113) Janecka, J.; Vyas, H. M.; Fujimoto, M. ESR Spectra of Methyl Radicals in γ-Irradiated Crystals of Sodium Acetate 3D2 O. J. Chem. Phys. 1971, 54, 3229–3230.

52

ACS Paragon Plus Environment

Page 53 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

(114) Tanimoto, M.; Saito, S.; Hirota, E. Microwave spectrum of the boron monoxide radical, BO J. Chem. Phys. 1986, 84, 1210-1214. (115) Suter, H. U.; Engels, B. Theoretical study of electron spin resonance parameters: H2 CN and H2 CO+ . J. Chem. Phys. 1994, 100, 2936-2942. (116) Müller, H. S. P.; Löblein, K.; Hübner, H.; Hüttner, W.; Brown, J. M. The rotational spectrum of the NF2 free radical: Determination of molecular structure. J. Mol. Spectrosc. 2008, 251, 185-193. (117) Bird, G. R.; Baird, J. C.; Jache, A. W.; Hodgeson, J. A.; Curl, R. F.; Kunkle, A. C.; Bransford, J. W.; Rastrup-Andersen, J.; Rosenthal, J. Microwave Spectrum of NO2 : Fine Structure and Magnetic Coupling.J. Chem. Phys. 1964 40, 3378–3390. (118) Müller, H. S. P.; Kobayashi, K.; Takahashi, K.; Tomaru, K.; Matsushima, F. Terahertz spectroscopy of N18 O and isotopic invariant fit of several nitric oxide isotopologs. J. Mol. Spectrosc. 2015, 310, 92-98. (119) Kawaguchi, K.; Saito, S.; Hirota, E. Microwave spectroscopy of the NCO radical in the ν2 = 0 2 Π, ν2 = 1 2 ∆, and ν2 = 2 2 Φ vibronic states. Mol. Phys. 1985, 55, 341-350. (120) Allen, M. D.; Evenson, K. M.; Gillet, D. A.; Brown, J. M. Far-Infrared Laser Magnetic Resonance ˜ 2 Πr State. J. Mol. Spectroscopic Study of the ν2 Bending Fundamental of the CCN Radical in Its X Spectrosc. 2000, 201, 18-29. (121) Tanimoto M.; Saito, S. Microwave spectroscopic study of the SiF3 radical: Spin-rotation interaction and molecular structure. J. Chem. Phys. 1999, 111, 9242-9247. (122) Kaupp, M.; Arbuznikov, A. V.; Hesselmann, A.; Görling, A. Hyperfine coupling constants of the nitrogen and phosphorus atoms: A challenge for exact-exchange density-functional and post-HartreeFock methods. J. Chem. Phys. 2010, 132, 184107. (123) Nguyen, M. T.; Creve, S.; Eriksson, L. A.; Vanquickenborne, L. G. Calculation of the hyperfine constants of phosphorus-containing radicals. Mol. Phys. 1997, 91, 537-550. (124) Hermosilla, L.; Calle, P.; García de la Vega, J. M.; Sieiro, C. Theoretical isotropic hyperfine coupling constants of third-row nuclei (29 Si, 31 P, and 33 S). J. Phys. Chem. A 2005, 109, 7626-7635. (125) Tanoura, M.; Chiba, K.; Tanaka, K.; Tanaka, T. Microwave spectroscopy of chlorine dioxide. J. Mol. Spectrosc. 1982, 95, 157-181. (126) McCarthy, M. C.; Fuchs, G. W.; Kucera, J.; Winnewisser, G.; Thaddeus, P. Rotational spectra of C4 N, C6 N, and the isotopic species of C3 N. J. Chem. Phys. 2003, 118, 3549-3557. (127) Gottlieb, C. A.; Gottlieb, E. W.; Thaddeus, P.; Kawamura, H. Laboratory detection of the C3 N an C4 H free radicals. Astrophys. J. 1983, 275, 916-921. (128) McCarthy, M. C.; Gottlieb, C. A.; Thaddeus, P.; Horn, M.; Botschwina, P. Structure of the CCCN and CCCCH radicals: Isotopic substitution and ab initio theory. J. Chem. Phys. 1995, 103, 7820-7827. (129) Fessenden, R. W. Determination of the relative signs of ESR hyperfine constants from higher-order effects. J. Mag. Res. 1969, 1, 277-290.

53

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 54 of 61

Table 1: Comparison of the HFC tensor components (in MHz) for CH (2 Π) radical obtained via COS-CC method with benchmark FCI and icMRCC results from Ref. 85 The CC results are indicated as deviations from the FCI values. The COS-CC and the spin-orbital based CC calculations include orbital-relaxation effects. 13

Method FCI COS-CCSD COS-CCSD(T) Tˆ3 -only COS-CCSD(T) Tˆ3e -only COS-CCSD(T) Tˆ3x -only COS-CCSD(T) (all) icMRCCSD(1e, 1o) icMRCCSDpsT(1e, 1o) ROHF-CCSD ROHF-CCSD(T) UHF-CCSD UHF-CCSD(T)

aF 34.05 -15.07 -12.89 -14.07 -9.34 -6.17 -14.31 -5.43 -3.63 -0.44 +1.46 +0.09

C Taa -84.27 -1.18 -1.25 -1.13 -0.95 -0.98 -0.57 -0.25 -0.33 -0.15 -0.29 -0.15

1

aF -61.81 -3.74 -1.03 -3.76 -5.15 -2.45 +2.61 +0.99 -1.11 -0.41 -1.53 -0.29

H Taa 40.32 -0.31 -0.37 -0.31 -0.12 -0.18 -0.35 -0.07 +0.16 +0.04 +0.26 +0.04

Table 2: Proton isotropic hyperfine coupling constants (in MHz) in radicals obtained in correlated electronic structure calculations. All electrons were correlated. UHFUHF-CC COS-CC 2 ˆ State hS iUHF MP2 CCSD CCSD(T) CCSD CCSD(T) Expt. 2 OH Π 0.757 -68.02 -75.39 -72.85 -70.58 -70.74 -73.25 101 2 NH2 B1 0.759 -61.94 -68.96 -66.49 -69.06 -66.77 -67.18 102 2 00 CH3 A2 0.762 -68.08 -72.96 -69.83 -77.15 -74.89 -70.06a 103 2 0 HCO A 0.766 380.55 359.89 361.53 354.62 363.36 388.9, 105 383.9 106 2 PH2 B1 0.771 -58.12 -48.85 -48.53 -50.41 -47.71 -48.87 107 2 0 HSiS A 0.780 317.84 312.17 314.66 307.05 317.27 335.74 108 2 00 CH3 SO A 0.788 54.39 40.95 43.09 41.13 51.53 48.48b 98 2 CH2 N B2 0.949 159.26 213.05 216.66 199.38 211.29 233.15 109 2 CH2 P B2 1.161 79.06 98.95 96.27 95.35 106.92 104.94 110 a The reported value corresponds to the planar (D3h ) equilibrium structure without out-of-plane vibrations as derived by Chipman 103 from the EPR spectroscopic results of Fessenden; 104 b The EPR spectrum was recorded at a low temperature (-185◦ C), at which the CH3 group acquires a fixed conformation (see Supporting Information). The reported value corresponds to the two equivalent out-of-plane protons.

54

ACS Paragon Plus Environment

Page 55 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 3: Boron, carbon, and nitrogen isotropic hyperfine coupling constants (in MHz) in radicals obtained in correlated electronic structure calculations. All electrons were correlated. Nucleus 11 B 13

14

a

C

BO CH3 HCO FCO NCO CH2 N CO+ CH2 P

State 2 Σ

hSˆ2 iUHF 0.800

UHFMP2 998.76

UHF-CC CCSD CCSD(T) 1043.85 1007.81

COS-CC CCSD CCSD(T) 1054.54 1032.67

A002 A0 2 0 A 2 Π 2 B2 2 Σ 2 B2

0.762 0.766 0.773 0.843 0.949 1.077 1.161

58.16 337.55 731.32 47.70 10.79 1321.01 -27.29

76.65 381.56 767.72 -63.91 -78.99 1557.67 -41.60

62.78 368.71 769.36 -53.57 -79.98 1602.51 -40.15

2

2

72.40 371.95 750.24 -51.82 -73.74 1469.88 -35.21

72.85 368.47 759.38 -52.87 -82.83 1566.28 -41.58

Expt. 1033(1), 90 1027.40 114 75.67a 103 377.5 91 802.5 91 -80.9b 1573(3), 1506±15 112 92

2 NH2 B1 0.759 24.54 28.60 27.51 22.10 26.67 28.06 102 2 NF2 B1 0.764 44.73 51.40 46.68 40.05 46.52 46.61 116 2 97 NO2 A1 0.771 147.61 145.67 144.95 142.69 147.31 152, 146.53±0.6 117 2 NO Π 0.798 -57.28 24.06 18.78 17.53 21.11 22.27 118 2 0 F2 NO A 0.816 321.49 240.15 249.08 246.56 272.02 261.24, 93 262.87 95 2 NCO Π 0.843 -1.56 15.49 13.93 10.57 16.41 15.13 119 2 CH2 N B2 0.949 1.48 26.64 24.59 19.55 26.26 25.92 109 2 CCN Π 1.128 -29.54 11.58 8.38 7.02 8.65 9.13 120 The reported value corresponds to the planar (D3h ) equilibrium structure without out-of-plane vibrations as derived by Chipman 103 from the EPR spectroscopic results of Fessenden; 104 b Taken from Ref. 115

N

Table 4: Oxygen and fluorine isotropic hyperfine coupling constants (in MHz) in radicals obtained in correlated electronic structure calculations. All electrons were correlated. Nucleus 17 O

19

F

OH HCO NO2 FCO ClO2 NO BO F2 NO NCO CO+ SiF3 NF2 FCO F2 NO

State 2 Π 2 0 A 2 A1 2 0 A 2 B1 2 Π 2 Σ 2 0 A 2 Π 2 Σ

hSˆ2 iUHF 0.757 0.766 0.771 0.773 0.789 0.798 0.800 0.816 0.843 1.077

UHFMP2 -45.57 -59.98 -64.75 -89.56 -6.70 -309.52 -80.98 -8.14 17.14 -94.32

UHF-CC CCSD CCSD(T) -52.08 -51.25 -44.26 -44.04 -58.79 -64.12 -56.97 -58.91 -35.29 -31.38 -38.41 -45.82 -12.49 -15.11 -45.95 -40.67 -32.01 -27.49 32.39 26.75

COS-CC CCSD CCSD(T) -38.96 -50.72 -38.37 -52.69 -57.91 -76.50 -56.55 -68.89 -27.69 -40.69 -26.07 -39.43 -19.39 -26.57 -40.93 -55.57 -20.66 -29.67 17.16 12.22

A0 B1 2 0 A 2 0 A

0.752 0.764 0.773 0.816

350.04 181.61 884.47 438.93

365.25 170.25 895.71 335.67

366.64 105.89 920.17 359.12

2

2

55

366.81 160.77 933.19 390.15

ACS Paragon Plus Environment

447.42 206.64 1039.80 485.41

Expt. -51.17 101 -42.32 94 -62.2 97

-35.46 118 -19(3) 90

18.5(5) 92 386.83 121 164.53 116 1016.(7), 99 935.9 100 401.69, 93 402.11 95

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 56 of 61

Table 5: The isotropic hyperfine coupling constant for the fluorine atom (2 P ) obtained from correlated electronic structure calculations using the uncontracted aug-cc-pCVQZ basis set. All electrons were correlated. Results in MHz.

a

Taken from Ref.; 78

b

From Ref. 78

Method COS-CCSD COS-CCSD(T) Experimenta UHF-CCSD UHF-CCSD(T) MRSD-CIb Result obtained with

aF 204.23 301.38 301.7 308.55 311.08 285.3 an uncontracted (23s,12p,10d,4f,2g) basis set.

Table 6: Isotropic hyperfine coupling constants (in MHz) for third-row elements obtained in correlated electronic structure calculations. All electrons were correlated. Nucleus 29 Si 31

P

33

S

35

Cl

SiF3 HSiS PH2 CH2 P HSiS CH3 SO ClO2

State 2 0 A 2 0 A 2 B1 2 B2 2 0 A 2 00 A 2 B1

UHFUHF-CC COS-CC hSˆ2 iUHF MP2 CCSD CCSD(T) CCSD CCSD(T) 0.752 -1342.12 -1332.34 -1314.62 -1348.35 -1332.82 0.780 -589.26 -568.59 -566.51 -611.81 -612.80 0.771 146.91 213.07 209.70 155.86 209.30 1.161 52.88 180.91 173.78 97.45 161.76 0.780 13.43 12.18 15.05 7.16 18.53 0.788 49.63 29.87 30.38 16.99 27.55 0.789 18.73 51.90 46.99 35.71 49.79 a The EPR spectrum was recorded at -185◦ C.

56

ACS Paragon Plus Environment

Expt. -1394.34 96 207.31 107 177.26 110 22.42a 98 46.12 125

Page 57 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Table 7: Components of the spin-dipolar hyperfine coupling tensor (in MHz) for radicals obtained in correlated electronic structure calculations. All electrons were correlated. UHFUHF-CC COS-CC MP2 CCSD CCSD(T) CCSD CCSD(T) Expt. OH Taa 87.85 89.86 89.30 87.64 88.52 87.09 101 Tbb − Tcc 58.96 58.37 58.28 58.65 58.88 56.68 NH2 Taa 18.23 19.22 18.99 18.66 18.79 18.22 102 Tbb − Tcc -10.88 -9.16 -9.29 -10.15 -10.10 -7.88 CH3 Tcc 3.94 3.06 3.12 3.30 3.31 2.2a 113 HCO Tcc -16.28 -16.89 -16.69 -16.67 -16.42 -15.5b , -11.8b 105 Taa 13.15 13.17 12.88 13.44 13.12 11.6b , 10.2b PH2 Taa -1.22 -0.99 -1.07 -1.11 -1.10 -1.00 107 Tbb − Tcc -10.87 -10.63 -10.58 -11.14 -11.20 -9.91 CH2 N Taa 0.80 9.48 8.37 9.01 9.33 8.29 109 Tbb − Tcc 10.00 2.62 3.69 2.46 2.00 4.15 CH2 P Taa 3.35 4.11 2.43 2.57 2.97 3.68 110 Tbb − Tcc -15.14 -13.65 -11.64 -10.78 -11.16 -11.70 11 B BO Adip 33.05 25.89 27.19 26.66 26.41 25, 90 27.13c 114 13 C CH3 Tcc 155.79 155.23 154.86 155.72 155.44 121.6a 113 HCO Tcc -53.89 -40.97 -41.29 -41.12 -39.43 -38.6 d 91 + 92 CO Adip 88.22 44.82 50.73 48.53 47.65 46(1), 48.2±0.7c 112 14 N NH2 Taa -43.77 -43.67 -43.49 -43.97 -43.80 -43.04 102 Tbb − Tcc -130.31 -129.41 -129.07 -130.06 -129.54 -132.29 NF2 Taa -49.28 -48.49 -47.53 -48.94 -47.77 -47.71 116 Tbb − Tcc -154.55 -151.59 -147.97 -151.44 -148.56 -148.66 NO2 Taa -18.64 -23.95 -21.76 -22.03 -20.59 -19.77 117 Tbb − Tcc 51.62 62.98 58.29 60.32 55.05 57.23 NO Taa -35.58 -39.99 -35.84 -38.38 -37.97 -39.25 118 Tbb − Tcc 276.55 100.79 118.04 111.60 110.08 112.59 NCO Taa -29.30 -33.58 -32.34 -32.11 -31.72 -32.13 119 Tbb − Tcc 109.95 91.29 90.62 91.85 89.17 88.49 CH2 N Taa -35.75 -47.27 -45.10 -45.48 -45.55 -45.14 109 Tbb − Tcc 158.49 116.10 118.77 121.99 118.66 115.59 CCN Taa -5.43 -19.36 -15.66 -14.57 -15.20 -20.61 120 Tbb − Tcc 78.34 45.62 48.91 44.06 46.65 46.77 17 O OH Taa 142.30 140.55 140.35 141.43 140.83 147.20 101 Tbb − Tcc -426.40 -424.21 -423.49 -427.33 -425.32 -429.49 NO Taa 330.98 60.58 71.58 63.72 64.45 61.91 118 Tbb − Tcc 365.35 -207.19 -182.38 -207.68 -206.39 -206.12 BO Adip -6.79 -21.72 -19.82 -17.72 -18.74 -12c 90 + CO Adip -4.85 -43.10 -37.12 -33.27 -34.37 -33.0(5)c 92 19 F SiF3 Tcc 108.60 114.12 115.33 108.12 112.95 115.39 121 NF2 Taa -221.76 -237.09 -241.58 -221.04 -233.61 -241.97 116 Tbb − Tcc -665.86 -665.95 -692.31 -622.58 -655.14 -694.64 31 P PH2 Taa -308.48 -314.08 -309.64 -314.84 -314.01 -300.14 107 Tbb − Tcc -902.17 -902.86 -896.45 -906.87 -900.54 -944.02 CH2 P Taa -362.39 -359.80 -330.68 -328.87 -329.85 -341.77 110 Tbb − Tcc 767.15 783.14 844.63 878.82 849.64 873.91 35 Cl ClO2 Taa -75.18 -78.99 -76.34 -75.50 -72.64 -77.74 125 Tbb − Tcc -227.79 -248.50 -237.62 -234.52 -226.28 -243.95 a The components of the 13 C and 1 H spin-dipolar HFC tensors perpendicular to the molecular plane; b The Tcc component is perpendicular to the molecular plane and Taa is the larger among the two in-plane components. The values Tcc =-11.8 MHz and Taa =10.2 MHz are from Ref., 106 which were recalculated in the principal inertial axis system in Ref.; 105 c The AN dip tensor takes the form (-Adip , -Adip , 2Adip ) for a 57 nucleus N located on an n-fold symmetry axis with n ≥ 3; d Direct comparison with experiment is possible ACS Paragon Plus Environment only for the out-of-plane component due to different choice of axis systems. Nucleus 1 H

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 58 of 61

Table 8: The HFC tensors for the ground 2 Σ state of C3 N radical obtained in spin-adapted COS-CC and spin-orbital based CC calculations using the uncontracted aug-cc-pCVTZ basis set. All electrons were correlated. Results in MHz. aF Taa (=2Adip ) 13 13 13 14 13 13 13 14 Method hSˆ2 i C1 a C2 C3 N C1 a C2 C3 N COS-CCSD 0.750 1105.72 232.59 28.05 -1.68 99.29 38.72 0.86 2.24 COS-CCSD(T) 0.750 1100.20 235.31 33.72 -1.05 94.52 40.67 0.95 2.57 Experimentb 973(2) 188.6(2) 23.55(2) -1.20(3)c 93.(2) 35.2(6) 1.44(2) 1.89(6)c Experimentd 999(3) 199(1) UHF-SCF 1.623 UHF-CCSD 0.786 1113.08 232.30 32.68 -1.18 85.52 53.19 4.30 5.58 UHF-CCSD(T) 0.738 1073.31 251.71 13.89 1.75 95.63 43.32 3.72 1.60 a The C atoms in the linear C3 N molecule were numbered as C1 C2 C3 N starting from the terminal C atom; b From Ref.; 126 c From Ref. 127 Ref. 126 reported the 14 N iHFCC to be -1.26(6), -1.234(6), -1.182(8) MHz for the species 13 CCCN, C13 CCN, and CC13 CN, respectively, and the Taa value to be 2.26(6), 1.88(2), and 1.92(1) MHz, respectively; d From Ref. 128

58

ACS Paragon Plus Environment

Page 59 of 61

0.40

Spin density at nucleus (au-3)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

0.30

COS-CCSD COS-CCSD(T) MRD-CI

0.20

0.18

0.40 0.37

0.19

0.10

0.05

0.02 Total

0.10 0.09 0.05

-0.05 -0.03 1s

2s

N

-0.10 -0.20

2s

0.08

0.07

0.00 1s

0.22

Total

F

1s

2s

3s

Total

Cl

-0.17

-0.30 -0.33 -0.33 -0.40

Figure 1: Contributions from the spin polarization of the 1s, 2s, and 3s shells to the spin density at the nuclei. The 1s shell contribution was computed by keeping the 2s (as well as 3s in Cl) atomic orbitals doubly occupied (i.e., frozen) in the electron correlation treatment and so on. The results for the 14 N nucleus are taken from Ref. 80

59

ACS Paragon Plus Environment

Journal of Chemical Theory and Computation

1.20 UHF-MP2 UHF-CCSD UHF-CCSD(T) COS-CCSD COS-CCSD(T)

1.00

0.80 Mean Absolute Error

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 60 of 61

0.60

0.40

0.20

0.00

H

1 Proton

11 Boron B

C

13 Carbon

N

14 Nitrogen

O

17 Oxygen

F

19 Fluorine

Third-row elements

Third-row elements

Figure 2: Mean absolute errors of the computed iHFCCs for various nuclei relative to experiment. The errors have been scaled by the element-specific prefactor PN (MHz/au3 ).

60

ACS Paragon Plus Environment

Page 61 of 61 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

For Table of Contents Only 207x138mm (300 x 300 DPI)

ACS Paragon Plus Environment