New Physics Constraints from Atomic Parity Violation in $^{133}$Cs

Our improved calculation of the nuclear spin-independent parity violating electric dipole transition amplitude ($E1_{PV}$) for $6s ~ ^2S_{1/2} - 7s ~ ^2S_{1/2}$ in $^{133}$Cs in combination with the most accurate (0.3\%) measurement of this quantity yields a new value for the nuclear weak charge $Q_W=-73.71(26)_{ex} (23)_{th}$ against the Standard Model (SM) prediction $Q_W^{\text{SM}}=-73.23(1)$. The advances in our calculation of $E1_{PV}$ have been achieved by using a variant of the perturbed relativistic coupled-cluster theory which treats the contributions of the core, valence and excited states to $E1_{PV}$ on the same footing unlike the previous high precision calculations. Furthermore, this approach resolves the controversy regarding the sign of the core correlation effects. We discuss the implications of the deviation of our result for $Q_W$ from the SM value by considering different scenarios of new physics.

The neutral current weak interactions due to the exchange of a Z 0 boson between the electrons and the nucleus in an atom leads to parity violation [1]. This phenomenon is referred to as atomic parity violation (APV). The nuclear spin-independent (NSI) APV has been measured to an accuracy of 0.35% in the 6s 2 S 1/2 − 7s 2 S 1/2 transition in 133 Cs [2]. This is the most accurate APV measurement to date, but two recent proposals [3,4] have the potential to surpass this accuracy. Thus, the stage is now clearly set to take the APV calculations in Cs to the next level. This indeed provides the motivation for our present work. The principal quantity of interest in the APV studies is the nuclear weak charge (NWC), which is a linear combination of the NSI neutral current weak interaction coupling coefficients between electrons and upand down-quarks in an atom [1,5]. The difference in the model independent value of NWC obtained from APV and that obtained from the Standard Model (SM) could shed light on new physics beyond the SM (BSM).
Following the conventions using 4-fermion operators defined in the PDG [6], the parity-violating leptonhadron interactions at low energies can be described by eγ µ γ 5 e q g eq AV qγµq + eγ µ e q g eq V A qγµγ 5 q , (1) where G F is the Fermi constant and the sum over q includes the interaction of electrons with up (u), down (d) and strange (s) quarks. Note that in earlier editions of the PDG [7] a different notation was used, writing C 1q instead of g eq AV and C 2q instead of g eq V A . The coupling constants g eq AV and g eq V A are defined in the static limit and are universal. The temporal component of the quark currents determines the NSI weak interaction Hamiltonian which is used in the atomic calculations. Assuming that the nucleons can be treated non-relativistically and point-like, Eq. (1) allows us to define parity-violating couplings for protons and neutrons at vanishing momentum transfer, g ep AV = 2g eu AV + g ed AV and g en AV = g eu AV + 2g ed AV . (2) These add up coherently across the nucleus giving rise to the NWC of the nucleus. With Z protons and N neutrons, the NWC at leading order is given by There are corrections due to the fact that the matrix elements of the electromagnetic and axial-vector neutral current operators between electronic states involved in the atomic transition carry a spatial dependence on the electric and weak charge distributions in the nucleus [8]. The dominating part of this correction will be taken into account in the calculation by using a non-pointlike charge distribution inside the nucleus and assuming that the proton and the neutron distributions are the same. However, we will have to add a small correction due to the difference between proton and neutron distributions. It turns out that this correction is dominated by the difference of radii of the proton and neutron distributions, i.e. the neutron skin (NSKIN) effect [9]. Therefore, it can be treated as an additive contribution where q p and q n are determined from the electronic axial form factor weighted by the proton and neutron distributions.
In the SM, the effective low-energy couplings are determined by the weak neutral-current gauge couplings of the Z-boson to quarks and leptons. They are fixed by the charge and isospin quantum numbers and the weak mixing angle, sin 2 θ W . We have g eu AV = − 1 2 + 4 3 sin 2 θ W , g ed AV = 1 2 − 2 3 sin 2 θ W , g eu V A = −g ed V A = 1 2 − 2 sin 2 θ W and therefore, g ep AV = − 1 2 + 2 sin 2 θ W and g en AV = 1 2 . At leading-order, the SM predicts relations between the weak mixing angle, the normalization of the effective 4-fermion operators given by the Fermi constant and the weak boson masses. The Fermi constant is given by G F = πα/( √ 2m 2 W sin 2 θ W ) with the fine structure constant α = e 2 /(4π) and the weak mixing angle is related to the gauge boson masses by m W = m Z cos θ W . High-precision measurements require to take into account higher-order electroweak radiative corrections, both in these parameter relations as well as in predictions for observables. It is convenient to choose the MS renormalization scheme where the weak mixing angle becomes a scale-dependent running coupling, usually denoted by sin 2θ W (µ). Effects due to Feynman diagrams with loops (γZ mixing, vertex corrections, box graphs) contributing to the observable have been calculated in [10] and can be absorbed into corrected effective couplings, as described in Ref. [11]. Numerically, these corrections can be taken into account by replacing Eq. (3) with The scale µ has to be chosen equal to the typical momentum scale of the experiment. We follow Ref. [12] and set µ = 2.4 MeV for Cs, but the precise value is not important since sin 2θ W (µ) depends only very weakly on µ at low scales.
In a specific model, the coefficients of 4-fermion operators are predictions. They would be related to parameters of an underlying theory. For example, 4-fermion contact interactions can originate from Feynman diagrams describing the exchange of a heavy particle at tree-level. Models with extra heavy Z ′ -bosons are well-motivated, for example in string-inspired grand unified models with an E 6 gauge group [13]. Spontaneous symmetry breaking generates two extra U(1) factors whose Z ′ -bosons mix with each other in general. The lighter of them, with mass M Z θ , contributes to the weak charge of the nucleon where z, A and B are parameters which depend on the mixing angle of the two extra Z-bosons, their gauge coupling constant g θ and on the mass [13]. The extra Z ′ is denoted Z χ for the special case of no mixing, and one finds [13,14] ∆Q Z,N W (Z χ ) = (Z + 2N ) APV measurements can therefore set a limit on the mass of such an extra heavy Z ′ boson.
A different type of heavy new particles without direct couplings to the ordinary fermions can enter at the loop level through the W and Z self energies. Examples would be SUSY models at high mass scales, or technicolor models. Just three parameters are needed to describe the corresponding effects on observables, usually called S, T and U [15]. We use the definition described in [6]. Such type of BSM physics can be absorbed in a modification of the neutral-current amplitudes by the factor ρ new = 1 + 0.00782T and by replacing the weak mixing angle with sin 2 θ W × (1 + 0.0157S − 0.0112T ), where the numerical coefficients in these relations are evaluated with the present world-average values of the weak mixing angle and the W -boson mass. This results in [14] Thus, APV is sensitive to S, the isospin-conserving parameter, while the dependence on the isospin-violating parameter T is very small.
Finally, we discuss the case of BSM physics at low mass scale which has caused considerable interest recently and is motivated by the search for a dark matter particle. A light vector boson associated with a U(1) gauge symmetry in the dark sector could couple to ordinary SM matter via kinetic mixing with the photon and mass mixing with the SM Z-boson [16,17]. Such a new boson is known as dark-Z boson and could be an additional source of parity violation. Its effect can be described by a modification of the running of the weak mixing angle in the intermediate to low mass range without visible effects at high-energy Z-pole measurements. A possible realization of such a scenario can be found with a 2-Higgs doublet model where mixing is generated through loop diagrams. The effective weak mixing angle seen at the energy scale µ would be shifted by [17] ∆ sin 2 where M Z d is the mass of the dark-Z, and ǫ and δ are model parameters, depending for example on the charged Higgs-boson mass if the model is realized with two Higgs doublets. For APV we can assume that µ ≪ M Z d which leaves us with Several experiments have narrowed down the parameter space for a dark-Z recently, but there is still room for a significant modification of sin 2 θ W that can be tested with APV. The NSI neutral current weak interaction Hamiltonian in an atom is given by [1]  Energy values This work 31357 (50) 20243(20) 12861(15) 9641(10) 5697(10)  where ρ nuc (r e ) is the electron density within the nucleus The charge in the nucleus is described by a Fermi distribution. The atomic wave function (|Ψ v ) of a state in the Cs atom is calculated by splitting the total Hamiltonian into two parts, where H em represents the dominant electromagnetic interactions in an atom and H N SI We have considered the Dirac-Coulomb-Breit interaction Hamiltonian along with lower-order QED corrections due to the self-energy and vacuum polarization effects as H em in our calculations (for details, see [18]). Since the strength of H N SI AP V is much weaker than that of H em in an atomic system, the wave function |Ψ v represents a state corresponding to the total Hamiltonian H = H em + λH w and its energy (say, E v ) can be expressed as (14) where the superscripts 0 and 1 stand for the zeroth-order and first-order contributions due to H w , respectively. The electric dipole transition amplitude corresponding to two same nominal parity states |Ψ i and |Ψ f in the presence of H N SI AP V can be written as [1,5] In the sum-over-states approach, the firstorder wave function is expanded as |Ψ , where I denotes all possible intermediate states, that can be divided into core states (contributions from these states are designated as "Core"), low-lying bound states (contributions from these states are given as "Main"), and the remaining high-lying states including continuum (whose contributions are mentioned as "Tail") for computational simplicity.  [42] 0.3464(10) 8P 1/2 ↔ 7S 0.9623 (20) 0.2296 (05) The latest two high precision calculations, reported in Refs. [19,20], are carried out by estimating the "Core", "Main" and "Tail" contributions by applying mixed many-body methods. The calculations in Ref. [19] included the valence triple excitation effects to "Main" by employing the relativistic coupled-cluster (RCC) theory, and it was found that these effects to the atomic properties of 133 Cs were relatively important in reducing the uncertainty in the E1 P V amplitude to 0.27% [19]. This result was in good agreement with the SM, however the calculation on which it is based had used a sum-overstates approach in which the 'Main' contributions were estimated only from the excited states up to the principal quantum number n = 9. Later Dzuba et al. reported another result in Ref. [20] with 0.5% accuracy by using the "Main" contribution from Ref. [19], but with different "Core" (opposite sign than [19]) and "Tail" contributions by taking into account certain sub-classes of correlation effects. They found substantial differences in these contributions from Ref. [19]; especially the "Core" contribution differed by about 200% (due to opposite sign). This resulted in 0.8% difference between the final results of Porsev et al. [19] and Dzuba et al. [20]. In addition, both the above works did not include double core-polarization (DCP) effects [21], and contributions from the Breit and QED effects were taken from the earlier works. To include all these neglected contributions and to treat all the electron correlation effects on an equal footing, we solve the inhomogeneous equation for the first-order wave function where E (1) v = 0 in the present case owing to the oddparity nature of H w . This is achieved by expressing the unperturbed and the first-order wave function of the Cs  (27) Ref. [20] 0.0018 (8)  a Contains additional contribution from the 9p 2 P 1/2 state. b Taken from previous calculation [43].
atom in the RCC theory framework as [22][23][24] and where |Φ v is obtained by determining the Dirac-Hartree-Fock (DHF) wave function of the closed-core (|Φ 0 ) and then, appending the corresponding valence orbital v to it are the core and the valence excitation operators. The superscript 0 represents the absence of any external perturbation. Similarly, T (1) and S (1) v are the core and the valence excitation operators with the superscript 1 representing the order of perturbation in H w . In our previous calculations, we had successfully employed this approach based on RCC theory with singles and doubles approximation (RCCSD method) for the evaluation of E1 P V amplitudes in Ba + [22], Ra + [23] and Yb + [24] and had achieved results within 1% accuracy. In the present work, we have implemented additional triple excitations beyond the RCCSD method (RCCSDT method) to achieve sub-one percent accurate E1 P V in 133 Cs as there is a renewed interest in the inclusion of the neglected correlation effects in this atom (e.g. see discussions in [25,26]). It is worth mentioning here that we excite all the electrons in the RCCSD method to account for the electron correlation effects, but correlate all the electrons except from the 1 − 3s, 2 − 3p, and 3d occupied orbitals and beyond n = 15 virtual orbitals for triple excitations due to limitations in the available computational resources. Also, we have considered active orbitals up to l = 5 in our RCC calculations and contributions from the orbitals belonging to higher angular momentum symmetries, quoted as 'Extra' hereafter, are estimated using low-order perturbative methods.
We first calculated the energies, electric dipole (E1) matrix elements and magnetic dipole hyperfine structure constants (A hyf ) of the states that give rise to dominant contributions to the determination of E1 P V in 133 Cs. By comparing these values with their corresponding experimental results [27][28][29][30][31][32][33][34][35][36][37][38][39][40][41][42], we have assessed the accuracies of the wave functions in the regions close to and far away from the nucleus. Our calculated values of these properties at different levels of approximations can be found in Ref. [18], however, the final values along with their experimental results are listed in Tables I and II. As can be seen from these two tables, comparison of our calculations of the above properties with their precisely known experimental data is very impressive and within sub-one percent accuracy. This strongly suggests that our atomic wave functions are very reliable and they can be used to determine the transition amplitude E1 P V accurately.
Keeping in mind our classification of the RCC terms, we find the 'Core' contribution to E1 P V and also the 'Main' contribution using our calculated energies, E1 matrix elements and amplitudes of H N SI AP V for the intermediate n(= 6, 7, 8)P 1/2 states that are quoted in Tables I  and II. After subtracting the "Core" and "Main" contributions from the final value, the remainder is taken as the "Tail" contribution. These contributions from the DHF, RCCSD and RCCSDT methods using the Dirac-Coulomb Hamiltonian are quoted in Table III. In addition, we give contributions from the Breit and QED interactions of our calculation using the RCCSDT method in the same table. The other neglected contributions due to 'Extra', possible neutral weak interactions among electrons (e − e), and the NSKIN effect that were not included in our RCC calculation, are also quoted in the above table. The small e − e contribution to E1 P V has been taken from Ref. [43].
In a seminal work, Fortson et al [44] had analyzed the NSKIN effect on APV. By adopting their analysis, the effect of NSKIN on E1 P V (δE1 N S P V ) can be estimated by where t is known as the neutron skin parameter which describes the relative difference of the r.m.s. radii of the neutron and proton distributions in the nucleus. This empirical formula was used in Refs. [19,20] for determining δE1 N S P V . Using t = 0.033(8) [45] and the uncorrected E1 P V value 0.8914, we get δE1 N S P V ≃ −0.0020 (5) in units of ×10 −11 i(−Q W /N )ea 0 . This is in good agreement with the recently estimated value by Brown et al. in Ref. [45]. However, Sil et al. have estimated ∆Q NS W by employing a more rigorous effective field theory framework [9]. We use the relation 20) and the numerical results from Ref. [9]. Interpreting the two model results considered there to define a central value and an uncertainty range, we find δE1 N S P V ≃ −0.00377(39)× 10 −11 i(−Q W /N )ea 0 by substituting Q W ≃ −73.23. This is a slightly larger correction than considered in the previous Refs. [19,20,46].
We compare individual contributions with the previously reported RCC results using the sum-over-states approach [19,46] and the latest reported result [20] of E1 P V in Cs. These calculations include the 9P 1/2 state in their "Main" contribution in the sum-over-states approach, whereas our "Tail" includes the contribution from this state. Our final E1 P V value is 0.8893 (27) in contrast to the results that have been reported previously as 0.8906(24) [19] and 0.8977 (40) [20] in units of ×10 −11 i(−Q W /N )ea 0 . The major difference between the results from Ref. [19] and ours is because of the fact that they account for different NSKIN effect. Large difference between the present calculation and that of Ref. [20] is mainly due to the 'Core' contributions, which have different signs in both the cases.
Nonetheless, one of the most important achievements of our calculation is that it resolves the ambiguity of the sign of the"Core" contribution from the calculations reported in Refs. [19] and [20]. We have adopted the same procedure as in Ref. [19] to estimate the uncertainty of E1 P V . This is also independently verified by analyzing major sources of uncertainties that can come from the neglected higher level excitations in the RCC theory and finite-size basis functions used in the calculation. If the difference between the RCCSD and RCCSDT values is assumed to be the maximum contribution from the neglected higher level excitations and considering 'Extra' as the maximum uncertainty due to incompleteness in the used basis functions, we also arrive at the same uncertainty of E1 P V .
It is necessary to combine our E1 P V value with the precisely measured Im(E1 P V /β) = 1.5935(56) mV/cm [2], where Im refers to the imaginary part and β is the vector polarizability of the 6s 2 S 1/2 − 7s 2 S 1/2 transition in 133 Cs to extract Q Z,N W . However, this also requires an accurate knowledge of β. We have determined this quantity by using the precisely available E1 matrix elements either from measurements or from our calculation as discussed in detail in Ref. [18] and obtain its value as β = 27.12 (4)  can be used to determine the weak mixing angle. We find sin 2θ W (2.4 MeV) = 0.2408(16) with a slightly smaller uncertainty and a significant shift of the central value compared with the previous determination [6].
The experimental value of Q Z,N W provides a constraint on the low-energy effective electron-quark couplings: 376g eu AV + 422g ed AV = 73.71 (35). Assuming the SM prediction for one of them, we find a value for the other: g eu AV = −0.1877(9) for g ed AV = 0.3419 and g ed AV = 0.3429(8) for g eu AV = −0.1888. We also find slightly improved limits on the BSM parameters described in the introduction. The isospin conserving oblique parameter S can be constrained to S ≃ 0.60(44) assuming T = 0 in Eq. (9). The central value of S is shifted to positive values, compared with the previous determination S ≃ −0.51(52) [6]. From Eq. (7), we obtain a limit on the mass of an extra Z ′ boson. Since the shift is always positive, we find a one-sided exclusion limit at 95 % confidence limit of M Zx > 2.36 TeV, using M W = 80.379 GeV [6], compared with recent limits from the ATLAS collaboration who found values ranging from 3.5 to 4.5 TeV [47]. Furthermore, using Eq. (11) we can constrain the dark-Z model parameter ǫδ MZ MZ d ≃ −0.0051 (37).
In conclusion, we used a perturbed version of the relativistic coupled-cluster theory to calculate the nuclear spin-independent parity violating electric dipole transition amplitude for the 6s 2 S 1/2 − 7s 2 S 1/2 transition in 133 Cs. The principal merit of this approach is that it treats the contributions of the core, valence and excited states to the above parity violating transition amplitude on the same footing, thereby overcoming the limitations of the previous high precision calculations of this quantity. Our work resolved the ambiguity in the sign difference for the contribution from core states. In addition, we estimated the uncertainty in our calculation of the parity violating transition amplitude in 133 Cs more rigorously than those in previous calculations. The salient implications of the deviation of the nuclear weak charge from the standard model, that is obtained in the present work, for probing possible new physics have been discussed. Our result, in combination with measurements from proposed new high-precision experiments, has the potential to improve the constraints on beyond the standard model physics in the future.
We acknowledge support by the Mainz Institute for Theoretical Physics (MITP). This work was initiated during the MITP Virtual Workshop "Parity Violation and Related Topics". Computations reported in this work were performed using the PRL Vikram-100 HPC cluster.