Flexible nuclear screening approximation to the two-electron spin-orbit coupling based on ab initio parameterization.

THE JOURNAL OF CHEMICAL PHYSICS 139, 204106 (2013)

Flexible nuclear screening approximation to the two-electron spin–orbit coupling based on ab initio parameterization Jakub Chalupskýa) and Takeshi Yanai Department of Theoretical and Computational Molecular Science, Institute for Molecular Science, Okazaki, Aichi 444-8585, Japan

(Received 5 September 2013; accepted 8 November 2013; published online 27 November 2013) The derivation, implementation, and validation of a new approximation to the two-electron spin– orbit coupling (SOC) terms is reported. The approximation, referred to as flexible nuclear screening spin–orbit, is based on the effective one-electron spin–orbit operator and accounts for two-electron SOC effects by screening nuclear charges. A highly flexible scheme for the nuclear screening is developed, mainly using parameterization based on ab initio atomic SOC calculations. Tabulated screening parameters are provided for contracted and primitive Gaussian-type basis functions of the ANO-RCC basis set for elements from H to Cm. The strategy for their adaptation to any other Gaussian basis set is presented and validated. A model to correct for the effect of splitting of transition metal d orbitals on their SOC matrix elements is introduced. The method is applied to a representative set of molecules, and compared to exact treatment and other approximative approaches at the same level of relativistic theory. The calculated SOC matrix elements are in very good agreement with their “exact” values; deviation below 1% is observed on average. The presented approximation is considered to be generally applicable, simple to implement, highly efficient, and accurate. © 2013 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4832737] I. INTRODUCTION

Spin–orbit coupling (SOC) is believed to be the most important spin-dependent relativistic effect.1–3 SOC gives rise or significantly contributes to many interesting physical and chemical phenomena, such as intersystem crossing (and spinforbidden chemical reactions), phosphorescence, zero-field splitting, and electron paramagnetic and nuclear magnetic resonance spectra. It is thus evident that SOC should often be included in quantum chemical calculations to correctly describe the problems under examination. The rigorous treatment of SOC can be formulated on the basis of the many-electron Dirac– Coulomb theory with the fully relativistic four-component scheme4–10 or the quasirelativistic two-component scheme based on the exact decoupling of positive- and negative-state components.11–18 However, the complexity of the quantum chemical implementation of these relativistic schemes is in a strict sense quite high and thus is considered to limit their applicability. Although recent progress in the development of two- and four-component methods has provided significant improvements to their computational efficiency, facilitating rigorous calculations of the relativistic effects for larger systems,10, 18–20 there is a continuing demand for highly simplified methods that can determine SOC. Such methods are useful not only in the two-component framework, but also in the approach that incorporates SOC as a correction to the onecomponent non-relativistic, or scalar relativistic, wave functions. This approach is widely used because it can be easily implemented within the existing non-relativistic quantum a) Electronic mail: [email protected]

0021-9606/2013/139(20)/204106/14/$30.00

chemical codes. The post-treatment of SOC is typically conducted by employing perturbation theory21–23 or configuration interaction.24–26 While the scalar relativistic effects are well characterized using one-electron terms alone, a correct description of the SOC effects requires consideration of two-electron SOC, which is known to be generally non-negligible.27, 28 However, evaluation of the two-electron SOC terms involves complex computer implementation and is computationally demanding. Therefore, the use of approximations is required to dramatically simplify the inclusion of SOC effects. The most widely used approximation consists of the so-called mean-field spin– orbit (SO) operator proposed by Heß et al.,29 and inspired by the original ideas of Blume and Watson.27 In this approximation, the full one- plus two-electron SO operator is truncated into an effective one-electron-only operator. This simplifies the evaluation of SOC matrix elements, especially on the side of wave function, because it neglects the effect of mutually doubly excited configurations. However, dependence on the two-electron SOC integrals in the mean-field SO operator remains, so that no reduction arises in the required types of SOC integrals. Although the evaluation of one-electron integrals does not constitute a problem, the same is not true for two-electron SOC integrals, the evaluation of which is one of the large practical tasks in quasirelativistic quantum chemistry. This requires substantially, approximately by an order of magnitude, more computational effort than the evaluation of electron-electron repulsion integrals; therefore, it often becomes the most time and resource consuming part of the entire electronic-structure calculation, and thus limits the size of practically tractable molecules and basis sets.

139, 204106-1

© 2013 AIP Publishing LLC

204106-2

J. Chalupský and T. Yanai

From this perspective, it becomes necessary to search for efficient, although still accurate, approximations to the two-electron SOC terms. The first type of already established approximation is based on neglecting part of the two-electron SOC integrals. The integrals that are typically neglected are due to the short-range SO interaction of multicenter type. This forms the basis of the so-called one-center approximation, which is, in combination with the mean-field approach, implemented, for example, in the favorite code AMFI of Schimmelpfennig.30 The second type of approximation avoids explicit evaluation of some (or all) two-electron SOC integrals by using simplified physical models that mimic two-electron SOC effects. This type covers a broad range of approximations, such as methods that use effective nuclear charges,31 relativistic effective core potentials (ECP),32 or effective potential SO operator.33 In this paper, we present a new type of mean-field-based SOC approximation, which is inspired by and closely related to the screened nuclear spin–orbit (SNSO) approximation proposed by Boettger.34 The important features of this approach, which has been named the flexible nuclear screening spin–orbit (FNSSO) approximation, can be summarized as follows. 1. SOC matrix elements are evaluated using only oneelectron SOC integrals along with a feasible parameterization to account for the two-electron contributions, so that implementation is fairly simple and highly efficient. A similar simplification is employed in some other methods, such as effective nuclear charge and SNSO approximations, whereas it is the feature that contrasts with the complexity of general mean-field approach or AMFI that employ the two-electron SOC integrals in part. 2. Our approach is based on a newly proposed parameterization scheme for the nuclear screening model. It has excellent flexibility, which allows the accurate and systematic adjustment of the screening effect towards the rigorous mean-field description that can be obtained by expensive ab initio calculations. In this work, the screening parameters are predetermined by atomic SOC calculations that involve the evaluation of two-electron SOC integrals, and compiled into tabulated data. 3. This study provides the ready-for-use parameter sets for nuclear screening in atom-centered Gaussian-type basis representations. These have been determined using the first-order Douglas–Kroll–Heß (DKH1) SO operator and Roos’ ANO-RCC basis set in contracted and uncontracted (primitive) forms to cover almost all the elements, i.e., from H to Cm where ANO-RCC is available. An inexpensive scheme for adapting the screening parameters to any other Gaussian basis sets is also provided. We demonstrate the accuracy of the FNSSO approximation on the values of SOC matrix elements for several molecules composed of elements from the lightest to the heaviest. This paper is organized as follows. The details of the present theory are described in Sec. II and the computer implementation in Sec. III. The computational details for the

J. Chem. Phys. 139, 204106 (2013)

determination of screening parameters and benchmark calculations are presented in Sec. IV. Results and discussion are given in Sec. V and conclusions are drawn in Sec. VI. II. THEORY A. Spin–orbit coupling and mean-field approach

