# Influence of the atomic modeling on the electron capture process

X. Mougeot<EMAIL_ADDRESS>Université Paris-Saclay, CEA, List,
Laboratoire National Henri Becquerel (LNE-LNHB), F-91120 Palaiseau, France  L.
Mouawad  A. Andoche  P.-A. Hervieux<EMAIL_ADDRESS>Université de
Strasbourg, CNRS, Institut de Physique et Chimie des Matériaux de Strasbourg,
UMR 7504, F-67000 Strasbourg, France

(August 27, 2024)

###### Abstract

Ongoing experimental efforts to measure with unprecedented precision electron
capture probabilities will challenge the current theoretical models. The short
range of the weak interaction necessitates an accurate description of the
atomic structure down to the nucleus region. In this study, an atomic modeling
has been developed in the framework of the relativistic density functional
theory within which several exchange-correlation functionals and self-
interaction-corrected models have been explored. A recent electron capture
model has been modified in order to only test the influence of the atomic
modeling on the calculated decay probabilities. In particular, two important
possibilities have been tested, the “no vacancy” and the “frozen orbital”
approximations. Results have been compared to experimental values evaluated
from available published data. New high-precision measurements are strongly
encouraged as the accuracy of the current experimental values does not allow
to draw a firm conclusion on the theoretical modeling to be employed.

††preprint: APS/123-QED

##  I Introduction

Electron capture is a low-energy weak interaction process in which a proton
embedded in a nucleus absorbs an atomic electron, leading to a nucleus of a
lower atomic number and the emission of a neutrino with a well-defined energy.
If the daughter nucleus is not in its ground state after the transition, its
deexcitation goes through a cascade of $\gamma$ transitions that leads to the
emission of either $\gamma$-rays or internal conversion electrons, the latter
creating atomic vacancies. In addition, the daughter atom just after the
transition is overall neutral but in an excited state, the capture process
being seen in its electronic configuration as the promotion of the captured
electron to the valence shell. The subsequent relaxation of the daughter atom
then leads to the emission of X-rays and Auger electrons all along the
propagation of the vacancies, until their final disappearance.

Ground-state to ground-state transitions are very difficult to measure,
especially for low-mass nuclei, because only X-rays and Auger electrons of
very low energies can be detected. Standardization of such pure electron
capture emitters in radionuclide metrology is challenging and directly depends
on the knowledge of the capture probabilities, which most often come from
theory. Recently, a new international reference system, called ESIR [Coulon
_et al._ , 2020], has been established by the International Bureau of Weights
and Measures (BIPM) based on activity measurements carried out with the liquid
scintillation counting (LSC) technique. The computed LSC efficiency strongly
decreases at very low energy and a precise modeling of the atomic
rearrangement is necessary, which depends on the number of initial vacancies
[Broda  _et al._ , 2007].

What is problematic for radionuclide metrology can become an advantage in
nuclear medicine for cancer treatment when the radioactive nuclei are
correctly vectorized into the human body. The vast majority of Auger electrons
are emitted with a few keV kinetic energy, deposited in a range from
nanometers to micrometers. This property makes them highly promising for an
accurate, well-controlled internal radiotherapy at the cell or even the DNA
level. Many short-lived nuclei decaying by electron capture have been studied
in this perspective [Bavelaar  _et al._ , 2018,Ku  _et al._ , 2019].

Long-lived isotopes can also be useful in other contexts. The decommissioning
of legacy nuclear sites necessitates the radioactivity inventory of various
materials. When non-destructive techniques such as $\gamma$ spectrometry are
impossible, precise measurements of samples are used together with scaling
methods to estimate the total amount of radioactivity [Hou, 2005,Lee  _et al._
, 2021]. Some long-lived pure electron capture emitters have been used to
investigate the age of the Solar System [Huss  _et al._ , 2009,Jörg  _et al._
, 2012] and the irradiation history of meteorites and lunar samples [Kubik
_et al._ , 1986].

Every application refers to evaluated nuclear decay data, which can be taken
either in the ENSDF (Evaluated Nuclear Structure Data File) database [Bhat,
1992] or in the DDEP (Decay Data Evaluation Project) database [Kellett, Mark
A. and Bersillon, Olivier, 2017], the latter being recommended by the BIPM for
metrology purposes. In ENSDF evaluations, the properties of $\beta$
transitions and electron captures are still calculated with the LogFT code
[Gove and Martin, 1971], which is based on an approximate theoretical model
limited by the computing power available in the 1970s. Recently, the BetaShape
code has been developed by one of the authors in order to improve the
theoretical predictions for $\beta$ transitions and electron captures
[Mougeot, 2015,Mougeot, 2019]. Its modeling takes into account several
additional physical phenomena and provides much more detailed information
required by different communities. The DDEP collaboration has already adopted
this new code for its nuclear data evaluations.

Capture probabilities from BetaShape have been compared to a selection of
measurements available in the literature, concluding to the need of new high-
precision measurements to validate and constrain the theoretical models
[Mougeot, 2018]. This is in part at the origin of the European metrology
project MetroMMC [Ranitzsch  _et al._ , 2020], which aims at improving the
knowledge of electron capture decay and the subsequent atomic relaxation
processes. The high-precision measurements from this project will challenge
the theoretical predictions, for which the accuracy of the atomic modeling is
pivotal.

Indeed, as the electron capture process takes place inside the nucleus, it is
mandatory to develop models that are as precise as possible and capable of
describing the electronic properties of atoms in this space region. In
addition, we can also wonder about the role played by electron correlations in
the decay process. Trying to answer these two issues constitute the main goal
of this work.

We first describe briefly the electron capture modeling and the determination
of the atomic wave functions in the BetaShape code. Next, a model developed in
the framework of relativistic density functional theory within the local
density approximation (RLDA) is described in details. Within the latter, most
popular exchange-correlation functionals and self-interaction-corrected models
have been explored. Binding energies, wave functions and capture probabilities
are finally compared between each other, and to experimental values when
available.

##  II Theoretical methods

Our modeling of electron capture decays is based on the formalism from
[Bambynek  _et al._ , 1977,Behrens and Bühring, 1982] and has already been
described in detail in [Mougeot, 2018,Mougeot, 2019]. In this section, we
focus on the dependency of this modeling on the atomic wave functions, and on
the determination of the wave functions themselves.

###  II.1 Electron capture

The low rest mass of the beta particle compared to the transition energies
necessitates a fully relativistic formalism of electron capture decay, usually
expressed in relativistic units ($\hbar=m_{e}=c=1$). Considering spherical
symmetry, radial and angular parts of the wave functions can be separated. An
atomic orbital is then characterized by its quantum numbers ($n$,$\kappa$),
the latter being the eigenvalue of the operator
$K=\beta(\vec{\sigma}\cdot\vec{L}+1)$ defined from the Dirac ($\beta$) and
Pauli ($\vec{\sigma}$) matrices and the angular momentum operator ($\vec{L}$).

Atomic levels being degenerated, one has to introduce their relative
occupation number $n_{n\kappa}$. In addition, each wave function is
characterized by its Coulomb amplitude $\beta_{n\kappa}$, as defined in
[Bambynek  _et al._ , 1977,Behrens and Bühring, 1982]. This quantity has to be
computed and resumes to the value of the wave function at the origin ($r=0$)
for $\kappa=-1$, i.e. in the case of $s_{1/2}$ orbitals (see Section III.2).

The transition probability per unit of time can be derived by applying first-
order time-dependent perturbation theory and results in the sum of the capture
probability of each subshell:

$\lambda_{\varepsilon}=\dfrac{G^{2}_{\beta}}{2\pi^{3}}\sum\limits_{n,\kappa}{\frac{\pi}{2}n_{n\kappa}C_{n\kappa}q^{2}_{n\kappa}\beta^{2}_{n\kappa}B_{n\kappa}\left(1+\sum\limits_{m,\kappa}{P_{m\kappa}}\right)}$
(1)

with $G_{\beta}$ the Fermi coupling constant and $q_{n\kappa}$ the neutrino
momentum. We do not consider in the present work any electron capture decay
for which the available energy would be sufficient for a competition with a
$\beta^{+}$ transition.

The quantity $C_{n\kappa}$ couples the lepton and nucleon wave functions. For
allowed and forbidden unique transitions, the nuclear dependency acts as a
constant factor. Therefore, determining ratios of relative capture
probabilities prevents from introducing any nuclear structure. In such cases,

$C_{n\kappa}\propto\dfrac{p_{n\kappa}^{2(k-1)}\leavevmode\nobreak\
q_{n\kappa}^{2(L-k)}}{(2k-1)![2(L-k)+1]!}$ (2)

with $k=|\kappa|$, $p_{n\kappa}$ the electron momentum and
$L=\left|J_{i}-J_{f}\right|$ the difference between the total angular momenta
of initial and final nuclear levels.

The quantity $B_{n\kappa}$ is a correction that accounts for overlap and
exchange effects. The former is induced by the capture process that changes
the nucleus charge, leading to an imperfect overlap of the initial and final
wave functions of the spectator electrons. The latter arises from electrons
being identical particles. The two following methods, available in the
literature, have been considered in the BetaShape code for calculating
$B_{n\kappa}$. Bahcall’s approach [Bahcall, 1965] only considers the first
three $s_{1/2}$ orbitals, assumes a complete set of states for the other
orbitals and makes use of the closure property to sum over the continuum
states. Vatai’s approach [Vatai, 1970] considers in addition the exchange term
of the $4s_{1/2}$ orbital and the overlap correction of each subshell, but
does not use the closure property. The generalization of theses approaches
[Mougeot, 2018] leads to
$B_{n\kappa}=\left|b_{n\kappa}/\beta_{n\kappa}\right|^{2}$ with

$b_{n\kappa}=t_{n\kappa}\left[\prod\limits_{m\neq
n}{<(m,\kappa)^{\prime}|(m,\kappa)>}\right]\times\left[\beta_{n\kappa}-\sum\limits_{m\neq
n}{\beta_{m\kappa}\dfrac{<(m,\kappa)^{\prime}|(n,\kappa)>}{<(m,\kappa)^{\prime}|(m,\kappa)>}}\right]\;,$
(3)

where $<(m,\kappa)^{\prime}|(m,\kappa)>$ is the overlap of the atomic wave
functions between the initial ($m$,$\kappa$) state and the final
($m$,$\kappa$)’ state, and correspondingly for the other overlaps. Bahcall’s
model is given by $t_{n\kappa}=1$ and Vatai’s model by

$\begin{array}[]{rl}t_{n\kappa}=<(n,\kappa)^{\prime}|(n,\kappa)>^{n_{n\kappa}-1/2|\kappa|}&\times\left[\prod\limits_{m\neq
n}{<(m,\kappa)^{\prime}|(m,\kappa)>^{n_{m\kappa}-1}}\right]\\\ &\\\
&\times\left[\prod\limits_{\tiny\begin{array}[]{c}m,\mu\\\ \mu\neq\kappa\\\
\end{array}\normalsize}{<(m,\mu)^{\prime}|(m,\mu)>^{n_{m\mu}}}\right]\;.\end{array}$
(4)

