Relativistic effects on electron-nucleus hyperfine ... - ACS Publications

Subscriber access provided by UNIV OF REGINA

Article

Relativistic effects on electron-nucleus hyperfine coupling studied with an `exact 2-component' (X2C) Hamiltonian Jochen Autschbach J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b01014 • Publication Date (Web): 14 Dec 2016 Downloaded from http://pubs.acs.org on December 16, 2016

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 free 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 accessible to all readers and 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.

Journal of Chemical Theory and Computation 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 23

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

Relativistic effects on electron-nucleus hyperfine coupling studied with an ‘exact 2-component’ (X2C) Hamiltonian

Jochen Autschbach Department of Chemistry University at Buffalo State University of New York Buffalo, NY 14260-3000, USA email: [email protected] Abstract: An ‘exact 2-component’ (X2C) transformation of the one-electron Hamiltonian is used to transform nuclear hyperfine magnetic field operators from the 4-component Dirac picture to 2-component form. Numerical applications are concerned with hyperfine coupling constants of one-electron and many-electron atoms, as well as the HgH radical, using spin-unrestricted scalar X2C Hartree-Fock and Kohn-Sham theory. Reference data for 2-component generalized-collinear X2C calculations, including spin-orbit coupling, are also provided for selected cases. Calculations for one-electron atomic 𝑛𝑠 states with 𝑛 = 1 – 3 show that the X2C transformed hyperfine operators give accurate hyperfine coupling constants. Kohn-Sham one-electron self-interaction errors for these states are small. The performance of the X2C transformed hyperfine operator for manyelectron systems is also promising. The method is straightforward to implement in codes using spin-unrestricted (‘1-component) or 2-component spinor orbitals.

1 ACS Paragon Plus Environment


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

1 Introduction For many decades already, there has been tremendous interest in atomic and molecular relativistic calculations using a 2-component (2-spinor) rather than the 4-component Dirac (4-spinor) formalism.1–10 As far as the one-electron part of the Hamiltonian is concerned, the transformation of the Hamiltonian furnished by the ‘exact 2-component’ (X2C) approach is an efficient and practical solution to this problem.3, 4, 11–22 Two-electron contributions require more effort, but relatively simple approximations, such as augmenting the nuclear potential in the X2C transformation by a model potential18, 19, 23 or carrying out the transformation based on a converged 4-component mean-field potential,20 have also been shown to be effective and accurate. Molecular property calculations related to static magnetic field interactions, chiefly NMR chemical shifts and indirect nuclear spin-spin coupling (𝐽 -coupling), as well as EPR 𝑔-tensors and electron-nucleus hyperfine interactions, additionally require the transformation of odd 4-component perturbation operators to 2-component form. This is in principle straightforward,24 but the procedure can become somewhat cumbersome when the first-order perturbed X2C transformation is required and when, for external fields, distributed gauge-origin atomic orbital (AO) methods and magnetic balance25–27 are applied. There are indications in the literature, based on an analysis of the normalized elimination of the small component (NESC, which is in principle numerically equivalent to X2C) and corresponding scalar NESC calculations, that a reasonably accurate 2-component hyperfine perturbation operator can be obtained by using only the field-free transformation of the one-electron Hamiltonian to 2-component form.11, 28 For electric-field perturbations, the contributions from the first-order terms in the transformation were also shown to be quite small29 when considering, for instance, typical magnitudes of picture-change or electron correlation effects on these properties. In the X2C method, the transformation matrices for the Hamiltonian are explicitly available in an atomic orbital (AO) basis. Given an AO matrix representation of a 4-component operator for a molecular property, it is then particularly straightforward to calculate the same property within the 2-component scheme simply by using the field-free transformation matrices, in order to furnish the all-important picture-change. Further justification of this approach is provided in Section 2. For the present work, we have implemented such an approach for electron-nucleus hyperfine coupling in the recently developed X2C module30 of the NWChem package,31 for scalar X2C spin-unrestricted (‘1-component’) self-consistent field (SCF) methods, i.e. Hartree-Fock (HF) and Kohn-Sham calculations, as well as for generalized-collinear 2-component X2C SCF calculations including spin-orbit (SO) coupling variationally. Numerical X2C calculations of the hyperfine integrals for 𝑛𝑠 states of one-electron atoms, with 𝑛 = 1, 2, 3 and nuclear charges 𝑍 = 54, 80, are compared to analytic results for Dirac spinors and indicate that the field-free X2C-transformed 2 ACS Paragon Plus Environment

Page 2 of 23

Page 3 of 23

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


hyperfine operator is very accurate for the highly important relativistic analog of the Fermi contact (FC) mechanism. Kohn-Sham calculations for the same one-electron states with a variety of commonly used approximate functionals introduce relative one-electron self-interaction errors of, at most, 7 parts per thousand. Calculations for Cs, Fr, and Hg+ as well as the HgH radical demonstrate an easy and practical route to relativistic hyperfine coupling calculations for many-electron atoms and molecules based on X2C ground-state electronic structure calculations. Section 2 provides the formalism underlying the calculations. Further technical and computational details are provided in Section 3. One-electron (hydrogen-like) atoms are discussed in Section 4, while results for many-electron systems are provided in Section 5. Conclusions and a brief outlook can be found in Section 6.

2 Decoupling transformations applied to magnetic field operators ̀ 𝝁 is used for vectors, vector operators, and their matrix represenBold-italic notation such as 𝒓, 𝑺, tations in a scalar (one component) atomic orbital (AO) basis set. Arrangements of AO matrices in a 2 × 2 super-matrix structure are indicated by upright-bold as in 𝐫, 𝐒, 𝛍. The 3-vector of the Pauli matrices is 𝝈 = (𝛔𝑥 , 𝛔𝑦 , 𝛔𝑧 ). A matrix representation of an operator in the Dirac picture in 4 × 4 block structure is denoted by notation such as 𝕙, 𝕌. For 2-component relativistic calculations, one seeks a completely de-coupled electrons-only operator.32, 33 For the Dirac one-electron Hamiltonian 𝕙0𝐷 matrix without magnetic fields, the transformation is [ ] 𝐡 𝟎 𝕌† 𝕙0𝐷 𝕌 = X2C (1) 𝟎 Ȃ with

[ ] 𝐔𝑈 𝑈 𝐔𝑈𝐿 𝕌= 𝐔𝐿𝑈 𝐔𝐿𝐿

(2)

using split notation for the 2 × 2 sub-blocks (𝟎 means that a block is filled with zeros, and 𝟏 (used further below) is a 2 × 2 unit matrix, and Ȃ means that we are not concerned with this block). In Equation (2), 𝑈 is an index for the ‘upper’ or ‘large’ components of the Dirac 4-spinor, and 𝐿 is an index for the ‘lower’ or ‘small’ components. (Indices 𝐿 and 𝑆 for ‘large’ and ‘small’ are more commonly used, but in the context of this work the different notation should not lead to confusion.) The 𝑈 𝑈 component of the transformed Hamiltonian is the desired operator; the 𝐿𝐿 block is not of interest here. The procedure is referred to as X2C if the decoupling transformation is achieved by matrix algebra.3 The specific form of the blocks 𝐔𝑈 𝑈 and 𝐔𝐿𝑈 needed for the present work has been presented in References 21 and 30, among others. A generalized matrix self-consistent-field (SCF) eigenvalue equation using 𝐡X2C for the one-electron part and the upper 3 ACS Paragon Plus Environment


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 4 of 23

component basis set overlap matrix is then solved in HF and KS X2C calculations, either as the target electronic structure method or as a starting point for correlated wavefunction calculations, possibly with some form of mean-field or model potentials two-electron contributions included in the transformation. X2C, like other relativistic Hamiltonians, can be implemented in a scalar relativistic (SR) approximation, i.e. ignoring spin-orbit (SO) effects, or with SO coupling included from the onset. The unperturbed Hamiltonian to be de-coupled is that of Dyall’s modified Dirac Equation11, 34 (mDE), in which the Dirac small / lower component is replaced by a pseudo-large / -upper component. The transformation is carried out on matrix representations of the operators, using the same atomic orbital (AO) basis set {𝜒𝜇 } for the upper and pseudo-upper component. The mDE and the unmodified Dirac equation in a basis set representation are equivalent if in the latter a restricted kinetic balance35, 36 (RKB) basis is adopted for the lower component. RKB means that if the upper component basis is {𝜒𝜇 }, then the lower components are represented in the basis {𝜉𝜇 }, with {𝜉𝜇 } = {

1 ̀ 𝜇} (𝝈 Ȃ 𝒑)𝜒 2𝑚𝑒 𝑐

(3)

Here, 𝒑̀ = −𝑖ℏ𝛁 is the linear momentum operator. The blocks of the transformation 𝕌 are then likewise represented in this basis set. The ‘minimal substitution’ 𝒑̀ Ă 𝒑̀ + 𝑒𝑨 for an electron is commonly used to incorporate magnetic effects into the quantum theoretical treatment of molecules. Here, 𝑨 is the vector potential and 𝑒 the unit charge. The magnetic-field dependent one-electron Hamiltonian matrix in the 4-component Dirac picture then reads mag

(4)

𝕙𝐷 = 𝕙0𝐷 + 𝕙𝐷 with

[

mag

𝕙𝐷

𝟎 𝝈Ȃ𝑨 = 𝑒𝑐 𝝈Ȃ𝑨 𝟎

] (5)

In Equation (4), 𝕙0𝐷 is the usual field-free Dirac one-electron Hamiltonian. As a reminder, here and in the following the 2 × 2 blocks of operators are represented in the RKB basis set. A general notation for the vector potential, 𝑨(𝒓) = 𝒃 × 𝑭 (𝒓)

(6)

can be used for the vector potential for common application scenarios, where 𝒃 is a constant vector representing the perturbation, and 𝒓 is the electron position vector. For instance, for an external 4 ACS Paragon Plus Environment

Page 5 of 23

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


static homogeneous magnetic field 𝑩, 1 1 𝑨(𝒓) = 𝑩 × (𝒓 − 𝑮) i.e. 𝒃 = 𝑩 , 𝑭 (𝒓) = (𝒓 − 𝑮) 2 2

(7)

Here, 𝑮 is the chosen gauge origin for the vector potential. For external magnetic fields, typical choices for 𝑮 are the laboratory coordinate origin, or the center of mass of the molecule, or its center of nuclear charges. When Gaussian units are used instead of SI – which remains common for magnetic properties – then the vector potential contains an additional factor of 𝑐 −1 . For the hyperfine magnetic field generated by the nuclear spin magnetic moment 𝒎𝐾 of a nucleus no. 𝐾 located at 𝑹𝐾 , the vector potential is 𝑨(𝒓) = 𝜘𝒎𝐾 × with

𝒓𝐾 𝑟3𝐾

=

𝒓𝐾 𝑟3𝐾

i.e. 𝒃 = 𝒎𝐾 , 𝑭 (𝒓) = 𝜘

𝜇 𝒓 − 𝑹𝐾 , 𝜘= 0 3 4𝜋 |𝒓 − 𝑹𝐾 |

𝒓𝐾 𝑟3𝐾

(8)

(9)

These expressions assume a point magnetic dipole. It is possible to ‘smear out’ the magnetic moment with a spherical Gaussian, for instance, in order to treat finite nucleus effects37–39 on the same approximate footing as in a commonly used model for finite nuclear charge distributions.40 The gauge origin for the hyperfine vector potential is usually chosen to coincide with the position of the nucleus, as in Equation (9). In Hartree atomic units, 𝜘 = 𝜇0 Ȃ(4𝜋) has the same numerical value as 𝑐 −2 , with 𝑐 Ȃ 137.036 = 𝛼 −1 atomic units. 𝛼 is the dimensionless fine structure constant. Both the external and the hyperfine field are given here in the Coulomb gauge, which implies 𝛁 Ȃ 𝑨 = 0. In principle, the minimal substitution to include the vector potential should also be included in the definition of the kinetic balance basis, leading to magnetic balance,25 in particular when vector potential dependent basis functions (e.g. gauge-including atomic orbitals (GIAOs)41, 42 for external magnetic fields) are already used for the upper / large components in order to distribute the gauge origin of an external field over the molecule. However, for properties that involve only nuclear magnetic fields, i.e. electron-nucleus hyperfine coupling and 𝐽 -coupling, field dependent basis sets are not commonly used because the position of the nucleus is a natural choice for the gauge origin. There is also the issue of diamagnetic contributions in second-order properties, such as the magnetizability, NMR shielding, and 𝐽 -coupling, which enter only implicitly – via coupling in the linear response function with negative-energy states – in calculations based on a magnetic operator as in Equation (5). Through transformations of the 4-component operators, magnetic balance, or a combination thereof, diamagnetic operators or operator matrices that are bi-linear in 5 ACS Paragon Plus Environment


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 6 of 23

the vector potential can be brought about in the 4-component formalism in a variety of ways.43–45 However, for 𝐽 -coupling it is known that the diamagnetic contribution is usually negligible.46, 47 Furthermore, for second-order properties such as NMR shielding and 𝐽 -coupling, the Sternheim approximation25, 48 has been applied, in which a diamagnetic contribution is obtained from the expectation value of the nonrelativistic operator 𝑒2 Ȃ(2𝑚𝑒 )𝑨 Ȃ 𝑨. A first-order one-electron property can be defined via an energy derivative with respect to a suitable chosen perturbation parameter 𝜆. For the 4-component case – if the basis set is not dependent on the perturbation – the Hellman-Feynman theorem can be used, which gives 𝑑𝐸 || = ‫܂ 𝐷𝜓| )𝜆(𝐷 ̀𝐻| 𝐷𝜓܂‬ 𝑑𝜆 ||𝜆=0

(10)

with 𝐻̀ 𝐷(𝜆) = (𝜕 𝐻̀ 𝐷 )Ȃ(𝜕𝜆)|𝜆=0 , using here for convenience a notation with the unperturbed 4component wavefunction 𝜓𝐷 and operators, not their representations in the basis. Assuming that the 4-component and X2C energies and their derivatives are equal, the unitary transformation operator can simply be inserted twice in the expectation value, to give 𝑑𝐸 || = ‫𝜓܂‬X2C |𝑈̀ † 𝐻̀ 𝐷(𝜆) 𝑈̀ |𝜓X2C ‫܂‬ | 𝑑𝜆 |𝜆=0

(11)

with 𝜓X2C = 𝑈̀ † 𝜓𝐷 . For a one-electron system, the X2C transformation matrix matching a 4component RKB calculation with a given basis set can be constructed exactly, and therefore the transformation of the 4-component perturbation operator matrix with the unperturbed 𝕌 should produce the correct 2-component perturbation operator. For many-electron systems, due to commonly applied approximations, the X2C transformation does not produce exactly the same energy as a corresponding 4-component calculation. Moreover, without taking care of perturbationdependent terms entering the X2C transformation, the expectation value of 𝑈̀ † 𝐻̀ 𝐷(𝜆) 𝑈̀ is not exactly the derivative of the X2C energy, as already emphasized in Reference 29, and if the basis set depends on 𝜆 there are additional terms. However, as mentioned in the Introduction, for properties such as hyperfine coupling and electric field gradients the errors from neglecting a perturbationdependence of a transformation to a 2-component picture appear to be quite small, in particular when considering other approximations that affect these properties. The discussion in the preceding two paragraphs suggests that a hyperfine operator matrix as in Equation (5), transformed to 2-component form without further adjustments, may constitute a mag reasonable choice for first-order properties. Consequently, we are seeking an AO matrix 𝐡X2C for the magnetic perturbation operator, obtained from the same transformation 𝕌 used to transform 𝕙0𝐷 . I.e., mag (12) 𝐡X2C = 𝐔𝐿𝑈 † 𝐡𝐿𝑈 𝐔𝑈 𝑈 + 𝐔𝑈 𝑈 † 𝐡𝑈𝐿 𝐔𝐿𝑈 6 ACS Paragon Plus Environment

Page 7 of 23

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


Here, 𝐡𝐿𝑈 is the AO matrix representation of 𝝈 Ȃ 𝑨 in the mixed lower / upper component basis, and 𝐡𝑈𝐿 is the AO matrix representation of 𝝈 Ȃ 𝑨 in the mixed upper / lower component basis. In the 𝑈𝐿 basis, 𝝈 Ȃ 𝒑̀ acts on the ket. In the 𝐿𝑈 basis, 𝝈 Ȃ 𝒑̀ acts on the bra. The factor of 𝑐 −1 in Equation (3) enters the 𝑈𝐿 and 𝐿𝑈 AO matrices and cancels the factor of 𝑐 in the magnetic perturbation operator of Equation (5). For calculations of magnetic properties, one needs in first order the operator derivatives with respect to the components 𝑏𝑢 of the vector 𝒃. With 𝝈Ȃ𝑨=𝒃Ȃ𝑭 ×𝝈 one obtains in the 4-component framework (the RKB basis is implicit) mag,𝑢

𝕙𝐷

[ ] ] [ mag 𝑈𝐿,𝑢 𝜕𝕙𝐷 | 𝟎 [𝑭 × 𝝈] 𝟎 𝐡 𝑢 | = 𝑒𝑐 = = 𝑒𝑐 𝐿𝑈 ,𝑢 | 𝜕𝑏𝑢 |𝒃=0 𝐡 𝟎 [𝑭 × 𝝈]𝑢 𝟎

(13)

with [𝑭 × 𝝈]𝑢 being the component 𝑢 of the vector 𝑭 × 𝝈. The matrix elements of 𝐡𝐿𝑈 ,𝑢 and 𝐡𝑈𝐿,𝑢 are given by 𝑒 ̀ ‫ 𝝈( 𝑢]𝝈 × 𝑭[|𝜇܂‬Ȃ 𝒑)|𝜈‫܂‬ 2𝑚𝑒 𝑒 ̀ = ‫ 𝝈(܂‬Ȃ 𝒑)𝜇|[𝑭 × 𝝈]𝑢 |𝜈‫܂‬ 2𝑚𝑒

Ă𝑈𝐿,𝑢 = 𝜇𝜈

(14a)

,𝑢 Ă𝐿𝑈 𝜇𝜈

(14b)

assuming that the basis set does not depend on 𝑏𝑢 . It is convenient to separate the AO matrix elements needed for the X2C-transformed magnetic operator into spin-dependent and spinindependent contributions. In units of the Bohr magneton 𝜇𝐵 = 𝑒ℏȂ(2𝑚𝑒 ), the matrix elements read Ă𝑈𝐿,𝑢 = −𝑖𝟏‫ ܂𝜈| 𝑢)𝛁 × 𝑭(|𝜇܂‬− 𝛔𝑢 ‫ 𝑭|𝜇܂‬Ȃ 𝛁|𝜈‫ ܂‬+ 𝜇𝜈 ,𝑢 Ă𝐿𝑈 = 𝑖𝟏‫ ܂𝜈|𝜇 𝑢)𝛁 × 𝑭(܂‬− 𝛔𝑢 ‫ 𝑭܂‬Ȃ 𝛁𝜇|𝜈‫ ܂‬+ 𝜇𝜈

∑

∑

𝛔𝑣 ‫܂𝜈 𝑢𝜕| 𝑣𝐹|𝜇܂‬

(15a)

𝑣

𝛔𝑣 ‫܂𝜈| 𝑣𝐹|𝜇 𝑢𝜕܂‬

(15b)

𝑣

with 𝑢, 𝑣, Ȃ {𝑥, 𝑦, 𝑧}. Hence, the X2C magnetic perturbation operator can be calculated easily from AO matrices with elements ‫ ܂𝜈 𝑣𝜕| 𝑢𝐹|𝜇܂‬for the 𝑈𝐿 block, and its transpose is used to calculate the 𝐿𝑈 block (assuming real basis functions). Subsequent application of the transformation in Equation (12) gives the desired 2-component magnetic perturbation operator matrix. A comment on Equations (12), (15) and the nonrelativistic limit is in order. We have 𝐡𝑈𝐿,† = 𝐡𝐿𝑈 . The matrices are generally not Hermitian. In the nonrelativistic limit, both 𝐔𝑈 𝑈 and 𝐔𝐿𝑈 7 ACS Paragon Plus Environment


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

become unit matrices.49 In this case, in Equation (12) 𝐡𝐿𝑈 and 𝐡𝑈𝐿 are simply added to create a Hermitian operator matrix in the AO basis. AO integrals entering these matrix elements combine 𝐹𝑢 and the derivative 𝜕𝑣 applied in the bra in one term and in the ket in the other term. The AO integrals entering the electron spin independent terms are the same as those needed in nonrelativistic calculations for the paramagnetic nuclear spin – electron orbital (PSO) hyperfine interaction (or the orbital-Zeeman (OZ) interaction with an external magnetic field, depending on whether (8) or (7) is chosen for the vector potential). For 𝑐 Ă Ȃ, the X2C transformation then simply produces the PSO operator matrix, but for finite 𝑐 there are important relativistic corrections from the transformation that occur mainly in the regions where the lower / small components ̀ times the upper / large components, that is, close to the atomic nuclei. differ from 1Ȃ(2𝑚𝑒 𝑐)(𝝈 Ȃ 𝒑) For comparison, in the quasi-relativistic zeroth-order regular approximation (ZORA), where the magnetic perturbation operators are easily derived at the operator level, the nonrelativistic PSO operator is multiplied with a relativistic ‘kinematic function’ that mainly differs from unity close to the nuclei.50 In the nonrelativistic limit for the X2C expression, after adding 𝐡𝐿𝑈 and 𝐡𝑈 𝐿 , the electron spin dependent terms afford AO integrals over 𝐹𝑢 𝜕𝑣 (𝜒𝜇 𝜒𝜈 ). With integration by parts, the derivative can be shifted to the 𝐹𝑢 term instead. This produces, up to constant factors, electric-field gradient AO integrals and, for 𝑣 = 𝑢, the familiar Dirac delta distribution ‘contact’ operator. Upon multiplication with the spin operators, this gives the well-known FC and SD contributions to the nonrelativistic hyperfine operator (or the spin-Zeeman operator, for an external field perturbation). This happens only when 𝐡𝐿𝑈 and 𝐡𝑈𝐿 can be added before the X2C transformation, meaning if ̀ times the upper the lower / small components of the orbitals are strictly proportional to (𝝈 Ȃ 𝒑) / large components. For comparison, in the quasi-relativistic ZORA approach the operator matrix elements representing the FC and SD interactions also afford AO integrals over 𝐹𝑢 𝜕𝑣 (𝜒𝜇 𝜒𝜈 ), but multiplied with the same relativistic kinematic function that enters the PSO-like terms, which suppresses a delta distribution for point nuclei as long as 𝑐 is finite.50 Even though the operators in the X2C AO integrals are not quite as singular as the nonrelativistic FC operator, the integrands nonetheless sample the regions very close to the nuclei. Indeed, it is known that relativistic effects, and even finite-nucleus effects, in particular on FC-like matrix elements can be very large.38, 51–53 In a SR framework, the spin independent terms (PSO, or OZ) do not contribute to the result unless the electronic state is orbitally degenerate. Furthermore, with a spin quantization axis in direction 𝑣 (usually chosen to be 𝑧), only terms dependent on 𝛔𝑣 contribute. In order to calculate the spin-dependent terms for the Zeeman or hyperfine interaction in a SR framework, the spin quantization directions are varied to be along 𝑣 = 𝑥, 𝑦, or 𝑧 (this can be done implicitly), the index 𝑢 for the perturbation is varied, and the AO integrals for 𝑢 are contracted with the spin density matrix, to generate a 3 × 3 matrix representing the interaction. Within a computational 8 ACS Paragon Plus Environment

Page 8 of 23

Page 9 of 23

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


framework where SO effects are included in the description of the electronic states, sizable PSO and OZ contributions to the hyperfine and Zeeman interaction, respectively, may be present even in the absence of orbital degeneracies,52, 55 and the distinction between electron spin-dependent and spin-independent terms blurs. In summary, the hyperfine coupling in the 2-component picture has contributions with similar characteristics as the nonrelativistic FC, SD, and PSO ones, and therefore the familiar terminology may be appropriate, keeping in mind that SO coupling blurs the distinction between spin-dependent and spin-independent mechanisms and that for point nuclei the operators are not singular like the nonrelativistic ones in their usual representation. As usual in quantum mechanics, transformations can be applied to the operators, and integration by parts can be used, such that the operator matrix elements remain the same but a different ‘mechanism’ is creating the observable. For example, the nonrelativistic limit of the hyperfine coupling can also be calculated with a continuum of global non-contact operators.54

3 Computational details The X2C calculations were carried out with a locally modified developer’s version of the NWChem quantum chemistry package.31, 56, 57 The implementation in NWChem for scalar relativistic and SO Hartree-Fock (HF) and Kohn-Sham (KS) X2C single point calculations as well as picture-change corrected electric field gradients in the same version of NWChem has been reported previously in Reference 30. The NWChem X2C module now includes new routines for the evaluation of the hyperfine operator matrix elements of Section 2, and corresponding routines to assemble hyperfine coupling tensors from spin-unrestricted scalar relativistic as well as from 2-component generalized-collinear SO X2C calculations. The integral code in NWChem was extended in order to produce integrals ‫ ܂𝜈 𝑣𝜕| 𝑢𝐹|𝜇܂‬for all combinations of 𝑢, 𝑣 Ȃ {𝑥, 𝑦, 𝑧} over Gaussian-type orbital (GTO) basis functions, utilizing the Rys polynomial algorithm for property integrals of Dupuis.58 The point-nucleus model was used for all calculations for the one-electron systems. For the many-electron systems, a spherical Gaussian nucleus model for the electronnucleus potential was used where indicated, for better comparison with data in the literature. For this work, the X2C transformation utilized the bare electron-nucleus potential. Tests with an available model potential in order to approximate two-electron contributions in the relativistic Hamiltonian will be left for future studies, along with a more detailed investigation of SO effects and the performance of 2-component calculations. For the present work, isotropic hyperfine mag,𝑢 coupling constants were calculated from the expectation values of the 𝑏𝑢 -derivatives of 𝐡X2C of Equation (12) taken at 𝒃 = 0. The hyperfine coupling in 2-component SO KS calculations was obtained from generalized-collinear SO density matrices and the integrals of Equations (15) as 9 ACS Paragon Plus Environment


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

described in Reference 53. Relativistic calculations utilized 𝑐 = 137.036 au exactly, whereas calculations with the X2C code simulating the nonrelativistic limit utilized a value of 𝑐 = 108 . In order to test the accuracy for one-electron Dirac atoms against known analytic 4-component results, an input was set up for a X2C Hartree-Fock (HF) calculation for the 1𝑠 state of Cs54+ , using an uncontracted set of 36 1𝑠 GTOs. All two-electron contributions in the Fock matrix were deleted in order to obtain the 𝑛𝑠 states with 𝑛 > 1 from the virtual orbitals. The basis set was comprised of the primitive 𝑠 functions of the Cs ANO-RCC basis,59 augmented by a geometric series of ten 1𝑠 GTOs deriving from the highest ANO exponent and a ratio of 5, resulting in a largest exponent of about 4.3Ȃ1014 . Test calculations showed that going much beyond this largest exponent produced increasing numerical noise, likely from the Rys quadrature routines for the hyperfine integrals. A similar procedure was used to create a basis set for Hg79+ with ten additional highexponent functions, resulting in a largest exponent of 5.1Ȃ1014 . One-electron self-interaction errors were also studied for these one-electron systems for the 𝑛𝑠 states with 𝑛 = 1, 2, 3, by switching to common approximate functionals and re-starting from the one-electron calculations to ensure rapid convergence to the intended state for 𝑛 > 1. The tested XC functionals are the Slater exchange local density approximation (LDA) functional, the VWN60 correlation LDA (incl. Slater exchange), three popular generalized gradient approximation (GGA) functionals (BLYP,61, 62 PBE,63 BP61, 64 ), the BLYP-based hybrid65 BHLYP with 50% exact exchange, and, for HgH, the PBE0 hybrid with 25% exact exchange.66 Calculations for many-electron systems were performed with uncontracted ANO-RCC basis sets59 for the heavy atoms, and – where indicated – augmented by additional 1𝑠 GTOs with high exponents using the same procedure as for the one-electron atoms. For hydrogen in the HgH calculations, the ‘pcJ4’ basis by Jensen67 was used. For the atomic calculations, polarization functions were removed, and for HgH basis functions beyond 𝓁 = 3 were removed. The HgH calculations used an inter-atomic distance of 1.7662 Å reported in Reference 68 as the equilibrium distance for the ground state, which is 0.02 Å less than a KS-optimized distance69 used by us in previous benchmarks52, 70 and 0.03 Å larger than an experimentally derived value from Reference 71. The X2C transformation should be carried out in an uncontracted AO basis. The NWChem implementation reported in Reference 30 required the use of a fully uncontracted basis. For the present work, additional functionality has been added to the code such that when a contracted basis set is used, then the X2C one-electron Hamiltonian and the transformation matrices used in Equation (1) are initially constructed in the corresponding uncontracted basis. The Hamiltonian matrix is subsequently transformed to the contracted basis, and the HF or KS calculation proceeds in the contracted basis. The construction of the hyperfine operator matrices follows the same recipe, i.e. the integrals ‫ ܂𝜈 𝑣𝜕| 𝑢𝐹|𝜇܂‬are obtained for the uncontracted basis, and the result 10 ACS Paragon Plus Environment

Page 10 of 23

Page 11 of 23

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


of Equation (12) is subsequently transformed to the contracted basis and used in conjunction with the self-consistent density matrix in order to produce elements of the hyperfine coupling tensor. It is noted that this procedure is not as efficient in the 2-component calculations as compared to the scalar relativistic ones, because the atomic spinors with the same angular momentum quantum number 𝓁 but different 𝑗 (e.g. 𝑝1Ȃ2 vs. 𝑝3Ȃ2 ) would ideally have separate contractions. The highly contracted ANO-RCC basis sets also caused problems with the module creating the initial guess in NWChem, and therefore results from calculations with fully uncontracted basis sets are reported herein for all calculations. For isotope-specific hyperfine coupling constants, given in the usual units of MHz, the following nuclear 𝑔-factors were used for the conversion (these are hard-coded in NWChem and agree with the isotope data of Reference 72): 𝑔( 1H) = 0.58569, 𝑔( 133Cs) = 0.737978, 𝑔( 199Hg) = 1.011770, 𝑔( 223Fr) = 0.780.

4 Hydrogen-like systems Consider a hydrogen-like (one-electron) atom with a point-nucleus of charge 𝑍. Hartree atomic units are used in this section, and the 𝜇0 Ȃ(4𝜋) factor from the vector potential is written as 𝜘 in order to avoid confusion with the different relativistic orders. Pyykkö et al.73–75 derived equations for the matrix elements of the hyperfine operator of Equation (5), with the point-nucleus vector potential of Equation (8), taken with the 4-component Dirac spinors |𝑛, 𝓁, 𝑗, 𝑚𝑗 ‫ ܂‬for a hydrogen-like atom. For the 𝑧 component of the nuclear spin magnetic moment, and an 𝑠 orbital spinor with principal quantum number 𝑛 and 𝑚𝑗 = ±1Ȃ2, the expectation value of the operator may be written as73 2 |𝑛, 0, 1Ȃ2, ±1Ȃ2‫ = ܂‬Ȃ𝜅 𝑎𝐷 3

(16)

𝑁 − 2(𝑛 − 1 + 𝛾) 𝛾(4𝛾 2 − 1)𝑁 4

(17)

mag,𝑧

‫𝑛܂‬, 0, 1Ȃ2, ±1Ȃ2|𝕙𝐷 with a radial integral

𝑎𝐷 = −2𝑍 3 and

Ȃ 1 − 𝑍 2 Ȃ𝑐 2 Ȃ 𝑁 = (𝑛 − 1)(𝑛 − 2𝛾) 𝛾=

A negative sign has been introduced in Equation (17) compared to the expression in Reference 73 in order to obtain a positive/negative sign for the hyperfine integrals of 𝑠 orbitals with 𝑚𝑗 equal to ±1Ȃ2. 11 ACS Paragon Plus Environment


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 12 of 23

An expansion of 𝑎𝐷 in terms of powers of 𝑍Ȃ𝑐 gives ( 2 ) ( ) 11𝑛 + 9𝑛 − 11 𝑍 5 203𝑛4 + 225𝑛3 − 134𝑛2 − 330𝑛 + 189 𝑍 7 2𝑍 3 −𝑎 = 3 + + 𝑛( 3𝑐 2 𝑛5 36𝑐 4 𝑛7 ) 6 5 4 3 2 1759𝑛 + 2205𝑛 − 721𝑛 − 3585𝑛 − 891𝑛 + 3969𝑛 − 1467 𝑍 9 + Ȏ(𝑐 −8 ) + 216𝑐 6 𝑛9 𝐷

(18)

The nonrelativistic limit of 𝑎𝐷 , i.e. the first term on the right-hand side, gives the nonrelativistic limit of the hyperfine integral as ±𝜘4𝑍 3 Ȃ(3𝑛), which is also easily derived from the nonrelativistic ‘Fermi contact’ (FC) operator 4𝜋 𝐡FC = 𝜘 𝛿(𝒓)𝛔𝑧 (19) 3 and the well-known expressions for the nonrelativistic hydrogen-like orbitals, evaluated at the position of the point-nucleus. The leading relativistic correction, of order 𝑐 −2 , has been known for a long time.73, 76 Table 1 collects for the Cs54+ and Hg79+ ions, i.e. the hydrogen-like atoms with 𝑍 = 55 and 𝑍 = 80, numerical values for the 𝑛𝑠 hyperfine integrals with 𝑛 = 1, 2, 3, not including the factor 𝜘 = 𝜇0 Ȃ(4𝜋) from the vector potential. The data listed as ‘exact’ are based on Equation (16), truncations thereof up to order 𝑐 −2 , 𝑐 −8 , and 𝑐 −12 , and the exact nonrelativistic limit, respectively. The X2C data were obtained with the NWChem implementation described in Section 3. Within the displayed number of significant figures, the exact relativistic hyperfine integrals for a speed of light of 𝑐 = 108 are identical to the exact nonrelativistic ones and therefore not listed separately. The small deviations between ‘exact nonrel. and ‘X2C nrel.’ in Table 1 are therefore indicative of the basis set incompleteness. For the relativistic calculations, the basis set incompleteness is likely more pronounced. The exact 𝑛𝑠1Ȃ2 Dirac spinors for finite 𝑐 have weak singularities at a point nucleus, which require very high exponents in the non-singular GTO basis set in order to obtain converged hyperfine integrals. Nonetheless, the data in Table 1 show very good agreement between the finite-basis numerical X2C calculations and the analytic Dirac results, with X2C in most cases being close to the Ȏ(𝑐 −12 ) truncations of the exact data. As already noted in Section 3, going to even higher exponents than the ones used for the results in Table 1 introduced noticeable numerical noise, likely from the Rys quadrature for the hyperfine AO integrals. The X2C data in Table 1 may also be slightly impacted by numerical noise, but not to a degree that should cause concern for using the code in production calculations. In this context it is also important to note that production calculations are likely to employ finite nucleus models, in which case the convergence of the hyperfine integrals with respect to the high-exponent set of basis functions is faster and extremely high exponents are not required. Further, in practical applications other approximations introduced via the basis set and the electronic structure model are going to be much more of a concern. 12 ACS Paragon Plus Environment

Page 13 of 23

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 relativistic effects are seen to be extremely important, which has long been known, increasing the 𝑛𝑠 hyperfine integrals for 𝑍 = 55 by factors of 1.31, 1.46, and 1.46 for 𝑛 = 1 to 3, respectively. For 𝑍 = 80, the increases due to relativistic effects are by factors of 1.97, 2.53, and 2.54, respectively. Indeed, relativistic hyperfine integrals for truly heavy elements can be several times as large as the nonrelativistic limit.73 The data in Table 1 emphasize that the leading relativistic order, 𝑐 −2 , is not accurate for the 𝑠-state hyperfine integrals of Cs54+ and a very poor approximation for Hg79+ . KS calculations for one-electron systems with commonly used approximate functionals suffer from the one-electron self-interaction error (SIE), which constitutes one facet of the delocalization error (DE).77 The SIE appears because the approximate exchange in the XC potential does not fully mag,𝑧

Table 1: Hyperfine integrals ‫𝑛܂‬, 0, 1Ȃ2, 1Ȃ2|𝕙𝐷 hydrogen-like systems with 𝑍 = 55 and 80 𝑎 𝑍=

|𝑛, 0, 1Ȃ2, 1Ȃ2‫ ܂‬for 𝑛𝑠 states with 𝑛 = 1 to 3 of

55

1𝑠

80

exact rel. 291 156 1 347 874 −12 exact Ȏ(𝑐 ) 291 148 1 341 293 −8 exact Ȏ(𝑐 ) 290 989 1 316 443 exact Ȏ(𝑐 −2 ) 275 434 1 031 655 X2C rel. 290 944 1 335 787 exact nrel. 221 883 682 667 X2C nrel. 221 830 682 663 2𝑠 exact rel. 40 538.7 216 227 −12 exact Ȏ(𝑐 ) 40 536.7 214 547 exact Ȏ(𝑐 −8 ) 40 497.8 208 437 −2 exact Ȏ(𝑐 ) 37 221.0 147 133 X2C rel. 40 508.2 214 158 exact nrel. 27 729.2 85 333.3 X2C nrel. 27 728.8 85 332.9 3𝑠 exact rel. 12 020.1 64 158.9 exact Ȏ(𝑐 −12 ) 12 019.6 63 660.3 −8 exact Ȏ(𝑐 ) 12 008.0 61 845.9 exact Ȏ(𝑐 −2 ) 11 024.6 43 634.9 X2C rel. 12 011.4 63 538.8 exact nrel. 8 216.05 25 284.0 X2C nrel. 8 216.37 25 284.4 𝑎 Point nuclei. Exact results from Equation (16) versus numerical results from the X2C implementation in NWChem described in Section 3. SI-based atomic units, excluding the 𝜘 = 𝜇0 Ȃ(4𝜋) factor in the vector potential. ‘rel.’ implies 𝑐 = 137.026 au, ‘nrel.’ implies 𝑐 = 108 au. 13 ACS Paragon Plus Environment


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

eliminate the spurious self-repulsion of the electron in the Coulomb part of the KS potential. We take the opportunity here to investigate briefly SIEs for the 𝑛𝑠 hyperfine integrals with 𝑛 = 1, 2, 3 of the two hydrogen-like systems of Table 1. The results are collected in Tables S1 and S2 in the Supporting Information. Generalized gradient approximation (GGA) functionals tend to produce larger SIEs overall than the local density approximation (LDA) functionals. In hybrid functionals, the exact exchange component serves to compensate for the SIE but the performance for the hyperfine integrals is in many cases worse than that of the LDA functionals. Further, the SIEs become amplified by relativistic effects. Overall, however, the deviations do not exceed 7 parts per thousand for the test cases and are therefore small compared to typical electron correlation effects and many-electron SIEs that can be observed in hyperfine calculations of many-electron systems.

5 Selected many-electron systems As a further test of the X2C hyperfine coupling approach, computational results for well-studied benchmark systems are provided in this section. X2C Hartree-Fock (HF) test calculations that were additionally performed for Cs and Hg+ with point nuclei and a high speed of light (108 au) produced numerically the same hyperfine coupling constants with the X2C and with nonrelativistic FC and SD integrals, within a relative accuracy of 10−5 or better. Likewise, the results agreed with nonrelativistic HF calculations using nonrelativistic hyperfine operators. Table 2 collects hyperfine constants for selected atoms in 2 𝑆 states in the usual units of MHz (isotope data are given at the end of Section 3), with reference data from the literature given in the footnote. A comparison is made between calculations using a point nucleus (PN), and calculations where, in the ground state calculation and in the X2C transformation of the one-electron Hamiltonian, the potential corresponding to a Gaussian finite nucleus (FN) model was applied. The hyperfine AO integrals used 𝑭 (𝒓) of Equation (8), corresponding to a point magnetic dipole, but previous work51 has shown that adopting a FN model for these integrals as well constitutes a secondary finite nuclear volume effect when compared to the larger impact via the potential and the ground state orbitals. 4-component KS calculations for Cs and Fr with the non-hybrid functional BP were previously reported in Reference 51, and calculations for the two alkali atoms at various non-hybrid and hybrid KS as well as HF and correlated wavefunction levels were reported in References 79 and 80 using a scalar quasi-relativistic Hamiltonian (IORA-mm) and nonrelativistic calculations. FN results for Hg+ obtained from a scalar NESC formalism can be found in Reference 28. Results from the scalar NESC and scalar X2C calculations should be numerically equivalent, all else being equal. Reference 28 demonstrated for selected Hg systems that errors from neglecting the dependence of the transformation to 2-component form on the perturbation 14 ACS Paragon Plus Environment

Page 14 of 23

Page 15 of 23

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


Table 2: Calculated hyperfine coupling constants (MHz) for many-electron 2 𝑆 atoms.𝑎 Cs Fr Hg+ PN FN PN FN PN FN HF 1847 1803 7533 6533 44644 40616 Slater 2130 2080 7764 6741 44431 40434 VWN 2541 2481 9202 7989 46199 42045 BLYP 2455 2397 8912 7728 45414 41301 PBE 2296 2242 8695 7543 45811 41679 BP 2424 2367 9016 7818 46109 41935 BHLYP 2368 2312 9061 7858 45907 41758 Scalar X2C calculations with 𝑐 = 137.036 au. PN = point nucleus, FN = finite nucleus, for the electron-nucleus potential. Uncontracted ANO-RCC basis sets. For the PN calculations, 10 high-exponent 1s GTOs were added. For further details and isotope information see Section 3. Selected data from the literature: Hg+ : 38312 (scalar NESC, HF, uncontr. SARC basis) and 43400 (contr. SARC basis),28 42366 (4-component numerical restricted active space).78 Cs: 2235 (4-component KS, BP),51 2410 (scalar IORA-mm QCISD) and 2052 (scalar IORA-mm HF),79 2571 (IORA-mm, KS, BPW91 functional),80 2298 (expt.)81 Fr: 9165 (PN) and 8096 (FN, 4-component KS, BP),51 7628 (expt., from 211 Fr value reported in Ref. 82 and converted to 223Fr with the ratio 0.780/0.889 of the respective nuclear 𝑔-factors for the 223/211 isotope. 𝑎

were small, as already pointed out. Correlation effects as furnished via the KS functionals – including the associated SIEs – increase the hyperfine coupling significantly for Cs and Fr, whereas the effects on the Hg+ ion are less pronounced on a relative scale (the Slater exchange data are slightly below HF). In agreement with 4-component data from Reference 51, the finite nucleus effects on the Cs hyperfine coupling are small compared to correlation / SIE effects. However, there is a reduction of the Fr hyperfine coupling by over 1000 MHz in the KS calculations due to the finite nuclear volume, which approaches the magnitude of the correlation effects. The overall agreement with the FN and PN 4-component KS calculations for Fr from Reference 51 is reasonable when considering that contracted basis sets were used in the latter. (Basis set differences also affect comparisons with other sets of available calculations; with the uncontracted ANO-RCC basis presumably providing the greatest flexibility among the tested sets.) The X2C calculations for Cs with PBE and BHLYP are close to the experiment. Likewise, the Fr KS results for the FN model bracket the experimental hyperfine coupling of 7628 MHz,82 whereas the absence of correlation in the HF calculation leads to a severe underestimation. For Hg+ , the FN results for KS calculations that include a correlation functional (i.e. all but those with pure Slater exchange) are in reasonable agreement with, if slightly below, a 42366 MHz reference value from a fully numerical 4-component multi-configurational wavefunction calculation.78 Results from scalar X2C calculations for the HgH doublet radical are collected in Table 3. 15 ACS Paragon Plus Environment


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

Table 3: HgH radical: Calculated isotropic hyperfine resapectively.𝑎 Hg PN FN HF 8357 7623 Slater 5515 5027 VWN 7305 6656 BLYP 7248 6600 PBE 7266 6619 BP 7395 6734 PBE0 7755 7067 BHLYP 8170 7445

coupling constants (MHz) for Hg and H, H PN 679 703 685 797 724 721 703 753

FN 678 702 684 796 723 720 702 752

Scalar X2C calculations with 𝑐 = 137.036 au. PN = point nucleus, FN = finite nucleus, for the electron-nucleus potential. Uncontracted ANO-RCC basis sets. For the PN calculations, 10 high-exponent 1s GTOs were added. For further details and isotope information see Section 3. 2-component X2C KS calculations (this work, FN) give 6427 (PBE) and 6934 (PBE0) for Hg. Selected data (FN) from the literature: Hg: 8238 (scalar NESC, HF, contr. SARC basis),28 6244 (4-component KS, BP),51 6615 (unscaled scalar ZORA with SO corrections, KS, PBE) and 7047 (PBE0),52 6368 (SO ZORA, ‘MA’ KS, PBE) and 6813 (PBE0),53 7002 (expt.)83 H: 745 (scalar ZORA with SO corrections, KS, PBE) and 727 (PBE0),52 668 (SO ZORA, ‘MA’ KS, PBE) and 653 (PBE0),53 710 (expt.)83 𝑎

Computational results for the Hg hyperfine coupling in HgH from 2-component X2C calculations, i.e. including SO coupling, for the PBE and PBE0 functional and otherwise identical computational settings are provided in the footnote of Table 3 along with selected reference data from the literature. Finite nuclear volume effects are very pronounced for Hg and negligible for H, as expected. The scalar X2C data in the body of the table do not include any effects from the PSO mechanism, which probes the magnetization in the electronic system resulting from an orbital angular momentum. For orbitally non-degenerate electronic states, such as the ground state of HgH, SO coupling is the only source of an orbital angular momentum. Previous linear response calculations using the quasi-relativistic ZORA framework52 produced a first-order PSO-SO perturbation correction of -200 MHz for the Hg hyperfine coupling in this molecule. The present scalar X2C data for Hg are lower than the previous ZORA results excluding the PSO-SO correction by a comparable amount. This difference is not mainly due to ZORA vs. X2C or other technical aspects but because in the ZORA benchmarks in References 52 and 53 a slightly larger, optimized, interatomic distance of 1.7877Å was used instead of the 1.7662Å of the present work. For comparison, a scalar ZORA PBE0 calculation was carried out for this work with the 1.7662Å distance, but otherwise identical settings to those used in Reference 52. This calculation produced a Hg hyperfine coupling FC+SD contribution of 7129 MHz, which is not far from the FN X2C value of 7067


Page 16 of 23

Page 17 of 23

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


in Table 3, and a PSO-SO perturbation correction of -213 MHz. The corresponding data for the proton are 709 and -2 MHz. In calculations where SO coupling is included variationally in the ground state, the SO effects on the hyperfine coupling reflect a combination of SO effects on the FC+SD mechanism and the presence of a PSO contribution. With the PBE functional (see Table 3 footnote) the combination of these effects amounts to a -192 MHz correction for the Hg hyperfine coupling in HgH, relative to the scalar X2C FC+SD contribution. With the PBE0 functional, the correction becomes smaller, -133 MHz. The scalar and SO X2C calculations with the PBE0 functional give 7067 and 6934 MHz, respectively, for the FN model, and closely bracket the experimental isotropic Hg hyperfine coupling of 7002 MHz.83 Likewise, the experimental proton hyperfine coupling is well reproduced by the PBE0 and several other functionals. BLYP, and to a lesser degree BHLYP, noticeably overestimate this hyperfine coupling. It has been pointed out previously80 that the LYP correlation functional tends to be less suited for properties, such as hyperfine coupling, that require a balanced account of opposite-spin versus same-spin correlation.

6 Conclusions For the one-electron atomic 𝑛𝑠 states, the accuracy of the hyperfine operators calculated with the same X2C transformation as used for the field-free one-electron Hamiltonian is excellent. Since 𝑠 orbitals and the relativistic analog of the FC mechanism tend to dominate the hyperfine coupling in atoms and molecules (counter-examples are known), this finding points to an easy route toward accurate 2-component relativistic hyperfine calculations for molecular systems. One-electron SIEs in the hyperfine coupling, as introduced by approximate functionals in KS calculations, are very small for one-electron 𝑛𝑠 states when compared to the correlation effects on this property in manyelectron systems. In agreement with previous work, it is found that finite-nuclear volume effects on the hyperfine coupling can be rather large, on the order of -10% of the point-nucleus value, if the coupling is dominated by valence 𝑠 orbital contributions from a heavy element. The availability of suitably accurate 2-component basis set representations of electron-nucleus hyperfine operator matrices further appears to offer an easy route toward X2C based calculations of 𝐽 -coupling tensors, which are known to afford very large relativistic effects if they involve heavy nuclei.84

Acknowledgments The author acknowledges financial support from the U.S. Department of Energy, Office of Basic Energy Sciences, Heavy Element Chemistry program, under grant DE-SC0001136 (formerly DEFG02-09ER16066). Thanks to Dr. Michel Dupuis for helpful clarifications regarding Reference 58. The author is indebted to Dr. Stefan Knecht for helpful discussions on the topic of relativistic 17 ACS Paragon Plus Environment


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

hyperfine coupling calculations that started during a sabbatical spent at ETH Zürich in 2015, and he thanks Dr. Trond Saue for constructive comments on the manuscript.

Supporting Information Available Additional calculated data demonstrating self-interaction errors for the one-electron systems in KS calculations. This information is available free of charge via the Internet at http://pubs. acs.org/.

References (1) Pyykkö, P. Chem. Rev. 2012, 112, 371–384. (2) Kutzelnigg, W. Chem. Phys. 2012, 395, 16-34. (3) Saue, T. ChemPhysChem 2011, 12, 3077–3094. (4) Liu, W. Mol. Phys. 2010, 108, 1679-1706. (5) Autschbach, J. J. Chem. Phys. 2012, 136, 150902. (6) Pyykkö, P. Annu. Rev. Phys. Chem. 2012, 63, 45-64. (7) Moss, R. E. Advanced Molecular Quantum Mechanics; Chapman and Hall: London, 1973. (8) Grant, I. P. Relativistic Quantum Theory of Atoms and Molecules; Springer: New York, 2007. (9) Dyall, K. G.; Fægri, K. Jr. Relativistic Quantum Chemistry; Oxford University Press: New York, 2007. (10) Reiher, M.; Wolf, A. Relativistic quantum chemistry. The fundamental theory of molecular science; Wiley-VCH: Weinheim, 2009. (11) Dyall, K. G. J. Chem. Phys. 1997, 106, 9618-9626. (12) Jensen, H. J. A. “REHE 2005 conference, Mülheim, Germany”, April, 2005. (13) Kutzelnigg, W.; Liu, W. J. Chem. Phys. 2005, 123, 241102. (14) Kutzelnigg, W.; Liu, W. Mol. Phys. 2006, 104, 2225-2240.


Page 18 of 23

Page 19 of 23

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


(15) Liu, W.; Peng, D. J. Chem. Phys. 2006, 125, 044102. (16) Ilias, M.; Saue, T. J. Chem. Phys. 2007, 126, 064102. (17) Filatov, M.; Dyall, K. G. Theor. Chem. Acc. 2007, 117, 333-338. (18) Liu, W.; Kutzelnigg, W. J. Chem. Phys. 2007, 126, 114107. (19) Peng, D.; Liu, W.; Xiao, Y.; Cheng, L. J. Chem. Phys. 2007, 127, 104106. (20) Sikkema, J.; Visscher, L.; Saue, T.; Iliaš, M. J. Chem. Phys. 2009, 131, 124116. (21) Liu, W.; Peng, D. J. Chem. Phys. 2009, 131, 031104-4. (22) Peng, D.; Reiher, M. Theor. Chem. Acc. 2012, 131, 1081. (23) van Wüllen, C.; Michauk, C. J. Chem. Phys. 2005, 123, 204113. (24) Sun, Q.; Liu, W.; Xiao, Y.; Cheng, L. J. Chem. Phys. 2009, 131, 081101-4. (25) Aucar, G. A.; Saue, T.; Visscher, L.; Jensen, H. J. A. J. Chem. Phys. 1999, 110, 6208-6218. (26) Komorovsky, S.; Repisky, M.; Malkina, O. L.; Malkin, V. G.; Ondik, I. M.; Kaupp, M. J. Chem. Phys. 2008, 128, 104101. (27) Olejniczak, M.; Bast, R.; Saue, T.; Pecul, M. J. Chem. Phys. 2012, 136, 014108-13 Erratum: J. Chem. Phys. 136, 239902 (2012). (28) Filatov, M.; Zou, W.; Cremer, D. J. Phys. Chem. A 2012, 116, 3481-3486. (29) Cheng, L.; Gauss, J. J. Chem. Phys. 2011, 135, 084114–084114-7. (30) Autschbach, J.; Peng, D.; Reiher, M. J. Chem. Theory Comput. 2012, 8, 4239-4248. (31) NWChem quantum chemistry program package. URL http://www.nwchem-sw.org. Accessed 07/15. (32) Foldy, L. L.; Wouthuysen, S. A. Phys. Rev. 1950, 78, 29. (33) Heully, J. L.; Lindgren, I.; Lindroth, E.; Lundquist, S.; Mårtenson-Pendrill, A. M. J. Phys. B 1986, 19, 2799-2815. (34) Dyall, K. G. J. Chem. Phys. 1994, 100, 2118-2127. (35) Sun, Q.; Liu, W.; Kutzelnigg, W. Theor. Chem. Acc. 2011, 129, 423–436. 19 ACS Paragon Plus Environment


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

(36) Stanton, R. E.; Havriliak, S. J. Chem. Phys. 1984, 81, 1910–1918. (37) Hennum, A. C.; Klopper, W.; Helgaker, T. J. Chem. Phys. 2001, 115, 7356-7363. (38) Malkin, E.; Malkin, I.; Malkina, O. L.; Malkin, V. G.; Kaupp, M. Phys. Chem. Chem. Phys. 2006, 8, 4079 - 4085. (39) Autschbach, J. ChemPhysChem 2009, 10, 2274-2283. (40) Andrae, D. Phys. Rep. 2000, 336, 413-527. (41) Ditchfield, R. Mol. Phys. 1974, 27, 789-807. (42) Krykunov, M.; Autschbach, J. J. Chem. Phys. 2005, 123, 114103-10. (43) Xiao, Y.; Liu, W.; Cheng, L.; Peng, D. J. Chem. Phys. 2007, 126, 214101. (44) Kutzelnigg, W. Phys. Rev. A 2003, 67, 032109-12. (45) Visscher, L. Adv. Quantum Chem. 2005, 48, 369–381. (46) Kowalewski, J. Annu. Rep. NMR Spectrosc. 1982, 12, 81-176. (47) Autschbach, J.; Ziegler, T. J. Chem. Phys. 2000, 113, 9410-9418. (48) Sternheim, M. M. Phys. Rev. 1962, 128, 676-677. (49) Peng, D.; Reiher, M. J. Chem. Phys. 2012, 136, 244108. (50) Autschbach, J.; Ziegler, T. J. Chem. Phys. 2000, 113, 936-947. (51) Malkin, E.; Repisky, M.; Komorovsky, S.; Mach, P.; Malkina, O. L.; Malkin, V. G. J. Chem. Phys. 2011, 134, 044111-8. (52) Aquino, F.; Pritchard, B.; Autschbach, J. J. Chem. Theory Comput. 2012, 8, 598-609. (53) Verma, P.; Autschbach, J. J. Chem. Theory Comput. 2013, 9, 1932-1948. (54) Chipman, D. M.; Rassolov, V. A. J. Chem. Phys. 1997, 107, 5488-5495. (55) Sharkas, K.; Pritchard, B.; Autschbach, J. J. Chem. Theory Comput. 2015, 11, 538-549.


Page 20 of 23

Page 21 of 23

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


(56) Bylaska, E. J.; de Jong, W. A.; Govind, N.; Kowalski, K.; Straatsma, T. P.; Valiev, M.; van Dam, J. J.; Wang, D.; Apra, E.; Windus, T. L.; Hammond, J.; Autschbach, J.; Aquino, F.; Nichols, P.; Hirata, S.; Hackler, M. T.; Zhao, Y.; Fan, P.-D.; Harrison, R. J.; Dupuis, M.; Smith, D. M. A.; Glaesemann, K.; Nieplocha, J.; Tipparaju, V.; Krishnan, M.; VazquezMayagoitia, A.; Jensen, L.; Swart, M.; Wu, Q.; Van Voorhis, T.; Auer, A. A.; Nooijen, M.; Crosby, L. D.; Brown, E.; Cisneros, G.; Fann, G. I.; Fruchtl, H.; Garza, J.; Hirao, K.; Kendall, R.; Nichols, J. A.; Tsemekhman, K.; Wolinski, K.; Anchell, J.; Bernholdt, D.; Borowski, P.; Clark, T.; Clerc, D.; Dachsel, H.; Deegan, M.; Dyall, K.; Elwood, D.; Glendening, E.; Gutowski, M.; Hess, A.; Jaffe, J.; Johnson, B.; Ju, J.; Kobayashi, R.; Kutteh, R.; Lin, Z.; Littlefield, R.; Long, X.; Meng, B.; Nakajima, T.; Niu, S.; Pollack, L.; Rosing, M.; Sandrone, G.; Stave, M.; Taylor, H.; Thomas, G.; van Lenthe, J.; Wong, A.; Zhang, Z. “NWChem, A Computational Chemistry Package for Parallel Computers, Version 6 (2012 developer’s version)”, Pacific Northwest National Laboratory, Richland, Washington 99352-0999, USA., 2012. (57) Valiev, M.; Bylaska, E. J.; Govind, N.; Kowalski, K.; Straatsma, T. P.; Dam, H. J. J. V.; Wang, D.; Nieplocha, J.; Apra, E.; Windus, T. L.; de Jong, W. A. Comp. Phys. Comm. 2010, 181, 1477–1489. (58) Dupuis, M. Comp. Phys. Comm. 2001, 134, 150-166. (59) Roos, B. O.; Lindh, R.; Malmqvist, P.; Veryazov, V.; Widmark, P. J. Phys. Chem. A 2005, 109, 6575–6579. (60) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200-1211. (61) Becke, A. D. Phys. Rev. A 1988, 38, 3098-3100. (62) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B 1988, 37, 785–789. (63) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865-3868. (64) Perdew, J. P. Phys. Rev. B 1986, 33, 8822-8824. (65) Becke, A. D. J. Chem. Phys. 1993, 98, 5648–5652. (66) Adamo, C.; Barone, V. J. Chem. Phys. 1999, 110, 6158-6170. (67) Jensen, F. Theor. Chem. Acc. 2009, 126, 371–382.



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

(68) Huber, K.; Herzberg, G. Constants of Diatomic Molecules. In NIST Chemistry WebBook, NIST Standard Reference Database Number 69; Mallard, W. G.; Linstrom, P. J., Eds.; National Institute of Standards and Technology: Gaithersburg MD, 20899, 2005 data prepared by J.W. Gallagher and R.D. Johnson, III. URL: http://webbook.nist.gov. (69) Autschbach, J.; Pritchard, B. Theor. Chem. Acc. 2011, 129, 453-466. (70) Autschbach, J.; Patchkovskii, S.; Pritchard, B. J. Chem. Theory Comput. 2011, 7, 21752188. (71) Stwalley, W. C. J. Chem. Phys. 1975, 63, 3062-3080. (72) Stone, N. “Table of Nuclear Quadrupole Moments, International Atomic Energy Agency publication INDC(NDS)-650”, 2013 URL: http://www-nds.iaea.org/ publications/indc/indc-nds-0650.pdf, accessed 08/2016. (73) Pyykkö, P.; Pajanne, E.; Inokuti, M. Int. J. Quantum Chem. 1973, 7, 785-806. (74) Pyykkö, P. Chem. Phys. 1977, 22, 289-296. (75) Lipas, P. O.; Pyykkö, P.; Pajanne, E. J. Chem. Phys. 1973, 58, 3248-3254. (76) Breit, G. Phys. Rev. 1930, 35, 1447-1451. (77) Cohen, A. J.; Mori-Sánchez, P.; Yang, W. Chem. Rev. 2012, 112, 289-320. (78) Brage, T.; Proffitt, C.; Leckrone, D. S. Astrophys. J. 1999, 513, 524-534. (79) Filatov, M.; Cremer, D. J. Chem. Phys. 2004, 121, 5618-5622. (80) Filatov, M.; Cremer, D. J. Chem. Phys. 2005, 123, 124101. (81) Fuller, G. H.; Cohen, V. W. Nucl. Data Tabs. 1969, A5, 433-612. (82) Liberman, S.; Pinard, J.; Duong, H. T.; Juncar, P.; Pillet, P.; Vialle, J.-L.; Jacquinot, P.; Touchard, F.; Büttgenbach, S.; Thibault, C.; de Saint-Simon, M.; Klapisch, R.; Pesnelle, A.; Huber, G. Phys. Rev. A 1980, 22, 2732–2737. (83) Knight, L. B. Jr.; Weltner, W. Jr. J. Chem. Phys. 1971, 55, 2061–2070. (84) Autschbach, J.; Zheng, S. Annu. Rep. NMR Spectrosc. 2009, 67, 1-95.


Page 22 of 23

Page 23 of 23

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


For Table of Contents µe

µN


Relativistic effects on electron-nucleus hyperfine ... - ACS Publications

Recommend Documents