The term spin–orbit coupling represents the interaction of an electron spin with the magnetic field caused by its movement relative to other charged particles of the system. In the quasirelativistic theory at the Breit–Pauli (BP) level,35 the SO operator can be written as a simple sum of one- and twoelectron contributions and reads ˆ +1 ˆ j ), Hˆ SO = Hˆ 1el-SO + Hˆ 2el-SO = g(i, h(i) 2 i j =i i (1) where the one-electron operator is (in SI units) defined by ˆ = hˆ i · sî = h(i) hˆ A (i) = hˆ iA · sî A

=

A

e2 −3 ZA riA (riA × pˆ i ) · sî , 2 2 8π ε0 m c A

(2)

and the two-electron operator by ˆ j ) = gˆ ij · (ˆsi + 2ˆsj ) g(i, =−

e2 r −3 (rij × pˆ i ) · (ˆsi + 2ˆsj ). 4π ε0 m2 c2 ij

(3)

Here, pˆ i and sî (pˆ j and sˆj ) are the momentum and spin operators of the ith (jth) electron, ZA is the charge of the Ath nucleus, and riA of magnitude riA and rij of magnitude rij are the positions of the ith electron relative to nucleus A and electron j, respectively. Boldface printed quantities are vectorial with three spatial components. The magnetic field, with which the electronic spin interacts, is thus caused by its relative movement with respect to the nuclei for the one-electron term, and with respect to other electrons for the two-electron term. Moreover, the one- and two-electron operators have opposite signs, and one-electron term is usually the leading term; therefore, the interactions with other electrons may be seen as a screening of the interactions with the nuclei. This perspective forms the basis of some approximations to the twoelectron SOC terms (see, for example, Refs. 31, 32, and 34). The two-electron term is typically divided into two contributions, the so-called spin-same-orbit and spin-other-orbit interactions, which is apparent from Eq. (3). Although this type of SO operator is not variationally stable and can be used only for molecules composed of relatively light elements, the simplicity and explicit form make it both attractive and useful for conceptual and qualitative discussions of SOC effects. It is in this spirit that we refer to the Breit–Pauli SO operator here. For the practical applications presented in this work, the rather more advanced Douglas– Kroll–Heß SO operator36–38 has been used. However, these two types of SO operators have the same structure, and the difference between them is covered only by partially different

204106-3


J. Chem. Phys. 139, 204106 (2013)

schemes for the evaluation of SOC integrals. Thus, all the formulas and arguments given below apply equally well to both. The evaluation of SOC matrix elements using the full SO operator is often too expensive to be practical. Heß et al.29 have introduced an effective one-electron SO operator that substantially simplifies the evaluation. Although this operator is slightly more complicated, because the two-electron SOC has lower symmetry with respect to particle interchange than the electron-electron repulsion, it is a Fock operator very similar to that used in Hartree–Fock theory.39 The SOC matrix element of a pair of Slater determinants that differ by a single spin orbital excitation ψ i → ψ j is given by 1 ˆ ˆ 2)|ψj ψk (ψi ψk |g(1, HijSO = ψi |h(1)|ψ j + 2 k

mean-field SOC matrix elements are then evaluated using Eq. (5). Note that the first two two-electron integrals are coulomb type, while the other two are exchange type. Although this mean-field SO operator is generally applicable to any type of wave function, it does have certain limitations. Most importantly, only an effective one-electron operator is used for construction of all SOC matrix elements, so that the mean-field approach does not account for the effect of mutually doubly excited configurations. Moreover, it is generally necessary to apply spin-averaging to arrive at Eq. (7). However, these limitations are of minor importance to a wide range of problems relevant to chemistry and physics, and the mean-field SO approach represents a powerful tool for quasirelativistic calculations, as well as a good starting point for further simplifications.

ˆ 2)|ψk ψj − ψi ψk |g(1, ˆ 2)|ψk ψj + ψk ψi |g(1, ˆ 2)|ψj ψk ), − ψk ψi |g(1,

(4)

where the summation runs over occupied spin orbitals ψ k . If we now consider an effective one-electron operator ˆf, which is used to approximate the total SO operator by ˆfi · sî , Hˆ SO ≈ Hˆ MF = (5)

B. Flexible nuclear screening spin–orbit approximation

To develop the approximations presented in this work, we begin with only a simple reorganization of the terms from Eq. (7). If we formally rewrite the Fock matrix elements as 2el fij = f1el ij + fij ,

i

and integrate out those parts of the two-electron SOC integrals that depend on the spin of the second electron, then more practically convenient equations involving only integrals over spatial parts of the orbitals may be obtained. For the given pair of determinants, matrix elements of the Fock-like operator ˆf that correspond to its three Cartesian or spherical spatial components are expressed as ˆ + 1 (2ik|ˆg|j k + 4ki|ˆg|kj fij = i|h|j 2 k − 3ik|ˆg|kj − 3ki|ˆg|j k),

(6)

where the summation runs over occupied spatial orbitals k.40 Equation (6) holds exactly only if all spatial orbitals k are occupied by two electrons, which is a relatively special case. In more general cases, such as when a multiconfigurational approach is used, this condition is not fulfilled. However, only a simple modification based on the introduction of the occupation numbers to Eq. (6) is necessary, if we assume equal occupation of spatial orbitals k by both α and β electrons, even for the partially occupied orbitals (this approximation is usually referred to as spin-averaging). Furthermore, by averaging these occupation numbers over more bra and ket determinants, the Fock matrix elements of the mean-field SOC approach may be introduced as ˆ + 1 Nk (2ik|ˆg|j k + 4ki|ˆg|kj fij = i|h|j 4 k − 3ik|ˆg|kj − 3ki|ˆg|j k),

(7)

where Nk ∈ 0, 2 are fixed occupation numbers of spatial orbitals k common to bra and ket wave functions, which reflect the average effect of the two-electron contributions from Eq. (6), i.e., two-electron contribution averaged, and also spin-averaged, over all considered determinant pairs. The

(8)

then we may define the quantity 1el(μ)

(μ)

Qij =

fij

1el(μ)

fij

(μ)

− fij

2el(μ)

=

−fij

1el(μ)

fij

,

(9)

where μ is one of the three spatial components. This quantity represents the ratio of the two- and one-electron parts of the Fock matrix elements, i.e., the degree of screening of the one-electron SO interactions caused by the corresponding two-electron interactions. The Qij quantities are therefore referred to as screening parameters. We may now define a new effective one-electron SO Fock operator based on the screening parameters in two different ways. First, the general two-dimensional screening parameters can be used to define the operator with the relation FNSSO2D (μ) 1el(μ) (μ) 1 − Qij , (10) = fij fij or, second, the screening effect can be approximated using one-dimensional screening parameters, and the operator is defined as FNSSO1D (μ) 1el(μ) (μ) (μ) fij 1 − Qii Qjj . (11) = fij While the former is clearly equivalent to the mean-field SO Fock operator, which is apparent from Eqs. (8)–(10), the latter is a parallel to the operator used in the SNSO approximation (see Ref. 34). Now the most important question is how to obtain the screening parameters without the need to evaluate all of the two-electron SOC integrals at every calculation. The short answer is that the screening parameters can be calculated once, tabulated, and then used continually. However, this is not possible when using only the general expressions given, because the values of the screening parameters are dependent on the molecular geometry, the bra and ket electronic states,

204106-4


the computational method, and the basis set. Thus, some additional approximations must be made. Let us consider a molecule with two atoms, A and B, and examine the interactions of an electron occupying orbital i (we will denote this electron ei here), which belongs to atom A. When ei interacts with nucleus A, which is the one-center type of one-electron interaction, it can be expected that this interaction is screened mainly by other electrons of atom A, especially by inner electrons with respect to ei . On the other hand, for the multi-center type of one-electron interaction, when ei interacts with nucleus B, it can be expected that the dominant contributions to screening come from the electrons of atom B (and less importantly, also from outer electrons of A). This can be described theoretically by formal decomposition of the mean-field SOC Fock matrix elements from Eq. (7) to the form 1 i|hˆ A |j + NkA (2ikA |ˆg|j kA fij = 4 k A A + 4kA i|ˆg|kA j −3ikA |ˆg|kA j −3kA i|ˆg|j kA ) , (12) where kA represent orbitals belonging to atom A. Here, for example, one-center one-electron SOC integrals, for which both i and j belong to A, can be screened “only” by corresponding one-center two-electron integrals, i.e., summation for k runs only over orbitals of A (kA ). These presumptions lead to the first simplification to be performed. The screening of one-center and multi-center one-electron SO interactions are treated separately, because different strategies are applied for each. 1. Screening of one-center spin–orbit interactions

It is apparent from Eq. (7) that Nk occupation numbers of the mean-field approach mediate the screening effect of orbitals k, and thus play a crucial role in this approximation. However, as previously suggested in the original work (see Ref. 29), many reasonable choices of their values lead to fairly accurate SOC results. In a typical multiconfigurational calculation involving two electronic states, Nk would be set to 2 for all core orbitals, and the average of the occupation in the bra and ket wave functions would be used for the valence orbitals. However, the occupation numbers of the valence orbitals may also be set simply to N/M, without losing much of the accuracy, if the valence space comprises N electrons in M orbitals, and this occupation is used for many pairs of electronic states. The robustness of the mean-field approach to the choice of Nk occupation numbers is very convenient, because we may employ the one-center approximation, neglect the influence of chemical bonding on the occupation numbers, and arrive at a situation of independent atoms, as conducted in the AMFI30 approximation. In such a case, the one-center screening parameters are basically dependent only on the basis set and can be determined from atomic calculations. The number of unique screening parameters of an atom is substantially reduced due to its spherical symmetry and the properties of the SO operator. If we assume that the Fock operator in Eq. (7) is defined in atomic orbitals, we can use a

J. Chem. Phys. 139, 204106 (2013)

more specific definition of the one-center screening parameters Q(μ) ij = Qij = Q(l, ni , nj ) =

2el −fi(n i ,l,mi )j (nj ,l,mj ) 1el fi(n i ,l,mi )j (nj ,l,mj )

, (13)

where only one of the (2l + 1)2 combinations of magnetic quantum numbers is necessary for the definition, as reflected by the dependence of Q on only the angular momentum l and shells (principal quantum numbers) n of atomic orbitals i and j. Here, any combination with a non-zero one-electron SOC integral (a necessary but not sufficient condition is mi = mj ) can be selected, for example, the screening parameter of the 2p shell of an atom can be defined by the z component of one- and two-electron SOC contributions between the 2px and 2py atomic orbitals. Thus, for any given combination of shells with the same angular momentum, there is only one unique screening parameter, although there are generally more nonzero and unique values of one- and two-electron SOC integrals. Such screening parameters can already be very simply tabulated for any element and basis set, and then used according to Eqs. (10) and (11) for the screening of one-electron SOC integrals in atomic-orbital basis. Such a table comprises separate blocks for different values of angular momentum, and each block is a square symmetric matrix with dimension equal to the number of basis functions of the particular angular momentum. The difference between the schemes of two-dimensional (Eq. (10)) and one-dimensional (Eq. (11)) screening thus lies in using the entire blocks (Qij ) or only their diagonal elements (Qii ), respectively. Efficiency and accuracy are not the only concerns, but also the general applicability of this approach; therefore, an additional approximation is now introduced that removes the dependence on the choice of basis set. If we consider a stationary electronic state of an atom, its electron density has a constant spherically symmetric distribution. Consequently, for a given electron at a given distance from the nucleus, the ratio of its one-electron and mean-field two-electron SOC contributions, which are used to define the degree of screening, is a constant number. Therefore, a representative set of primitive functions may be selected, in our case Gaussian-type orbitals (GTOs), of which the exponents reflect the electron-nucleus distance, and the values of screening parameters for these are then tabulated. The screening parameters for the primitive GTOs of any basis set can then be determined from this tabulated dependence, for example, by employing simple interpolations. The typical shape of such dependence is shown for diagonal (one-dimensional) screening parameters in the upper panel of Fig. 1. As expected, the screening effect vanishes close to the nucleus, i.e., for compact GTOs with large exponents. On the other hand, electrons at larger distance from the nucleus are basically screened fully, because the screening parameter values approach 1 for small exponents. Moreover, it is interesting that the shape of these dependencies does not differ significantly for most atoms, as shown in the lower panel of Fig. 1. Therefore, it could be considered that one (fitted) distribution per angular momentum could be used to determine the screening parameters of primitive GTOs for any

204106-5


J. Chem. Phys. 139, 204106 (2013)

1

preted as a screening caused by all electrons of neutral atom A. For strongly charged atoms, slightly different screening parameter values may be more adequate. These may thus be set, for example, to pA /ZA , where pA can be either the atomic population obtained from a population analysis, or the number of electrons of neutral atom A lowered by its formal oxidation state. However, different choices of these values lead only to fairly small changes in the SOC matrix elements in most cases.

p d f

0.8

Qii

0.6 0.4 0.2 0 -5

0

5 ln(αi)

10

1

15

p

0.8

Qii

0.6 0.4 0.2 0 -5

0

5 ln(αi)

10

15

FIG. 1. Dependence of the one-dimensional screening parameter values on the logarithm of the exponent of primitive GTOs for p, d, and f angular momenta of a uranium atom (upper panel), and for the p angular momentum of all atoms from hydrogen to curium (lower panel).

element. However, to obtain maximal accuracy, we do not adopt such an approach throughout this work. 2. Screening of multi-center spin–orbit interactions

The situation for screening of multi-center SOC integrals is much more theoretically complicated, because the dependence of the screening parameters on the geometry and composition of the molecule under investigation cannot be avoided. However, a very simple and relatively accurate scheme for such screening can be based on few practical observations. As shown in Fig. 1, electrons become fully screened at larger distances from the nucleus. Thus, all the multi-center integrals should have a very small effect on the SOC values, unless they involve orbitals (and nuclei) of two atoms lying relatively close to each other, e.g., atoms forming a bond. However, even if it is so, analysis of the integrals indicates that the screening is very strong (almost complete for most orbitals up to the valence shell) if neither i nor j belongs to the nucleus involved in the one-electron SOC integral. Based on these findings, we propose to simply set the one-dimensional screening parameters that correspond to multi-center interactions either to 0, if the multi-center integrals should be left unscreened, or to 1 to account for their screening. The latter means that the screening of electron interaction in orbital i with nucleus A is complete, and it is inter-

3. Scheme for the screening of all one-electron spin–orbit integrals

The one-center and multi-center interactions have been treated separately; therefore, the results for these two types of interactions must be combined to arrive at the screening for any type of one-electron SOC integral. Generally, none of the one-electron integrals are neglected. However, some of the integrals may be screened completely, depending on the choice of screening parameters. For the calculations presented in this work, we have employed a model for the screening of multi-center interactions caused by all electrons of neutral atoms. In this case, the fully multi-center one-electron integrals that typically have the smallest effect are screened down to zero. The screening of all other integrals, including the semi-multi-center integrals, is described in terms of the tabulated one-center screening parameters, as presented in Table I, which summarizes the scheme for screening available in this model. Thus, in any calculation, all the screening factors are constructed using only the database of onedimensional (and possibly also two-dimensional) atomic screening parameters, and are then used for the multiplication of one-electron SOC integrals. Note that the methods employing only one-dimensional, and a combination of one- and twodimensional screening parameters differ only in the screening of one-center integrals involving atomic orbitals with the same angular momentum. In other cases, the same type of screening is used. Although we refer to these two types of methods as one-dimensional and two-dimensional methods of screening throughout this work, we keep in mind that a significant part of the scheme for screening is common to both. At the end of this section, we emphasize the advantages of the approach presented in this work, some of which are relatively significant. First, this approximation to the twoelectron SOC terms is highly efficient, because only the oneelectron SOC integrals need to be evaluated explicitly, and may thus be used even for very large molecules. Despite the TABLE I. Scheme for the screening of one-electron SOC integrals i|hˆ A |j , where i and j are atom-centered GTOs. Label A refers to the set of GTOs centered on atom A. Multi-center interactions are screened by all electrons of neutral atoms (Qii = 1 for i ∈ / A). Integral type

Conditions

Screening available

Screening factor

One-center

i ∈ A, j ∈ A

Semi multi-center Fully multi-center

i ∈ A, j ∈ /A i∈ / A, j ∈ /A

Two-dimensional One-dimensional One-dimensional One-dimensional

(1 − Qij )

(1 − Qii Qjj ) √ (1 − Qii ) 0

204106-6


J. Chem. Phys. 139, 204106 (2013)

simplicity of this approximation, it can provide highly accurate SOC results, because it is based on a flexible scheme of the screening and respects the relative importance of different types of SOC integrals. The one-center approximation is not generally employed, which thus allows to also account for the effect of multi-center SOC integrals. Moreover, coulomb- and exchange-type two-electron SOC integrals are treated on the same footing, and the same applies to the spin-same-orbit and spin-other-orbit contributions. Finally, but definitely not least, the approximation can be used without further modifications not only in any method based on wave functions, but also in methods based on density functional theory (DFT), because only the one-electron SO operator is used for the construction of SOC matrix elements. C. Comparison to other approximative methods

An alternative strategy for the evaluation of SOC, which employs the effective potential SO operator, has been developed in the framework of DFT. Here, formally only oneˆ must be electron SOC integrals of the form i|∇VKS × p|j calculated using the Kohn–Sham potential ZA ρ(r ) + dr

VKS (r) = −

| |r − r | |r − r A A +

δEX [ρ] δEC [ρ] + , δρ(r) δρ(r)

TABLE II. Range of one-dimensional screening parameter values that correspond to valence atomic orbitals with various angular momenta determined using C, Fe, Os, Pb, and Cm atoms (only values of occupied orbitals have been included). Qcoulomb was calculated using the numerical integration of the coulomb term and the same wave functions as that for QFNSSO . Improved values of SNSO screening parameters (ISNSO) proposed in this work are reported as well. l

QSNSO

Qcoulomb

ZA · QFNSSO

QISNSO

p d f

2 10 28

1.7–2.3 9.8–12.0 28.4–31.0

2.8–3.4 11.6–13.5 30.0–32.6

2.85 11.7 31.1

of the SNSO approach, we cannot generally expect to achieve the accuracy of FNSSO, due to the much smaller flexibility of the screening (the same argument can also be used for an approach based on effective nuclear charges). Despite the incompleteness of the description based only on coulomb potential, it does provide a reasonable interpretation of the screening parameters. These represent the number of (inner) electrons, which screen the corresponding one-electron SO interaction (see Ref. 34). The screening parameters provided here are thus tabulated in the same form, i.e., the one-center FNSSO screening parameters defined in Eq. (13) are multiplied by the nuclear charge of particular atoms, as is shown in Table II.

(14)

where the four terms represent electron-nuclear attraction, coulomb, exchange, and correlation potentials, respectively. Although this form of SO operator is computationally more attractive than the general mean-field operator, its use is restricted to DFT, and the evaluation of SOC integrals requires substantially more computational effort than that required for the methods presented in this work. Moreover, it has been noted by Neese41 that the exchange contribution is of the wrong sign and size. The importance of this deficiency can be easily understood, when we compare the screening parameters obtained by evaluating the coulomb term alone, and all two-electron integrals (see Table II). The coulomb contribution, or more precisely, the contribution of spin-same-orbit coulomb type integrals, encompasses approximately 60%– 70%, 85%–90%, and 95% of the total two-electron SOC corresponding to the p, d, and f orbitals of atoms, respectively. Thus, this approach can lead to a 20%–30% error in the overall SOC,41, 42 especially in the case of light atoms from the p block, unless some empirical corrections are introduced (Neese41 suggested to scale the exchange contribution by a factor of −2). A very similar picture can also be found for the SNSO approximation. The screening parameters used in this approach (see Ref. 34 and Table II) should approximate the effect of only the coulomb potential, so that it tends to overestimate SOC. Therefore, we propose to use higher values of SNSO screening parameters that would cover the screening effect of all the terms. These could be 2.85 for p, 11.7 for d, and 31.1 for f, which represent the average FNSSO values for the valence orbitals of the first rows of particular blocks in the periodic table. Although these significantly improve the accuracy

III. IMPLEMENTATION

The implementation of the present approach is relatively straightforward and requires only fairly simple modifications of the existing code for the evaluation of one-electron SOC integrals. Because only BP and DKH types of SO operators are used in our quantum chemistry package named ORZ, we restrict ourselves to these here. Van Wüllen and Michauk43 noted that it is more general, and also convenient, to apply screening to the matrix elements ˆ × pˆ (in a one-component framework, V is the nuclear of pV attraction potential) before these are passed to the subroutines for DKH transformation. In this case, the same scheme for the screening can be principally used for BP and arbitrary-order DKH SO operators. However, a question remains in what level of relativistic theory should be used to determine the screening parameters. In this work, we use DKH1 one- and two-electron SOC integrals for the tabulation, and their values thus represent screening at the DKH1 level. These screening parameters can also be used for the screening of BP SOC integrals in molecules composed of lighter atoms. It is not yet known whether their values are also sufficiently suited for the “direct” screening of higher-order DKH one-electron SOC integrals in systems that contain the heaviest atoms. However, this possible problem can be simply overcome by accounting for the two-electron SOC by screening at the DKH1 level and adding this approximate first-order DKH two-electron contribution to the arbitrary-order DKH one-electron contribution. This approach is based on quite a general scheme of screening, and because the summation over nuclei is typiˆ cally conducted in the subroutines for evaluation of the i|pV ˆ matrix elements (as also in ORZ), we apply the × p|j

204106-7


screening already to the potential V , i.e., the nuclear charges are scaled using the corresponding factors from Table I beˆ fore passing them to the subroutines for evaluation of i|pV ˆ . Moreover, the matrix elements of pV ˆ × pˆ are at first × p|j calculated in the basis of atom-centered primitive GTOs, even for the BP operator, and general contraction of the basis set is not used. This strategy allows us to distinguish between the screening of one-center and multi-center one-electron SOC integrals, and to use different screening parameter values for any pair of not only primitive, but also contracted Gaussiantype basis functions. Basically, the same procedure is used in our own implementation of the SNSO approximation, except the definition of the screening parameters. However, we account for the effect of multi-center SOC integrals even in this approach, unlike the original formulation suggests. The algorithm for the evaluation of “screened” oneelectron SOC integrals, which is implemented in ORZ and covers all SNSO and FNSSO models used in this work, can be concisely summarized as follows. prepare all necessary screening parameters Qij for primitive GTOs i do for primitive GTOs j do scale nuclear charges by screening factors from Table I ˆ × p|j ˆ integrals calculate submatrix of i|pV end for end for perform DKH SOC calculation transform integrals to the basis of contracted GTOs

In the case of FNSSO, screening parameters that correspond to one-center SO interactions are, either for contracted or primitive GTOs, read from the database. Those for multicenter interactions are set to 1 in the present model (see prior text for other alternatives). If we employ FNSSO based on screening parameters for primitive GTOs and a basis set other than that used for their tabulation, then the screening parameters for the primitive GTOs are determined from their tabulated dependence on the two-dimensional space of GTO exponents in the logarithmic scale by linear interpolations between the nearest exponents (a sequence of three one-dimensional linear interpolations is used in our implementation). IV. COMPUTATIONAL DETAILS A. Ab initio tabulation of screening parameters

Calculations were performed on atoms and the two different previously suggested strategies were used to tabulate the screening parameters Qij that correspond to one-center SO interactions. The ANO-RCC basis set44 was applied for both strategies, because it is well suited to these types of calculations and is available for all elements up to curium. While ANO-RCC based on contracted GTOs with the quadruple-ζ quality plus polarization functions (QZP) has been used for the first strategy, which leads to basis-set-dependent screening parameters, all primitive GTOs of ANO-RCC have been used as the representative set of primitive basis functions for

J. Chem. Phys. 139, 204106 (2013)

the second strategy, which allows for basis-set-independent screening. Except for the difference in basis set, the same computational protocol has been applied to both of these strategies. Tabulated screening parameters for all atoms up to Cm have been determined from calculations based on complete active space self-consistent field (CASSCF)45 wave functions. CASSCF active spaces were composed of the most important atomic orbitals, i.e., valence s and p for main group elements (1s and 2s for H and He), valence s and d for transition metals (and also p for Zn, Cd, and Hg), and valence s, d, and f for lanthanoids and actinoids (only s and d for La, Lu, Ac, and Th). For atoms with closed shell configuration, the lowest excited states of triplet multiplicity have been used, while ground states, varying from doublet to nonet, have been calculated for all other atoms. In these calculations, all of the degenerate electronic states (from one to several) that form the corresponding atomic term were included and the orbitals were averaged over all of them. Scalar relativistic effects have been included via the second-order DKH2 spinless Hamiltonian.36, 46 All atomic CASSCF calculations have been performed using the Molcas quantum chemistry package.47 The CASSCF wave functions obtained were then used for the construction of the mean-field operator from Eq. (7), where the averaged occupation numbers of core and active orbitals were utilized for the definition of Nk . Finally, the values of screening parameters were determined from the one-electron and mean-field two-electron SOC integrals in atomic-orbital basis and tabulated. For this purpose, the firstorder DKH1 SO operator was employed. All SOC calculations used for the tabulation have been performed using the SOSS program of Kývala.48 Both sets of tabulated screening parameters, which correspond to contracted and primitive GTOs of ANO-RCC, are provided in the supplementary material.49 B. Benchmark calculations

For the benchmark calculations, a fairly representative set of molecules was selected, which contain atoms from the lightest to the heaviest. Furthermore, the set of molecules comprises systems with dominant contributions to SOC coming from p, d, and f electrons; the set contains a series of biradicals XH2 (X = C, Si, Ge, Sn, and Pb), H2 X2 (X = O, S, Se, Te, and Po), CO+ , NO, octahedral hexaaqua complexes X(H2 O)6 (X = Fe3+ , Fe2+ , Co2+ , Ru3+ , and Os3+ ), the tetrahedral [NiCl4 ]2− complex, and NdO2 and UO2 oxides. The geometries of all molecules, except the transition metal complexes, have been obtained from DFT optimizations employing the B3LYP50 functional and def2-TZVP51 basis set with the corresponding ECP for heavier atoms (defTZVP52 with the ECP for Nd and U). Ground electronic states have been used for these optimizations, i.e., doublet for CO+ and NO, and triplet for CH2 , NdO2 , and UO2 (the triplet is expected to be the ground state for the latter two, although it is not certain due to the high density of states), and singlets for all others. The geometries of the octahedral hexaaqua complexes of Ru3+ and Os3+ have been taken from previous

204106-8


work,53 where the strategy used has also been applied to optimizations of other complexes presented here, i.e., the PBE54 functional, def2-SVP51 basis set, and RI-J approximation52, 55 have been used. These calculations have also been performed for the ground electronic state, which is singlet in the case of [Fe(H2 O)6 ]2+ , triplet for [NiCl4 ]2− , and doublet in other cases. The Turbomole program was utilized for all geometry optimizations.56 Molecular geometries used for the SOC calculations presented in this work can be found in the supplementary material.49 CASSCF wave functions, which were employed for subsequent evaluation of the SOC matrix elements, were obtained using the Molcas and ORZ programs. While for molecules comprising only main group elements, all valence electrons and orbitals have been included in the active space, more limited active spaces have been used for other molecules. Active space was constructed from 3 valence d orbitals and the corresponding electrons of the t2g subshell for complexes of Fe3+ , Ni2+ , Ru3+ , and Os3+ , from the valence eg subshell for complexes of Co2+ , and from the valence d shell for complexes of Fe2+ . In the case of NdO2 and UO2 molecules, 4 and 3 valence orbitals and their electrons were used, respectively; the 4f orbitals of Nd, and the 7s and 5f orbitals of U. The choice of low-lying electronic states, for which the wave functions and SOC matrix elements were calculated, is apparent from the tables provided in Sec. V. For all considered states of a particular molecule, one common set of molecular orbitals was used to provide directly comparable SOC values.57 These have been obtained from state-averaged CASSCF calculations (averaging over all included states) for CO+ , NO, transition metal complexes except that of Fe2+ , and NdO2 and UO2 molecules, while state-specific molecular orbitals of the lowest triplet have been used for biradicals, and those of the ground state singlet for [Fe(H2 O)6 ]2+ , hydrogen peroxide, and its heavier congeners. To provide a deeper insight into the accuracy of the presented methods, several strategies have been used for the evaluation of SOC matrix elements, all employing DKH1 type of SO Hamiltonian. The following SO operators have been applied: (1) full (labeled as SOfull in the tables provided later); (2) one-electron (SO1el ); (3) mean-field (SOMF ); (4) full using one-center approximation (SO1c ); (5) AMFI; (6) SNSO; (7) SNSO with improved values of screening parameters (ISNSO); (8) FNSSO using one-dimensional screening parameters of contracted GTOs (FNSSOC1D ); (9) FNSSO using two-dimensional screening parameters of primitive GTOs (FNSSOP2D ). Our choice of the two methods based on the scheme of screening presented here will be explained later. SOC matrix elements have been evaluated using the SOSS (operators 1–4), Molcas (5), and ORZ (2 and 6–9) programs. For the major part of the benchmark calculations, the ANO-RCC basis set contracted to triple-ζ quality plus polarization functions (TZP) has been used for all atoms, except ligands of transition metal complexes, for which doubleζ quality plus polarization functions (DZP) contraction has been employed. Furthermore, the accuracy of the presented methods has been tested on selected molecules using various basis sets, including STO-3G,58 3-21G,59 and 6-31G,60 cc-pVXZ (X = D, T, and Q),61 and ANO-RCC with DZ and

J. Chem. Phys. 139, 204106 (2013)

DZP contraction. In all calculations involving the ANO-RCC basis set, scalar relativistic effects have been included via the DKH2 spinless Hamiltonian.

V. RESULTS AND DISCUSSION A. Calculations employing ANO-RCC basis set

The main data set, given in Table III, covers the wide range of SO interaction strength. According to results obtained by exact treatment of SOC at the DKH1 level (SOfull ), the magnitude of the SOC matrix elements is varied from about 10 to almost 6000 cm−1 . Although there is a general increase with the nuclear charge, it is also strongly dependent on the character of the electronic states, especially on the angular momentum of electrons involved in the interaction. We have already seen that the relative importance of two-electron SOC, which is represented in this work by the screening parameters, increases with the angular momentum (see Table II or the supplementary material49 ). The SOC in the molecules examined is predominantly from electrons of just one angular momentum (p up to NO, d for transition metal complexes, and f for NdO2 and UO2 ), and thus the very same effect on the values of molecular SOC can be observed here. While for molecules composed of light main group elements, the size of the twoelectron SOC contribution of p electrons is approximately 30%–50% of the size of their one-electron contribution, this decreases to only approximately 4% for their heaviest congeners. However, that contribution is still approximately 18% for the d electrons of [Os(H2 O)6 ]3+ , and almost 36% for the f electrons of UO2 . Quite surprisingly, the largest two-electron contribution to SOC is found here for NdO2 , which is more than 50% of that for the one-electron SOC. The importance of an appropriate description of the two-electron SOC is reflected by the relative error of the one-electron approximation (SO1el ), which can sometimes be even over 100%, i.e., the error of the prediction can be larger than the actual value to predict. On average, an error of approximately 45% is observed for the set, if the two-electron SOC integrals are neglected. It has been already reported (see, for example, review of Heß and Marian1 ) that the mean-field approximation is usually very accurate. When the average of the bra and ket CASSCF occupation numbers is used for the definition of Nk from Eq. (7), the error is far below 1% for most of the SOC matrix elements presented here (see SOMF column of Table III). Maximal errors of −1.4% and 1.1% are found for [NiCl4 ]2− and NdO2 , respectively, and the mean absolute error is only 0.3%. Thus, the role of doubly excited configurations is fairly negligible for SOC between the studied electronic states. This is of course very convenient, because the present approach is based on further approximations to the mean-field SO operator. In many cases, similar performance is also observed for the one-center approximation, either exact (SO1c ) or combined with the mean-field approach and predetermined occupation numbers of atomic orbitals (AMFI). However, the size of the error for the one-center approximation can in problematic cases, where multi-center SOC integrals are relatively large, and one- and two-electron multicenter contributions do not compensate each other, exceed

H2 S2 H2 Se2 H2 Te2 H2 Po2 CO+ NO [Fe(H2 O)6 ]3+

[Fe(H2 O)6 ]2+

[Co(H2 O)6 ]2+ [NiCl4 ]2−

[Ru(H2 O)6 ]3+

NdO2 UO2 MPE MAPE

SO1el

PE

SOMF

PE

SO1c

PE

AMFI

PE

SNSO

PE

ISNSO

PE

FNSSOC1D

PE

FNSSOP2D

PE

S0 /T1 S0 /T1 S0 /T1 S0 /T1 S0 /T1 S0 /T1 S0 /T2 S0 /T1 S0 /T2 S0 /T1 S0 /T2 S0 /T1 S0 /T2 S0 /T1 S0 /T2 D1 /D2 D2 /D3 D1 /D2 D1 /D2 D1 /D3 D2 /D3 S0 /T1 S0 /T2 S0 /T3 D1 /D2 T1 /T2 T1 /T3 T2 /T3 D1 /D2 D1 /D3 D2 /D3 D1 /D2 D1 /D3 D2 /D3 T1 /T2 T1 /T2

10.9 57.5 391.9 936.6 3310.4 41.1 52.8 162.6 184.9 784.6 885.4 1698.9 1968.8 5225.9 5994.4 18.1 60.3 62.6 235.7 230.1 224.0 533.6 543.9 543.9 15.0 322.2 321.7 321.7 553.8 533.9 511.0 1805.1 1714.3 1610.8 2442.3 3005.1

21.6 74.6 439.7 1006.1 3453.0 61.4 79.2 198.9 226.2 867.0 978.2 1813.0 2100.8 5433.3 6231.2 25.5 97.7 101.8 410.6 400.4 389.4 973.4 989.7 989.7 26.8 553.8 553.0 552.9 770.3 742.5 710.6 2190.5 2080.3 1954.7 5052.7 4670.6

97.6 29.7 12.2 7.4 4.3 49.4 49.8 22.4 22.3 10.5 10.5 6.7 6.7 4.0 4.0 40.9 61.9 62.5 74.2 74.0 73.8 82.4 82.0 82.0 78.9 71.9 71.9 71.9 39.1 39.1 39.0 21.4 21.4 21.3 106.9 55.4 44.7 44.7

11.0 57.6 392.0 936.6 3310.4 40.9 52.6 162.5 184.9 784.6 885.4 1698.8 1968.8 5225.9 5994.4 18.2 60.1 62.9 235.4 229.8 223.7 534.9 545.1 545.1 15.1 317.7 317.2 317.2 553.7 533.8 510.9 1805.0 1714.3 1610.7 2468.1 3006.0

0.2 0.1 0.0 0.0 0.0 − 0.5 − 0.4 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.5 − 0.4 0.4 − 0.1 − 0.1 − 0.1 0.3 0.2 0.2 0.4 − 1.4 − 1.4 − 1.4 0.0 0.0 0.0 0.0 0.0 0.0 1.1 0.0 − 0.1 0.3

9.1 57.9 391.7 937.1 3310.5 42.9 54.6 161.8 184.1 785.2 886.2 1698.5 1968.3 5227.3 5995.7 17.9 56.7 66.4 235.7 230.1 224.0 533.4 543.8 543.8 14.9 322.3 321.8 321.8 553.7 533.9 511.0 1805.2 1714.4 1610.8 2443.9 3005.4

− 16.7 0.7 0.0 0.1 0.0 4.3 3.3 − 0.5 − 0.5 0.1 0.1 0.0 0.0 0.0 0.0 − 0.9 − 6.1 6.0 0.0 0.0 0.0 0.0 0.0 0.0 − 0.6 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.1 0.0 − 0.3 1.1

11.0 57.8 391.6 937.2 3309.9 42.0 53.7 161.6 183.8 785.3 886.1 1699.1 1968.7 5228.7 5994.2 17.0 58.0 66.4 228.5 222.9 216.8 539.7 549.6 549.6 15.5 321.7 321.2 321.1 550.9 531.0 508.2 1799.4 1708.9 1605.6 2438.5 2977.3

0.3 0.5 − 0.1 0.1 0.0 2.1 1.6 − 0.6 − 0.6 0.1 0.1 0.0 0.0 0.1 0.0 − 6.3 − 3.8 6.0 − 3.1 − 3.1 − 3.2 1.1 1.0 1.0 3.1 − 0.2 − 0.2 − 0.2 − 0.5 − 0.5 − 0.6 − 0.3 − 0.3 − 0.3 − 0.2 − 0.9 − 0.2 1.2

14.0 63.8 411.8 965.8 3368.7 47.2 60.6 175.2 199.2 818.4 923.5 1745.7 2023.0 5309.0 6089.4 19.9 69.8 74.7 252.6 246.4 239.8 599.1 609.7 609.7 16.9 354.0 353.4 353.4 595.0 573.6 549.0 1901.9 1806.3 1697.3 2691.7 3246.9

28.1 10.8 5.1 3.1 1.8 14.7 14.7 7.8 7.7 4.3 4.3 2.8 2.8 1.6 1.6 9.7 15.7 19.3 7.2 7.1 7.0 12.3 12.1 12.1 12.9 9.9 9.9 9.9 7.4 7.4 7.4 5.4 5.4 5.4 10.2 8.0 8.7 8.7

11.0 59.3 400.1 948.7 3332.8 40.6 52.1 164.6 187.2 796.9 899.2 1716.2 1988.8 5254.3 6026.8 17.5 59.5 63.4 225.7 220.2 214.3 535.6 545.0 545.0 15.2 320.2 319.7 319.6 565.1 544.8 521.5 1852.8 1759.6 1653.5 2430.9 3089.7

0.8 3.0 2.1 1.3 0.7 − 1.3 − 1.4 1.2 1.2 1.6 1.6 1.0 1.0 0.5 0.5 − 3.2 − 1.4 1.2 − 4.2 − 4.3 − 4.3 0.4 0.2 0.2 1.6 − 0.6 − 0.6 − 0.6 2.0 2.0 2.0 2.6 2.6 2.7 − 0.5 2.8 0.4 1.7

11.1 57.8 392.3 937.3 3311.8 40.5 52.2 161.0 183.1 782.9 883.5 1697.3 1966.8 5224.4 5992.6 17.6 59.8 63.6 227.6 222.0 216.0 538.4 548.2 548.1 15.3 321.5 321.0 320.9 551.1 531.2 508.5 1802.8 1712.1 1608.8 2448.2 2998.8

1.8 0.4 0.1 0.1 0.0 − 1.4 − 1.3 − 1.0 − 1.0 − 0.2 − 0.2 − 0.1 − 0.1 0.0 0.0 − 3.0 − 0.9 1.5 − 3.5 − 3.5 − 3.6 0.9 0.8 0.8 2.3 − 0.2 − 0.2 − 0.2 − 0.5 − 0.5 − 0.5 − 0.1 − 0.1 − 0.1 0.2 − 0.2 − 0.4 0.9

11.0 57.7 391.7 936.1 3316.5 40.9 52.6 162.2 184.4 784.2 885.0 1697.8 1967.5 5235.2 6004.8 17.8 59.7 64.1 228.3 222.7 216.6 538.9 548.8 548.8 15.4 321.4 320.9 320.9 551.5 531.7 508.9 1802.2 1711.6 1608.2 2459.0 3001.9

1.0 0.2 0.0 0.0 0.2 − 0.5 − 0.5 − 0.3 − 0.3 0.0 0.0 − 0.1 − 0.1 0.2 0.2 − 1.7 − 1.0 2.3 − 3.2 − 3.2 − 3.3 1.0 0.9 0.9 2.7 − 0.3 − 0.3 − 0.3 − 0.4 − 0.4 − 0.4 − 0.2 − 0.2 − 0.2 0.7 − 0.1 − 0.2 0.8

J. Chem. Phys. 139, 204106 (2013)

[Os(H2 O)6 ]3+

SOfull


CH2 SiH2 GeH2 SnH2 PbH2 H2 O2

States

204106-9

TABLE III. Magnitude of SOC matrix elements (cm−1 ) calculated at DKH1 level with exact (SOfull ) and various approximate treatments of two-electron SOC (see text), using CASSCF wave functions and the ANO-RCC basis set, and the percentage error (PE) of the approximative approaches. Mean (MPE) and mean absolute (MAPE) percentage errors are also reported.

204106-10


5%. Such behavior can be found in molecules composed of light atoms (and with delocalized orbitals), for example, in CH2 , CO+ , and NO. Table III shows that the AMFI approach seems to be slightly more accurate than the exact one-center approximation in these cases (the size of the maximal error is approximately 6%). This is true especially for CH2 , where almost the exact result is observed with AMFI, while the exact one-center approximation underestimates the S0 /T1 SOC matrix element by almost 17%. The screening parameter values used in the SNSO approximation are systematically lower than those obtained by the mean-field SO atomic calculations (compare SNSO and FNSSO values in Table II). Consequently, the SNSO approximation tends to overestimate the SOC matrix elements due to insufficient screening of the one-electron terms. Although the error of the one-electron approximation is significantly reduced, the overestimation for molecules composed of light atoms can be as large as 15%–30%. However, this is not so surprising, because this approximation was originally meant only for SOC in heavy-element-containing systems. On the other hand, Table III shows that the error can be 5%–10%, even for the heaviest elements (9% on average for the set). It is noteworthy that we can generally expect only slightly better performance of methods based on the effective potential SO operator if we do not apply any corrections to the exchange contribution, or if only the coulomb term is used for the evaluation of two-electron SOC effects. Despite the somewhat lower accuracy of the original formulation, the concept of the SNSO approach is attractive due to its simplicity, both theoretically and practically (concerning its implementation and efficiency), and the possibility of simple systematic improvements. This possibility can be demonstrated by using the improved values of SNSO screening parameters suggested earlier in this work (ISNSO column). The mean absolute error is decreased to less than 2%, and errors of 1%–2% are observed for most of the SOC matrix elements reported here (maximal error size is only approximately 4%), which is already comparable to some of the more sophisticated methods. Moreover, our choice of the “improved” SNSO screening parameter values is based on the FNSSO screening parameters for lighter atoms of particular blocks in the periodic table. This choice is quite reasonable, because those are the elements with the largest (relatively) two-electron SOC contributions. However, rather even larger screening parameter values may be chosen when higher accuracy is desirable for heavier elements (see the supplementary material49 for more details on screening parameters). With regard to the improvement of accuracy achieved solely by a simple modification of the SNSO screening parameters, it is not so surprising that the accuracy of the meanfield approach can be almost attained with a more flexible model of screening. In the SNSO approximation, only a few different screening parameters are used, each of which correspond to a particular angular momentum. On the other hand, each shell (more precisely GTO), or combination of shells (GTOs) in the case of two-dimensional screening, of each angular momentum of each element has its own screening parameter in the present FNSSO approach. This flexibility of screening allows the description of any SOC matrix element

J. Chem. Phys. 139, 204106 (2013)

with relatively high accuracy. The maximal error size of the FNSSO method presented here is approximately 3%–4%, and the mean absolute error is below 1% (see Table III). Most of the matrix elements are predicted to be within 0.5% deviation from their exact values, especially for molecules containing heavier atoms. Another positive indication is that FNSSO seems to conserve the trends of the mean-field approach from which it originates. Slightly higher than average errors of 1%–3% and 1%– 2% for FNSSOC1D and FNSSOP2D , respectively, are observed for molecules with significant multi-center contribution to SOC, e.g., for CH2 , CO+ , and NO. That is quite logical because we use a fairly simple model for the screening of multi-center one-electron SOC integrals, so that their lower accuracy can be generally expected. However, the error of the one-center approximation seems to be significantly reduced, even with this simple model, and we thus consider it quite satisfactory. As a demonstration, the multi-center one-electron SOC integrals can be left unscreened. In that case, the error exceeds 8% for S0 /T1 SOC in CH2 and for D2 /D3 SOC in CO+ , which have significant contributions from the twoelectron multi-center integrals. The complexes of the first row transition metals are special cases, for which we observe maximal errors of the FNSSO approach (approximately between −3% and +3%). However, these deviations have a principal rather than numerical origin and reflect the fact that we use screening parameters determined for isolated atoms, even when these atoms were placed into a chemical environment without spherical symmetry (the same is also true for the AMFI approach due to its predetermined atomic-orbital occupation numbers – compare the AMFI and FNSSO results for these molecules in Table III). In the environment with lower symmetry, the screening parameter value that corresponds to the valence d shell splits (as does its energy), because its occupation changes. Thus, there should, in principal, be additional screening parameters that correspond to the d subshells of transition metal atoms (ions) in various chemical environments (for various types of d splitting), which is somewhat unfortunate. Nevertheless, according to the present findings, a simple model can be used to correct the d screening parameter values based on the ligand field stabilization energy (LFSE), which describes the extent of the d occupation non-uniformity (see standard textbooks of inorganic chemistry, for example, Ref. 62). In the case of the simplest d splitting that corresponds to octahedral and tetrahedral ligand fields, the original one-dimensional (similarly for two-dimensional) screening parameter of the d shell can be split into screening parameters of t2g and eg subshells according to LFSE t2g d Qii = Qii 1 ∓ 2 , 100 LFSE e Qiig = Qdii 1 ± 3 , 100 where the upper and lower signs represent the octahedral and tetrahedral fields, respectively, and the LFSE in relative scale (in the units of d splitting) corresponds to certain d configurations and ligand-field types. If we consider, for example, the

204106-11


J. Chem. Phys. 139, 204106 (2013)

TABLE IV. Magnitude of SOC matrix elements (cm−1 ) calculated at the DKH1 level using the FNSSO approach with screening parameters corrected for splitting of d orbitals, CASSCF wave functions and the ANO-RCC basis set for first row transition metal complexes. Percentage (PE), mean percentage (MPE), and mean absolute percentage (MAPE) errors are also reported.

[Fe(H2 O)6 ]3+

[Fe(H2 O)6 ]2+

[Co(H2 O)6 ]2+ [NiCl4 ]2−

MPEa MAPEa

States

FNSSOC1D

PE

FNSSOP2D

PE

D1 /D2 D1 /D3 D2 /D3 S0 /T1 S0 /T2 S0 /T3 D1 /D2 T1 /T2 T1 /T3 T2 /T3

234.5 229.1 223.0 534.1 544.2 544.4 15.2 317.7 317.2 317.2

− 0.5 − 0.4 − 0.4 0.1 0.1 0.1 1.6 − 1.4 − 1.4 − 1.4 − 0.3 0.6

235.2 229.8 223.7 534.6 544.9 545.1 15.3 317.6 317.1 317.1

− 0.2 − 0.1 − 0.2 0.2 0.2 0.2 2.0 − 1.4 − 1.4 − 1.4 − 0.1 0.5

a

Calculated using corrected values presented here together with the values for the rest of the molecules from Table III.

d5 configuration of a transition metal ion in a strong octahedral ligand field, for which LFSE is equal to 2.0, then screening of the t2g and eg subshells is reduced by 4% and enhanced by 6%, respectively. This describes the electronic situation in the [Fe(H2 O)6 ]3+ complex; therefore, it can be expected that SOC between the three near degenerate (strictly degenerate in ideal symmetry) lowest doublets, which differ in occupation of the orbitals of the t2g subshell, should increase by approximately 4%. Table IV, which summarizes the SOC results for first row transition metal complexes calculated using FNSSO with screening parameters corrected for splitting of the transition metal d orbitals, shows that the error for [Fe(H2 O)6 ]3+ is reduced by 3% in the actual calculation and becomes 0.5% or less. Similar improvement is also achieved in the case of the [Fe(H2 O)6 ]2+ complex with d6 configuration, where SOC between the states examined is predominantly from configurations differing in t2g → eg single-excitations and a small enhancement of the screening (reduction of SOC) can be expected. Somewhat worse results are observed for the Co2+ complex, which is probably due to the higher sensitivity of its weak D1 /D2 SOC matrix element, because this matrix element is in the ideal Th geometry symmetrically forbidden. Finally, almost perfectly corrected results are obtained for the tetrahedral [NiCl4 ]2− complex, if we realize that the target of the FNSSO approximation is by its nature the meanfield SO operator, rather than the exact SO operator. Table IV indicates that the mean absolute error for the set is decreased to 0.5% for FNSSOP2D (0.3% with respect to SOMF ), which is already fairly close to the error of the mean-field approach. Although it has been shown that the error arising from the splitting of transition metal d orbitals can in principle be corrected, it is not expected that this correction would have to be generally applied. We consider the accuracy of the uncorrected FNSSO approach, which leads to errors of up to 3%–4% for the first row transition metal complexes (less than 1% for heavier transition metals), to be quite satisfac-

tory. Moreover, the error observed for metals with d5 configuration placed in a strong octahedral field (for example, the [Fe(H2 O)6 ]3+ complex) should be the maximal expected error, because further splitting of the d subshells does not significantly change the discrepancy between their occupation. Last but not least, this type of error is systematic and dependent on the d configuration and type of ligand (crystal) field, and is thus predictable.

B. Calculations that employ various basis sets

We have used only the ANO-RCC basis set for the tabulation of all screening parameters used in this work; therefore, it is important to also test the accuracy of this approach for other basis sets. For the purpose of this testing, four molecules composed of lighter atoms were selected, for which twoelectron SOC contributions play a significant role, namely, CH2 , H2 O2 , NO, and [Fe(H2 O)6 ]3+ . Table V shows that the use of very small basis sets, such as STO-3G and 3-21G, can lead to errors in SOC as large as 40% (results obtained with the largest ANO-RCC used as a reference). From the SOC matrix elements studied, it seems that smaller basis sets tend to systematically underestimate SOC. Moreover, this underestimation is further enhanced when the AMFI approach is used, and the total error with respect to a reasonably accurate prediction may exceed 50%. However, this is not the case for only the AMFI approach. For all approximative SO operators used here, the strongest inconsistency of results is obtained for the smallest basis sets. This is due to the use of atomic quantities (occupation numbers or screening parameters) obtained from much more accurate calculations, which are, in case of very small basis sets, inadequate due to different descriptions of electron density. Moreover, we do not obtain very stable relative errors for the ISNSO and FNSSOC1D operators, even with larger basis sets. This is not surprising, because ISNSO uses a less flexible scheme of screening, and the accuracy of FNSSOC1D can be conserved only when the structure of the basis set in use is highly similar to that of the contracted ANO-RCC basis set. A much better situation is realized for the other two approximations. The relative error with AMFI seems to be stable from the cc-pVDZ basis set, and that of FNSSOP2D is already stable from the 3-21G basis. Although the STO-3G basis seems to be the only problematic one for FNSSOP2D , the relative error of the SOC predictions is still not significantly high, but varies from approximately 0.5% to 8%, and as shown in Table V, is approximately one order of magnitude smaller than the error of the basis set in all four cases. For other basis sets, the relative errors for particular molecules span the range within 1%, and the accuracy of these predictions is thus well consistent with that for ANO-RCC. That leads us to the conclusion that the FNSSOP2D approach may be considered to be basis-setindependent. Furthermore, FNSSOP2D seems to provide the least biased results from the approximative operators used here, which indicates that the screening parameters of primitive GTOs are, in terms of their transferability, well chosen quantities (perhaps preferable to the occupation numbers of atomic orbitals). Inspection of the results provided in

204106-12


J. Chem. Phys. 139, 204106 (2013)

TABLE V. Magnitude of SOC matrix elements (cm−1 ) calculated at the DKH1 level with exact (SOfull ) and approximate treatments of two-electron SOC (see text), using CASSCF wave functions and various basis sets. Matrix elements between the S0 and T1 electronic states are reported for CH2 and H2 O2 , and between the D1 and D2 for NO and [Fe(H2 O)6 ]3+ . Percentage (PE), mean percentage (MPE), and mean absolute percentage (MAPE) errors of approximative operators are also reported.

CH2

H2 O2

NO

[Fe(H2 O)6 ]3+

Basis set

SOfull

Baserra

AMFI

PE

ISNSO

PE

FNSSOC1D

PE

FNSSOP2D

PE

STO-3G 3-21G 6-31G cc-pVDZ cc-pVTZ cc-pVQZ ANO-RCC DZ ANO-RCC DZP ANO-RCC TZP STO-3G 3-21G 6-31G cc-pVDZ cc-pVTZ cc-pVQZ ANO-RCC DZ ANO-RCC DZP ANO-RCC TZP STO-3G 3-21G 6-31G cc-pVDZ cc-pVTZ cc-pVQZ ANO-RCC DZ ANO-RCC DZP ANO-RCC TZP STO-3G 3-21G 6-31G cc-pVDZ cc-pVTZ ANO-RCC DZ ANO-RCC DZP ANO-RCC TZP/DZP

10.3 9.0 10.0 10.5 10.7 10.9 11.1 11.1 10.9 30.7 34.1 36.9 38.9 40.0 40.4 41.0 41.8 41.1 41.3 51.0 57.7 57.9 60.4 61.6 64.8 63.5 62.6 132.4 203.7 223.6 239.8 237.9 239.6 238.3 235.7

− 5.9 − 17.9 − 8.5 − 3.7 − 2.1 − 0.7 1.8 1.9 0.0 − 25.3 − 17.0 − 10.4 − 5.5 − 2.7 − 1.7 − 0.2 1.7 0.0 − 34.1 − 18.6 − 7.9 − 7.5 − 3.6 − 1.7 3.4 1.3 0.0 − 43.8 − 13.6 − 5.1 1.7 0.9 1.7 1.1 0.0

8.9 7.1 8.5 10.5 10.7 10.8 11.1 11.1 11.0 29.3 29.6 34.0 39.6 40.7 41.3 41.9 42.6 42.0 40.9 45.7 55.2 61.9 64.5 65.6 68.7 67.4 66.4 113.4 174.0 201.0 233.1 231.2 232.6 231.1 228.5

− 13.1 − 20.6 − 14.7 − 0.3 − 0.2 − 1.1 − 0.1 − 0.3 0.3 − 4.7 − 13.3 − 7.8 2.0 1.8 2.1 2.1 1.8 2.1 − 0.9 − 10.3 − 4.3 6.9 6.7 6.6 6.1 6.2 6.0 − 14.3 − 14.6 − 10.1 − 2.8 − 2.8 − 2.9 − 3.0 − 3.1 − 1.9 5.2

10.9 9.6 10.4 10.9 10.9 10.9 11.2 11.2 11.0 31.8 35.1 37.1 39.0 39.7 40.0 40.3 41.1 40.6 44.5 54.3 60.0 60.1 61.9 62.7 65.4 64.2 63.4 143.8 204.2 218.6 228.8 227.3 229.4 227.8 225.7

5.5 7.2 4.3 3.3 1.8 0.6 0.1 0.2 0.8 3.5 2.7 0.7 0.3 − 0.8 − 1.1 − 1.7 − 1.7 − 1.3 7.7 6.5 4.0 3.7 2.4 1.8 1.0 1.2 1.2 8.6 0.2 − 2.2 − 4.6 − 4.4 − 4.3 − 4.4 − 4.2 1.1 2.8

11.0 9.9 10.6 11.1 11.3 11.6 11.3 11.3 11.1 31.8 35.3 37.2 39.0 40.1 40.9 40.3 41.1 40.5 44.7 54.9 60.5 60.6 63.0 64.7 65.6 64.4 63.6 144.3 206.9 220.3 230.8 229.3 231.3 229.7 227.6

6.7 9.6 6.3 5.3 5.5 6.7 1.0 1.3 1.8 3.4 3.3 0.9 0.5 0.2 1.1 − 1.7 − 1.7 − 1.4 8.1 7.6 4.8 4.6 4.2 5.1 1.2 1.4 1.5 9.0 1.6 − 1.5 − 3.8 − 3.6 − 3.5 − 3.6 − 3.5 2.3 3.6

10.3 9.1 10.1 10.7 10.8 11.0 11.2 11.2 11.0 29.5 33.9 36.5 38.6 39.8 40.2 40.7 41.5 40.9 41.5 52.4 59.1 59.5 62.0 63.2 66.2 64.9 64.1 121.7 196.5 217.2 232.3 230.3 232.3 230.8 228.3

− 0.5 0.8 0.9 1.3 1.0 1.2 0.3 0.5 1.0 − 3.9 − 0.7 − 1.0 − 0.7 − 0.5 − 0.5 − 0.7 − 0.7 − 0.5 0.4 2.8 2.4 2.7 2.7 2.6 2.2 2.3 2.3 − 8.1 − 3.6 − 2.8 − 3.1 − 3.2 − 3.0 − 3.1 − 3.2 − 0.2 1.9

MPE MAPE a

Error (%) of the basis set (with respect to ANO-RCC with the most accurate contraction used here).

Tables III–V indicates why FNSSOP2D with its very promising accuracy and general applicability is our method of choice. Here, we briefly discuss the two FNSSO models used throughout this work, because this is closely related to the structure of the GTO basis sets. If we consider the screening parameters for the contracted GTOs of ANO-RCC (see supplementary material49 ), we note that they usually span a relatively narrow range of values, especially up to the “valence” shell. This means that the ANO-RCC basis with QZP contraction does not contain GTOs that are too compact or diffuse, at least from the point of SOC screening. In such a case, a one-dimensional scheme that uses geometric mean of diagonal screening parameters according to Eq. (11) seems to sufficiently describe the screening (only small differences from the

two-dimensional screening are observed for our benchmark set). On the other hand, exponents of the primitive GTOs span a wide range of values, as with their screening parameters (see also Fig. 1). In this case, the geometric-mean approximation is not so accurate, although it has qualitatively correct behavior in most parts of the dependencies of two-dimensional screening parameters on primitive GTO exponents, and onedimensional screening scheme leads to slightly higher errors (deviations over 4% from the two-dimensional version can be observed). Thus, when a basis set with fairly diffuse or compact GTOs, either primitive or contracted, is employed, then a two-dimensional scheme of screening is preferable. Moreover, the set of screening parameters in use should be adapted to these diffuse or compact functions to avoid strong deviations of their SOC contributions, most probably

204106-13


underestimation for compact and overestimation for diffuse GTOs. C. General remarks on the accuracy of SOC predictions

We have shown that the magnitude of SOC matrix elements is quite sensitive to the quality of the basis set used, and errors as large as 40% can be obtained with very small basis sets. Although this is quite an extreme case, we should expect that our SOC predictions can easily suffer from a deviation of 1%–5% when standard basis sets of double-ζ or triple-ζ , or sometimes even higher quality are used. We can also suppose that the quality of the wave function matters in the same way and this can also lead to percentage scale errors. For example, when state-specific CASSCF wave functions are used for both S0 and T1 states of CH2 , the corresponding SOC matrix element is (with the same basis set as in Table III) approximately 12.4 cm−1 , i.e., larger by almost 13%. This is also a reason why the SOC matrix element values presented in this work, even the exact ones, should not be taken as highly accurate predictions of their “real” values. The key role plays the comparison of performance of exact and approximate (DKH1) SO operators here. Moreover, it is not uncommon to employ a Breit–Pauli SO operator for calculations on molecules containing atoms up to the fourth period. Although this is relatively safe for first row transition metals, for which it usually leads to errors below 1%, deviations of more than 7% can be observed for the corresponding main group elements, and still approximately 1% for lighter atoms of the second row of the p block, such as silicon and sulfur. Similarly, the DKH1 SO operator may lead to non-negligible errors for the heaviest atoms. If all of these limitations are considered, then in a typical calculation (of which the primary goal is not providing SOC values as accurate as possible), it is probably not worth performing the expensive evaluation of all two-electron SOC integrals, but rather a highly efficient and reasonably accurate approximation. As demonstrated in this work, one such alternative is FNSSO. VI. CONCLUSIONS

A new approximation to the two-electron SOC terms has been presented, which is referred to as FNSSO. This approximation, inspired by the SNSO approach and formed on the basis of the mean-field SO operator, uses an effective one-electron SO operator and accounts for the effect of two-electron SOC by screening the nuclear charge involved in the corresponding one-electron SOC integrals. This screening is based on a highly flexible scheme, where basically any pair of atomic orbitals has its own screening parameter. While screening parameters that correspond to one-center SO interactions have been parameterized by ab initio atomic calculations, fairly simple assumptions have been applied to the screening of multi-center interactions. The high flexibility of the screening used in FNSSO allows mostly spectroscopic accuracy (error within 1 cm−1 ) to be achieved for SOC matrix elements in molecules composed of light atoms, and results were obtained with deviations of typically a few wavenum-

J. Chem. Phys. 139, 204106 (2013)

bers for heavy-atom containing systems with respect to the exact treatment of two-electron SOC at the same level of relativistic theory (DKH1 SO operators were used throughout this work). It has been shown that the version based on two-dimensional screening parameters for primitive GTOs (FNSSOP2D ) is independent of basis set selection, and thus, only one simple database of screening parameters is used for any calculation. Moreover, we have introduced a model to correct the effect of transition metal d orbital splitting on the SOC matrix elements, which arises from the atomic nature of the screening parameters. In addition to high accuracy, the FNSSO approach is highly efficient, because only one-electron SOC integrals have to be explicitly evaluated. Last but not least, this approach is easily implemented, because only fairly simple modifications of existing codes for one-electron SOC calculations are required. Both sets of screening parameters that correspond to contracted and primitive Gaussian-type basis functions of the ANO-RCC basis set are available for all elements up to curium via the supplementary material.49 ACKNOWLEDGMENTS

The authors would like to acknowledge our colleague Mojmír Kývala for many inspiring discussions regarding SOC effects, and also for allowing us to use his program SOSS. One of us (T.Y.) was supported by a Grant-in-Aid from the Core Research for Scientific Research (B) (Grant No. 25288013) of the Ministry of Education, Culture, Sports, Science and Technology of Japan (MEXT). 1 B.

A. Heß and C. M. Marian, in Computational Molecular Spectroscopy, edited by P. Jensen, and P. R. Bunker (John Wiley & Sons, New York, 2000), pp. 169–219. 2 K. G. Dyall and K. Fægri, Introduction to Relativistic Quantum Chemistry (Oxford University Press, New York, 2007). 3 M. Reiher and A. Wolf, Relativistic Quantum Chemistry (Wiley, Weinheim, 2009). 4 L. Visscher, O. Visser, P. Aerts, H. Merenga, and W. Nieuwpoort, Comput. Phys. Commun. 81, 120 (1994). 5 T. Saue, K. Fægri, T. Helgaker, and O. Gropen, Mol. Phys. 91, 937 (1997). 6 H. Quiney, H. Skaane, and I. Grant, Adv. Quantum Chem. 32, 1 (1998). 7 W. Liu, G. Hong, D. Dai, L. Li, and M. Dolg, Theor. Chem. Acc. 96, 75 (1997). 8 T. Yanai, H. Iikura, T. Nakajima, Y. Ishikawa, and K. Hirao, J. Chem. Phys. 115, 8267 (2001). 9 T. Nakajima, T. Yanai, and K. Hirao, J. Comput. Chem. 23, 847 (2002). 10 M. S. Kelley and T. Shiozaki, J. Chem. Phys. 138, 204113 (2013). 11 T. Saue, ChemPhysChem 12, 3077 (2011). 12 D. Peng and M. Reiher, Theor. Chem. Acc. 131, 1081 (2012). 13 W. Liu, Mol. Phys. 108, 1679 (2010). 14 J. Seino and M. Hada, Chem. Phys. Lett. 461, 327 (2008). 15 W. Kutzelnigg and W. Liu, J. Chem. Phys. 123, 241102 (2005). 16 M. Barysz and A. J. Sadlej, J. Chem. Phys. 116, 2696 (2002). 17 J. Sikkema, L. Visscher, T. Saue, and M. Iliaš, J. Chem. Phys. 131, 124116 (2009). 18 J. Seino and H. Nakai, J. Chem. Phys. 137, 144101 (2012). 19 L. Storchi, L. Belpassi, F. Tarantelli, A. Sgamellotti, and H. M. Quiney, J. Chem. Theory Comput. 6, 384 (2010). 20 D. Peng, N. Middendorf, F. Weigend, and M. Reiher, J. Chem. Phys. 138, 184105 (2013). 21 S. J. Havriliak and D. R. Yarkony, J. Chem. Phys. 83, 1168 (1985). 22 B. A. Heß, C. M. Marian, and S. D. Peyerimhoff, in Modern Electronic Structure Theory, edited by D. Yarkony (World Scientific, Singapore, 1995).

204106-14 23 P.-Å.


Malmqvist, B. O. Roos, and B. Schimmelpfennig, Chem. Phys. Lett. 357, 230 (2002). 24 B. A. Hess, R. J. Buenker, C. M. Marian, and S. D. Peyerimhoff, Chem. Phys. Lett. 89, 459 (1982). 25 C. Teichteil, M. Pelissier, and F. Spiegelmann, Chem. Phys. 81, 273 (1983). 26 R. M. Pitzer and N. W. Winter, J. Phys. Chem. 92, 3061 (1988). 27 M. Blume and R. E. Watson, Proc. R. Soc. London, Ser. A 270, 127 (1962); 271, 565 (1963). 28 M. Dolg, H. Stoll, H.-J. Flad, and H. Preuss, J. Chem. Phys. 97, 1162 (1992). 29 B. A. Heß, C. M. Marian, U. Wahlgren, and O. Gropen, Chem. Phys. Lett. 251, 365 (1996). 30 B. Schimmelpfennig, AMFI: An Atomic Mean-Field Integral Program (University of Stockholm, Stockholm, Sweden, 1996). 31 S. Koseki, M. W. Schmidt, and M. S. Gordon, J. Phys. Chem. 96, 10768 (1992). 32 Y. S. Lee, W. C. Ermler, and K. S. Pitzer, J. Chem. Phys. 67, 5861 (1977). 33 G. Schreckenbach and T. Ziegler, J. Phys. Chem. A 101, 3388 (1997). 34 J. C. Boettger, Phys. Rev. B 62, 7809 (2000). 35 W. Pauli, Z. Phys. 43, 601 (1927); G. Breit, Phys. Rev. 34, 553 (1929). 36 M. Douglas and N. M. Kroll, Ann. Phys. (N.Y.) 82, 89 (1974). 37 B. A. Hess, Phys. Rev. A 32, 756 (1985); J. Almlöf, K. Fægri, Jr., and H. H. Grelland, Chem. Phys. Lett. 114, 53 (1985). 38 R. Samzow, B. A. Heß, and G. Jansen, J. Chem. Phys. 96, 1227 (1992); U. Wahlgren, M. Sjøvoll, H. Fagerli, O. Gropen, and B. Schimmelpfennig, Theor. Chem. Acc. 97, 324 (1997). 39 D. R. Hartree, Proc. Cambridge Philos. Soc. 24, 89 (1928); V. Fock, Z. Phys. 61, 126 (1930); 62, 795 (1930). 40 A similar expression can also be derived for a pair of equivalent Slater determinants; however, their SOC contribution vanishes if only one set of real orthonormal orbitals common to bra and ket wave functions is used, as is conducted in many implementations. If real orbital basis is used, the second two-electron term in Eqs. (4) and (6) vanishes as well. 41 F. Neese, J. Chem. Phys. 122, 034107 (2005). 42 O. L. Malkina, J. Vaara, B. Schimmelpfennig, M. Munzarová, V. Malkin, and M. Kaupp, J. Am. Chem. Soc. 122, 9206 (2000). 43 C. van Wüllen and C. Michauk, J. Chem. Phys. 123, 204113 (2005). 44 P.-O. Widmark, P.-Å. Malmqvist, and B. O. Roos, Theor. Chim. Acta 77, 291 (1990); B. O. Roos, V. Veryazov, and P.-O. Widmark, Theor. Chem. Acc. 111, 345 (2003); B. O. Roos, R. Lindh, P.-Å. Malmqvist, V. Veryazov, and P.-O. Widmark, J. Phys. Chem. A 108, 2851 (2004); 109, 6575 (2005); Chem. Phys. Lett. 409, 295 (2005); J. Phys. Chem. A 112, 11431 (2008). 45 B. O. Roos, P. R. Taylor, and P. E. M. Siegbahn, Chem. Phys. 48, 157 (1980); B. O. Roos, Int. J. Quantum Chem., Symp. 17, 175 (1980); P. E. M. Siegbahn, A. Heiberg, B. O. Roos, and B. Levy, Phys. Scr. 21, 323 (1980); P. E. M. Siegbahn, J. Almlöf, A. Heiberg, and B. O. Roos, J. Chem. Phys. 74, 2384 (1981).

J. Chem. Phys. 139, 204106 (2013) 46 B.

A. Hess, Phys. Rev. A 33, 3742 (1986); G. Jansen and B. A. Hess, ibid. 39, 6016 (1989). 47 G. Karlström, R. Lindh, P.-Å. Malmqvist, B. O. Roos, U. Ryde, V. Veryazov, P.-O. Widmark, M. Cossi, B. Schimmelpfennig, P. Neogrády, and L. Seijo, Comput. Mater. Sci. 28, 222 (2003). 48 M. Kývala, SOSS: A Program for ab initio Quasirelativistic Calculations (Institute of Organic Chemistry and Biochemistry, v.v.i., Academy of Sciences of the Czech Republic, Prague, Czech Republic, 2002–2012. 49 See supplementary material at http://dx.doi.org/10.1063/1.4832737 for values of tabulated screening parameters and molecular geometries. 50 A. D. Becke, Phys. Rev. A 38, 3098 (1988); C. T. Lee, W. T. Yang, and R. G. Parr, Phys. Rev. B 37, 785 (1988); A. D. Becke, J. Chem. Phys. 98, 5648 (1993); P. J. Stephens, F. J. Devlin, C. F. Chabalowski, and M. J. Frisch, J. Phys. Chem. 98, 11623 (1994). 51 A. Schäfer, H. Horn, and R. Ahlrichs, J. Chem. Phys. 97, 2571 (1992); F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys. 7, 3297 (2005). 52 K. Eichkorn, F. Weigend, O. Treutler, and R. Ahlrichs, Theor. Chem. Acc. 97, 119 (1997). 53 M. Srnec, J. Chalupský, M. Fojta, L. Zendlová, L. Havran, M. Hocek, M. Kývala, and L. Rulíšek, J. Am. Chem. Soc. 130, 10947 (2008). 54 J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). 55 K. Eichkorn, O. Treutler, H. Öhm, M. Häser, and R. Ahlrichs, Chem. Phys. Lett. 240, 283 (1995). 56 TURBOMOLE V6.4 2012, a development of University of Karlsruhe and Forschungszentrum Karlsruhe GmbH, 1989–2007, TURBOMOLE GmbH, since 2007; available from http://www.turbomole.com. 57 This was necessary due to slightly different strategies for the evaluation of SOC matrix elements in the programs used. While Molcas and SOSS can calculate transition matrix elements between states with mutually nonorthogonal molecular orbitals, only one set of molecular orbitals (averaged, for example, even over states with different multiplicity) has to be used in ORZ to make evaluation of the transition matrix elements also possible with CASSCF wave functions based on density matrix renormalization group theory (DMRG). 58 W. J. Hehre, R. F. Stewart, and J. A. Pople, J. Chem. Phys. 51, 2657 (1969); W. J. Pietro and W. J. Hehre, J. Comput. Chem. 4, 241 (1983). 59 J. S. Binkley, J. A. Pople, and W. J. Hehre, J. Am. Chem. Soc. 102, 939 (1980); K. D. Dobbs and W. J. Hehre, J. Comput. Chem. 8, 861 (1987). 60 W. J. Hehre, R. Ditchfield, and J. A. Pople, J. Chem. Phys. 56, 2257 (1972); V. A. Rassolov, J. A. Pople, M. A. Ratner, and T. L. Windus, ibid. 109, 1223 (1998). 61 T. H. Dunning, Jr., J. Chem. Phys. 90, 1007 (1989); N. B. Balabanov and K. A. Peterson, ibid. 123, 064107 (2005). 62 F. A. Cotton, G. Wilkinson, and P. L. Gaus, Basic Inorganic Chemistry, 2nd ed. (John Wiley & Sons, New York, 1987).

The Journal of Chemical Physics is copyrighted by the American Institute of Physics (AIP). Redistribution of journal material is subject to the AIP online journal license and/or AIP copyright. For more information, see http://ojps.aip.org/jcpo/jcpcr/jsp

Flexible Force Field Parameterization through Fitting on the Ab Initio-Derived Elastic Tensor.

AUTOMATED FORCE FIELD PARAMETERIZATION FOR NON-POLARIZABLE AND POLARIZABLE ATOMIC MODELS BASED ON AB INITIO TARGET DATA.

Ab Initio Screening of CO2-philic Groups.

Ab initio water pair potential with flexible monomers.

Understanding Anharmonicity in fcc Materials: From its Origin to ab initio Strategies beyond the Quasiharmonic Approximation.

An ab initio study of magneto-electric coupling of YMnO3.

High-throughput ab initio screening for two-dimensional electride materials.

Ab initio RNA folding.

Ab initio study on the low-lying excited states of gas-phase PH(+) cation including spin-orbit coupling.

Vibronic bandshape of the absorption spectra of dibenzoylmethanatoboron difluoride derivatives: analysis based on ab initio calculations.

An ab initio reinvestigation of PC and PH coupling constants in phosphates.

Implantable flexible pressure measurement system based on inductive coupling.

Ab initio density matrix renormalization group study of magnetic coupling in dinuclear iron and chromium complexes.

Ab Initio Screening Approach for the Discovery of Lignin Polymer Breaking Pathways.

PaFlexPepDock: parallel ab-initio docking of peptides onto their receptors with full flexibility based on Rosetta.

Benchmarking ab initio binding energies of hydrogen-bonded molecular clusters based on FTIR spectroscopy.

The hydrogen abstraction reaction O(3P) + CH4: a new analytical potential energy surface based on fit to ab initio calculations.

IPred - integrating ab initio and evidence based gene predictions to improve prediction accuracy.

Sphinx: merging knowledge-based and ab initio approaches to improve protein loop prediction.

Revealing Missing Human Protein Isoforms Based on Ab Initio Prediction, RNA-seq and Proteomics.

AIDA: ab initio domain assembly server.

Ab-initio SCF investigation ofγ-aminobutyric acid.

Ab initio and all-atom modeling of detergent organization around Aquaporin-0 based on SAXS data.

Langevin spin dynamics based on ab initio calculations: numerical schemes and applications.