The difference between Bahcall and Vatai models in the resulting capture
probabilities is used in BetaShape as a theoretical uncertainty component. The
other component comes from input data and corresponds to the propagation of
the uncertainty on the transition energy.

The last term in Equation (1) corrects for the shaking effects, i.e. internal
excitation (shake-up) and internal ionization (shake-off). These effects take
into account additional final states in which other atomic vacancies are
present with the initial vacancy due to capture. The modeling implemented in
BetaShape is limited to a single additional vacancy due to shaking. The
separate calculation of each process requires the knowledge of all unoccupied
bound states up to the Fermi level for shake-up, and summing over the infinite
continuum states for shake-off. However, the determination of capture
probabilities only requires the total probability $P_{m\kappa}$ of creating a
secondary vacancy in an orbital ($m$,$\kappa$) consecutive to the capture of
an electron in an orbital ($n$,$\kappa$). Such a probability is simpler to
calculate when considering the non-shaking probability, and can be expressed
as [Crasemann  _et al._ , 1979]

$P_{m\kappa}=1-<(m,\kappa)^{\prime}|(m,\kappa)>^{2n_{m\kappa}}-\sum\limits_{l\neq
m}{<(l,\kappa)^{\prime}|(m,\kappa)>^{2}}\;.$ (5)

These effects are already taken into account in Bahcall modeling of overlap
and exchange with the use of the closure property. On the contrary, they are
not included in Vatai’s approach. The shaking correction is therefore only
applied along with Vatai modeling.

Finally, the possibility of radiative capture process is taken into account
with radiative corrections, for which calculation depends on the atomic
energies. However, they are only significant when a $\beta^{+}$ transition
competes, i.e. for the determination of capture-to-positron ratios.

###  II.2 BetaShape wave functions

The relativistic electron wave functions in BetaShape are determined following
essentially the method of Behrens and Bühring [Behrens and Bühring, 1982], as
described in detail in [Mougeot, 2019,Mougeot and Bisch, 2014]. A neutral atom
is considered, which orbitals are filled in following the Madelung rule. The
Coulomb potential of the nucleus is modeled from a uniformly charged sphere
and a screened potential parameterized on Dirac-Hartree-Fock-Slater self-
consistent results from [Salvat  _et al._ , 1987]. For the bound states, an
exchange potential is added that includes a constant prefactor. This prefactor
is adjusted in order to force the convergence procedure to tabulated atomic
energies.

As underlined by Vatai in [Vatai, 1970], the vacancy created in a subshell by
the capture process has a significant influence on all the orbitals. This hole
effect is corrected employing first-order time-independent perturbation theory
[Mougeot, 2018], as proposed by Vatai.

The tabulated atomic energies from Desclaux [Desclaux, 1973] for relativistic
single electrons were used in a first study [Mougeot, 2018]. However,
comparison with a selection of experimental values showed that the accuracy of
the atomic wave functions was not sufficient to distinguish between Bahcall
and Vatai models of overlap and exchange correction. In a second study
[Mougeot, 2019], atomic energies from Kotochigova et al. [Kotochigova  _et
al._ , 1997] recommended on the NIST website [Kotochigova  _et al._ , 2021]
were considered. They had been determined for a point charge nucleus using the
relativistic local-density approximation, including some electron
correlations. Good agreement was found with experimental values and as
expected, Vatai model of overlap and exchange, corrected for shaking, was
found to be more accurate than Bahcall model.

It is clear that this method for determining the atomic wave functions is
inconsistent: first, because the nucleus is a spherical charge distribution
while Kotochigova’s modeling is for a point charge nucleus; second, because no
correlations are present in the adjusted Coulomb potential. On the contrary,
the model described below is totally consistent. A specific version of the
BetaShape code has been developed for using its wave functions as input, read
from external files together with the energies and electronic configuration.
The approximate correction of the hole effect has been removed and this effect
is exactly taken into account by determining the wave functions for an
electronic configuration with a vacancy due to capture.

###  II.3 A realistic atomic modeling

In order to describe the electron capture process, one has to work in the
framework of relativistic quantum mechanics for describing the electronic part
of the atom. In this context, it is also mandatory to include correlation
effects – beyond mean-field – and nuclear finite size effects by using a
precise description of the nucleus charge distribution – beyond the usual
point charge approximation. Indeed, the penetration of the electron wave
function inside the nucleus can be important for heavy elements.

One of the most popular and powerful method for including correlation effects
at a relativistic level is the multiconfiguration Dirac-Fock (MCDF) model,
which uses a linear combination of Slater determinants to approximate the
electron many-body wave function [Desclaux, 1975]. However, due to the
exponential increase with the atomic number of configurations involved in the
computation, the use of MCDF is restricted to atoms having a relatively low
atomic number ($Z<30$).

An alternative approach, although a priori less accurate, is provided by the
relativistic density functional theory (RDFT) within the local density
approximation (RLDA) [Dreizler and Engel, 1998]. The latter is much less time
consuming, more easy to use and able to cover a large number of radionuclides.
In the present work, we have developed a specific calculation code based on
this method. The main ingredients of RLDA are reminded below and we refer to
[Rajagopal, 1978,MacDonald and Vosko, 1979,Engel, 2002] for more details. In
the following, a spherical symmetry is assumed and atomic units are used if
not otherwise stated.

In RDFT, one must solve the one electron Dirac equation

$\left[T+\beta/\alpha^{2}+V_{\mathrm{eff}}[\rho](r)\right]\Psi_{i}(\vec{r})=\varepsilon_{i}\Psi_{i}(\vec{r})\;,$
(6)

which can be expressed as

$\left[\vec{\alpha}.\vec{p}/\alpha+\beta/\alpha^{2}+V_{\mathrm{eff}}[\rho](r)\right]\Psi_{i}(\vec{r})=\varepsilon_{i}\Psi_{i}(\vec{r})\;,$
(7)

where $\vec{\alpha}$ and $\beta$ are the Dirac matrices [Strange, 2008],
$\alpha$ is the fine structure constant, $\Psi_{i}$ is a bi-spinor and
$V_{\mathrm{eff}}[\rho]$ is an effective spherical one-particle potential. $T$
is the kinetic energy that introduces a coupling between the small and large
components of the Dirac wave function and $\beta/\alpha^{2}$ is the electron
rest energy. The electron density is given by

$\rho(r)=\sum_{i=1}^{N_{\mathrm{orb}}}\Psi_{i}^{\dagger}(\vec{r})\Psi_{i}(\vec{r})=\sum_{i=1}^{N_{\mathrm{orb}}}\rho_{i}(r)\;,$
(8)

where $N_{\mathrm{orb}}$ is the number of orbitals.

Within the framework of RLDA, the effective potential can be expressed as

$V_{\mathrm{eff}}[\rho](r)=V_{\mathrm{H}}[\rho](r)+V_{\mathrm{xc}}[\rho](r)+V_{\mathrm{ext}}(r)\;,$
(9)

where
$V_{\mathrm{H}}[\rho](r)=\int\frac{\rho(r^{\prime})}{\left|\vec{r}-\vec{r}^{\prime}\right|}d\vec{r}^{\prime}$
is the Hartree potential, $V_{\mathrm{xc}}[\rho]$ is the relativistic exchange
and correlation potential, and $V_{\mathrm{ext}}$ is the external potential
due to the interaction of the electrons with the nucleus.

The nuclear charge distribution $\rho_{\mathrm{nuc}}(r)$ is modeled by using a
two-parameter Fermi-type distribution [Pálffy, 2010]. The value of the root
mean square charge (RMS) radius of the nucleus, which is needed as an input
parameter, is taken from [Angeli and Marinova, 2013] where existing
experimental data are quoted. Otherwise, the empirical formula from [Johnson
and Soff, 1985] is used. The external potential $V_{\mathrm{ext}}$ is then
computed thanks to

$V_{\mathrm{ext}}(r)=-\int\dfrac{\rho_{\mathrm{nuc}}(r^{\prime})}{\left|\vec{r}-\vec{r}^{\prime}\right|}d\vec{r}^{\prime}\;.$
(10)

For $V_{\mathrm{xc}}[\rho]$, we have used different exchange-correlation
functionals: Vosko, Wilk and Nusair (VWN) [Vosko  _et al._ , 1980], Gunnarsson
and Lundqvist (GL) [Gunnarsson and Lundqvist, 1976] and Perdew and Wang (PW)
[Perdew and Wang, 1992]. We have also implemented the gradient-dependent
exchange-correlation functional of van Leeuwen and Baerends (LB) [van Leeuwen
and Baerends, 1994]. This functional is non-local and depends not only on the
electron density but also on the electron density gradient.

We have elaborated a computer code for generating electron energies and wave
functions. It is based on the numerical schemes of [Čertík  _et al._ , 2013]
and is complemented by the incorporation of a realistic description of the
nucleus and by the implementation of the optimized effective potential method
(see Section II.4). We have used an exponential grid well adapted for radius
values within the area of the nucleus, which are essential for the modeling of
the electron capture process. The possibility to create atomic vacancies, as
needed for the modeling of electron capture, has been also implemented.

###  II.4 Optimization of the effective potential

We have first tested and benchmarked our code by reproducing the atomic
energies from Kotochigova et al. [Kotochigova  _et al._ , 1997] for atoms from
Hydrogen to Uranium. They have been obtained within the framework of RLDA
using the Vosko, Wilk and Nusair exchange-correlation functional and a point
charge model for describing the charge distribution of the nucleus.

In a second step, we have optimized our modeling in regards to the effective
potential. For a neutral atom, the asymptotic behavior of the effective
potential is given by the exchange contribution to $V_{\mathrm{xc}}$, which
behaves at large distance as $\rho^{1/3}$. As a consequence, the effective
potential $V_{\mathrm{eff}}$ decreases exponentially to zero and does not
reproduce the expected $1/r$ asymptotic behavior. This problem does not appear
in the Dirac-Fock theory because the Dirac-Fock exchange potential exactly
compensates the self-interaction term contained in the Hartree potential (see
Figure 1).

Generally, the use of this correction leads to a much better description of
the outermost shells. In particular, the experimental ionization potential can
be generally well reproduced [van Leeuwen and Baerends, 1994]. In addition,
although it has not been studied in the present work, this correction is
crucial for a correct description of unoccupied and continuum states that are
essential for shake-up and shake-off processes. However, the drawback of this
correction and its variants (ADSIC [Ciofini  _et al._ , 2005], SIC [Perdew and
Zunger, 1981], GAM [Politis  _et al._ , 1998]) is an inaccurate description of
the inner shells, which are mainly involved in the electron capture process
(i.e. $1s_{1/2}$, $2s_{1/2}$, $2p_{1/2}$). A way to overcome this
inconvenience is to use the optimized effective potential method, which was
originally described within a non-relativistic framework by J. B. Krieger et
al. [Krieger  _et al._ , 1992,Li  _et al._ , 1993] and extended to the
relativistic case by Xiao-Min Tong and Shih-I Chu [Tong and Chu, 1998].

