${\cal{P,T}}$-Violating and Magnetic Hyperfine Interactions in Atomic Thallium

We present state-of-the-art configuration interaction and coupled cluster calculations of the electron electric dipole moment, the nucleon-electron scalar-pseudoscalar, and the magnetic hyperfine interaction constants ($\alpha_{d_e}, \alpha_{C_S}, A_{||}$, respectively) for the thallium atomic ground state $^2P_{1/2}$. Our present best values are $\alpha_{d_e} = -559 \pm 28$, $\alpha_{C_S} = 6.77 \pm 0.34$ $[10^{-18} e$ cm$]$, and $A_{||} = 21172 \pm 1059$ [MHz]. These findings lead to a significant reduction of the theoretical uncertainties for ${\cal{P,T}}$-odd interaction constants but not to stronger constraints on the electron electric dipole moment, $d_e$, or the nucleon-electron scalar-pseudoscalar coupling constant, $C_S$.

In most BSM models [11] the dominant CP-odd sources in open-shell systems are the electron EDM, d e , and the nucleon-electron scalar-pseudoscalar (Ne-SPS) coupling, C S . The panel [55] shows the constraints (yellow surface) on d e and C S using the combined information from measurements [3,[12][13][14] and calculations [4,5,[15][16][17][18][19][20][21][22][23], including the associated experimental and theoretical uncertainties, on the open-shell systems ThO (green), YbF (red), HfF + (orange) and Tl (blue) through a global fit in the d e / C S plane. Results from a single system do, therefore, not constrain d e or C S individually at all in this multiple-source interpretation [24], but lead to a fanshaped surface of allowed combinations. The width of this surface is a function of the experimental and theoretical uncertainties.
This means that a substantial reduction of an uncertainty for an individual system could lead to more stringent constraints on the unknown CP-violating parameters. The main reason for this is that the surfaces for different systems are not fully aligned, which is due to the different dependency of electron EDM and Ne-SPS atomic interactions on the electric charge of the respective heavy nuclei [25,26].
A substantial part of the width of the surface for the Tl atom is due to the great spread of theoretical values for the electron EDM atomic enhancement, R, calculated in the past by various groups using different electronic-structure approaches [1,[19][20][21][22][23]27]. Strikingly, Nataraj et al. [22] used a high-level many-body approach, the Coupled Cluster (CC) method, and produced a value for R that strongly disagrees with the results from all other groups, on the order of 20%.
The purpose of this paper is twofold: 1. We use state-of-the-art Configuration Interaction (CI) and Coupled-Cluster approaches for largescale applications to determine the mentioned atomic interaction constants. We put particular emphasis on the electron EDM enhancement R and a conclusive resolution of the major discrepancy between literature values. Claims about physical effects that purportedly underlie these discrepancies are scutinized.
2. We investigate whether a reduced uncertainty for R(Tl) impacts the above-described constraints on d e and C S .
The paper is structured as follows. In section II we lay out the theory underlying the atomic electron EDM, Ne-SPS, and magnetic hyperfine interactions constants. The following section III contains technical details about our calculations, results, and a discussion of these results in comparison with literature values. The final section IV concludes on our findings.

II. THEORY
An atomic electric dipole moment (EDM) is defined [28] (p. 16) as where ∆ε is a P, T -odd energy shift and E ext is an external electric field. In atoms with nuclear spin I ≤ 1 2 [29] and in an electronic state with unpaired electrons, this energy shift is dominated by and originates from either the electron EDM, d e , or a P, T -odd nucleon-electron (Ne) interaction, or a combination of the two [10,11]. The two cases are presented separately.

A. Atomic EDM due to electron EDM
The Hamiltonian for the interaction of the electron electric dipole moment, d e , is for an atomic system where γ 0 is a Dirac matrix, Σ = σ 0 0 σ is a vector of spin matrices in Dirac representation, j is an electron index, E(r j ) the electric field at position r j and the bare fermion's electric dipole moment is expressed as d = d e γ 0 Σ, necessarily linearly dependent on the particle's spin vector Σ [1,30]. Supposing a non-zero electron EDM d e , the resulting energy shift can be evaluated as where ψ(E ext ) is the field-dependent atomic wavefunction of the state in question. The expectation value in Eq. (3) has the physical dimension of electric field and can be regarded as the mean interaction of each electron EDM with this field in the respective state. Following stratagem II of Lindroth et al. [31] the expectation value is recast in electronic momentum form as an effective one-body operator where the approximation lies in assuming that ψ is an exact eigenfunction of the field-dependent Hamiltonian of the system. This momentum-form EDM operator has already been used as early as in 1986, by Johnson et al. [27]. In the present work the field-dependent Hamiltonian is the Dirac-Coulomb (DC) Hamiltonian (in a.u. , e = m e = = 1) with E ext weak and homogeneous, the indices j, k run over n electrons, Z the proton number with the nucleus K placed at the origin, and α are standard Dirac matrices. E ext is not treated as a perturbation but included a priori in the variational optimization of the atomic wavefunction. Furthermore, the final results reported in this work include high excitation ranks in the correlation expansion of ψ. For these reasons, the approximation in Eq. (4) is considered very good in the present case. Within the so-defined picture and using Eqs. (1), (3), and (4) the atomic EDM becomes The (dimensionless) atomic EDM enhancement factor is defined as R : for the sake of simplicity, the enhancement factor is The external field used in the experiment on Tl [14] was E ext = 1.23 × 10 7 V m ≈ 0.2392 × 10 −4 a.u. In the present work E ext = 0.24 × 10 −4 a.u. is used. This is a very small field which is well within the linear regime considering the derivative in Eq. (7). The enhancement factor may under these circumstances be written as a function of two field points We set E ext (1) := 0 from which it follows that E eff (1) = 0, and so E eff is calculated as described in reference [32]. ψ is an approximate configuration interaction (CI) eigenfunction of the Dirac-Coulomb Hamiltonian including E ext . Alternatively, E eff can be calculated within the finite-field approach [33,34]. The latter has been used in coupled cluster calculations. The electron EDM enhancement factor R is in the particle physics literature often denoted as the atomic-scale interaction constant of the electron EDM.

B. Nucleon-Electron Scalar-Pseudoscalar Interaction
The effective Hamiltonian for a P, T -odd nucleon-electron scalar-pseudoscalar interaction is written as and the resulting atomic energy shift is accordingly where A is the nucleon number, C S is the S-PS nucleon-electron coupling constant, G F is the Fermi constant[56] and ρ N (r j ) is the nucleon density at the position of electron j. Note that in the present work we define γ 5 := ıγ 0 γ 1 γ 2 γ 3 , whereas Flambaum and co-workers [21,25,26] define γ 5 := −ıγ 0 γ 1 γ 2 γ 3 which explains the sign difference between the present Ne-SPS atomic interaction constants and those of Flambaum and co-workers. Next, we define (see also reference [36]) in analogy with Eq. (7) an Ne-SPS ratio [57] and so one can write, using Eq. (1), and in the linear regime The initial implementation of this expectation value in the latter expression has been described in reference [37]. The independent implementation of the matrix elements of the Hamiltonian (11) has been developed in ref. [4]. For comparison with literature results we also define the S-PS nucleon-electron interaction constant

C. Magnetic Hyperfine Interaction
Minimal substitution according to p −→ p − q c A in the Dirac equation and representing the vector potential in magnetic dipole approximation as A D (r) = m×r r 3 with m the nuclear magnetic dipole moment leads to the magnetic hyperfine Hamiltonian for a single point charge q at position r outside the finite nucleus. Given the nuclear magnetic dipole moment vector as m = µ I µ N I = g I µ N I where µ is the magnetic moment in nuclear magnetons (µ N ), g I is the nuclear g-factor, and I is the nuclear spin, Eq. (17) for a single electron is written aŝ Based on Eq. (17) we now define the magnetic hyperfine interaction constant for n electrons in the field of nucleus K (in a.u. ) where 1 2cmp is the nuclear magneton in a.u. and m p is the proton rest mass. The term 1 M J in the prefactor of Eq. (18) is explained as follows.
The vector operator αi×r iK which depends linearly on the total electronic angular momentum projection quantum number M J . However, the magnetic hyperfine energy has to be independent of M J which is assured by the above prefactor 1 M J . Magnetic hyperfine interaction matrix elements have been calculated based on the implementations in references [4,39] which do not make direct use of the Wigner-Eckart theorem and reduced matrix elements.

A. Technical Details
Gaussian atomic basis sets of double-, triple-, and quadruple-ζ quality [40][41][42] (including correlating functions for 4f and 5d shells in the case of CI and cvDZ/CC) [43] have been used in the present work.
The atomic spinor basis is obtained in Dirac-Coulomb (DC) Hartree-Fock (HF) approximation where the Fock operator is defined by averaging over 6p 1 j=1/2 and 6p 1 j=3/2 open-shell electronic configurations. A locally modified version of the DIRAC program package [44] has been used for all electronic-structure calculations. Interelectron correlation effects are taken into account through Configuration Interaction (CI) theory as implemented in the KRCI module [45] of DIRAC. Coupled cluster (CC) calculations have been carried out within the mrcc code [46][47][48].
We use the experimental value [49] for the nuclear magnetic moment of 205 Tl with nuclear spin I = 1 2 , µ = 1.63821[µ N ], in calculations of the magnetic hyperfine interaction constant.

B. Results for Atomic Interaction Constants
The results from the systematic study of many-body effects on atomic EDM enhancement (R), Ne-SPS interaction ratio (S) and magnetic hyperfine interaction constant (A) are compiled in Table I. The general strategy is to first qualitatively investigate the relative importance of various many-body effects on the properties using a rather small atomic basis set. Then, in a second step, accurate models are developed that include all important many-body effects using the insight from the first step and larger atomic basis sets. Since EDM enhancement and Ne-SPS interaction ratio are analytically related [25,26] it is sufficient to discuss the trends for R only. C.
Step 1: Many-Body effects in cvDZ basis

Valence electron correlation
The result of R = −388 for CAS1in3 which is a singles CI expansion for the electronic ground state can be regarded as close to a DC Hartree-Fock result. The Full-CI (FCI) result including only the three valence electrons (CAS3in4 SDT3/60au) of R = −487 shows that valence correlation effects lead to a considerable change by more than 25% (in the large cvQZ basis by more than 35%). The valence FCI enhancement in cvQZ basis of R = −587 is, therefore, a benchmark. This value is closely reproduced using the universal basis set of reference [22]. Further effects can be considered as modifications of this benchmark result and will be studied one by one.

Subvalence electron correlation
Subvalence electrons of the Tl atoms are those occuping the 5s, 5p, and 5d shells. All other electrons will be considered core electrons. Correlations among the 5d electrons and in particular of the 5d and the valence electrons lead to a strong decrease of R, on the absolute, on the order of 10%. Corresponding contributions from the 5s and 5p electrons are significantly smaller.   Subsets of important CI models based on the findings of the previous subsection have been repeated using the larger atomic basis sets, cvTZ and cvQZ. The single best values from these calculations are given by the model SD18 CAS3in4 SDT21/35au. These latter values V are then corrected by a "correction shift", calculated as follows: The final best CI values are obtained by adding the above sum of individual corrections to the value from the model SD18 CAS3in4 SDT21/35au. Table II gives values of R, S and A || ( 205 Tl) constants obtained within the all-electron coupled cluster with single, double and non-iterative triple cluster amplitudes, CCSD(T), method employing several basis sets. One can see a good convergence of the results in the series of the Dyall's DZ, TZ and QZ basis sets: values of R obtained within the QZ and TZ basis sets differ by about 2%. Table II also gives values of the constants obtained within the Nataraj's universal basis set [22]. Note that the latter basis set is the even-tempered basis set (geometry progression). One can see a good agreement of the results obtained within the QZ basis set and Nataraj's universal basis set. Table III gives values of R calculated with different number of correlated electrons. As can be seen contributions from subvalence and outer-core electrons are close to those obtained within the CI approach above.

E. Accurate CC results
To check the convergence with respect to electron correlation effects we performed a series of successive 21-electron coupled cluster calculations within the TZ basis set (see Table IV). In these calculations two sets of atomic bispinors were used. The first one was obtained within the DCHF approximation where the Fock operator is defined by averaging over 6p 1 j=1/2 and 6p 1 j=3/2 open-shell electronic configurations as in the CI case above. The second one was obtained within the closed-shell DCHF method for the Tl + cation. One can see that CC values gives almost identical result for each set at any level. Besides, the Calculations were performed using DCHF spinors for the neutral Tl atom (V N potential) and for the Tl + cation (V N −1 potential) cases. In both cases the external field perturbs both the spinor coefficients and the CC amplitudes.
-418 -402 CCSD -531 -530 CCSD(T) -521 -522 CCSDT -523 -523 CCSDT(Q) -522 -522  [51] 21310.8 ± 0.0 contribution of correlation effects beyond the CCSD(T) model is almost negligible in the considered case. We considered models up to coupled cluster with Single, Double, Triple and perturbative Quadruple cluster amplitudes, CCSDT(Q). Contribution of the effect of the Breit interaction on R has been estimated in reference [23] as 0.36%. Based on the uncertainties discussed above we conservatively estimate the uncertainty of our final CC value for R to be less than 5%.

F. Discussion in Comparison with Literature Results
Our present best results are shown in Table V in comparison with previous work. The earlier controversy between different groups over results for R(Tl) can be condensed into three main points which we address one by one.

Basis sets
From the results in Tables I and II it is evident that a large atomic basis set, at least of quadruple-zeta quality, must be used for obtaining very accurate interaction constants. The results in VI and II obtained with our correlation methods demonstrate that the basis set used by Nataraj et al. in ref. [22] fulfills this requirement, yielding interaction constants that are very close to those obtained with Dyall's cvQZ basis set and the same correlation expansion. The earlier suggestion of Porsev et al. about an inadequate basis set used in ref. [22] can, therefore, be excluded as a possible reason for the outlier result in ref. [22].

Treatment of correlation effects by the many-body method
It is claimed in reference [22] that the treatment of electron correlation effects were more complete than in references [20] and [21]. We have therefore first attempted to reproduce the electron EDM enhancement calculated by Nataraj et al. by using the same many-body Hamiltonian and EDM operator, the same atomic basis set ("Nataraj universal") and the the same method, CCSD(T). A persisting difference with the approach of Nataraj et al. is the use of CC amplitudes for the closed shells of neutral Tl (our case) or the closed shells of the singly-ionized Tl + (Nataraj case). These results are shown in Table II Table III, becomes R = −555) and in good agreement with the best results of Liu et al. [20], Dzuba et al. [21], and Porsev et al. [23], see Table V. The correct evaluation of the electron EDM enhancement in our codes has been assured by comparative tests of the independent implementations of present CI and CC, as well as with the DIRRCI module [32,52] in the DIRAC program package. All three independent implementations produce the same values of R for small test cases using Full CI / Full CC expansions. These findings strongly suggest that the CC wavefunctions used by us and by Nataraj et al. are almost identical, but that the evaluation of R in reference [22] is flawed.
Since correlation effects have been treated at a very similar (but physically more accurate) level in the present work as in ref. [22] and the result is very different, the claim of correlation effects being responsible for the large difference between previous results is untenable.
First, given a fixed atomic basis set and a fixed many-body Hamiltonian [58], the Full CI expansion delivers the exact solution in the N -particle sector of Fock space [53], independent of the orbital/spinor basis used for this Full CI expansion. This implies that a many-body expansion that closely approximates the Full CI expansion, such as CCSDT or CCSDT(Q), must also be nearly independent of the employed Dirac-Fock potential.
Our results in Table IV clearly confirm this conjecture and demonstrate that even in the more approximate CCSD expansion the electron EDM enhancement factor R is almost independent (0.2% difference) of the underlying spinor set. As the many-body expansion becomes more approximate, such as in the CI model SD18 CAS3in4 SD21 (see Table I) basic theory leads us to expect that the difference in R should increase which is indeed the case (roughly 2% difference). Adding external Triple excitations to the CI expansion, model SD18 CAS3in4 SDT21, quenches the difference to a mere 1.2%, again in accord with expectation. Even the use of a V N −3 potential (i.e., spinors optimized for the Tl 3+ system) changes R by only about 3% relative to spinors for the neutral atom in the SD18 CAS3in4 SDT21 model. This difference is expected to be even smaller in CC models.
Despite of the unimportance of the employed spinor set in highly-correlated calculations, we have used the physically most accurate spinors for the neutral Tl atom in obtaining our best final results. The ratios of our calculated P, T -odd interaction constants are

IV. CONCLUSIONS
We conclude from our findings that the result of Nataraj et al. in ref. [22] is unreliable and should be excluded from the dataset used to constrain the CP-odd parameters d e and C S . Likewise, the result by Sahoo et al. [50] (see Table V) -presumably obtained with a similar code as R(Tl) by Nataraj et al. -is also significantly too small and should be excluded from the reliable dataset. A further proof of this conclusion is the ratio obtained using the values from refs. [22] and [50] which amounts to The panel is the updated version of the one shown in the introduction, using the dataset of reliable calculations of α de and α C S for the Tl atom. The strongly reduced uncertainty of atomic interaction constants for Tl leads to a discernable shrinking of the associated parameter surface (blue), but does not lead to modified constraints. The essential reason for this is the extremely high sensitivity of the experiments on ThO (green) and HfF + (orange). However, tighter constraints on d e and C S can be obtained by including experimental and theoretical results for closed-shell atomic systems as discussed in ref. [54].