In order to compute the optimized effective potential
$V_{\mathrm{eff}}^{\mathrm{KLI}}$ – giving what we call the KLI model in the
following, one has to solved self-consistently the equation

$\begin{array}[]{rl}V_{\mathrm{eff}}^{\mathrm{KLI}}(r)=&V_{\mathrm{eff}}[\rho](r)+\dfrac{fo(N_{\mathrm{orb}})\rho_{N_{\mathrm{orb}}}(r)}{\rho(r)}\times(-V_{\mathrm{H}}[\rho_{N_{\mathrm{orb}}}](r)-V_{\mathrm{xc}}[\rho_{N_{\mathrm{orb}}}](r))\\\
&\\\
&+\sum\limits_{i=1}^{N_{\mathrm{orb}}-1}{\dfrac{fo(i)\rho_{i}(r)}{\rho(r)}}\times\\{-V_{\mathrm{H}}[\rho_{i}](r)-V_{\mathrm{xc}}[\rho_{i}](r)\\\
&\\\
&-\int_{0}^{\infty}\rho_{i}(r)(-V_{\mathrm{xc}}[\rho_{i}](r)-V_{\mathrm{H}}[\rho_{i}](r))r^{2}dr+\int_{0}^{\infty}\rho_{i}(r)V_{\mathrm{eff}}^{\mathrm{KLI}}(r)r^{2}dr\\}\;,\end{array}$
(11)

where $fo(i)$ is the occupation number of the orbital $i$. The term
$V_{\mathrm{eff}}^{\mathrm{KLI}}(r)$ that appears on both sides of the above
equation will be quickly adjusted after a few iterations. All the details
concerning the derivation of this equation are given in [Tong and Chu, 1998].

Figure 1:  (Color online) KLI (solid) and Kotochigova (dashed) effective
potentials for 109Cd.

##  III Results

We have carried out many simulations employing different exchange-correlation
functionals (VWN, GL, PW and LB) and self-interaction-corrected models (ADSIC,
SIC, GAM and KLI). It turns out that KLI is the best compromise for obtaining:
i) electron binding energies of inner subshells close to the experimental
values (see Table 1); and ii) a long-range Coulomb behavior of the effective
potential (see Figure 1). Moreover, we have found that within the KLI model,
the results depend only very slightly on the exchange-correlation functional
utilized. Therefore, we have employed the VWN functional for the comparison of
the results obtained with both Kotochigova and KLI models.

Table 1: Binding energies (in atomic units) of the innermost orbitals for Kotochigova and KLI models. Experimental values have been taken in [Sevier, 1979] for 7Be and all values with uncertainties, and in [Bearden and Burr, 1967] in other cases. |  | $1s$ |  |  | $2s$ |  |  | $2p_{1/2}$ |
---|---|---|---|---|---|---|---|---|---
Element | Koto. | KLI | Exp. | Koto. | KLI | Exp. | Koto. | KLI | Exp.
7Be | -3.856 | -4.815 | -4.384 | -0.206 | -0.377 | -0.343 | – | – | –
41Ca | -144.40 | -148.93 | -148.57 (7) | -15.16 | -16.14 | -16.26 (7) | -12.375 | -13.394 | -13.230 (4)
54Mn | -235.08 | -240.81 | -240.70 (6) | -27.22 | -28.38 | -28.72 (6) | -23.368 | -24.606 | -24.357 (37)
55Fe | -255.90 | -261.88 | -261.80 (4) | -29.99 | -31.20 | -31.49 (7) | -25.92 | -27.21 | -26.94 (4)
109Cd | -968.69 | -981.42 | -981.60 | -143.97 | -146.64 | -147.66 | -23.37 | -24.61 | -24.36 (4)
125I | -1204.13 | -1218.71 | -1218.93 | -186.24 | -189.33 | -190.65 | -174.98 | -178.32 | -178.31
138La | -1414.24 | -1430.37 | -1430.46 | -225.36 | -228.79 | -230.27 | -212.81 | -216.51 | -216.49

###  III.1 Binding energies

Binding energies of the inner subshells of the atoms considered in this work
are given in Table 1 for Kotochigova and KLI models, and are compared to
experimental values taken in [Sevier, 1979,Bearden and Burr, 1967]. It is
clear that KLI energies are much more accurate, which validates the approach
developed to correct for the asymptotic behavior of the effective potential
and to obtain accurate binding energies for the inner subshells.

###  III.2 Wave functions

The use of the KLI wave functions in the electron capture modeling described
in Section II.1 requires, of course, to go from the atomic to the relativistic
units. However, one also has to determine the Coulomb amplitudes
$\beta_{n\kappa}$ for each wave function. In Behrens and Bühring formalism
[Behrens and Bühring, 1982], the relativistic electron wave function is
defined as

$\psi_{n\kappa}(\vec{r})=\left(\begin{array}[]{c}S_{\kappa}f_{n\kappa}(r)\chi^{\mu}_{-\kappa}\\\
g_{n\kappa}(r)\chi^{\mu}_{\kappa}\\\ \end{array}\right)\;,$ (12)

where $f_{n\kappa}(r)$ and $g_{n\kappa}(r)$ are respectively the small and
large radial components, and $S_{\kappa}$ is the sign of $\kappa$. The usual
spin-angular functions $\chi_{\kappa}^{\mu}$ are defined on an orthonormal
basis of spherical harmonics and are not of direct interest for the present
work.

To compute these wave functions, Behrens and Bühring express the effective
potential as a Taylor series

$V_{\text{eff}}=\sum\limits_{m=0}^{\infty}{v_{m}r^{m}}$ (13)

and introduce a local power series expansion that gives close to the origin

$\begin{Bmatrix}f_{n\kappa}(r)\\\
g_{n\kappa}(r)\end{Bmatrix}=\beta_{n\kappa}\dfrac{(p_{n\kappa}r)^{k-1}}{(2k-1)!!}\leavevmode\nobreak\
\leavevmode\nobreak\ \sum\limits_{j=0}^{\infty}{\begin{Bmatrix}a_{j}\\\
b_{j}\end{Bmatrix}r^{j}}\;.$ (14)

The series coefficients are defined by the following recurrence relationships:

$\begin{array}[]{lll}a_{0}=1\text{, }&b_{0}=0&\text{if }\kappa>0\\\
a_{0}=0\text{, }&b_{0}=1&\text{if }\kappa<0\\\ \end{array}$ (15)

and

$\begin{array}[]{l}(j+k-\kappa)a_{j}=-(W-1-v_{0})b_{j-1}+\sum\limits_{m=1}^{j-1}{v_{m}b_{n-m-1}}\\\
(j+k+\kappa)b_{j}=(W+1-v_{0})a_{j-1}-\sum\limits_{m=1}^{j-1}{v_{m}a_{n-m-1}}\\\
\end{array}$ (16)

with $W$ the total energy of the electron.

The effective KLI potential is locally fitted, after convergence, with a third
order polynomial. The determination of the Coulomb amplitude $\beta_{n\kappa}$
is then straightforward knowing the radial wave functions. In Table 2, we
compare some Coulomb amplitudes for different atoms determined with BetaShape,
Kotochigova and KLI models. As explained below, the differences are due to the
binding energies and the nucleus modeling.

We present in Figure 2 a comparison of the $1s$ and $2s$ wave functions in
neutral 7Be ($Z=4$) for the three different atomic models considered in this
work. As expected for low-$Z$ atoms, relativistic effects are very weak,
leading to $f\ll g$. We can also expect that the finite size of the nucleus
should play a minor role. It is indeed the case since the large component $g$
of the wave function exhibits the same behavior in the three models.
Concerning the small component $f$ (which is approximately 100 times smaller),
a different behavior can be seen between the wave function from Kotochigova
model and those from KLI and BetaShape. It shows up within the area of the
nucleus. This difference is due to the modeling of the nucleus charge
distribution in the three models: point charge for Kotochigova, sphere with
uniform charge density for BetaShape and Fermi distribution for KLI.

Table 2: Coulomb amplitudes $\beta_{n\kappa}$ for different atoms considered in the present work. The electronic configuration corresponds to neutral, fully relaxed atoms. |  | $\beta_{1s}$ |  |  | $\beta_{2s}$ |  |  | $\beta_{2p_{1/2}}$ |
---|---|---|---|---|---|---|---|---|---
Element | BS | Koto. | KLI | BS | Koto. | KLI | BS | Koto. | KLI
7Be | 0.008763 | 0.009073 | 0.009295 | 0.001584 | 0.001862 | 0.002110 | – | – | –
41Ca | 0.1174 | 0.1186 | 0.1174 | 0.03324 | 0.03449 | 0.03417 | 0.001823 | 0.001879 | 0.001929
54Mn | 0.1690 | 0.1742 | 0.1707 | 0.05263 | 0.05299 | 0.05200 | 0.003756 | 0.003743 | 0.003775
55Fe | 0.1806 | 0.1867 | 0.1825 | 0.05665 | 0.05721 | 0.05601 | 0.004225 | 0.004227 | 0.004249
109Cd | 0.5825 | 0.6414 | 0.5830 | 0.1955 | 0.2180 | 0.1983 | 0.02898 | 0.03214 | 0.02971
125I | 0.7224 | 0.8161 | 0.7248 | 0.2507 | 0.2829 | 0.2514 | 0.04140 | 0.04668 | 0.04213
138La | 0.8576 | 0.9864 | 0.8585 | 0.2985 | 0.3474 | 0.3025 | 0.05382 | 0.06231 | 0.05509

Similarly, we present in Figure 3 the same comparison for neutral 138La
($Z=57$). As expected for medium-$Z$ atoms, relativistic effects start to
become important leading to a non-negligible small component $f<g$. With
respect to the case of 7Be, the difference between the binding energies ($1s$
and $2s$) of Kotochigova (also used in BetaShape) and KLI is relatively much
smaller (see Table 1). Consequently, one also expect a smaller difference
between the wave functions, what can clearly be seen on Figure 3 for distances
greater than the nuclear radius ($1.6\times 10^{-2}$ rel. u.). For smaller
distances, the wave functions from Kotochigova model and those from KLI and
BetaShape differ, here again, due to the different modeling of the nuclear
charge distribution.

In order to see more differences between BetaShape and KLI models, let us now
investigate the higher $2p_{1/2}$ and $5s$ orbitals, for which wave functions
are presented in Figure 4. The analysis of $2p_{1/2}$ is similar to the one
given for $1s$ or $2s$. For the $5s$ orbital, relatively large differences can
be observed. This tendency is certainly due to the large difference (13%)
between the binding energies: 38.924 eV for Kotochigova and BetaShape, and
44.115 eV for KLI. The choice of $5s$ was motivated by the fact that this
outer $s$ orbital has a non-negligible probability of presence inside the
nucleus.

Figure 2:  (Color online) Large ($g$) and small components ($f$) of the radial
wave functions for the $1s$ and $2s$ orbitals of 7Be. Dash-dotted, dashed and
solid curves are from KLI, Kotochigova and BetaShape models, respectively.
Figure 3:  (Color online) Large ($g$) and small components ($f$) of the radial
wave functions for the $1s$ and $2s$ orbitals of 138La. Dash-dotted, dashed
and solid curves are from KLI, Kotochigova and BetaShape models, respectively.

###  III.3 Capture probabilities

Internal conversion theory has been recently tested comparing precise
measurements with different predictions [Hardy  _et al._ , 2018]. Two
hypotheses were considered: i) the “no vacancy” approximation, which assumes a
final state corresponding to the daughter atom after full atomic relaxation;
and ii) the “frozen orbital” (FO) approximation, which assumes that the final
atomic orbitals in the internal conversion process correspond to those of the
parent atom with a vacancy. It was found that FO approximation gives the most
accurate and consistent results over a wide range of atomic numbers, as
expected from the time scales of the physical processes [Hardy  _et al._ ,
2014].

The situation is comparable for electron capture. The weak interaction
transforms the parent nucleus with atomic number $Z$ into the daughter nucleus
with atomic number ($Z-1$). This interaction is mediated by massive bosons and
is therefore of very short range. The time scale of the electron capture decay
might thus be seen as instantaneous compared to the atomic time scale, i.e.
the vacancy lifetime. The BetaShape modeling assumes the FO approximation, and
we have considered with the KLI modeling both the FO and the “no vacancy”
approximations. We underline that in the latter case, the hole effect is
obviously not taken into account.

The dominant electron capture transition in each of the radionuclides of
interest for the European project MetroMMC have been studied. Calculation of
capture probabilities have been performed using the recommended Q-values
established in the latest Atomic Mass Evaluation AME2020 [Wang  _et al._ ,
2021]. Nuclear level energies have been taken from the latest ENSDF
evaluations for the following decays: 54Mn [Dong and Junde, 2014], 109Cd
[Kumar  _et al._ , 2016], 125I [Katakura, 2011] and 138La [Chen, 2017]. The
results are given in Table 3 together with experimental values when available
and are discussed below.

Figure 4:  (Color online) Large ($g$) and small components ($f$) of the radial
wave functions for the $2p_{1/2}$ and $5s$ orbital of 138La. Dash-dotted,
dashed and solid curves are from KLI, Kotochigova and BetaShape models,
respectively.

7Be. This decay has been measured a lot for decades. Its chemical form has
been demonstrated to have a significant influence on the half-life, which
seems reasonable with such a small number of atomic orbitals. The
$P_{L}/P_{K}$ ratio is thus expected to be influenced by environmental
effects. The experimental value is from a precise measurement in which 7Be was
implanted in a thin Ta metallic film [Fretwell  _et al._ , 2020]. It is clear
that none of the calculated values agrees and that the atomic modeling has a
very strong impact on the $P_{L}/P_{K}$ ratio, with a factor of 2 to 3 between
the predictions. It is noteworthy that such a difference appears only for this
radionuclide. The difference with the experimental value could be due in part
to the assumption of a decay in vacuum. In medium effects were estimated in
[Ray  _et al._ , 2002] and, if this simple modeling is correct, a strong
correction factor of 0.2986 is expected for metallic Ta. Applying this
correction to our predictions does not help to agree with the experimental
value.

41Ca. To our knowledge, no experimental value of any capture probability is
available in the literature for this radionuclide. We can simply observe that
the BetaShape value lies in between the “no vacancy” and the FO KLI results,
and that the FO approximation tends to decrease the $P_{L}/P_{K}$ ratio.

54Mn. The first experimental value comes from the weighted mean of two old
measurements with quite large uncertainties: $P_{L}/P_{K}=$ 0.098 (6) from
[Moler and Fink, 1963] and $P_{L}/P_{K}=$ 0.106 (3) from [Manduchi and
Zannoni, 1963]. Even if the “no vacancy” KLI result seems to be the best
prediction, it is difficult to be conclusive.

The experimental $P_{K}$ probability has been determined with measured values
listed in [Bambynek  _et al._ , 1977] and complemented with those in the DDEP
evaluation [Bé  _et al._ , 2004]. Only two measurements, consistent, do not
depend on the fluorescence yield $\omega_{K}$ and their weighted mean is
$P_{K}=$ 0.901 (6). Nine other measurements do depend on $\omega_{K}$. The
resulting dataset being discrepant, the uncertainty on their weighted mean has
been expanded in order to include the most precise value, leading to
$P_{K}\omega_{K}=$ 0.253 (6). It is interesting to note that the fluorescence
yield used in the DDEP evaluation from a semi-empirical fit from Bambynek
[Bambynek, 1984], $\omega_{K}=$ 0.289 (5), is significantly different to the
older experimental value from the same author [Bambynek, 1967], $\omega_{K}=$
0.2793 (17). With the former, one obtains a value of $P_{K}=$ 0.875 (26) while
with the latter, one has the higher value $P_{K}=$ 0.906 (22). Comparing with
the first result, it suggests that the correct fluorescence yield is the
experimental value, and we have chosen for our comparison the most precise
$P_{K}$ value. The FO KLI prediction is in very good agreement, which is
definitely not the case for the other two theoretical results.

55Fe. All three experimental capture probability ratios result from the
weighted means of two precise measurements [Pengra  _et al._ , 1972,Loidl  _et
al._ , 2018]. Surprisingly, the BetaShape predictions exhibit excellent
agreement for the $P_{L}/P_{K}$ and $P_{M}/P_{L}$ ratios while calculations
performed with the KLI atomic modeling not. It is noteworthy that the three
capture ratios are not consistent. Calculating $P_{M}/P_{K}$ from the weighted
means of the other two, one find a significantly higher value of
$P_{M}/P_{K}=$ 0.01813 (33), which is again in excellent agreement with the
BetaShape result but with none of the KLI predictions.

109Cd. The experimental values have been determined as weighted means of
several measurements listed in [Bé  _et al._ , 2016]. The $P_{LMNO}/P_{K}$
value has been reevaluated by revising one of the measured values as suggested
in [Bambynek  _et al._ , 1977]. Both BetaShape and FO KLI predictions are in
very good agreement with experiment, while the “no vacancy” KLI results are
significantly not in accordance.

125I. Seven measurements of the $P_{K}$ capture probability were performed in
the past, listed in [Bambynek  _et al._ , 1977,Bé  _et al._ , 2011]. The
experimental value in Table 3 results from their weighted mean. The value is
largely dominated by a single old measurement from [Leutz and Ziegler, 1964]
and only a slight dependency in the $\omega_{K}$ fluorescence yield is
expected. For this radionuclide, the “no vacancy” KLI result is the only one
being in good agreement, the other two models predicting a higher K capture
probability.

138La. This transition is the only one considered in the present study that is
not allowed, being second forbidden unique. In this case, $C_{n\kappa}$ is not
equal to unity and must be taken into account. All the measured values come
from a single recent experiment [Quarati  _et al._ , 2016]. The situation is
quite unclear as it depends on the capture probability ratio that is
considered. For $P_{L}/P_{K}$, the BetaShape result is in remarkable
accordance with the experimental value. It is the contrary for $P_{M}/P_{K}$,
KLI results being in excellent agreement with a slight advantage for the FO
approximation. For $P_{M}/P_{L}$, no theoretical prediction agrees with the
measurement.

However, as underlined in [Quarati  _et al._ , 2016] 138La decay seems to be
sensitive to its particular nuclear structure and additional terms in
$C_{n\kappa}$ might be necessary due to accidental cancellation of the
dominant nuclear matrix element. In [Mougeot, 2019], the nuclear levels were
described quite simply as single particle states and the BetaShape results
were found to be in better agreement with the measured values.

Table 3: Comparison of calculated and measured capture probabilities for the different isotopes considered in the present work. The three models and the experimental values are described in the text. Isotope | Energy (keV) | Quantity | BetaShape | KLI | KLI | Experimental
---|---|---|---|---|---|---
|  |  |  | “no vacancy” | “frozen orbital” |
7Be | 861.89 (7) | $P_{L}/P_{K}$ | 0.105 (8) | 0.1606 (41) | 0.0509 (20) | 0.070 (7)
41Ca | 421.64 (14) | $P_{L}/P_{K}$ | 0.09800 (40) | 0.10415 (16) | 0.09078 (16) | –
54Mn | 542.2 (10) | $P_{L}/P_{K}$ | 0.11219 (31) | 0.10785 (8) | 0.09590 (19) | 0.1044 (27)
|  | $P_{K}$ | 0.88419 (34) | 0.88623 (10) | 0.90005 (21) | 0.901 (6)
55Fe | 231.12 (18) | $P_{L}/P_{K}$ | 0.11629 (31) | 0.11236 (8) | 0.10073 (20) | 0.1165 (9)
|  | $P_{M}/P_{K}$ | 0.01824 (12) | 0.019390 (32) | 0.014824 (45) | 0.01786 (29)
|  | $P_{M}/P_{L}$ | 0.1568 (11) | 0.17257 (31) | 0.14716 (49) | 0.1556 (26)
109Cd | 127.1 (18) | $P_{K}$ | 0.8148 (14) | 0.8097 (11) | 0.8164 (12) | 0.815 (2)
|  | $P_{LMNO}/P_{K}$ | 0.2274 (12) | 0.2350 (11) | 0.2250 (12) | 0.2279 (21)
125I | 150.28 (6) | $P_{K}$ | 0.79927 (41) | 0.79798 (7) | 0.80376 (23) | 0.7971 (14)
138La | 312.60 (30) | $P_{L}/P_{K}$ | 0.3913 (25) | 0.4077 (15) | 0.4242 (49) | 0.391 (3)
|  | $P_{M}/P_{K}$ | 0.0965 (9) | 0.09908 (41) | 0.1002 (11) | 0.102 (3)
|  | $P_{M}/P_{L}$ | 0.2465 (20) | 0.2430 (22) | 0.2362 (24) | 0.261 (9)

##  IV Conclusions and perspectives

For the past few years, the DDEP international collaboration has employed the
BetaShape code in their evaluated nuclear decay data in order to improve the
accuracy of the beta and electron capture properties. In the present work, we
have performed an exhaustive study of the role played by the modeling of the
atomic electrons on the electron capture process.

A realistic model has been developed within the framework of RDFT, in which we
have implemented few exchange-correlation functionals and self-interaction-
corrected models, among the most popular in the community. From our analysis,
we found that the best choice which fulfills some physical constraints
(binding energies of inner subshells, asymptotic behavior of the effective
potential) is provided by the optimized effective potential method originally
developed by Krieger, Li and Iafrate (KLI). It turns out that, within this
model, correlation effects seem to play a minor role. However, a more precise
study would be required in order to obtain a more decisive conclusion on that
issue.

The electron capture model has been adapted in order to consistently use the
binding energies, wave functions and electronic configurations from the KLI
code. It allowed us to precisely include the hole effect on the other
orbitals, due to the vacancy created by the capture process, while this effect
is accounted for approximately with the BetaShape atomic modeling. We have
demonstrated a posteriori why the latter modeling can lead to surprisingly
good results. When Kotochigova and KLI binding energies are relatively close
to each other, i.e. for the innermost shells, introducing a nucleus spatial
extension instead of considering a point charge gives the dominant distortion
on the radial wave functions close to the origin. However, it also seems that
some inaccuracies (hole and shaking effects, which depend on the different
possible overlaps) somehow compensate each other.

Several transitions of interest have been calculated and the theoretical
capture probabilities have been compared to experimental values. Such a test
covers a wide range of atomic numbers, $3\leq Z\leq 57$. As in the early stage
of the study of internal conversion [Hardy  _et al._ , 2014], the situation is
unclear and unsatisfactory and we face to similar difficulties. The KLI model
with the “frozen orbital” approximation should in principle provide the most
accurate predictions. Our analysis of the published experimental values has
shown that new high-precision measurements are definitely required. It is also
necessary to have a more complete set of measured capture probabilities for
each transition because one or two are not sufficient to be conclusive,
especially for the outer shells. One can expect that these new accurate
experimental values and further improvements of the theoretical models will
clarify the situation, as it did with internal conversion [Hardy  _et al._ ,
2018].

In the near future, the electron capture model should benefit from ongoing
developments on beta decay that aim at including a realistic nuclear structure
for the calculation of forbidden non-unique transitions. Such an improvement
will necessitate a revision of the overlap and exchange correction, with a
precise coupling of lepton and nucleon wave functions.

Within the framework of the KLI atomic model, efforts will be devoted to the
separate computation of shake-up and shake-off processes. Indeed, good
description of the unoccupied levels and the continuum states is mandatory for
such processes, shake-off being in particular conditioned by the asymptotic
behavior of the effective potential.

A more precise study devoted to the role played by the correlation effects on
electron capture probabilities is also required. A comparison between the
present results and those obtained with the MCDF method for light atoms (e.g.
7Be and 41Ca) will be investigated, what could lead to a revised
parameterization of the correlation potential.

Finally, another important issue concerns the role played by the environment
on the electron capture process, as it has been seen e.g. in 7Be decay. Such a
complicated problem deserves in itself a detailed specific study.

###### Acknowledgements.

This work was performed within the European project MetroMMC, “Measurement of
fundamental nuclear decay data using metallic magnetic calorimeters”. This
project 17FUN02 MetroMMC has received funding from the EMPIR programme co-
financed by the Participating States and from the European Union’s Horizon
2020 research and innovation programme.

## References

  * Coulon  _et al._ (2020) R. Coulon, R. Broda, P. Cassette, S. Courte, S. Jerome, S. Judge, K. Kossert, H. Liu, C. Michotte, and M. Nonis, Metrologia 57, 035009 (2020).
  * Broda  _et al._ (2007) R. Broda, P. Cassette, and K. Kossert, Metrologia 44, S36 (2007).
  * Bavelaar  _et al._ (2018) B. M. Bavelaar, B. Q. Lee, M. R. Gill, N. Falzone, and K. A. Vallis, Frontiers in Pharmacology 9, 996 (2018).
  * Ku  _et al._ (2019) A. Ku, V. J. Facca, Z. Cai, and R. M. Reilly, EJNMMI Radiopharmacy and Chemistry 4, 27 (2019).
  * Hou (2005) X. Hou, Radiochimica Acta 93, 611 (2005).
  * Lee  _et al._ (2021) Y.-J. Lee, J.-M. Lim, J.-H. Lee, S.-B. Hong, and H. Kim, Nuclear Engineering and Technology 53, 1210 (2021).
  * Huss  _et al._ (2009) G. R. Huss, B. S. Meyer, G. Srinivasan, J. N. Goswami, and S. Sahijpal, Geochimica et Cosmochimica Acta 73, 4922 (2009).
  * Jörg  _et al._ (2012) G. Jörg, Y. Amelin, K. Kossert, and C. Lierse v. Gostomski, Geochimica et Cosmochimica Acta 88, 51 (2012).
  * Kubik  _et al._ (1986) P. W. Kubik, D. Elmore, N. J. Conard, K. Nishiizumi, and J. R. Arnold, Nature 319, 568 (1986).
  * Bhat (1992) M. R. Bhat, in  _Nuclear Data for Science and Technology_, edited by S. M. Qaim (Springer, Berlin, Heidelberg, 1992) pp. 817–821.
  * Kellett, Mark A. and Bersillon, Olivier (2017) Kellett, Mark A. and Bersillon, Olivier, EPJ Web Conf. 146, 02009 (2017).
  * Gove and Martin (1971) N. Gove and M. Martin, Atomic Data and Nuclear Data Tables 10, 205 (1971).
  * Mougeot (2015) X. Mougeot, Phys. Rev. C 91, 055504 (2015).
  * Mougeot (2019) X. Mougeot, Appl. Radiat. Isot. 154, 108884 (2019).
  * Mougeot (2018) X. Mougeot, Appl. Radiat. Isot. 134, 225 (2018).
  * Ranitzsch  _et al._ (2020) P. C.-O. Ranitzsch, D. Arnold, J. Beyer, L. Bockhorn, J. J. Bonaparte, C. Enss, K. Kossert, S. Kempf, M. Loidl, R. Mariam, O. J. Nähle, M. Paulsen, M. Rodrigues, and M. Wegner, Journal of Low Temperature Physics 199, 441 (2020).
  * Bambynek  _et al._ (1977) W. Bambynek, H. Behrens, M. H. Chen, B. Crasemann, M. L. Fitzpatrick, K. W. D. Ledingham, H. Genz, M. Mutterer, and R. L. Intemann, Rev. Mod. Phys. 49, 77 (1977).
  * Behrens and Bühring (1982) H. Behrens and W. Bühring,  _Electron Radial Wave Functions and Nuclear Beta Decay_ (Clarendon, Oxford, 1982).
  * Bahcall (1965) J. N. Bahcall, Nuclear Physics 71, 267 (1965).
  * Vatai (1970) E. Vatai, Nuclear Physics A 156, 541 (1970).
  * Crasemann  _et al._ (1979) B. Crasemann, M. H. Chen, J. P. Briand, P. Chevallier, A. Chetioui, and M. Tavernier, Phys. Rev. C 19, 1042 (1979).
  * Mougeot and Bisch (2014) X. Mougeot and C. Bisch, Phys. Rev. A 90, 012501 (2014).
  * Salvat  _et al._ (1987) F. Salvat, J. D. Martínez, R. Mayol, and J. Parellada, Phys. Rev. A 36, 467 (1987).
  * Desclaux (1973) J. Desclaux, Atomic Data and Nuclear Data Tables 12, 311 (1973).
  * Kotochigova  _et al._ (1997) S. Kotochigova, Z. H. Levine, E. L. Shirley, M. D. Stiles, and C. W. Clark, Phys. Rev. A 55, 191 (1997).
  * Kotochigova  _et al._ (2021) S. Kotochigova, Z. H. Levine, E. L. Shirley, M. D. Stiles, and C. W. Clark,  _NIST Chemistry WebBook, NIST Standard Reference Database 141_, edited by P. Linstrom and W. Mallard (National Institute of Standards and Technology, Gaithersburg MD, 20899, 2021) Chap. Atomic Reference Data for Electronic Structure Calculations.
  * Desclaux (1975) J. Desclaux, Computer Physics Communications 9, 31 (1975).
  * Dreizler and Engel (1998) R. M. Dreizler and E. Engel, “Relativistic density functional theory,” in  _Density Functionals: Theory and Applications: Proceedings of the Tenth Chris Engelbrecht Summer School in Theoretical Physics Held at Meerensee, near Cape Town South Africa, 19–29 January 1997_, edited by D. Joubert (Springer, Berlin, Heidelberg, 1998) pp. 147–189.
  * Rajagopal (1978) A. K. Rajagopal, Journal of Physics C: Solid State Physics 11, L943 (1978).
  * MacDonald and Vosko (1979) A. H. MacDonald and S. H. Vosko, Journal of Physics C: Solid State Physics 12, 2977 (1979).
  * Engel (2002) E. Engel,  _Relativistic Electronic Structure Theory, Part 1. Fundamentals_, edited by P. Schwerdtfeger (Elsevier, Amsterdam, 2002) pp. 524–624.
  * Strange (2008) P. Strange,  _Relativistic Quantum Mechanics: With Applications in Condensed Matter and Atomic Physics_ (Cambridge University Press, 2008).
  * Pálffy (2010) A. Pálffy, Contemporary Physics 51, 471 (2010).
  * Angeli and Marinova (2013) I. Angeli and K. Marinova, Atomic Data and Nuclear Data Tables 99, 69 (2013).
  * Johnson and Soff (1985) W. Johnson and G. Soff, Atomic Data and Nuclear Data Tables 33, 405 (1985).
  * Vosko  _et al._ (1980) S. H. Vosko, L. Wilk, and M. Nusair, Canadian Journal of Physics 58, 1200 (1980).
  * Gunnarsson and Lundqvist (1976) O. Gunnarsson and B. I. Lundqvist, Phys. Rev. B 13, 4274 (1976).
  * Perdew and Wang (1992) J. P. Perdew and Y. Wang, Phys. Rev. B 45, 13244 (1992).
  * van Leeuwen and Baerends (1994) R. van Leeuwen and E. J. Baerends, Phys. Rev. A 49, 2421 (1994).
  * Čertík  _et al._ (2013) O. Čertík, J. E. Pask, and J. Vackář, Computer Physics Communications 184, 1777 (2013).
  * Ciofini  _et al._ (2005) I. Ciofini, C. Adamo, and H. Chermette, Chemical Physics 309, 67 (2005), electronic Structure Calculations for Understanding Surfaces and Molecules. In Honor of Notker Roesch.
  * Perdew and Zunger (1981) J. P. Perdew and A. Zunger, Phys. Rev. B 23, 5048 (1981).
  * Politis  _et al._ (1998) M. F. Politis, P. A. Hervieux, J. Hanssen, M. E. Madjet, and F. Martín, Phys. Rev. A 58, 367 (1998).
  * Krieger  _et al._ (1992) J. B. Krieger, Y. Li, and G. J. Iafrate, Phys. Rev. A 46, 5453 (1992).
  * Li  _et al._ (1993) Y. Li, J. B. Krieger, and G. J. Iafrate, Phys. Rev. A 47, 165 (1993).
  * Tong and Chu (1998) X.-M. Tong and S.-I. Chu, Phys. Rev. A 57, 855 (1998).
  * Sevier (1979) K. D. Sevier, Atomic Data and Nuclear Data Tables 24, 323 (1979).
  * Bearden and Burr (1967) J. A. Bearden and A. F. Burr, Rev. Mod. Phys. 39, 125 (1967).
  * Hardy  _et al._ (2018) J. Hardy, N. Nica, V. Iacob, and M. Trzhaskovskaya, Applied Radiation and Isotopes 134, 406 (2018).
  * Hardy  _et al._ (2014) J. Hardy, N. Nica, V. Iacob, S. Miller, M. Maguire, and M. Trzhaskovskaya, Applied Radiation and Isotopes 87, 87 (2014).
  * Wang  _et al._ (2021) M. Wang, W. Huang, F. Kondev, G. Audi, and S. Naimi, Chinese Physics C 45, 030003 (2021).
  * Dong and Junde (2014) Y. Dong and H. Junde, Nuclear Data Sheets 121, 1 (2014).
  * Kumar  _et al._ (2016) S. Kumar, J. Chen, and F. Kondev, Nuclear Data Sheets 137, 1 (2016).
  * Katakura (2011) J. Katakura, Nuclear Data Sheets 112, 495 (2011).
  * Chen (2017) J. Chen, Nuclear Data Sheets 146, 1 (2017).
  * Fretwell  _et al._ (2020) S. Fretwell, K. G. Leach, C. Bray, G. B. Kim, J. Dilling, A. Lennarz, X. Mougeot, F. Ponce, C. Ruiz, J. Stackhouse, and S. Friedrich, Phys. Rev. Lett. 125, 032701 (2020).
  * Ray  _et al._ (2002) A. Ray, P. Das, S. K. Saha, S. K. Das, and A. Mookerjee, Phys. Rev. C 66, 012501 (2002).
  * Moler and Fink (1963) R. B. Moler and R. W. Fink, Phys. Rev. 131, 821 (1963).
  * Manduchi and Zannoni (1963) G. Manduchi and G. Zannoni, Il Nuovo Cimento 27, 251 (1963).
  * Bé  _et al._ (2004) M.-M. Bé, V. Chisté, C. Dulieu, E. Browne, V. Chechev, N. Kuzmenko, R. Helmer, A. Nichols, E. Schönfeld, and R. Dersch,  _Table of Radionuclides_, Monographie BIPM-5, Vol. 1 (Bureau International des Poids et Mesures, Pavillon de Breteuil, F-92310 Sèvres, France, 2004) ; Evaluation update in 2014\.
  * Bambynek (1984) W. Bambynek, in  _X-84 Proc. X-ray and Inner-Shell Processes in Atoms, Molecules and Solids, Leipzig Aug. 20-23, 1984_ , edited by A. Meisel (VEB Druckerei, Thomas Münzer, Langensalza, 1984) post-deadline paper P-1.
  * Bambynek (1967) W. Bambynek, Zeitschrift für Physik 206, 66 (1967).
  * Pengra  _et al._ (1972) J. G. Pengra, H. Genz, J. P. Renier, and R. W. Fink, Phys. Rev. C 5, 2007 (1972).
  * Loidl  _et al._ (2018) M. Loidl, M. Rodrigues, and R. Mariam, Applied Radiation and Isotopes 134, 395 (2018).
  * Bé  _et al._ (2016) M.-M. Bé, V. Chisté, C. Dulieu, M. Kellett, X. Mougeot, A. Arinc, V. Chechev, N. Kuzmenko, T. Kibédi, A. Luca, and A. Nichols,  _Table of Radionuclides_, Monographie BIPM-5, Vol. 8 (Bureau International des Poids et Mesures, Pavillon de Breteuil, F-92310 Sèvres, France, 2016).
  * Bé  _et al._ (2011) M.-M. Bé, V. Chisté, C. Dulieu, X. Mougeot, V. Chechev, N. Kuzmenko, F. Kondev, A. Luca, M. Galán, A. Nichols, A. Arinc, A. Pearce, X. Huang, and B. Wang,  _Table of Radionuclides_, Monographie BIPM-5, Vol. 6 (Bureau International des Poids et Mesures, Pavillon de Breteuil, F-92310 Sèvres, France, 2011).
  * Leutz and Ziegler (1964) H. Leutz and K. Ziegler, Nuclear Physics 50, 648 (1964).
  * Quarati  _et al._ (2016) F. Quarati, P. Dorenbos, and X. Mougeot, Applied Radiation and Isotopes 109, 172 (2016).


11institutetext: Institute for Theoretical Physics, Technische UniversitÃt
Wien, 1040 Vienna, Austria 22institutetext: Department of Physics, P.O. Box
35, 40014 University of Jyväskylä, Finland 33institutetext: Helsinki Institute
of Physics, P.O. Box 64, 00014 University of Helsinki, Finland
44institutetext: FakultÃt fÃr Physik, UniversitÃt Bielefeld, D-33615
Bielefeld, Germany

# Fermion and gluon spectral functions far from equilibrium

Kirill Boguslavski Speaker, .
K.B. would like to acknowledge support by the Austrian Science Fund (FWF)
under project P 34455-N.
The authors are grateful to A. Kurkela, M. Mace and J. Peuron for
collaboration and discussions. The authors wish to acknowledge the National
Energy Research Scientific Computing Center, a U.S. Department of Energy
Office of Science User Facility supported under Contract No. DE-
AC02-05CH11231, the CSC - IT Center for Science, Finland, and the Vienna
Scientific Cluster (VSC, Austria) under the project number 71444 for
computational resources. 11<EMAIL_ADDRESS>Tuomas Lappi  T.L.
is supported by the Academy of Finland, project 321840 and under the European
Unionâs Horizon 2020 research and innovation programme by the European
Research Council (ERC, grant agreement No. ERC-2015-CoG-681707) and by the
STRONG-2020 project (grant agreement No 824093). The content of this article
does not reflect the official opinion of the European Union and responsibility
for the information and views expressed therein lies entirely with the
authors. 2233 Sören Schlichting  S.S. is supported by the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) through the CRC-TR
211 “Strong-interaction matter under extreme conditions” Project number
315477589. 44

###### Abstract

Motivated by the quark-gluon plasma, we develop a simulation method to obtain
the spectral function of (Wilson) fermions non-perturbatively in a non-Abelian
gauge theory with large gluon occupation numbers Boguslavski:2021kdd . We
apply our method to a non-Abelian plasma close to its non-thermal fixed point,
i.e., in a far-from-equilibrium self-similar regime, and find mostly very good
agreement with perturbative hard loop (HTL) calculations. For the first time,
we extract the full momentum dependence of the damping rate of fermionic
collective excitations and compare our results to recent non-perturbative
extractions of gluonic spectral functions in two and three spatial dimensions
Boguslavski:2021buh ; Boguslavski:2018beu .

##  1 Introduction

Nonperturbatively strong bosonic fields can be encountered in a variety of
systems, including heavy-ion collisions at early times Schlichting:2019abc ;
Berges:2020fwq , cosmological reheating Allahverdi:2010xz , or the infrared
sector of gluonic and scalar fields in thermal equilibrium. Understanding
their interactions with fermions is important for phenomenological
applications including jets (energy loss, jet quenching) or electromagnetic
observables, to mention a few in heavy-ion collisions.

For a microscopic description of the dynamics, the calculation of spectral
functions $\rho(\omega,p)$ is very useful, since they encode medium
interactions and the spectrum of collective excitations.

In this talk we will present a new method to compute the spectral function of
Wilson fermions non-perturbatively in a weakly-coupled non-Abelian gauge
theory with large gluon occupation numbers $f\gg 1$ (Sec. 2). In this far-
from-equilibrium situation, we can use classical-statistical lattice
simulations to describe the underlying dynamics Berges:2004yj ; Berges:2013lsa
. Coupling them to the Dirac equation allows us to extract the fermion
spectral function nonperturbatively. Our results are then discussed and
compared to perturbative hard-loop (HTL) computations (Sec. 3) and to gluonic
spectral functions in different dimensions (Sec. 4, for the comparison see
Sec. 5). We conclude in Sec. 5. More details on our method and results for the
fermion spectral function can be found in Ref. Boguslavski:2021kdd . Our
discussion of the gluon spectral function results is based on Refs.
Boguslavski:2021buh ; Boguslavski:2018beu .

##  2 Method & setup

###  2.1 Theory

We consider $\mathrm{SU}(N_{c})$ gauge theory on a cubic lattice of size
$N_{s}^{3}$ and lattice spacing $a_{s}$. The Yang-Mills and Dirac Hamiltonians
read in $A_{0}=0$ gauge

$\displaystyle H_{\text{Y}M}=$
$\displaystyle\frac{1}{g^{2}a_{s}}\sum_{\mathbf{x},i}\text{Tr}[E_{i}(t^{\prime},\mathbf{x})^{2}]+\frac{1}{2}\sum_{j}\text{ReTr}[1-U_{ij}(t^{\prime},\mathbf{x})]$
$\displaystyle\hat{H}_{W}=$
$\displaystyle\frac{1}{2}\sum_{\mathbf{x}}\left[\hat{\psi}^{\dagger}(t^{\prime},\mathbf{x}),\gamma^{0}\left(-i\not{D}_{s}[U]+m\right)\hat{\psi}(t^{\prime},\mathbf{x})\right],$
(1)

with the link fields $U_{j}\approx\exp\left(ig\,a_{s}A_{j}\right)$ replacing
gauge fields, plaquettes $U_{ij}$ computed from link fields, and
$-i\not{D}_{s}[U]\hat{\psi}$ being the tree-level improved Wilson Dirac
operator Mace:2016shq ; Mace:2019cqo . We neglect the backreaction of fermions
on the evolution of gauge fields since it is subleading when gluonic
occupation numbers are large, as in our case. We will consider $N_{c}=2$ in
the following to reduce the computational complexity of the calculations.

The fermionic fields are formulated in terms of a mode expansion Aarts:1998td

$\hat{\psi}(t^{\prime},\mathbf{x})=\frac{1}{\sqrt{V}}\sum_{\lambda,\mathbf{p}}\hat{b}_{\lambda,\mathbf{p}}(t)\,\phi_{\lambda,\mathbf{p}}^{u}(t^{\prime},\mathbf{x})+\hat{d}^{\dagger}_{\lambda,\mathbf{p}}(t)\,\phi_{\lambda,\mathbf{p}}^{v}(t^{\prime},\mathbf{x}),$
(2)

where $\lambda=1,\cdots,2N_{c}$ labels spin and color indices. At the
reference time $t^{\prime}=t$, where $t$ will be the time when we extract the
spectral function, the creation and annihilation operators $\hat{b}$,
$\hat{d}$ satisfy the usual anti-commutation relations at $t$ and the wave
functions are set to plane waves

$\left.\phi_{\lambda,\mathbf{p}}^{u}(t^{\prime},\mathbf{x})\right|_{t^{\prime}{=}t}=u_{\lambda}(\mathbf{p})e^{+i\mathbf{p}\cdot\mathbf{x}}\,,\quad\left.\phi_{\lambda,\mathbf{p}}^{v}(t^{\prime},\mathbf{x})\right|_{t^{\prime}{=}t}=v_{\lambda}(\mathbf{p})e^{-i\mathbf{p}\cdot\mathbf{x}}\,.$
(3)

###  2.2 Fermion spectral function

The spectral function of fermions is defined as

$\rho^{\alpha\beta}(x,y)=\left\langle\left\\{\hat{\psi}^{\alpha}(t^{\prime},\mathbf{x}),\hat{\bar{\psi}}^{\beta}(t,\mathbf{y})\right\\}\right\rangle,$
(4)

where we have written the Dirac indices $\alpha$, $\beta$ explicitly. In our
classical-statistical framework, the spectral function can be computed
efficiently in momentum space. Plugging in the mode expansion (2) into (4),
Fourier transforming with respect to $\mathbf{x}-\mathbf{y}$, and using (3)
for time $t$, one arrives at

$\displaystyle\rho^{\alpha\beta}(t^{\prime},t,\mathbf{p})$ (5)
$\displaystyle=\frac{1}{V}\sum_{\lambda,\mathbf{q}}\left\langle\tilde{\phi}^{u,\alpha}_{\lambda,\mathbf{q}}(t^{\prime},\mathbf{p})\left(\tilde{\phi}^{u,\gamma}_{\lambda,\mathbf{q}}(t,\mathbf{p})\right)^{*}+\tilde{\phi}^{v,\alpha}_{\lambda,\mathbf{q}}(t^{\prime},\mathbf{p})\left(\tilde{\phi}^{v,\gamma}_{\lambda,\mathbf{q}}(t,\mathbf{p})\right)^{*}\right\rangle\gamma_{0}^{\gamma\beta}$
$\displaystyle=\frac{1}{V}\sum_{\lambda}\left\langle\tilde{\phi}^{u,\alpha}_{\lambda,\mathbf{p}}(t^{\prime},\mathbf{p})\,u_{\lambda}^{\dagger,\gamma}(\mathbf{p})+\tilde{\phi}^{v,\alpha}_{\lambda,-\mathbf{p}}(t^{\prime},\mathbf{p})\,v_{\lambda}^{\dagger,\gamma}(-\mathbf{p})\right\rangle\gamma_{0}^{\gamma\beta}\,,$

where $\langle.\rangle$ denotes a classical-statistical average over gluonic
configurations. The initialization and the specific form of the initial
conditions in Eq. (3), where $\phi_{\lambda,\mathbf{q}}^{u/v}(t,\mathbf{x})$
is given by a plane wave, allows to simplify the expression in the last line
above.

###  2.3 Simulation algorithm

The simulation algorithm can be summarized as follows:

  1. 1.

Set initial conditions for gluons at $t^{\prime}=0$, generating a
configuration with $\langle
E_{T}^{*}(t^{\prime}{=}0,\mathbf{p})E_{T}(t^{\prime}{=}0,\mathbf{q})\rangle\propto
p\,f(t^{\prime}{=}0,p)(2\pi)^{3}\delta(\mathbf{p}-\mathbf{q})$.

  2. 2.

Solve classical equations of motion (EOMs) for gauge fields for $0\leq
t^{\prime}\leq t$, set Coulomb-type gauge
$\left.\partial^{j}A_{j}\right|_{t}=0$ at $t^{\prime}=t$, i.e., at the time
when the spectral function is measured.

  3. 3.

For each momentum mode $\mathbf{p}$ initialize
$\phi_{\lambda,\mathbf{p}}^{u/v}$ at $t^{\prime}=t$ using Eq. (3).

  4. 4.

Evolve gauge and fermionic fields for $t^{\prime}>t$ using a leap-frog scheme
to solve the classical EOMs for gauge and the Dirac equation for fermionic
fields.

  5. 5.

Calculate the fermionic spectral function $\rho(t^{\prime},t,\mathbf{p})$ for
each momentum mode $\mathbf{p}$ according to Eq. (5).

##  3 Nonperturbatively computed fermion spectral functions

###  3.1 Non-equilibrium state: self-similar turbulent attractor

We now consider spectral functions in a state of high gluonic occupancy
$f_{g}(t{=}0,p\lesssim Q)\sim\frac{1}{g^{2}}\gg 1$. Such a system approaches a
self-similar attractor with

$\displaystyle f(t,p)=(Qt)^{\alpha}f_{s}\left((Qt)^{\beta}p\right),$ (6)

after some short transient time. The scaling exponents and scaling function
then become insensitive to details of initial conditions and their universal
values in Minkowski space-time are given by $\beta=-1/7$ and $\alpha=4\beta$
for a 3+1 dimensional Kurkela:2012hp ; Schlichting:2012es ; Berges:2013fga
and by $\beta=-1/5$ and $\alpha=3\beta$ for a 2+1D plasma Boguslavski:2019fsb
. Since we know the time-evolution of characteristic momentum scales at those
attractor states, like the hard scale $\Lambda\sim Q(Qt)^{-\beta}$ and the
soft mass scale $m\sim Q(Qt)^{\beta}$, these states are useful assets to study
spectral functions in more detail.

Figure 1: Left: Components of the spectral function $\rho(t{+}\Delta t,t,p)$
at different momenta. Right: The spectral function
$\rho_{+}(t,\omega,p)\equiv\rho_{V}^{0}+\rho_{V}$ in frequency space for the
same momenta. The values at negative frequency stem from
$\rho_{+}(t,-\omega,p)=\rho_{-}(t,\omega,p)\equiv\rho_{V}^{0}-\rho_{V}$.
Figures taken from Ref. Boguslavski:2021kdd .

###  3.2 Fermion $\rho$ in 3+1D at an attractor

We run our simulations on a $256^{3}$ lattice with $Qa_{s}=0.75$ with almost
vanishing fermion mass $m=0.003125\,Q$ to avoid ambiguities of the Wilson-
Dirac equation. The fermion spectral function is extracted at the reference
time $Qt=1500$, i.e., sufficiently late to make sure that we are close to the
self-similar regime of Eq. (6). The (nonvanishing) vector components of the
spectral function

$\rho^{0}_{V}=\frac{1}{4}\text{Tr}(\rho\gamma^{0})\,,\quad\rho_{V}=-\frac{E_{\mathbf{p}}\,p^{j}}{4\,p^{2}}\,\text{Tr}(\rho\gamma^{j})\,,$
(7)

are shown in the left panel of Fig. 1 as functions of $\Delta t\equiv
t^{\prime}-t$ for different momenta. One observes that they exhibit damped
oscillations with a well defined frequency (dispersion relation) and damping
rate.

This can be made more quantitative after Fourier transforming it with respect
to $\Delta t$ into the frequency domain, as shown in the right panel of Fig.
1. One finds two quasiparticle peaks at $\omega=\pm\omega_{\pm}(p)$ with
respective width (damping rates) $\gamma_{\pm}(p)$. To extract these values,
we recall that the spectral function can be computed perturbatively in HTL as
Blaizot:2001nr

$\displaystyle\rho_{+}^{\rm HTL}(\omega,p)=2\pi\,\beta_{+}(\omega/p,p)$
$\displaystyle\quad+2\pi\left[Z_{+}(p)\delta(\omega-\omega_{+}(p))+Z_{-}(p)\delta(\omega+\omega_{-}(p))\right],$
(8)

with the Landau damping term $\beta_{+}$. We extend this expression to allow
for a finite peak width, denoting it HTL+$\gamma$,

$\displaystyle\rho_{+}^{{\rm
HTL}+\gamma}(\omega,p)=2\pi\,\beta_{+}(\omega/p,p)$ (9)
$\displaystyle\quad+\frac{2Z_{+}(p)\gamma_{+}(p)}{(\omega-\omega_{+}(p))^{2}+\gamma_{+}^{2}(p)}+\frac{2Z_{-}(p)\gamma_{-}(p)}{(\omega+\omega_{-}(p))^{2}+\gamma_{-}^{2}(p)},$

and show it in the right panel of Fig. 1, together with its Landau damping and
quasiparticle contributions separately. One finds very good agreement with our
data. The good agreement with our data thus allows us to extract
$\omega_{\pm}(t,p)$, $Z_{\pm}(t,p)$ and $\gamma_{\pm}(t,p)$ by fitting to (9).

Figure 2: Extracted dispersion relations $\omega_{\pm}(t,p)$, residues
$Z_{\pm}(t,p)$ and damping rates $\gamma_{\pm}(t,p)$ as functions of $p$ by
fitting to (9). Figures taken from Ref. Boguslavski:2021kdd .

These extracted quantities are shown in Fig. 2, together with results from a
smaller and coarser $64^{3}$ lattice with $Qa_{s}=1$ for comparison. The HTL
predictions for the dispersion relations $\omega_{\pm}^{\mathrm{HTL}}(p)$ and
residues $Z_{\pm}^{\mathrm{HTL}}(p)$ are included as dashed lines. Note that
they are completely determined after calculating the HTL fermion mass

$m_{F}=\left[C_{F}\int\dfrac{\mathrm{d}^{3}p}{(2\pi)^{3}}\,\dfrac{g^{2}f_{g}(p)}{p}\right]^{1/2}$
(10)

whose value is indicated by red arrows. Hence, the nice agreement between our
simulation results and leading order HTL expressions shows that these
quantities are dominated by perturbative physics.

In addition, our simulations allow for a first-principles insight into the
momentum dependence of damping rates $\gamma_{\pm}(p)$, which is hard within a
perturbative treatment. In general, we find that $\gamma_{+}(p)$ decrease with
momentum and that its values are much smaller than the respective dispersion
$\omega_{+}(p)$, validating a quasiparticle picture for (almost) all momentum
modes.

Figure 3: Left: Fermion spectral function $\rho_{+}(t,\omega,p{=}0)$ at
different times in terms of a) the constant scale $Q$ and b) the time-
dependent fermion mass $m_{F}(t)$. Right: A summary of our results: $\rho_{+}$
as a function of frequency and momentum. Figure taken from Ref.
Boguslavski:2021kdd .

From HTL Braaten:1992gd  one would expect the damping rate to scale with the
effective temperature $\gamma^{\mathrm{HTL}}(t,p{=}0)\propto g^{2}T_{*}(t)\sim
Q(Qt)^{-3/7}$. Instead, we find that it rather scales with the fermion mass
$\gamma(t,p{=}0)\sim m_{F}(t)\sim Q(Qt)^{-1/7}$, as we show in the left panel
of Fig. 3. By rescaling all dimensionful quantities with
$m_{F}(t)\equiv\omega_{\pm}(t,p{=}0)$, all curves corresponding to different
(reference) times $t$ collapse on a single curve. Note that this also confirms
a self-similar evolution of the spectral function at the attractor.

##  4 Gluonic spectral functions

We will now discuss the gluon spectral function
$\rho(x^{\prime},x)=i\left\langle\left[\hat{A}(x^{\prime}),\hat{A}(x)\right]\right\rangle$
in different dimensions at corresponding self-similar attractors. Its
measurement follows a very similar algorithm to the one for fermions and
combines classical-statistical simulations with linear response theory. We
refer to Refs. Boguslavski:2018beu ; Kurkela:2016mhu ; Boguslavski:2021buh
for more details. Closely related algorithms were also applied to scalar
theories at self-similar attractors in PineiroOrioli:2018hst ;
Boguslavski:2019ecc .

We will discuss the color averaged transversely polarized (dotted) gluonic
spectral functions
$\dot{\rho}_{T}(t^{\prime},t)=\partial_{t^{\prime}}\rho_{T}(t^{\prime},t)$.
One obtains $\dot{\rho}_{T}(t,\omega)\approx\omega\rho_{T}$ by Fourier
transforming in $\Delta t\equiv t^{\prime}-t$.

Apart from the spectral function one can also measure the statistical
correlation function $\langle
E(t^{\prime})E(t)\rangle\equiv\ddot{F}(t^{\prime},t)$. While in thermal
equilibrium it is related to the spectral function via the fluctuation-
dissipation relation (FDR), far from equilibrium it was found numerically in
Refs. Boguslavski:2018beu ; Boguslavski:2021buh  that statistical and spectral
correlation functions are related via a generalized FDR (also seen in the left
panels of Figs. 4 and 5).

Figure 4: Gluonic simulation results in a 3+1D plasma. Left: The spectral
function at different momentum modes as a function of $\omega$. In contrast to
the right panel of Fig. 1, we show only the positive frequency part for the
gluon spectral function, which agrees with the values at negative frequency
due to $\dot{\rho}_{T}(t,-\omega,p)=\dot{\rho}_{T}(t,\omega,p)$. Right: The
extracted damping rate at different times as a function of momentum, rescaled
by powers of time (main plot) and without rescaling (inset). Figures taken
from Ref. Boguslavski:2018beu .

Figure 5: Gluonic simulation results in 2+1D plasmas. Left: The transversely
polarized spectral function at different momentum modes as a function of
$\omega$ for a 2+1D (left) and a 2+1D simulation coupled to an adjoint scalar
field (right). Right: The extracted damping rate of the genuinely 2+1D plasma
at different times as a function of momentum, rescaled by the time-dependent
plasmon mass $\omega_{\mathrm{pl}}(t)$ (main plot) and without rescaling
(inset). Figures taken from Ref. Boguslavski:2021buh .

###  4.1 Gluon spectral function in isotropic 3+1D plasmas

For the 3+1D simulations, the same lattice parameters and initial condition
are used as in Sec. 3.2 (see Ref. Boguslavski:2018beu  for details on the
gluonic simulation results). The gluonic spectral function
$\dot{\rho}_{T}\approx\omega\rho_{T}$ is shown in the left panel of Fig. 4,
extracted at the same reference time $Qt=Qt_{\text{ref}}=1500$ as we did for
the fermion $\rho_{+}$. One finds that $\dot{\rho}_{T}$ is dominated by narrow
Lorentzian quasiparticle peaks. The damping rate (peak width) is hence seen to
be much smaller than the dispersion (peak position)
$\gamma_{T}(t,p)\ll\omega_{T}(t,p)$. It is further shown in Ref.
Boguslavski:2018beu  that HTL provides an accurate description, with the
dispersion relation $\omega_{T}(p)$, the quasiparticle residue $Z_{T}(p)$ and
the Landau damping contribution for $|\omega|\leq p$ (here suppressed by one
power of $\omega$) agreeing well with HTL predictions.

Similarly to the fermionic case in Eq. (9), the damping rate $\gamma_{T}(t,p)$
is extracted from the gluon spectral function by fitting Lorentzian shapes to
the data. The results are shown in the right panel of Fig. 4. The full
momentum dependence is obtained and plotted at different times in the inset.
Rescaling it as $\gamma_{T}(t,p)/g^{2}T_{*}$, with the effective temperature
$g^{2}T_{*}\sim Q(t/t_{\text{ref}})^{-3/7}$, and plotting it as a function of
$p/\omega_{\mathrm{pl}}(t)$, with the plasmon mass
$\omega_{\mathrm{pl}}(t)\equiv\omega_{T}(t,p{=}0)\sim m_{F}(t)\sim
Q(t/t_{\text{ref}})^{-1/7}$, all the curves are seen to fall on top of each
other. Hence, as time grows the ratio
$\frac{\gamma_{T}(t,p)}{\omega(t,p{=}0)}\sim(Qt)^{-2/7}$ decreases. This shows
that the peaks become even more narrow with time and quasiparticle
approximations remain valid.

###  4.2 Gluon spectral function in 2+1D plasmas

In Sec. 4.1 we considered a 3+1D plasma with a gluonic distribution $f(t,p)$
that is isotropic in momentum space. If we now allow for momentum anisotropy
with the extreme case being $f(t,p_{T},p_{z}{=}0)$, the system becomes
effectively 2+1 dimensional. Here we can distinguish between a genuinely 2+1D
plasma consisting of only gluonic excitations and one that is coupled to an
adjoint scalar field, which we refer to as Glasma-like. The latter corresponds
to the case of extreme momentum anisotropy mentioned above.

We consider these 2+1D systems far from equilibrium and discuss the spectral
functions close to their self-similar attractors (see Ref. Boguslavski:2021buh
for details). The numerical results are summarized in Fig. 5. In the left
panel $\dot{\rho}_{T}$ is shown in the frequency domain for different momenta.
Different from gluonic and fermionic spectral functions in 3+1D, it exhibits
broad non-Lorentzian peaks of approximately Gaussian shape and differs
considerably from perturbative HTL predictions that are displayed as green
dashed lines. The failure of HTL perturbation theory is not surprising in
2+1D. There the HTL approximation fails because soft momentum modes provide
important contributions Boguslavski:2019fsb .

The damping rate $\gamma_{T}(t,p)$ is extracted by performing Gaussian fits to
the data and is shown at different times in the right panel of Fig. 5.
Rescaled by the time-dependent plasmon mass $\omega_{\mathrm{pl}}(t)$, the
damping rate becomes stationary, showing
$\gamma(t,p)\sim\omega_{\mathrm{pl}}(t)$ even quantitatively. Hence, the
damping rate remains of the order of $\omega_{\mathrm{pl}}(t)$. Since the
dispersion is $\omega_{T}(t,p)\approx\omega_{\mathrm{pl}}(t)$ for modes
$p\lesssim\omega_{\mathrm{pl}}(t)$, quasiparticle approximations break down
for these soft momenta.

##  5 Comparison and conclusion

We have presented a new tool to extract the fermion spectral function in
highly occupied plasmas nonperturbatively and we have revisited a similar
method for gluonic spectral functions. They enable a first-principles
determination of damping rates and were applied to classical self-similar
attractors in different dimensions.

The extracted fermion and gluon spectral functions share many properties in
3+1D. They exhibit Lorentzian quasiparticle peaks and are well described by
perturbative HTL expressions. Their damping rates are in general much smaller
than the dispersions, confirming a valid quasiparticle picture. However, there
are also interesting differences in the momentum dependence of the damping
rates. While the gluonic damping rate $\gamma_{T}$ increases with momentum,
the fermionic one $\gamma_{+}$ decreases. Moreover, the gluonic damping rate
decreases in time faster than the mass, while the fermionic one scales with
the mass for $p=0$.

Interestingly, the scaling $\gamma(t,p{=}0)\sim\omega_{\mathrm{pl}}(t)$ has
also been observed for the gluonic spectral function in 2+1D systems (see Sec.
4.2). However, there the peaks in the spectral function are generally much
broader than in 3+1D and a perturbative HTL description breaks down.

So far, we have applied the methods to systems in Minkowski space-time. In the
future, we want to apply our framework to Bjorken expanding (anisotropic)
systems and to heavy quarks, which are of phenomenological relevance to heavy-
ion collisions. Moreover, nonperturbative properties of the spectral or
statistical correlators out of equilibrium can have a sizable impact on
transport coefficients, as was shown in Boguslavski:2020tqz . Thus, it will be
interesting to study how observables like the photon production rate are
affected by our numerical results.

## References

  * (1) K. Boguslavski, T. Lappi, M. Mace, S. Schlichting (2021), 2106.11319
  * (2) K. Boguslavski, A. Kurkela, T. Lappi, J. Peuron, JHEP 05, 225 (2021), 2101.02715
  * (3) K. Boguslavski, A. Kurkela, T. Lappi, J. Peuron, Phys. Rev. D98, 014006 (2018), 1804.01966
  * (4) S. Schlichting, D. Teaney, Ann. Rev. Nucl. Part. Sci. 69, 447 (2019), 1908.02113
  * (5) J. Berges, M.P. Heller, A. Mazeliauskas, R. Venugopalan (2020), 2005.12299
  * (6) R. Allahverdi, R. Brandenberger, F.Y. Cyr-Racine, A. Mazumdar, Ann. Rev. Nucl. Part. Sci. 60, 27 (2010), 1001.2600
  * (7) J. Berges, AIP Conf. Proc. 739, 3 (2004), hep-ph/0409233
  * (8) J. Berges, K. Boguslavski, S. Schlichting, R. Venugopalan, JHEP 05, 054 (2014), 1312.5216
  * (9) M. Mace, N. Mueller, S. Schlichting, S. Sharma, Phys. Rev. D 95, 036023 (2017), 1612.02477
  * (10) M. Mace, N. Mueller, S. Schlichting, S. Sharma, Phys. Rev. Lett. 124, 191604 (2020), 1910.01654
  * (11) G. Aarts, J. Smit, Nucl. Phys. B 555, 355 (1999), hep-ph/9812413
  * (12) A. Kurkela, G.D. Moore, Phys. Rev. D 86, 056008 (2012), 1207.1663
  * (13) S. Schlichting, Phys. Rev. D 86, 065008 (2012), 1207.1450
  * (14) J. Berges, K. Boguslavski, S. Schlichting, R. Venugopalan, Phys. Rev. D 89, 114007 (2014), 1311.3005
  * (15) K. Boguslavski, A. Kurkela, T. Lappi, J. Peuron, Phys. Rev. D 100, 094022 (2019), 1907.05892
  * (16) J.P. Blaizot, E. Iancu, Phys. Rept. 359, 355 (2002), hep-ph/0101103
  * (17) E. Braaten, R.D. Pisarski, Phys. Rev. D 46, 1829 (1992)
  * (18) A. Kurkela, T. Lappi, J. Peuron, Eur. Phys. J. C 76, 688 (2016), 1610.01355
  * (19) A. Piñeiro Orioli, J. Berges, Phys. Rev. Lett. 122, 150401 (2019), 1810.12392
  * (20) K. Boguslavski, A. Piñeiro Orioli, Phys. Rev. D 101, 091902 (2020), 1911.04506
  * (21) K. Boguslavski, A. Kurkela, T. Lappi, J. Peuron, JHEP 09, 077 (2020), 2005.02418