Cosmological decay of Higgs-like scalars into a fermion channel

TL;DR

This paper investigates the decay of Higgs-like scalars into fermions in a cosmological setting, deriving non-perturbative survival probabilities and highlighting deviations from Minkowski decay laws due to cosmological effects.

Contribution

It introduces a non-perturbative method to compute decay probabilities of scalar particles in expanding universes, revealing memory effects and providing bounds on particle survival based on cosmological parameters.

Findings

Derived explicit survival probability formulas accounting for cosmological redshift.

Identified that phenomenological Minkowski decay laws underestimate particle lifetimes.

Established an upper bound on survival probability as a function of redshift.

Abstract

We study the decay of a Higgs-like scalar Yukawa coupled to massless fermions in post-inflationary cosmology, combining a non-perturbative method with an adiabatic expansion. The renormalized survival probability $\mathcal{P}_\Phi(t)$ of a (quasi) particle ``born'' at time $t_b$ and decaying at rest in the comoving frame, $\mathcal{P}_\Phi(t) = \Big[\frac{t}{t_b}\Big]^{-\frac{Y^2}{8\pi^2}}~ e^{ \frac{Y^2}{4\pi^2}\,\big(t/t_b\big)^{1/4} } \,e^{-\Gamma_0\,(t-t_b)}~ \mathcal{P}_\Phi(t_b) $, with $\Gamma_0$ the decay rate at rest in Minkowski space-time. For an ultrarelativistic particle we find $\mathcal{P}_\Phi(t) = e^{-\frac{2}{3}\Gamma_0\,t_{nr}\,(t/t_{nr})^{3/2}}~ \mathcal{P}_\Phi(t_b)$ before it becomes non-relativistic at a time $t_{nr}$ as a consequence of the cosmological redshift. For $t\gg t_{nr}$ we find $\mathcal{P}_\Phi(t) = \Big[\frac{t}{t_{nr}}\Big]^{-\frac{Y^2}{8\pi^2}}~ e^{ \frac{Y^2}{4\pi^2}\,\big(t/t_{nr}\big)^{1/4} }~\Big[\frac{t}{t_{nr}}\Big]^{\Gamma_0 t_{nr}/2} \,e^{-\Gamma_0\,(t-t_{nr})}~ \mathcal{P}_\Phi(t_{nr})$. The extra power is a consequence of the memory on the past history of the decay process. We compare these results to an S-matrix inspired phenomenological Minkowski-like decay law modified by an instantaneous Lorentz factor to account for cosmological redshift. Such phenomenological description \emph{under estimates the lifetime of the particle}. For very long lived, very weakly coupled particles, we obtain an \emph{upper bound} for the survival probability as a function of redshift $z$ valid throughout the expansion history $\mathcal{P}_\Phi(z) \gtrsim e^{-\frac{\Gamma_0}{H_0}\,\Upsilon(z,z_b)}\,\mathcal{P}_\Phi(z_b)$, where $\Upsilon(z,z_b)$ only depends on cosmological parameters and $t_{nr}$.

Equations

$$S=\int d^{3}x\;dt\;\sqrt{-g}\,\Bigg{\{}\frac{1}{2}{\dot{\phi}^{2}}-\frac{(\nabla\phi)^{2}}{2a^{2}}-\frac{1}{2}\,\big{[}M^{2}+\xi\,R\big{]}\phi^{2}+\overline{\Psi}\Big{[}i\,\gamma^{\mu}\;\mathcal{D}_{\mu}-m_{f}-Y\phi\Big{]}\Psi\Bigg{\}}\,,$$

$$R=6\Big{[}\frac{\ddot{a}}{a}+\Big{(}\frac{\dot{a}}{a}\Big{)}^{2}\Big{]}\,,$$

$$g^{\mu \nu }(x)=e_a^{\mu }(x)e_b^{\nu }(x){\eta }^{ab},$$

$${\gamma }^{\mu }(x)={\gamma }^a{e}_a^{\mu }(x),\{{\gamma }^{\mu }(x),{\gamma }^{\nu }(x)\}=2g^{\mu \nu }(x),$$

$${\mathcal{D}}_{\mu }={\partial }_{\mu }+\frac{1}{8}[{\gamma }^c,{\gamma }^d]e_c^{\nu }\left({\partial }_{\mu }{e}_{d\nu }-{\Gamma }_{\mu \nu }^{\lambda }{e}_{d\lambda }\right),$$

$$g_{\mu \nu }=C^2(\eta ){\eta }_{\mu \nu },C(\eta )\equiv a(t(\eta )),$$

$$e_a^{\mu }=C^{-1}(\eta ){\delta }_a^{\mu }\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}{e}_{\mu }^a=C(\eta ){\delta }_{\mu }^a.$$

$$S_{f}=\int d^{3}x\,d\eta\,C^{4}(\eta)\,\overline{\Psi}(\vec{x},\eta)\,\Bigg{[}i\frac{\gamma^{0}}{C(\eta)}\,\Big{(}\frac{d}{d\eta}+3\frac{C^{{}^{\prime}}(\eta)}{2C(\eta)}\Big{)}+i\,\frac{\gamma^{i}}{C(\eta)}\nabla_{i}-m_{f}-Y\,\phi\Bigg{]}\Psi(\vec{x},\eta)\,.$$

$$\sqrt{-g}\;\overline{\Psi}\Big{(}i\;\gamma^{\mu}\;\mathcal{D}_{\mu}\Psi-m_{f}-Y\phi\Big{)}\Psi=\big{(}C^{3/2}(\eta)\,\overline{\Psi}(\vec{x},\eta)\big{)}\;\Big{[}i\;{\not\!{\partial}}-(m_{f}+Y\phi)\;C(\eta)\Big{]}\big{(}C^{3/2}(\eta)\,{\Psi}(\vec{x},\eta)\big{)}\,,$$

$$C(\eta )\varphi (\vec{x},t)=\chi (\vec{x},\eta )\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}{C}^{\frac{3}{2}}(\eta )\Psi (\vec{x},t)=\psi (\vec{x},\eta ),$$

$$S=\int d^{3}x\;d\eta\,\Big{\{}\mathcal{L}_{0}[\chi]+\mathcal{L}_{0}[\psi]+\mathcal{L}_{I}[\chi,\psi]\Big{\}}\;,$$

$${\mathcal{L}}_0[\chi ]$$

$${\mathcal{L}}_0[\psi ]$$

$${\mathcal{L}}_I[\chi ,\psi ]$$

$${\mathcal{M}}^2(\eta )=m_{\varphi }^2{C}^2(\eta )-\frac{C^{\prime \prime }(\eta )}{C(\eta )}(1-6\xi ),$$

$$M_f(\eta )=m_{f}C(\eta ).$$

$$\chi^{\prime\prime}_{\vec{k}}(\eta)+\Big{[}k^{2}+\mathcal{M}^{2}(\eta)\Big{]}\chi_{\vec{k}}(\eta)=0\,.$$

$$\chi(\vec{x},\eta)=\frac{1}{\sqrt{V}}\,\sum_{\vec{k}}\Big{[}a_{\vec{k}}\,g_{k}(\eta)\,e^{i\vec{k}\cdot\vec{x}}+a^{\dagger}_{\vec{k}}\,g^{*}_{k}(\eta)\,e^{-i\vec{k}\cdot\vec{x}}\Big{]}\,,$$

$$\Big{[}\frac{d^{2}}{d\eta^{2}}+k^{2}+\mathcal{M}^{2}(\eta)\Big{]}g_{k}(\eta)=0\,.$$

$$g_k^{{}^{\prime }}(\eta )g_k^{\ast }(\eta )-g_k^{\ast {}^{\prime }}(\eta )g_k(\eta )=-i,$$

$$\psi (\vec{x},\eta )=\frac{1}{\sqrt{V}}\sum _{\vec{k},\lambda =1,2}\left[b_{\vec{k},\lambda }{U}_{\lambda }(\vec{k},\eta )e^{i\vec{k}\cdot \vec{x}}+d_{\vec{k},\lambda }^{†}{V}_{\lambda }(\vec{k},\eta )e^{-i\vec{k}\cdot \vec{x}}\right],$$

$$\displaystyle\Bigg{[}i\;\gamma^{0}\;\partial_{\eta}-\vec{\gamma}\cdot\vec{k}-M_{f}(\eta)\Bigg{]}U_{\lambda}(\vec{k},\eta)$$

$$\displaystyle\Bigg{[}i\;\gamma^{0}\;\partial_{\eta}+\vec{\gamma}\cdot\vec{k}-M_{f}(\eta)\Bigg{]}V_{\lambda}(\vec{k},\eta)$$

$$U_{\lambda }(\vec{k},\eta )$$

$$V_{\lambda }(\vec{k},\eta )$$

$${\gamma }^0{\mathcal{U}}_{\lambda }={\mathcal{U}}_{\lambda },{\gamma }^0{\mathcal{V}}_{\lambda }=-{\mathcal{V}}_{\lambda }.$$

$$\left[\frac{d^2}{d{\eta }^2}+k^2+M_f^2(\eta )-i{M}_f^{\prime }(\eta )\right]f_k(\eta )$$

$$\left[\frac{d^2}{d{\eta }^2}+k^2+M_f^2(\eta )+i{M}_f^{\prime }(\eta )\right]h_k(\eta )$$

$$\frac{d}{d\eta }(U_{\lambda }^{†}(q,\eta )U_{\lambda }(q,\eta ))=0\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}\frac{d}{d\eta }(V_{\lambda }^{†}(q,\eta )V_{\lambda }(q,\eta ))=0.$$

$$U_{\lambda }^{†}(q,\eta )U_{{\lambda }^{\prime }}(q,\eta )=V_{\lambda }^{†}(q,\eta )V_{{\lambda }^{\prime }}(q,\eta )={\delta }_{\lambda ,{\lambda }^{\prime }},$$

$$\Big{(}\frac{\dot{a}}{a}\Big{)}^{2}=H^{2}(t)=H^{2}_{0}\,\Bigg{[}\frac{\Omega_{M}}{a^{3}(t)}+\frac{\Omega_{R}}{a^{4}(t)}+\Omega_{\Lambda}\Bigg{]}\,,$$

$$H_0=1.5\times 10^{-42}\mathrm{GeV}\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}{\Omega }_M=0.308\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}{\Omega }_R=5\times 10^{-5}\text{leavevmode}\text{leavevmode};\text{leavevmode}\text{leavevmode}{\Omega }_{\Lambda }=0.692.$$

$$\frac{d\,C(\eta)}{d\eta}=H_{0}\,\sqrt{\Omega_{M}}\,\Big{[}a_{eq}+C(\eta)+s\,C^{4}(\eta)\Big{]}^{1/2}\,,$$

Full text

I Introduction

The decay and scattering of particles are some of the most fundamental processes in particle physics, within and beyond the Standard Model, with profound impact in cosmology. These processes are ultimately responsible for establishing a state of local thermodynamic and chemical equilibrium and are fundamental ingredients in kinetic processes in the early universekolb ; linde ; bernstein . Particle decay is not only ubiquitous, but it plays an important role in big bang nucleosynthesis (BBN)kolb ; steigman ; sarkar ; fields ; serpicobbn ; pospelovbbn ; jeda , and the generation of the baryon and lepton asymmetriestrodden ; kolba ; lepgen ; lepgen2 . The decay of long-lived dark matter particles is constrained by various cosmological and astrophysical probeszentner ; dm1 ; dm2 ; dm3 ; dm4 , and recently it has been suggested that the two body decay of a long lived dark matter particle may relieve the tension between distance ladder and cosmic microwave background measurements of the Hubble constantloeb .

Most treatments of particle decay (and/or inverse decay) in cosmology implement the S-matrix quantum field theory approach as in Minkowski space-time. In this framework, the unstable decaying state is prepared at a time far in the past ($t\rightarrow-\infty$), and one obtains the transition amplitude to a given final state far in the future ($t\rightarrow\infty$). Taking the infinite time limit in the transition amplitude yields a total energy conserving delta function. Squaring this delta function to obtain the transition probability yields a total energy conservation delta function multiplied by the total elapsed time. Dividing by this large time and summing over all the final states for a given decay channel gives the total transition probability per unit time, namely a decay rate. Energy conservation, a consequence of the infinite time limit, yields kinematic constraints (thresholds) for decay and scattering processes.

In an expanding cosmology such an approach is at best approximate and at worst unreliable when the Hubble expansion rate is large even during a post-inflationary early stage of a radiation dominated cosmology, or if the lifetime of the particle is of the order of the Hubble time. In a spatially flat Friedmann-Robertson-Walker (FRW) cosmology there are three space-like Killing vectors associated with spatial translational invariance and spatial momentum conservation, however, as a consequence of cosmological expansion there is no global time-like Killing vector, therefore particle energy is not manifestly conserved in scattering or decay processes.

A consistent formulation of dynamic processes in an expanding cosmology requires implementing methods of quantum field theory in curved space timeparker ; zelstaro ; fulling ; birford ; bunch ; birrell ; fullbook ; mukhabook ; parkerbook . Early studies revealed a wealth of novel phenomena such as particle productionparker ; birford ; bunch and processes that are forbidden in Minkowski space time as consequence of strict energy conservation.

S-matrix theory was extended to simple cosmological space times to study the decay of a massive particle into two massless particles conformally coupled to gravity in ref.spangdecay . In referencesvilja ; viljafer these methods were adapted to calculate the decay of a massive bosonic particle at rest into two massless bosonic particles conformally coupled to gravity and into massless fermions Yukawa coupled to a scalar.

More recentlyherring the decay of bosonic particles into two other bosonic degrees of freedom during a radiation dominated era was studied by implementing a non-perturbative method. This method was adapted to quantum field theory from the study of linewidths in quantum opticsboyww ; boycosww , combined with a physically motivated adiabatic expansion. While the results of this reference agreed with those obtained in ref.vilja for a particle decaying at rest in the comoving frame in the long time limit, they revealed new phenomena for highly relativistic decaying particles as a consequence of the cosmological redshift, and the relaxation of kinematic thresholds as a consequence of energy uncertainties determined by the Hubble scale.

Our study in this article is a natural extension of that in ref.herring focusing on decay of a heavy bosonic particle into fermions, a more relevant case for standard model physics (and probably beyond) since most of the fermionic degrees of freedom in the standard model (with the possible exception of neutrinos) are Yukawa coupled to the Higgs boson.

Brief summary:

The study of fermionic degrees of freedom as decay products introduces several conceptually important distinctions with the bosonic case studied in refs.vilja ; herring that results in novel aspects of cosmological decay. First, fermionic degrees of freedom couple to the background gravitational field via the spin connectionweinbergbook ; parker ; duncan ; casta ; birrell ; parkerbook ; fercurved1 ; fercurved2 ; fercurved3 ; feradia ; feradia2 ; feradia3 ; boydVS ; baacke . Secondly, fermions Yukawa coupled to a bosonic degree of freedom yield a renormalizable theory. Recently the decay of a bosonic particle Yukawa coupled to fermions was studied within a non-perturbative real time framework in Minkowski space-timezenoqft . This study revealed novel transient dynamics associated with the dressing of the decaying particle by fermion-antifermion pairs into a quasiparticle state, which decays on a longer time scale. Such “dressing” leads to the necessity of an ultraviolet divergent renormalization of the decaying state and a detailed understanding of the various time scales to separate the many-particle dynamics of renormalization and dressing from that of the actual decay of the quasiparticle. Such dynamical effects cannot be addressed within an S-matrix framework since these effects are not secular in time and their contribution vanishes when the transition probability is divided by the total time in the infinite time limit. The dynamics of dressing and quasiparticle formation have been recently addressed in ref.collins for a consistent interpretation of the reduction formula in asymptotic quantum field theory.

We introduce a dynamical renormalization that absorbs the ultraviolet divergences associated with fermion pairs into a renormalized survival probability at a renormalization time scale $t_{b}$. The survival probability obeys a dynamical renormalization group equation with respect to $t_{b}$. The cosmological redshift encodes the memory of the transient dynamics of quasiparticle formation in the decay law not seen in Minkowski space-time. If the decaying particle is ultrarelativistic, the decay dynamics depends crucially on $t_{nr}$, the time scale at which it becomes non-relativistic as a consequence of the cosmological redshift. An S-matrix inspired, phenomenologically motivated, Minkowski-like decay law is shown to under estimate the lifetime of the decaying state. Section (II) introduces the model and the adiabatic approximation, section (III) summarizes the non-perturbative framework to obtain the time evolution of the survival probability. In section (IV) we obtain the decay function for massless fermions during radiation domination, section (V) describes the dynamical renormalization method, section (VI) analyzes the decay dynamics of the renormalized survival probability during radiation domination, compares the results to an S-matrix inspired decay function, and introduces an upper bound to the decay function for very long lived, very weakly coupled particles valid all throughout the expansion history. Section (VII) discusses the various results analyzing their regime of validity and highlighting several implications. Section (VIII) presents our conclusions summarizing the main results. Various appendices contain technical details, in particular appendix (B) derives the decay law in Minkowski space time, highlighting the renormalization aspects to compare to the curved space-time case.

II The model:

We consider a Higgs-like scalar field Yukawa coupled to one Dirac fermion in a spatially flat Friedmann-Robertson-Walker (FRW) cosmology with scale factor $a(t)$ in comoving time. Generalizing to include Majorana fermions and/or more fermionic species is straightforward.

In comoving coordinates, the action is given by

[TABLE]

where

[TABLE]

is the Ricci scalar, and $\xi$ is the coupling to gravity, with $\xi=0,1/6$ corresponding to minimal or conformal coupling, respectively. Introducing the vierbein field $e^{\mu}_{a}(x)$ defined as

[TABLE]

where $\eta_{ab}$ is the Minkowski space-time metric, the curved space time Dirac gamma- matrices $\gamma^{\mu}(x)$ and the fermionic covariant derivative $\mathcal{D}_{\mu}$ are given byweinbergbook ; birrell ; duncan ; casta

[TABLE]

where the $\gamma^{a}$ are the Minkowski space time Dirac matrices, chosen to be in the standard Dirac representation, and the covariant derivative $\mathcal{D}_{\mu}$ is given in terms of the spin connection by

[TABLE]

where $\Gamma^{\lambda}_{\mu\nu}$ are the usual Christoffel symbols.

For an (FRW) in conformal time $d\eta=dt/a(t)$, the metric becomes

[TABLE]

where $\eta_{\mu\nu}=\textrm{diag}(1,-1,-1,-1)$ is the flat Minkowski space-time metric, and the vierbeins $e^{\mu}_{a}$ are given by

[TABLE]

The fermionic part of the action in conformal coordinates now becomes

[TABLE]

The Dirac Lagrangian density in conformal time simplifies to

[TABLE]

where $i{\not\!{\partial}}=\gamma^{a}\partial_{a}$ is the usual Dirac differential operator in Minkowski space-time in terms of flat space time $\gamma^{a}$ matrices. Introducing the conformally rescaled fields

[TABLE]

and neglecting surface terms, the action becomes

[TABLE]

with

[TABLE]

The effective time dependent masses are given by

[TABLE]

and

[TABLE]

In the non-interacting case, $Y=0$, the Heisenberg equations of motion for the spatial Fourier modes with comoving wavevector $\vec{k}$ for the conformally rescaled scalar field are

[TABLE]

The Heisenberg fields are quantized in a comoving volume $V$, the real scalar field $\chi$ is expanded as

[TABLE]

where the mode functions $g_{k}(\eta)$ obey

[TABLE]

The mode functions are chosen to obey the Wronskian condition

[TABLE]

and $a,a^{\dagger}$ obey the usual canonical commutation relations.

For Dirac fermions the field $\psi({\vec{x}},\eta)$ is expanded as

[TABLE]

where the spinor mode functions $U,V$ obey the Dirac equationsfercurved1 ; fercurved2 ; fercurved3 ; feradia ; feradia2 ; feradia3 ; boydVS ; baacke

[TABLE]

These equations become simpler by writing

[TABLE]

with $\mathcal{U}_{\lambda};\mathcal{V}_{\lambda}$ being constant spinorsboydVS ; baacke obeying

[TABLE]

Inserting (II.23,II.24) into the Dirac equations (II.21,II.22) and using (II.25), it follows that the mode functions $f_{k}(\eta);h_{k}(\eta)$ obey the equations

[TABLE]

Multiplying the Dirac equations on the left by $\gamma^{0}$, it is straightforward to confirm that

[TABLE]

We choose the normalizations

[TABLE]

so that the operators $b,b^{\dagger},d,d^{\dagger}$ obey the canonical anticommutation relations. Furthermore, we will choose particle-antiparticle boundary conditions so that $h_{k}(\eta)=f^{*}_{k}(\eta)$ (see below). We note that for $m_{f}=0$ the conformally rescaled fermi fields obey the same equations as in Minkowski space-time but in terms of conformal time, whereas this only occurs for bosons if they are conformally coupled to gravity, namely with $\xi=1/6$, or for a radiation dominated cosmology (see below). The equivalence of massless fermions to those in Minkowski space-time will allow a direct comparison with the case of decay in flat space time studied in ref.zenoqft and summarized in appendix (B), and to interpret the differences with the curved space-time case.

II.1 Adiabatic approximation in post-inflationary cosmology:

The standard (post-inflation) cosmology is described by radiation (RD), matter (MD) and dark energy (DE) dominated stages, we take the latter to be described by a cosmological constant. Friedmann’s equation in comoving time is

[TABLE]

where the scale factor is normalized to $a_{0}=a(t_{0})=1$ today. We take as representative the following values of the parameters wmap ; spergel ; planck :

[TABLE]

Passing to conformal time $\eta$ with $d\eta=dt/a(t)$, where the metric is given by (II.5) and $C(\eta)\equiv a(t(\eta))$, it follows that

[TABLE]

with

[TABLE]

$a_{eq}$ is the scale factor at matter-radiation equality.

Hence the different stages of cosmological evolution, namely (RD), (MD), and (DE), are characterized by

[TABLE]

We will begin by studying the dynamics of particle decay during the (RD) dominated era, generalizing afterwards to the case of a very long lived, very weakly coupled particle. During (RD) and (MD) we find,

[TABLE]

where

[TABLE]

and conformal time in terms of the scale factor is given by

[TABLE]

During the (RD) stage

[TABLE]

and the relation between conformal and comoving time is given by

[TABLE]

a result that will prove useful in the study of the decay law during this stage.

Bosonic fields:

Solving the mode equations (II.18,II.26,II.27) with the cosmological scale factor (II.35) is obviously very challenging, instead we implement a physically motivated adiabatic expansion. To highlight the nature of the expansion let us consider first the bosonic mode equation (II.18). The term proportional to $C^{\prime\prime}/C$ in (II.18) vanishes identically in a radiation dominated cosmology or for conformally coupled bosonic fields for which $\xi=1/6$. We argue below that we can consistently neglect this term to leading order in the adiabatic expansion all throughout the cosmological evolution during (RD) and (MD) (see eqn. (II.54)). Neglecting this term, the mode equation (II.18) becomes

[TABLE]

We recognize that

[TABLE]

where

[TABLE]

is the local energy measured by a comoving observer, and $k_{p}(t)$ is the physical wavevector redshifting with the cosmological expansion.

Writing the solution of (II.40) in the WKB formbirrell ; fullbook ; mukhabook ; parkerbook ; birford

[TABLE]

and inserting this ansatz into (II.40) it follows that $W_{k}(\eta)$ must be a solution of the equationbirrell

[TABLE]

This equation can be solved in an adiabatic expansion

[TABLE]

We refer to terms that feature $n$-derivatives of $\omega_{k}(\eta)$ as of n-th adiabatic order. The nature and reliability of the adiabatic expansion is revealed by considering the term of first adiabatic order

[TABLE]

this is most easily recognized in comoving time $t$ in terms of the comoving local energy (II.41,II.42) and the Hubble expansion rate

[TABLE]

In terms of these variables, the first order adiabatic ratio (II.46) becomesherring

[TABLE]

where

[TABLE]

is the local Lorentz factor.

The adiabatic approximation relies on the smallness of the (time dependent) adiabatic ratio

[TABLE]

corresponding to the physical wavelength $\propto 1/k_{p}(t)$ and/or the Compton wavelength of the particle $1/m_{\phi}$ being much smaller than the size of the particle horizon $d_{H}(t)\propto 1/H(t)$ at a given time. During (RD) the particle horizon grows as $a^{2}(t)$ and during (MD) it grows as $a^{3/2}(t)$ whereas the physical wavelength grows as $a(t)$. Therefore, if at a given initial time the adiabatic approximation is valid and $H(t)\ll E_{k}(t)$ the reliability of the adiabatic expansion improves with the cosmological expansion.

To understand the origin of this approximation consider that the decaying particle is produced in the (RD) stage during which

[TABLE]

where $g_{\text{eff}}\lesssim 100$ is the number of ultrarelativistic degrees of freedom. Therefore,

[TABLE]

An upper bound on this ratio is obtained by considering that the decaying particle is produced at the scale of grand unification with $T\simeq 10^{15}\,\mathrm{GeV}$, assuming that this scale describes the onset of the (RD) era. Taking a typical comoving energy $E_{k}(t)\simeq T(t)$, one finds that $H(t)/E_{k}(t)\lesssim 10^{-3}$ and diminishes with cosmological expansion and diminishing temperature. This argument suggests that for typical particle physics processes the adiabatic ratio $H(t)/E_{k}(t)\ll 1$ throughout the post-inflation thermal history.

In terms of this adiabatic ratio, we find

[TABLE]

where $R(t)$ is the Ricci scalar (II.2). Furthermore, it is straightforward to find that

[TABLE]

therefore, this ratio is of second adiabatic order and can be safely neglected to the leading adiabatic order pursued in this study, justifying the simplification of the mode equations to (II.40) even for non-conformal coupling to gravity.

In this study we consider the zeroth-adiabatic order with the mode functions given by

[TABLE]

Since the decay function is $\propto Y^{2}$, keeping the zeroth adiabatic order yields the leading contribution to the decay law. Furthermore, as shown in detail in ref.herring particle production as a consequence of cosmological expansion is an effect of higher order in the adiabatic expansion, thus it can be safely neglected to leading order.

The phase of the mode function has an immediate interpretation in terms of comoving time and the local comoving energy (II.41,II.42), namely

[TABLE]

which is a natural generalization of the phase of positive frequency particle states in Minkowski space-time.

During the (RD) era with $C(\eta)$ given by (II.38) we find that the criterion (II.50) for the validity of the adiabatic approximation implies

[TABLE]

Fermi fields: The adiabatic expansion is straightforwardly applied to the fermionic case and has been discussed in the literaturefercurved2 ; fercurved3 ; feradia ; feradia2 ; feradia3 . Beginning with the mode equations (II.26,II.27) with $M^{\prime}_{f}(\eta)=m_{f}\,C^{\prime}(\eta)$ and, now with $\omega^{2}_{k}(\eta)=k^{2}+M^{2}_{f}(\eta)$, it follows that

[TABLE]

therefore the purely imaginary term in these mode equations are of first adiabatic order and will be neglected to leading (zeroth) adiabatic order. Hence, to leading order we find

[TABLE]

In what follows we will refer to $\omega^{2}_{k}(\eta)=k^{2}+M^{2}(\eta)$ for both bosonic and fermionic degrees of freedom with $M^{2}(\eta)=m^{2}C^{2}(\eta)$ for either case. To leading (zeroth) order in the adiabatic expansion the Dirac spinor solutions in the standard Dirac representation and with the normalization conditions (II.29) are found to be

[TABLE]

and

[TABLE]

To leading adiabatic order these spinors satisfy the completeness relations

[TABLE]

where the projector operators at different times $\Lambda^{+}_{k}(\eta,\eta^{\prime})\leavevmode\nobreak\ ;\leavevmode\nobreak\ \Lambda^{-}_{k}(\eta^{\prime},\eta)$ and their properties are given in appendix (A).

III Non-perturbative approach to the decay law:

In Minkowski space-time, the decay rate of a particle is typically computed via S-matrix theory by obtaining the transition probability per unit time from an in-state prepared in the infinite past to an out-state in the infinite future. Obviously, such an approach – taking the infinite time limit– is not suitable in a time dependent gravitational background. An alternative approach in Minkowski space-time considers the Dyson-resummed propagator in frequency space that includes radiative corrections through the self-energy. The imaginary part of the self-energy evaluated on the mass shell in frequency space is identified with the decay rate, and a Breit-Wigner approximation to the full propagator, namely approximating the self-energy near the (complex) pole yields the exponential decay law. Again such an approach is not available in an expanding cosmological background where the lack of a time-like Killing vector prevents Fourier transforms in time-frequency, and makes the self-energy explicitly dependent on two time arguments, not only on the difference.

Instead we implement a quantum field theory method that complements and extends the Wigner-Weisskopf theory of atomic linewidths, that is particularly suited to study time evolution in time dependent situations. This method is manifestly unitary and yields a non-perturbative description of transition amplitudes and probabilities directly in real time. We summarize below the main aspects of the method as it applies to this study, referring the reader toherring ; boyww ; boycosww for details. The total Hamiltonian in conformal time is given by $H_{0}+H_{I}$ where $H_{0}$ is the free field Hamiltonian and

[TABLE]

is the interaction Hamiltonian in the interaction picture. Passing to the interaction picture wherein a given state is expanded in the Fock states associated with the creation and annihilation operators $a,a^{\dagger},b,d$, etc. of the free theory, namely $|\Phi(\eta)\rangle_{I}=\sum_{n}\mathcal{C}_{n}(\eta)|n\rangle$, the amplitudes obey the coupled equations

[TABLE]

This is an infinite hierarchy of integro-differential equations for the coefficients $\mathcal{C}_{n}(\eta)$. Consider that initially the state is $|\Phi\rangle$ so that $\mathcal{C}_{\Phi}(\eta_{i})=1\,;\,\mathcal{C}_{\kappa}(\eta_{i})=0$ for $|\kappa\rangle\neq|\Phi\rangle$, and consider a first order transition process $|\Phi\rangle\rightarrow|\kappa\rangle$ to intermediate multiparticle states $|\kappa\rangle$ with transition matrix elements $\langle\kappa|H_{I}(\eta)|\Phi\rangle$. Obviously the state $|\kappa\rangle$ will be connected to other multiparticle states $|\kappa^{\prime}\rangle$ different from $|\Phi\rangle$ via $H_{I}(\eta)$. Hence for example up to second order in the interaction, the state $|\Phi\rangle\rightarrow|\kappa\rangle\rightarrow|\kappa^{\prime}\rangle$. Restricting the hierarchy to first order transitions from the initial state $|\Phi\rangle\leftrightarrow|\kappa\rangle$, and neglecting the contribution from vacuum diagrams which just yield a re-definition of the vacuum state111This is one of the main differences with the method used in referencesspangdecay ; vilja ; viljafer where a disconnected vacuum diagram is also included in the transition amplitude. (see discussion in ref.herring ) results in the following coupled equations

[TABLE]

Equation (III.4) with $C_{\kappa}(\eta_{i})=0$ is formally solved by

[TABLE]

and inserting this solution into equation (III.3) we find

[TABLE]

where we have introduced the self-energy

[TABLE]

We study the decay of a single particle bosonic state into a fermion-anti-fermion pair to leading order in the Yukawa coupling and the adiabatic approximation. Therefore the initial state is a single particle bosonic state with momentum $\vec{k}$, namely $|\Phi\rangle\equiv|1^{\chi}_{\vec{k}}\rangle$. The set of states $|\kappa\rangle$ with a non-vanishing matrix element of $H_{I}$ with this single particle state are $|\kappa\rangle\equiv|1^{f}_{\vec{p},\lambda},\overline{1^{f}}_{\vec{q},\lambda^{\prime}}\rangle$ where $\lambda,\lambda^{\prime}$ are the polarization of the fermion and antifermion states. The matrix elements entering in the evolution of the amplitudes are

[TABLE]

and the self-energy (III.7) to leading order in the adiabatic expansion becomes

[TABLE]

where $\vec{q}=\vec{k}-\vec{p}$. This is the fermionic one-loop self energy in curved space time to leading order in the adiabatic expansion.

Obviously the differential equation (III.6) cannot be solved exactly with the above self-energy. In Minkowski space time the self-energy is a function of the time difference allowing a solution via Laplace transformboyww ; boycosww . However, in a time dependent expanding cosmology such an approach is not available. This is a consequence of the lack of a global time-like Killing vector. Instead for weak coupling we resort to a Markov approximationherring . While details are available in ref.herring ; boycosww ; boyww to which the reader is referred, we summarize here the main aspects of this approximation.

We begin by introducing

[TABLE]

such that

[TABLE]

with the condition

[TABLE]

Then (III.6) can be written as

[TABLE]

which can be integrated by parts to yield

[TABLE]

Since $\mathcal{E}_{\Phi}\propto\mathcal{O}(Y^{2})$ the first term on the right hand side of (III.14) is of order $Y^{2}$, whereas the second is $\mathcal{O}(Y^{4})$ because $d\,\mathcal{C}_{\Phi}(\eta)/d\eta\propto Y^{2}$. Therefore up to $\mathcal{O}(Y^{2})$ the evolution equation for the amplitude $\mathcal{C}_{\Phi}$ becomes

[TABLE]

with solution

[TABLE]

This expression clearly highlights the non-perturbative nature of the Wigner-Weisskopf approximation. The imaginary part of the self energy $\Sigma_{\Phi}$ yields a renormalization of the adiabatic frequencies and will not be addressed hereboycosww ; boyww , whereas the real part determines the decay law

[TABLE]

where we introduced the survival probability $\mathcal{P}_{\Phi}(\eta)$ with $\mathcal{P}_{\Phi}(\eta_{i})=|C_{\Phi}(\eta_{i})|^{2}$. This final expression for the survival probability directly exhibits the non-perturbative nature of the method. The self-energy is given by (III.9) to leading order in Yukawa coupling.

In references (herring ,boyww ,boycosww ,zenoqft ) it has been established that this non-perturbative framework correctly describes the short, intermediate and long time dynamics in Minkowski space-time. It provides a real-time non-perturbative resummation of Feynman diagrams to a given order in the perturbative expansion and in Minkowski space-time is equivalent to the time evolution obtained from the inverse Fourier transform of the Dyson resummed propagator in momentum space. For a decaying particle, the propagator has a complex pole in the second (or higher) Riemann sheet, for weak coupling, in the narrow width approximation the long time behavior is completely determined by this pole. The Wigner-Weisskopf and Markov approximation yield exactly the same result including wave function renormalization and also describes correctly the early time behaviorzenoqft .

This equivalence has been discussed in greater detail in refs. (boyww ,boycosww ,zenoqft ). The expansion yielding the Markov approximation (III.14) can be systematically implemented (boyww ,boycosww ,zenoqft ) as an expansion in time derivatives of the amplitudes, which in turn is an expansion in powers of the coupling (squared). This expansion relies on a separation of time scales: the typical scale(s) in the self-energy kernel is the inverse mass of the decaying particle $1/M$, whereas the typical scale of time evolution of the amplitude is $\propto 1/Y^{2}M$, which determines the relaxation rate. This separation can be surmised from eqn. (III.15) which is tantamount to taking $\mathcal{C}_{\phi}$ outside the integral and evaluating it at $\eta^{\prime}=\eta$. Namely, the amplitude varies very slowly on the time scale of variation of the self-energy. For vanishing coupling, the amplitudes remain constant, thus vary slowly for weak coupling, as compared to the time variation of the self-energy.

In the adiabatic approximation in an expanding cosmology, the time scales in the self energy are completely determined by the adiabatic frequencies as explicitly shown by expression (III.9). Therefore, even with expansion the time scales in the self-energy are much shorter than the relaxation time scale $\propto 1/Y^{2}$ of the decaying state. This separation, even during an expanding cosmology, but in the adiabatic approximation validates, the Markov approximation.

IV Massless fermions:

Our goal in this article is to study the decay of a heavy Higgs-like scalar field into much lighter fermions, neglecting the fermion masses. This is a suitable scenario for the standard model where the Higgs scalar can decay into all the charged leptons and quarks but for the top, and the quark and lepton masses may be safely neglected. Such scenario also includes the possibility of decay into neutrinos in the case that neutrino masses originate in Yukawa couplings to a Higgs-like scalar beyond the standard model. We postpone the study of decay into heavier fermionic degrees of freedom to a companion article. Focusing on the case of massless fermions allows a direct comparison with results in Minkowski space time, which are summarized in appendix (B). Furthermore, understanding this simpler case provides a pathway towards the more general case of massive fermions to be studied elsewhere.

For massless fermions $\omega^{\psi}_{k}(\eta)=k$, in this case the projector operators $\Lambda^{\pm}$ in (III.9) are given by eqn. (A.12) in appendix (A), and the self-energy (III.9) can be written in dispersive form as

[TABLE]

where the spectral density is given by

[TABLE]

with the result

[TABLE]

We carry out the $k_{0}$ integral in (IV.1) by introducing an upper (comoving) ultraviolet cutoff $\Lambda$ and a short time convergence factor $\eta-\eta^{\prime}\rightarrow\eta-\eta^{\prime}-i\epsilon$ with $\epsilon\rightarrow 0^{+}$ and replacing $k^{2}_{0}\rightarrow-d^{2}/d\eta^{\,{}^{\prime}\,2}$ yielding the final result for the self-energy

[TABLE]

In our analysis we will keep $\Lambda$ fixed but large and take the limit $\epsilon\rightarrow 0^{+}$ first, clearly this is the correct limit when the theory is considered as an effective field theory valid below a cutoff $\Lambda$. We note that the flat space time limit is obtained by replacing $\eta\rightarrow t$, and the frequency $\omega^{\phi}_{k}$ to be time independent (see appendix (B)).

It remains to perform the time integrals to obtain $\Gamma_{\Phi}(\eta)$ and $\int^{\eta}_{\eta_{i}}\Gamma_{\Phi}(\eta^{\prime})\,d\eta^{\prime}$ given by eqn. (III.17). The total time derivative in (IV.4) is integrated by parts and consistently with keeping the leading order in the adiabatic expansion, terms of the form $\omega^{\prime}/\omega^{2}$ are neglected since these yield higher order adiabatic corrections. In the limit $\epsilon\rightarrow 0^{+}$ for fixed $\Lambda$ we find the decay function

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

In obvious notation the contributions $I_{2b}(k,\eta),I_{3b}(k,\eta)$ are the $\Lambda$ independent terms in $I_{2,3}$ respectively. These three contributions are studied separately below, analyzing their cutoff dependent and independent terms extracting the different physics of each term.

IV.1 Analysis of $I_{1,2,3}$:

In the following analysis we will take the cutoff $\Lambda$ to be the largest of all scales, in particular $\Lambda\gg\omega_{k}(\eta)$ at all times.

$I_{1}$: $I_{1}$ vanishes identically as $\eta\rightarrow\eta_{i}$ and the oscillatory terms become negligibly small for $\Lambda(\eta-\eta_{i})\gg 1$, therefore $I_{1}$ grows to its asymptotic value

[TABLE]

very rapidly, on a time scale $\eta-\eta_{i}\simeq 1/\Lambda$. This divergent contribution corresponds to a renormalization of the amplitude and is similar to a linearly divergent renormalization in Minkowski space timezenoqft (see appendix (B)).

$I_{2}$: The technical details of the analysis of $I_{2}$ are relegated to appendix (D). The main result is that for $\Lambda\,(\eta-\eta_{i})\gg 1$

[TABLE]

where $\gamma_{E}=0.577\cdots$ is Euler’s constant and $I_{2b}(k,\eta)$ is given by equation (D.5) in appendix (D) where this contribution is analyzed in detail. We discuss this contribution in further detail in sections (V,VI) below.

$I_{3}$: With $\Lambda\gg\omega_{k}$ the argument of the sine function in the first term in eqn. (IV.8), namely in $I_{3a}(\Lambda,k,\eta)$, simplifies to $\Lambda\,(\eta^{\prime}-\eta_{1})$, therefore

[TABLE]

Defining $\sigma=\Lambda\,(\eta^{\prime}-\eta_{1})\leavevmode\nobreak\ \leavevmode\nobreak\ ;\leavevmode\nobreak\ \leavevmode\nobreak\ \sigma_{f}=\Lambda\,(\eta^{\prime}-\eta_{i})$, and taking the limits $\Lambda\rightarrow\infty\leavevmode\nobreak\ \leavevmode\nobreak\ ;\leavevmode\nobreak\ \leavevmode\nobreak\ \sigma_{f}\rightarrow\infty$, the integral over $\eta_{1}$ in eqn. (IV.11) becomes

[TABLE]

therefore in this limit we find

[TABLE]

where we used $\omega^{\phi}_{k}(\eta)=C(\eta)E^{\phi}_{k}(t)=m_{\phi}\,C(\eta)\,\gamma_{k}(t)$ and $C(\eta)d\eta^{\prime}=dt^{\prime}$, with $\gamma_{k}(t)=\sqrt{1+k^{2}_{p}(t)/m^{2}}$ being the Lorentz factor whose time dependence is a consequence of the cosmological redshift.

In appendix (E) we provide the analysis for $I_{3b}$, gathering both terms we find that

[TABLE]

where $\mathcal{S}(k,\eta^{\prime})$ is given by (E.25) with asymptotic limit $\mathcal{S}(k,\eta^{\prime})\rightarrow 1$ for large $\eta^{\prime}$. Therefore $I_{3}=I_{3a}+I_{3b}$ does not depend on $\Lambda$ in the limit $\Lambda\rightarrow\infty$. This is similar to the case in Minkowski space time (see appendix (B) where the equivalent term is called $T_{3}(k,t)$, eqn. (B.7)).

V Renormalization: dynamics of “dressing”

The final result for the decay function in (IV.5), $I(\Lambda,k,\eta)$ is given by ,

[TABLE]

where $I_{fin}(k,\eta)$ is independent of the cutoff $\Lambda$ in the limit $\Lambda\rightarrow\infty$, and for $(\eta-\eta_{i})\gg 1/\Lambda$ it is given by

[TABLE]

The linear and logarithmic dependence on the cutoff $\Lambda$ are exactly the same as in Minkowski space timezenoqft , as obtained in the appendix (B). This similarity is expected as the cutoff dependence arises from the short distance behavior of the self-energy correction which should be insensitive to the curvature of space time. As discussed in ref.zenoqft the origin of this divergence is the “dressing” of the bare single particle state by a cloud of fermion-anti-fermion pairs into a renormalized quasiparticle state. In a renormalizable theory the growth of the density of states at high energy implies that this cloud of excitations contains high energy states. The dynamical build-up of the cloud of excitations occurs on a time scale $\eta-\eta_{i}\simeq 1/\Lambda$ at which the divergent contributions to $I_{1,2}$ saturate, see eqn. (IV.6) and the discussion in appendix (D).

The “dressing” of the bare into the physical renormalized quasiparticle state is accounted for by the wave-function renormalization of the amplitudezenoqft . For large cutoff scale $\Lambda$ and for a weakly coupled theory with $Y^{2}\ll 1$ there is a wide separation betweeen the time scales of formation of the dressed renormalized state $\eta-\eta_{i}\simeq 1/\Lambda$, the time scale of typical oscillations $\eta-\eta_{i}\simeq 1/\omega^{\phi}_{k}(\eta)$ and finally the decay time scale $\eta-\eta_{i}\propto 1/Y^{2}\omega^{\phi}_{k}(\eta)$, which for weak coupling is the longest scale. Therefore, we can evolve the initial state in time up to an intermediate time scale $\eta_{b}$ with $(\eta_{b}-\eta_{i})\gg 1/\Lambda$, but much smaller than the typical decay time scale $\propto 1/Y^{2}\omega^{\phi}_{k}(\eta_{i})$, so that the initial state had enough time to be “dressed” by fermion-antifermion pairs into the renormalized quasiparticle state, but did not have time to decay. For example, taking $\eta_{b}-\eta_{i}=1/\omega^{\phi}_{k}(\eta_{i})$ fulfills the conditions of time scale separation because $\omega^{\phi}_{k}\ll\Lambda$, and because for $Y^{2}\ll 1$ there will be many oscillations of the field before it decays. Taking this renormalization scale is tantamount to an “on-shell” renormalization scheme. We identify $\eta_{b}$ as the time of formation – or “birth” – of the “dressed” or quasiparticle statezenoqft , which after formation decays on a much longer time scale.

The time evolution of the “bare” single particle state until it is renormalized or “dressed” is implemented by the following procedure. Writing

[TABLE]

where, taking $(\eta_{b}-\eta_{i})\gg 1/\Lambda$, the subtracted quantity

[TABLE]

is independent of $\Lambda$ for $\eta>\eta_{b}$ and $\Lambda(\eta_{b}-\eta_{i})\gg 1$. The subtracted contributions $I_{2b}(k,\eta,\eta_{b})\leavevmode\nobreak\ ;\leavevmode\nobreak\ I_{3S}(k,\eta,\eta_{b})$ are defined as follows

[TABLE]

and are obtained explicitly in appendices (D,E) respectively. During (RD) we find (see appendix (D) for definitions and eqn. (D.12))

[TABLE]

with

[TABLE]

$J(\xi^{\prime})$ is given by eqn. (D.9) in appendix (D), and

[TABLE]

where $\mathcal{S}(\eta)$ is given by eqn. (E.25) in appendix (E). The contribution from $I(\Lambda,k,\eta_{b})$ is absorbed into wave-function renormalization $Z$ as follows. Writing equation (III.17) as

[TABLE]

where the renormalized probability is given by

[TABLE]

The exponent in the wave function renormalization $Z(\eta_{b})$ is given by

[TABLE]

yielding an ultraviolet divergent wave function renormalization. The renormalized probability obeys

[TABLE]

The decay function that describes the time evolution of the renormalized survival probability is given by

[TABLE]

it is finite and independent of $\Lambda$ in the large cutoff limit. The time scale $\eta_{b}$ acts as a renormalization scale, obviously the survival probability $\mathcal{P}_{\Phi,r}(\eta)$ is independent of this renormalization scale, hence it obeys a dynamical renormalization group equation, namely

[TABLE]

The solution of this equation is, obviously222Note the similarity with the usual renormalization group function associated with the running of the wave function renormalization that yields anomalous dimensions.,

[TABLE]

$\mathcal{P}_{\Phi,r}(\eta_{b})$ describes the probability of the renormalized quasiparticle state. This “dressed” state decays with the finite and cutoff independent decay function $\int^{\eta}_{\eta_{b}}\Gamma_{\Phi}(\eta^{\prime})d\eta^{\prime}$ on time scales much longer than the “dressing” or renormalization scale $\eta_{b}$.

In the following analysis we will drop the subscript $r$ from $\mathcal{P}_{\Phi,r}$ to simplify notation since we will be strictly dealing with the renormalized survival probability.

The decay function (V.13) depends explicitly on the initial time $\eta_{i}$ (see explicit expressions in appendix (D)). However, $\mathcal{P}_{\Phi,r}(\eta_{b})$ is defined at the renormalization scale $\eta_{b}$ and it is taken to be the initial probability of the fully renormalized state after all the short time transient dynamics that result in the “dressing” of the bare into the renormalized quasiparticle state have subsided. Therefore, the dependence of the contributions (V.6,V.8) on $\eta_{i}$ must be traded for a dependence on $\eta_{b}$.

Let us write

[TABLE]

with $\beta\gg 1$ so that the $\Lambda$ dependent terms in $I_{1,2}$ reached their asymptotic behavior. For example, the “on-shell” renormalization scheme corresponds to $\beta\equiv\Lambda/\omega_{k}(\eta_{i})$. Therefore, in terms of the Hubble rate and the physical cutoff $\Lambda_{ph}(\eta_{i})=\Lambda/C(\eta_{i})$ at the initial time $H(\eta_{i})$ we find in (RD)

[TABLE]

Since the cutoff scale $\Lambda$ is taken to be much larger than any of the energy scales and the adiabatic condition requires that $H(\eta)/E_{k}(\eta)\ll 1$ at all times, it follows that $H(\eta_{i})/\Lambda_{ph}(\eta_{i})\ll H(\eta)/E_{k}(\eta)\ll 1$. Furthermore, we find that

[TABLE]

the second term in the bracket is at most of first adiabatic order, this is the case for the “on-shell” renormalization scheme for which $\beta\,\omega_{k}(\eta_{b})/\Lambda=1$. Hence, to leading adiabatic order we can safely replace $\omega_{k}(\eta_{i})\rightarrow\omega_{k}(\eta_{b})$ in the expressions. Using the results of appendix (D) we find that similar arguments justify the replacement $\gamma_{k}(\eta_{i})\rightarrow\gamma_{k}(\eta_{b})$ along with $\eta_{i}\rightarrow\eta_{b}$ in all the quantities that enter in the decay function. In the limit of large cutoff $\Lambda$ the trade-off between the variables at the initial time and those at the renormalization scale $\eta_{b}$ does not depend on the cutoff as it must be for a consistent effective field theory description well below the cutoff scale. We note that the adiabatic approximation plays an important role in this separation and is a necessary ingredient because the frequencies depend on time unlike in Minkowski space-time. In particular for the “on-shell” renormalization scheme

[TABLE]

because the adiabatic condition (during (RD)) corresponds to $\omega_{k}(\eta_{i})\,\eta_{i}\gg 1$ (see equation (II.57)).

VI Dynamics of decay.

Once we have absorbed the ultraviolet divergences into a renormalization of the amplitude, we now proceed to analyze the main physical aspects of the decay dynamics leveraging the adiabatic approximation.

VI.1 Decay during radiation domination:

We assume that the decaying particle has been produced early during the (RD) stage by some (unspecified) particle physics process at a high energy/temperature scale, focusing first on the dynamics of decay during this era. The subtracted decay function $I_{S}(k,\eta,\eta_{b})$ (V.4) can be written in a compact manner amenable to a numerical study as

[TABLE]

with

[TABLE]

and

[TABLE]

where $\xi,W[\xi]$ are defined in eqn. (V.7), and we have introduced the following functions (see appendices D,E)

[TABLE]

where $J[\xi]$ is defined in eqn. (D.9) in appendix (D). To leading adiabatic order (see appendix (E))

[TABLE]

where $Si[x]$ is the sine-integral function (see equation (E.26) and discussion in appendix (E)).

We highlight that the contribution $I^{R}_{S}$ is a distinct feature of the renormalizable Yukawa interaction and of the fermionic density of states, whereas $I_{3S}$ in (V.8) is very similar to the decay function found in the scalar case studied in ref.herring .

As discussed above, to leading adiabatic order we set $\eta_{b}=\eta_{i}$ in $I_{3S}$ and obtain (see appendix (E))

[TABLE]

[TABLE]

where $\widetilde{I}_{3S}(k,\eta,\eta_{b})$ must be obtained numerically.

However, before we engage in a numerical study we analyze the different contributions to extract a physical picture of which terms dominate at different time scales. In order to analyze the behavior in the different regimes, we write the Lorentz factor both in terms of the variable $\xi=\frac{\eta}{\eta_{i}}-1$ (see appendix (C)) as well as in terms of comoving time with the equivalence $1+\xi\equiv\sqrt{t/t_{i}}$ (see also appendix (C)),

[TABLE]

where $t_{nr}$ is the comoving time scale at which the decaying particle becomes non-relativistic, given by

[TABLE]

Whence the limits

[TABLE]

Let us focus first on the contribution $I^{R}_{S}(k,\eta,\eta_{b})$ given by (VI.2). In Minkowski space-time the frequencies are time independent, therefore $W[\xi^{\prime}]=1$ and $J(\xi^{\prime})=(\omega_{k}-k)\eta_{i}\,\xi^{\prime}$. The analysis of appendix (B) shows that in Minkowski space time for $\xi\gg 1$ the second term in (VI.2), namely $F_{1}-F_{2}$, yields $2\ln[\xi/\xi_{b}]+\mathrm{constant}$, thereby cancelling the logarithmic time dependence of the first term (see appendix (B)). Such cancellation only occurs during a limited interval in time in the expanding cosmology as a consequence of the time dependence of the frequencies. This follows from the analysis of appendix (D) which shows that there are three distinct stages:

i) $\xi\lesssim\xi_{m}$: where $\xi_{m}$ given by (D.15,D.16) is the time scale at which $F_{2}[\xi,\xi_{b}]$ reaches a maximum. During this interval $F_{1}-F_{2}$ in (VI.2) is negligible and $I^{R}_{S}\simeq 2\ln[\xi/\xi_{b}]$.

ii) $\xi_{m}<\xi\lesssim\gamma_{i}$: during this interval the function $F_{1}[\xi,\xi_{b}]$ continues to rise monotonically whereas $F_{2}[\xi,\xi_{b}]$ oscillates around its constant asymptotic value $F_{2}[\xi,\xi_{b}]\simeq F_{2}[\xi_{m},\xi_{b}]\simeq 2\ln[\xi_{m}/\xi_{b}]$, a behavior summarized by figure (7) and equation (D.19) in appendix (D).

For $\omega_{i}\eta_{i}\gg 1$ the results (D.15,D.16) show that $\xi_{m}\ll\gamma_{i}$ for all values of $\gamma_{i}\geq 1$. Therefore, for $\gamma_{i}\gg 1$, during the interval $\xi_{m}\leq\xi<\gamma_{i}$ it follows that $W[\xi^{\prime}]\simeq 1$ and $F_{1}\simeq 2\,\ln[\xi/\xi_{b}]$ thereby (approximately) cancelling the logarithm from the first term in $I^{R}_{S}$, whereas $F_{2}\simeq 2\ln[\xi_{m}/\xi_{b}]$ remains constant, yielding a plateau in $I^{R}_{S}$. This approximate cancellation is effective during a time interval that increases for $\gamma_{i}\gg 1$ (see discussion in appendix (D)). According with eqn. (VI.9) and the limits (VI.11) during this interval, wherein $I^{R}_{S}$ is approximately constant, the decaying particle is in the ultrarelativistic regime. In this stage the constancy of $I^{R}_{S}$ is expected because in the ultrarelativistic regime the frequencies are nearly time independent since $\omega_{k}(\eta)\simeq k\simeq\omega_{i}$. Therefore $W[\xi]\simeq 1$ yielding $F_{1}\simeq 2\ln[\xi/\xi_{b}]$ thereby cancelling the logarithmic time dependence of the first term in (VI.2), similarly to Minkowski space-time.

If $\gamma_{i}\gg 1$ the decaying particle is “born” ultrarelativistically and there is a (long) time window $\xi_{m}<\xi<\gamma_{i}$ within which $\sqrt{W[\xi^{\prime}]}\simeq 1$ and $F_{1}[\xi,\xi_{b}]\simeq 2\ln[\xi/\xi_{b}]$ thereby approximately cancelling the first term in $I^{R}_{S}$ whereas $F_{2}[\xi,\xi_{b}]$ remains nearly constant. Therefore for $\gamma_{i}\gg 1$ it follows that $I^{R}_{S}(k,\eta,\eta_{b})$ rises rapidly on a time scale $\simeq\xi_{m}$ reaching a maximum and remaining nearly constant $I^{R}_{S}\simeq 2\,\ln[\xi_{m}/\xi_{b}]$ until $\xi\simeq\gamma_{i}$.

**iii) $\xi\gg\gamma_{i}$: ** The cosmological redshift eventually makes the decaying particle to become non-relativistic when $\xi\gg\gamma_{i}\gg 1$. During this stage the particle is non-relativistic as a consequence of the cosmological redshift. The time dependence of the frequency now yields $\sqrt{W[\xi^{\prime}]}+1/\sqrt{W[\xi^{\prime}]}\gg 2$, hence $F_{1}>2\,\ln[\xi]$. In this stage it follows that $W[\xi]\approx\xi/\gamma_{i}$, therefore for $\xi\gg\gamma_{i}\gg 1$ we find that $F_{1}[\xi]\simeq 2\sqrt{\xi/\gamma_{i}}$ and $F_{2}[\xi,\xi_{b}]\simeq 2\ln[\xi_{m}/\xi_{b}]$. For $\xi\gg\gamma_{i}$, the integral for $F_{1}[\xi,\xi_{b}]$ is estimated by splitting it into the stages $\xi_{b}\leq\xi\leq\gamma_{i}$ and $\xi>\gamma_{i}$. The first stage yields $2\ln[\gamma_{i}/\xi_{b}]$ since during this (ultrarelativistic) stage $W[\xi^{\prime}]\simeq 1$, and the second yields (approximately) $2\sqrt{\xi/\gamma_{i}}$ since during this (non-relativistic) stage $W[\xi]\approx\xi/\gamma_{i}$.

In summary, for a particle that is “born” ultrarelativistically, namely with $\gamma_{i}\gg 1$, the contribution $I^{R}_{S}$ rises rapidly up to a value $\simeq 2\,\ln[\xi_{m}/\xi_{b}]$ on a time scale $\xi_{m}\ll\gamma_{i}$ given by (D.16), remains nearly constant up to a time scale $\xi\simeq\gamma_{i}$ at which the particle becomes non-relativistic, and begins to fall-off as $-2\sqrt{\xi/\gamma_{i}}$ for $\xi\gg\gamma_{i}$.

In the opposite limit when $\gamma_{i}\simeq 1$ the decaying particle is non-relativistic already at the initial time and $\omega_{k}(\eta)\simeq m_{\phi}\,C(\eta)$. In this case $F_{2}[\eta,\eta_{b}]$ saturates rapidly, on a scale $\xi_{m}\simeq\pi/\omega_{i}\,\eta_{i}\ll 1$, and $F_{1}[\eta,\eta_{b}]$ grows faster than logarithmically, hence $F_{1}-F_{2}$ becomes larger than the logarithm in the first term of $I^{R}_{S}$ and negative. This behavior leads to an early suppression of decay.

This analysis is approximately summarized during the ultrarelativistic (UR) and non-relativistic (NR) regimes, by (see eqn. (D.19) in appendix (D)),

[TABLE]

The main aspects of this analysis are confirmed by a numerical study summarized in figures (1) and (2) for $\gamma_{i}=2,10$ respectively. Notice the different scales in the figures highlighting the emergence of the plateau and the crossover to a diminishing (negative) square root behavior at a scale $\xi\simeq\gamma_{i}$.

Decay at rest: for a very massive particle “born” and decaying at rest in the comoving frame, namely for $\gamma_{i}=1$, and $\omega_{i}\eta_{i}\gg 1$ we can provide an analytic form of the decay function for time scales $\xi\gg\xi_{b}\simeq 1/\omega_{i}\eta_{i}$ for on-shell renormalization. As discussed in appendices (D,E), $F_{2}[\eta,\eta_{b}]$ reaches its asymptotic limit on a time scale $\xi\simeq\pi/2\omega_{i}\eta_{i}\ll 1$ (see equation (D.15) in appendix (D)). Furthermore, the function $\mathcal{S}(\xi^{\prime})$ in (VI.3) reaches its asymptotic value $\mathcal{S}\simeq 1$ at a time scale $\xi^{\prime}\simeq\pi/\omega_{i}\eta_{i}\ll 1$. Therefore for $\xi^{\prime}\gg 1/\omega_{i}\eta_{i}$ we can neglect the contribution from $F_{2}$ and set $\mathcal{S}(\xi^{\prime})=1$ in (VI.3), hence $I_{3S}[k,\eta,\eta_{b}]$ is given by the first term in eqn. (VI.7) with $\gamma_{i}=1$ and multiplied by a factor $2$ to account for $\mathcal{S}=1$. Gathering all terms we find in this case ($\gamma_{i}=1\leavevmode\nobreak\ ;\leavevmode\nobreak\ \omega_{i}\eta_{i}\gg 1\leavevmode\nobreak\ ;\leavevmode\nobreak\ \xi\gg 1/\omega_{i}\eta_{i}$),

[TABLE]

where we have neglected a constant term of $\mathcal{O}(1)$. This expression displays all the features described above. Note that for $\xi\ll 1$ the logarithmic time dependence cancels out, but for $\xi\gg 1$ the first logarithm in (VI.13) continues to grow, however the negative square root eventually dominates the contribution of the first terms within brackets. These are precisely the terms arising form the renormalization and their time dependence is a consequence of the time dependence of the frequencies.

To compare to the decay law in Minkowski space time it is convenient to cast the result (VI.13) in terms of comoving time, using $1+\xi=\eta/\eta_{i}$, with $\eta=\sqrt{2t/H_{R}}$ (see eqn. (II.39) valid in (RD)), and the relation

[TABLE]

Setting $\eta_{i}=\eta_{b}$ to leading adiabatic order, we find for $\gamma_{i}=1$ and $t\gg t_{b}$

[TABLE]

leading to the survival probability for $t\gg t_{b}$

[TABLE]

This is one of the important results of this study. Remarkably $\Gamma_{0}$ is the same as the decay width at rest in Minkowski space time, however the power law with “anomalous dimension” $Y^{2}/8\pi^{2}$ and the stretched exponential with the power law $(t/t_{b})^{1/4}$ are a consequence of the renormalization and the time dependence of the frequencies, manifestly a consequence of the expanding cosmology. The combined effect of these two terms yields a slowing down of the decay as compared with the case of Minkowski space time with a concomitant enhancement of the lifetime of the decaying particle as compared to Minkowski space-time. This is a noteworthy result: as a consequence of the cosmological expansion the contribution from the renormalization and quasiparticle formation *slows down the decay * leading to an enhancement of the lifetime of the initial state.

Decay of particles with $\gamma_{i}\gg 1$: these are particles that are “born” ultrarelativistically. For $\gamma_{i}\gg 1$ the contribution from $I^{R}_{S}(k,\eta,\eta_{b})$ has been summarized by eqn (VI.12) and is displayed in fig. (2): a rapid rise on a time scale $\xi_{m}\ll\gamma_{i}$ given by (D.16) up to $I^{R}_{S}\simeq 2\,\ln[\xi_{m}/\xi_{b}]$ followed by a near plateau during the stage while $\xi\lesssim\gamma_{i}$. This contribution falls off slowly as $-\sqrt{\xi/\gamma_{i}}$ during the non-relativistic stage, $\xi\geq\gamma_{i}$ (see eqn. (VI.12)). While a quantitative analysis of $I_{3S}$ requires a numerical study, we can obtain a fairly accurate estimate as follows. The contribution from $\mathcal{S}$ to $I_{3S}$ (see equation (VI.3)) is discussed in appendix (E), and can be approximately summarized as: $\mathcal{S}\approx 0$ for $\xi<\xi_{s}$ and $\mathcal{S}(\eta)\approx 1$ for $\xi>\xi_{s}$ with $\xi_{s}$ given by (E.28,E.29).

With $\gamma_{i}\gg 1$, the ultrarelativistic stage corresponds to $\gamma_{i}\gg\xi$, during the stage $\gamma_{i}\gg\xi_{s}\gg\xi$ it follows that $\mathcal{S}\approx 0$, using $1+\xi=\sqrt{t/t_{i}}$, and eqns. (VI.10,VI.14), during this stage $I_{3S}$ is given in comoving time $t$ by

[TABLE]

where

[TABLE]

also describes the decay function in the case of a scalar field decaying into two massless scalarsherring . During this stage for $t\ll t_{nr}$ we find

[TABLE]

For $\gamma_{i}\gg\xi\gg\xi_{s}$ it follows that $\mathcal{S}\simeq 1$, therefore the above result is multiplied by a factor 2. Hence, during the ultrarelativistic stage with $\gamma(t)\gg 1$, or $t\ll t_{nr}$, and $\mathcal{S}=1$ in (VI.3), it follows that

[TABLE]

which when combined with the result (VI.12) yields in this ultrarelativistic regime, for $\gamma_{i}\gg\xi\gg\xi_{s},\xi_{m}$

[TABLE]

Neglecting the perturbatively small non-secular constant in the decay function from the first term in (VI.21)333Or absorbing it in a finite perturbatively small time independent wave function renormalization of $\mathcal{P}_{\Phi}$. , we find in the time interval for $t_{nr}\gg t\gg t_{b}$ during which the decaying particle is ultrarelativistic and the transient dynamics of quasiparticle formation has saturated

[TABLE]

We can now use the property (V.15) and write for $t>t_{nr}$

[TABLE]

where

[TABLE]

After the decaying particle becomes non-relativistic for $\xi\gg\gamma_{i}$ or $t\gg t_{nr}$ when $\gamma(t)\simeq 1$, the contribution $\mathcal{S}\simeq 1$ and $I_{3S}(\xi)-I_{3S}(\xi_{nr})$ becomes

[TABLE]

the dots in the above expression stand for terms of higher order in the ratio $t_{nr}/t$.

Finally, combining with the result given by eqn. (VI.12), the total decay function after the particle has become non relativistic $\xi\gg\gamma_{i}$ (or $t\gg t_{nr}\gg t_{b}$) is given in comoving time by

[TABLE]

where we have neglected a perturbatively small constant term and approximated $t_{i}\gamma^{2}_{i}\simeq t_{nr}$ for $\gamma_{i}\gg 1$. Hence for $t\gg t_{nr}\gg t_{b}$ we find

[TABLE]

It would be expected that after $t_{nr}$ when the particle has become non-relativistic as a consequence of the cosmological redshift, the time evolution of the survival probability would be similar to that of a particle born and decaying at rest. However, the result (VI.27) features an extra power law with exponent $\Gamma_{0}t_{nr}/2$ as compared to the decay function for the particle born at rest, eqn. (VI.16). This difference reflects the memory of the past evolution in the form of the integral (VI.24).

We can provide a measure of the impact of curved space time effects on the decay function by comparing the results above to a phenomenological, S-matrix inspired Minkowski decay law allowing for a local time dilation factor to account for the cosmological redshift, namely

[TABLE]

where $\Gamma_{0}=\frac{Y^{2}\,m_{\phi}}{8\pi}$ is the decay width at rest in Minkowski space time, and $\gamma(t)$ the local Lorentz factor (VI.9). The comparison to the cutoff independent subtracted decay function (VI.1) is facilitated by introducing

[TABLE]

so that the Minkowski-like decay function is given by

[TABLE]

where a factor is included in (VI.30) to ensure that $I_{M}(t_{i}=t_{b})=0$ consistently with the subtraction definining (VI.1). For $t\gg t_{i}$ this phenomenological decay function is interpreted as that of Minkowski space-time but with the instantaneous Lorentz time dilation factor. For $t\gg t_{i}$ it provides a “benchmark” to compare the results obtained above for the decay function to an S-matrix inspired instantaneous Minkowski decay law.

Before we engage in a numerical comparison, it is illuminating to analyze the cases discussed above.

Non-relativistic: $\gamma(t)=1$. For this case the $I_{S}(t)$ is given by (VI.15), the last term of which is precisely $I_{M}(t)$ for $\gamma(t)=1$. The first two terms in (VI.15) yield a negative contribution for $t\gg t_{b}=t_{i}$, therefore the cosmological decay function is smaller in this case than the phenomenological Minkowski function, leading to a longer lifetime.

Ultra-relativistic: $\gamma_{i}\gg 1$. During the ultrarelativistic regime $\gamma(t)\gg 1$ ($t\ll t_{nr}$), taking the time large enough so that the transient build-up of $\mathcal{S}$ in eqn. (VI.3) has saturated, the cosmological decay function is given by (VI.21) whereas I_{M}(t)\simeq\pi m_{\phi}t\,\big{(}t/t_{nr}\big{)}^{1/2}. The logarithmic term in (VI.21) could be fairly large for large $\gamma_{i}$ thereby yielding $I_{S}(t)>I_{M}(t)$ during a time interval. This can be understood from the following argument.

As discussed above and in appendix (D), for $\gamma_{i}\gg 1$ the contribution $I^{R}_{S}$ (see eqn. (VI.2)) rises on a time scale $\xi_{m}\simeq(3\pi\gamma^{2}_{i}/\omega_{i}\eta_{i})^{1/3}$ up to a maximum $\simeq 2\ln(\xi_{m}/\xi_{b})$ after which it remains nearly constant up to $\xi\simeq\gamma_{i}$ yielding the logarithmic term in (VI.21). For example, for $\gamma_{i}\simeq 200\leavevmode\nobreak\ ,\leavevmode\nobreak\ \omega_{i}\eta_{i}\simeq 100$ and “on-shell” renormalization with $\xi_{b}=1/\omega_{i}\eta_{i}$, the contribution from $I^{R}_{S}$ rises up to a value $\simeq\frac{4}{3}\ln[\sqrt{3\pi}\gamma_{i}\omega_{i}\eta_{i}]\simeq 14.7$ on a comoving time scale $t_{m}/t_{i}\approx 240$. Since the Hubble time scale $1/H(t)=2t$ during (RD), it follows that $I^{R}_{S}$ rises up to the plateau over $\simeq 240$ Hubble times, with the possibility that during this time $I_{S}(t)>I_{M}(t)$. However, after the particle has become non-relativistic, namely for $t\gg t_{nr}$, the cosmological decay function $I_{S}(t)$ is given by (VI.26) whereas

[TABLE]

showing that $I_{S}(t)\ll I_{M}(t)$ for $t\gg t_{nr}$. This suggests a crossover behavior for very large values of $\gamma_{i}$: there is an early time window during the ultrarelativistic stage wherein the cosmological decay function may be larger than the Minkowski one, however as the decaying particle eventually becomes non-relativistic the latter will ultimately dominate. This behavior is borne out by a detailed numerical study.

Figures (3,4,5) show a comparison between the phenomenological Minkowski decay function (VI.29), the total contribution $I_{S}$ (VI.1) along with $I_{3S}$ (VI.3) for on-shell renormalization with $\omega_{i}\eta_{i}=100$ and $\gamma_{i}=10,50,200$ respectively. For these values the transition time to the non-relativistic behavior is $t_{nr}/t_{i}\simeq 10^{2},2.5\times 10^{3},4\times 10^{4}$ respectively. For $\gamma_{i}=10,50$ figs. (3,4) show that $I_{S}$ and $I_{3S}$ are nearly indistinguishable, namely $I^{R}_{S}$ (VI.2) is subleading in these cases, and that the phenomenological $I_{M}$ is always larger than $I_{S}$.

However, for $\gamma_{i}=200$ fig. ( 5) shows that the contribution from $I^{R}_{S}$ dominates at early time, rising on a time scale $t/t_{i}\simeq 100$. In this case $I_{M}$ is smaller than $I_{S}$ during a substantial time window, $\approx 500$ Hubble times from the “birth” of the quasiparticle, before crossing over to becoming the largest decay function.

Therefore we conclude that in the ultrarelativistic case, for very large values of $\gamma_{i}$, the decay function is larger than the phenomenological Minkowski one within a substantial time interval but eventually becomes smaller at a time scale that depends on the various parameters. In either case, at long time the decaying particle lives longer than predicted by a Minkowski decay law extrapolated to the expanding cosmology. This is a generic result: after an intermediate time scale that depends on $\gamma_{i}$, the cosmological decay function is smaller than the phenomenological Minkowski-like one. Therefore the S-matrix inspired phenomenological Minkowski decay law under estimates the lifetime of the decaying particle.

VI.2 Long lived particles: decay during matter domination or beyond.

The discussion above focused on decay during the radiation dominated era that lasts until $C(\eta)=a_{eq}\simeq 10^{-4}$ corresponding to an ambient temperature $T\simeq\mathrm{eV}$ at a time $t_{eq}\approx 10^{12}\,\mathrm{secs}$. If the decaying particle is very long lived as would befit a dark matter candidate, it would continue to decay during the matter and perhaps dark energy dominated eras. This case corresponds to an extremely small Yukawa coupling, which allows to safely neglect early transient effects that saturate at early times. The general form of the decay function after renormalization is given by eqns. (V.13,V.4). Under the assumption of very weak Yukawa coupling we can neglect the contribution from the cosine term in $I_{2b}$, eqn. (V.6) (the contribution $F_{2}$ in eqn. (VI.5)) and we can set $\mathcal{S}=1$ in eqn. (V.8). This is because both terms saturate on short time scales therefore they yield perturbatively small corrections to the decay function for very weak Yukawa coupling as compared to the terms that continue to grow in time. Hence, neglecting these perturbatively small transient contributions for very weak Yukawa couplings, the decay function simplifies to

[TABLE]

where $\xi=(\eta-\eta_{i})/\eta_{i}$ and $F_{1}$ is given by (VI.4) and the dots in (VI.32) stand for constant terms that are of $\mathcal{O}(Y^{2})$.

For general scale factor $W[\xi]$ is given by

[TABLE]

Let us analyze each term separately in order to understand their behavior at long time during the (MD) era, taking as an upper bound $C(\eta)\simeq\mathcal{O}(1)$, or, upon using eqn. (II.37) $\eta\simeq\sqrt{a_{eq}}/H_{R}$. With “on-shell” renormalization ($\xi_{b}=1/\omega_{i}\eta_{i}$), we find

[TABLE]

Taking the initial time to correspond to an initial temperature $10^{15}\,\mathrm{GeV}$ yields $C(\eta_{i})\simeq 10^{-28}$ therefore the logarithm contribution to the decay function for $\eta\simeq\sqrt{a_{eq}}/H_{R}$ yields

[TABLE]

Obtaining the contribution from $F_{1}$ over the whole history from early (RD) into (MD) can be done numerically, although this is a rather challenging task because of the enormous dynamic range with the scale factor varying over twenty four orders of magnitude. However, we can provide a simple estimate of the remaining two terms of the decay function at long time during the (MD) era and/or beyond. If the particle remains ultrarelativistic, then as discussed in the previous sections the contribution from $F_{1}$ cancels the logarithmic time dependence of the first term, hence the combination of the first two terms saturates (this is the plateau in fig. (2)) and yields a perturbatively small time independent contribution to the decay function. Hence during this ultrarelativistic stage the last term in (VI.32) dominates the decay function.

After the particle has become non-relativistic then $W[\xi]\simeq C(\eta)/\gamma_{i}C(\eta_{i})\gg 1$ and

[TABLE]

during (MD) using eqn. (II.35) and taking as an upper bound $\eta\simeq\sqrt{a_{eq}}/H_{R}$ we find

[TABLE]

Finally, we can estimate the last term in (VI.32) during the stage when the particle is non-relativistic and (MD) dominated, taking $\gamma(\eta^{\prime})\simeq 1$ and taking $\eta\simeq\sqrt{a_{eq}}/H_{R}$, we find

[TABLE]

Since during the ultrarelativistic stage the time dependence of the first and second term cancel out and the last term dominates the decay dynamics, we conclude that the last term in (VI.32) dominates the decay dynamics of a very long-lived particle with very weak Yukawa coupling, all throughout the time evolution. Since the first (logarithmic) term is always subdominant, and the second term is negative, larger in magnitude than the logarithmic term but also subdominant at late time, the last term in (VI.32) yields an upper bound to the decay function throughout all the expansion history. It can be written as a function of the redshift by recalling that $C(\eta)\,d\eta=dt$, and using $dt=da/(a\,H(a))$ with $H(a)$ the Hubble expansion rate given by eqn. (II.30). Writing the local Lorentz factor as \gamma(a(t)))=\Big{[}\frac{a^{2}_{nr}}{a^{2}(t)}+1\Big{]}^{1/2}\leavevmode\nobreak\ ;\leavevmode\nobreak\ a_{nr}\equiv k/m_{\phi}, we find that the upper bound to the decay function at redshift $z$ is given by

[TABLE]

where $\Gamma_{0}=Y^{2}m_{\phi}/8\pi$ is the decay rate at rest in Minkowski space time, and

[TABLE]

depends solely on the cosmological parameters and $a_{nr}=k/m_{\phi}$ the scale factor at which the decaying particle transitions from ultrarelativistic to non-relativistic, and we have taken $z_{b}\gg 1$. In particular this result for the decay function is insensitive to the early transient dynamics.

The redshift evolution of the survival probability all throughout the expansion history is summarized concisely as

[TABLE]

The inequality in eqn. (VI.41) reflects that eqn. (VI.39) yields an upper bound to the decay function. For $a_{nr}=0$, namely when the decaying particle is “born” at rest, it follows that $\Upsilon(z,z_{b})=H_{0}(t-t_{b})$ independently of the cosmology, and we can compare the result (VI.39) for $a_{nr}=0$ to the case of the particle decaying at rest given by eqn. (VI.15) valid during the (RD) era. The discussion on dominant terms above clarifies that the last term in (VI.15) dominates the decay dynamics, whereas the first two terms combine into a negative contribution which becomes subleading at long time for very weak Yukawa couplings. Hence it is clear that for very weak Yukawa coupling and long time (VI.39) becomes the leading contribution and yields an upper bound to the decay function for long-lived particles decaying at rest. Furthermore, for $a\ll a_{nr}$, namely when the decaying particle is ultrarelativistic and taking this regime to be during the (RD) era with $a\propto t^{1/2}$ it follows that

[TABLE]

in agreement with the decay law (VI.22) during the ultrarelativistic regime in (RD). This analysis confirms the validity of the decay law (VI.41) with (VI.40) as an upper bound to describe the evolution of the survival probability for very weakly coupled, long lived particles all throughout the cosmological evolution, under the assumption that the fermionic decay products can be considered massless in the decay process.

VII Discussion

The final form of the renormalized decay function, eqn. (VI.1) describing the time evolution of the survival probability of the quasiparticle state is amenable to a straightforward numerical study. The analysis of section (VI) reveals a very rich dynamical evolution with various different time scales. The shortest time scales describe the build-up of the quasiparticle; this early transient dynamics is absorbed into a wave function renormalization of the quasiparticle survival probability at a time scale $t_{b}$. After this short time transient there remain the time scales over which $F_{2}$ (VI.5) saturates at a constant value and $\mathcal{S}$ (VI.6) rapidly approaches $\mathcal{S}\simeq 1$. The detailed dynamics over these scales was studied analytically and numerically in appendices (D) and (E)) respectively. The evolution of the survival probability on the intermediate and long time scales becomes simpler and can be summarized succinctly. Furthermore, because the short time transients saturate to constant values, for weak Yukawa coupling the largest contributions to the decay dynamics arises from terms that are secular (grow in time) over the intermediate and long time scales.

Decay at rest in the comoving frame ($\gamma_{i}=1$): The time evolution of the survival probability is given by

[TABLE]

where $\Gamma_{0}=\frac{Y^{2}}{8\pi}\,m_{\phi}\,$ is the decay width of a particle at rest in Minkowski space-time. The power law and stretched exponentials are both a remnant of the renormalization, or “dressing” of the bare into the quasiparticle state and a distinct consequence of the cosmological redshift. Indeed, in Minkowski space-time the terms that give rise to these contribution become time independent after the transient dynamics, whereas, in curved space-time, the origin of these contributions is the time dependence of the frequencies via the cosmological redshift.

The methods that we implemented in this study, a non-perturbative formulation combined with a physically motivated adiabatic expansion including a consistent treatment of renormalization, are very different from those implemented in ref.viljafer . The decay law of a particle at rest (VII.1) is also very different from that reported in ref.viljafer . The origin of the discrepancy is not clear to us. However, since the power law and stretched exponentials originate precisely from the contributions to the renormalization of the survival probability, we suspect that the discrepancy originates in the treatment of the ultraviolet divergences. These are of the same form as in Minkowski space-time (see appendix (B) and ref.zenoqft ) as expected since these are short distance divergences, but have not been discussed or addressed in ref.viljafer . As explained above, the time dependence of the frequency yields an unexpected contribution to the decay law on longer time scales that originates in the dynamics of quasiparticle formation.

Born ultrarelativistically: if the particle is “born” or produced ultrarelativistically, namely with $\gamma_{i}\gg 1$ during (RD), an important time scale is $t_{nr}=\frac{k^{2}}{2m^{2}_{\phi}H_{0}\sqrt{\Omega_{R}}}$, which determines when the particle transitions from being ultrarelativistic ($\gamma(t)\gg 1$ or $t\ll t_{nr}$) to non-relativistic ($\gamma(t)\simeq 1$ or $t\gg t_{nr}$) as a consequence of the cosmological redshift. The dynamical evolution of the survival probability is different in these stages. a) ultrarelativistic stage: ($\gamma(t)\gg 1$, or $t_{b}\ll t\ll t_{nr}$ )

[TABLE]

b) non-relativistic stage ($t\gg t_{nr}$ or $\gamma(t)\simeq 1$),

[TABLE]

Although for $t\gg t_{nr}$ the particle has become non-relativistic because of the cosmological redshift, as compared to the case of decay at rest (VII.1), this decay law features a new power with exponent $\Gamma_{0}t_{nr}/2$. Its origin is the memory of the decay function manifest in the form of the integral of the cosmological redshift in (V.8) over the whole history of the decay process. Therefore, even well after the decaying particle has become non-relativistic, the survival probability features an enhancement factor that “knows” about the past history when the particle was ultrarelativistic. The dynamics during the transition from the ultrarelativistic to the non-relativistic behavior must be studied numerically, and the previous section shows such study for several values of the parameters.

Massless fermions vs massless bosons: Ref.herring studied the decay of a scalar into two massless scalars, therefore we can now compare the results of that study to those obtained here for the case of scalar decay into massless fermions. The main difference is in the contribution $I^{R}_{S}$ in eqn.(VI.1) which is given by eqn. (VI.2). The contribution from $I_{3S}$ to the decay function is the same for fermions and bosons, for example the function $G[x]$ is the same that enters in scalar decayherring . The extra contribution, namely $I^{R}_{S}$ has the same origin as the ultraviolet divergent contributions that are absorbed in wave function renormalization. This is also the case in Minkowski space-timezenoqft as shown in appendix (B). Whereas in Minkowski space time this contribution becomes time independent after a short time transient and is absorbed into wave function renormalization, in a FRW cosmology, it is time dependent as a consequence of the cosmological redshift and becomes important for non-relativistic particles. Namely, $I^{R}_{S}$ is a remnant of the physical process of quasiparticle formation. There is no such contribution in the case of decay into two scalars because the theory in this case is superrenormalizable, hence there is no equivalent of the $I^{R}_{S}$ term. This contribution suppresses the decay function at long time, thereby enhancing the lifetime of the decaying particle. This behavior is yet another source of discrepancy with the results of ref.viljafer , which finds a larger rate in the fermionic case. The source of this discrepancy are precisely the “anomalous” power and stretched exponential which are a consequence of the quasiparticle formation and wave function renormalization. Although the decay probability requires an ultraviolet divergent wave function renormalization even in Minkowski space-time, this seems to be an aspect missing in the treatment of ref.viljafer . The cumulative effect of these differences results in that a meaningful comparison to our study has eluded us.

“Benchmarking” the decay law: The decay laws obtained above are very different from the usual exponential decay familiar in Minkowski space time, one of the reasons for the difference being the cosmological redshift. Thus a natural question arises: would an S-matrix inspired, phenomenologically motivated exponential decay law with a time dependent Lorentz factor to account for the cosmological redshift describe even approximately the decay of the particle?. This motivates the comparison of the previous results to the following Minkowski-like decay law (in (RD))

[TABLE]

For decay at rest $\gamma(t)=1$, this decay law misses the power with anomalous dimension and the stretched exponential, whose combination is negative. Therefore the Minkowski-like decay law over-estimates the suppression of the survival probability in the case of decay at rest. For a particle that is produced ultrarelativistically, during the stage wherein $\gamma(t)\gg 1$, namely $t\ll t_{nr}$ one finds

[TABLE]

which is smaller than (VII.2). For $t\gg t_{nr}$ when the decaying particle has become non-relativistic

[TABLE]

Comparing this result with (VII.3) clearly shows that the phenomenological Minkowski decay law including the instantaneous Lorentz factor over-estimates the suppression of the survival probability, namely under estimates the lifetime of the decaying state. The discrepancies with the cosmological decay law, both the factor $2/3$ in (VII.2) along with the powers and stretched exponential in (VII.3) are traced to i) the memory of quasiparticle formation, ii) the memory of the past evolution in the integral of the time dilation factor. None of these can be captured by a phenomenological Minkowski-like decay law including an instantaneous Lorentz factor as such description has no memory of the past evolution. We draw two important conclusions from this comparison: i) a phenomenological, S-matrix inspired Minkowski decay law under estimates the lifetime of the decaying particle since it over estimates the suppression of the survival probability, ii) describing particle decay in cosmology in terms of a decay rate, even one that includes the cosmological redshift in the time dilation factor, is not only not useful but is misleading insofar as missing important physical processes and yielding a substantial under estimate of the lifetime of the decaying particle.

On initial conditions: we have taken the initial state to correspond to a single particle state of a given momentum, to compare to the usual case in S-matrix theory. The calculation of a decay rate in S-matrix theory considers the transition amplitude from “in” single particle state (prepared at time $-\infty$) to an “out” multiparticle state at time $+\infty$. Our main point is that such calculation is not meaningful in an expanding cosmology, motivating the study of the previous sections. Thus the chosen initial condition allows us to directly compare to what would be expected from S-matrix theory, appended with an exponential decay law with the rate calculated with S-matrix.

Alternative initial conditions may be considered but all imply not only technical complexities, but also conceptual aspects: a single particle but spatially localized wave packet will not only decay via the decay of the different single particle components for different momenta, but its amplitude will also decay as a consequence of dispersion and spreading. Spatially narrow wave packets will decay the fastest and separating systematically the different physical processes is, in general, not only technically daunting but implies some ambiguity as to how to extract a “decay”. Another physically motivated initial condition would be to take the initial state to emerge from the decay of a heavier particle. Obviously, such choice would have inherent ambiguities from the choice of the parent particle and its decay kinematics. These aspects notwithstanding, the framework developed in the previous sections can be simply adapted to alternative initial conditions.

Modifications to BBN? Although the results obtained in this study do not apply directly to neutron decay, since we focused on scalar decay Yukawa coupled to massless fermions, and the small phase space available for three body neutron decay is a result of the small neutron-proton mass difference, let us explore the consequences of the results on this process, with all these caveats. First: the neutron is “born” after the QCD phase transition at $T_{QCD}\simeq 150\,\mathrm{MeV}$ at a time $t_{b}\simeq 10^{-5}\,\mathrm{secs}$, because the neutron mass $M_{N}\simeq\mathrm{GeV}\gg T_{QCD}$ it is “born” at rest in the plasma. Let us identify the dimensionless coupling $Y^{2}/8\pi\equiv\Gamma_{N}/M_{N}$ where $\Gamma_{N}\simeq 10^{-3}\,\mathrm{secs}^{-1}$ is the neutron’s lifetime. Hence $Y^{2}/8\pi\simeq 10^{-21}$, and taking $t/t_{b}\simeq 1/\Gamma_{N}t_{b}\simeq 10^{8}$ we see that the power law with “anomalous” dimension and the stretched exponential correction to the usual exponential decay law in eqn. (VI.16) for decay at rest are all but negligibly small and would not affect the dynamics of neutron decay during (BBN). Of course, there are the above mentioned caveats to this conclusion which should only be taken as an extrapolation and as a gross estimate of the effects. This analysis also suggests that the corrections to the decay law are more important for particles “born” very early during (RD) and very long lived a situation that befits most descriptions of a dark matter candidate.

Emergence of a local decay law with constant S-matrix decay rate. If a measurement of the time evolution of the survival probability is carried out during a * sufficiently short* time interval $\Delta t=t_{f}-ti$ and sufficiently long after the transient dynamics has subsided, we would expect that the decay law would be nearly exponential with a nearly constant decay rate. Namely we would expect that locally during such a short time interval the survival probability is given by

[TABLE]

where $\overline{\Gamma}_{\Phi}(k)$ is a constant, related to the S-matrix rate. This law cannot emerge during the transient stage dominated by the power laws in the survival probabilities (VI.22), (VI.27). However, in fact, it emerges naturally at longer time scales after the transient dynamics becomes negligible, from the expression (VI.39) when considered during time intervals $\Delta t=t_{f}-t_{i}\ll 1/H(t_{i})$, with $H(t_{i})$ the Hubble expansion rate at the time $t_{i}$. This is seen as follows: keeping only the last term in eqn. (VI.32) (neglecting transients), and passing to comoving time with $C(\eta^{\prime})d\eta^{\prime}=dt^{\prime}$, it follows that

[TABLE]

We now expand $\gamma(t^{\prime})$ around $t_{i}$, $\gamma(t^{\prime})=\gamma(t_{i})+\gamma^{\prime}(t_{i})(t^{\prime}-t_{i})+\cdots$ and integrate to obtain,

[TABLE]

therefore we clearly see that for time intervals $\Delta t\ll 1/H(t_{i})$ the decay law features small departures from the exponential S-matrix inspired one, with corrections of order $\Delta tH(t_{i})\ll 1$. This is expected on physical grounds: on very short time scales, much shorter than the Hubble time, a local Minkowski space-time approximation is warranted by the equivalence principle. For example, accelerator experiments today measure the lifetime of standard model particles, these experiments are obviously insensitive to the Hubble time today $\simeq 13.5\,\mathrm{Gyr}$. Therefore the S-matrix-inspired calculation in Minkowski space time is warranted as it describes the measurement of lifetimes much smaller than $1/H_{0}$. In the early stages of a radiation dominated cosmology during rapid expansion, or for lifetimes comparable to the Hubble time, such S-matrix inspired approximation does not describe correctly the dynamics of decay, as discussed in detail in the previous sections.

Caveats: We have focused on studying scalar decay into massless fermion pairs, a situation that approximates most of the fermionic decay channels of a Higgs scalar in the standard model. An important aspect of this decay process is that it does not feature thresholds. Including the mass for the decay products introduces kinematic thresholds, a consequence of strict energy-momentum conservation. In ref.herring it was argued that the Hubble rate of expansion introduces a natural energy uncertainty leading to a relaxation of the kinematic thresholds, thereby allowing processes that are forbidden in Minkowski space-time by energy conservation. Furthermore, ref.zenoqft has shown that energy uncertainties associated with transient non-equilibrium aspects of the decay allow decay into heavier particles during a time interval. In an expanding cosmology these effects may combine with the energy uncertainty from Hubble expansion to enhance the decay by opening up novel channels that would be otherwise forbidden by strict energy conservation. These aspects associated with the masses of the decay products will be the subject of further study.

The inclusion of masses for the decay products becomes a more pressing issue in the case of decay of very long lived particles studied in section (VI.2) where we have extended the results obtained for the (RD) era to provide an upper bound on the decay function all throughout the expansion history. Therefore, the decay law (VI.41) with the decay function (VI.40) must be understood within the context of decay of a heavy particle into massless or nearly massless fermionic channels with the caveat that such an approximation may be of limited validity during the (MD) or (DE) eras and should be interpreted as indicative of the decay dynamics.

In this study we have neglected finite temperature corrections to decay vertices and masses, their inclusion requires studying the time evolution of an initial density matrix. Furthermore, if the decay products thermalize with the medium, their population build-up will lead to Pauli blocking factors thereby suppressing the decay of the parent particle. These effects remain to be studied but are beyond the scope and goals of this article.

Possible implications: The time dependence of the decay function reveals non-equilibrium aspects that have not been previously recognized, not only from the transient build-up of the quasiparticle but also the memory effects that yield the unexpected power laws and stretched exponentials. These novel non-equilibrium effects may lead to interesting and perhaps important dynamics relevant to baryogenesis and leptogenesis. In particular, we envisage corrections to quantum kinetic processes for particle production and their inverse processes. Typically quantum kinetics inputs transition rates perhaps with finite temperature contributions but ultimately obtained from S-matrix theory. Namely, such transition rates are obtained in the infinite time limit and the forward and backward probabilities input strict energy conservation, and as a consequence, obey detailed balance . The richer time dependence of the decay function revealed by this study, with the hitherto unexplored novel non-equilibrium aspects, suggests that similar dynamical processes may enter in a modified quantum kinetic description in the early universe. We expect to report on these and other related issues in future studies.

VIII Summary, conclusions and further questions.

In this article we studied the decay of a bosonic particle into massless fermions via a Yukawa coupling in post-inflation cosmology. The approximation of massless fermions is warranted for a heavy Higgs-like scalar within or beyond the standard model decaying into most charged leptons or quarks (but for the top) of the standard model. We implemented a non-perturbative method that yields the time evolution of the survival probability $\mathcal{P}_{\Phi}(t)$ combined with a physically motivated adiabatic expansion. This expansion is justified when $H(t)/E_{k}(t)\ll 1$ where $H(t)$ is the Hubble rate and $E_{k}(t)$ the local energy of the particle as measured by a comoving observer. We have argued that this approximation is valid for typical particle physics processes during the radiation dominated era and beyond. In a standard cosmology the reliability of this approximation improves with the cosmological expansion, therefore if the adiabatic condition is fulfilled at the initial time when the decaying particle is produced, its reliability improves along the expansion history.

Particle decay into fermionic channels introduces novel phenomena associated with ultraviolet divergences requiring renormalization that result into two different physical processes: i) the build-up of a quasiparticle state out of the bare initial state by dressing with fermion-antifermion pairs, ii) the decay of this quasiparticle state via the emission of fermion pairs. These two different processes occur on widely separated time scales. We introduced a dynamical renormalization method that allows to separate the dynamics of formation of the quasiparticle from its decay on longer time scales. It relies on introducing a renormalization time scale $t_{b}$ to absorb the transient dynamics of formation into the wave function renormalization of the quasiparticle state. The survival probability obeys a dynamical renormalization group equation with respect to $t_{b}$. The decay function of this renormalized state is ultraviolet finite and cutoff independent.

We carried out a detailed analytic and numerical study of the decay function during the radiation dominated era. The dynamics of decay depends crucially on whether the particle is non-relativistic or relativistic. For a particle that is “born” at rest in the comoving frame during (RD) we find that after short time transients, the survival probability is given by

[TABLE]

where $Y$ is the Yukawa coupling and $\Gamma_{0}$ is the decay rate at rest in Minkowski space-time. The scale $t_{b}$ is an intermediate time scale that describes the build-up of the quasiparticle state, and $\mathcal{P}(t_{b})$ is the renormalized probability of such state. The power of $t/t_{b}$ with “anomalous” dimension and the stretched exponential with power $1/4$ are both a remnant of the formation of the quasiparticle on long time scales as a consequence of the cosmological redshift.

For the case in which the decaying particle is “born” ultrarelativistically the time evolution over the whole history during (RD) must be obtained numerically. Different regimes emerge depending on whether the particle is ultrarelativistic for $t\ll t_{nr}$ or non-relativistic for $t\gg t_{nr}$ where $t_{nr}=k^{2}/(2m^{2}_{\phi}H_{0}\sqrt{\Omega_{R}})$ is the time scale at which the decaying particle of mass $m_{\phi}$ that is born ultrarelativistically with comoving momentum $k$ transitions to being non-relativistic as a consequence of the cosmological redshift. During the ultrarelativistic regime ($t\ll t_{nr}$) we find for $t\gg t_{b}$ that the decay function is a stretched exponential

[TABLE]

whereas for $t\gg t_{b}$ and after the particle has become non-relativistic ($t\gg t_{nr}$) we find

[TABLE]

The extra power of $t/t_{nr}$ as compared to the case when the particle is born at rest (see eqn. (VIII.1)) is a consequence of the memory of the decay function on the past history during the ultrarelativistic stage.

The cosmological decay law is compared to a phenomenological Minkowski-like, S-matrix inspired decay law with an instantaneous Lorentz time dilation factor

[TABLE]

we found that this phenomenological law describes at long times a much faster decay thereby under estimating the lifetime of the decaying particle.

The decay dynamics revealed by this study during (RD) allows us to extrapolate to the case of very long lived, i.e. very weakly coupled particles. We obtain a decay function that yields an upper bound to the survival probability all throughout the expansion history under the assumption of two body decay into a massless fermions, it is given by

[TABLE]

where $\Upsilon(z,z_{b})$ is given by (VI.39) and depends only on the cosmological parameters and the scale factor at which the particle transitions from ultrarelativistic to non-relativistic.

One important conclusion from these results is that using a decay rate as measure of the decay dynamics is not a useful concept and misses the correct dynamical evolution. An S-matrix calculation of transition amplitudes or probabilities, where the time interval is taken to infinity not only it does not capture the various different dynamical scales and temporal behaviour of the survival probability, but substantially under estimates the lifetime of the decaying state.

An important corollary of this study is that the S-matrix approach to describe quantum decay in the cosmological setting is in general inadequate, while it may yield a good approximation for processes of decay at rest for weakly coupled particles late in the cosmological history, it misses important non-equilibrium dynamics. The non-equilibrium effects revealed by our study, from the transient dynamics of the formation to the quasiparticle, to the memory of the decay function on the past history of the decaying particle could be relevant in the quantum kinetics of processes in the very early universe. These could have potential impact in CP-violating non-equilibrium dynamics, baryogenesis and leptogenesis and merit further study.

Acknowledgements.

N.H. acknowledges discussions with A. Zentner. D. B. thanks John Collins for illuminating exchanges on aspects related to renormalization and quasiparticle formation.

Appendix A Projectors

Introducing the notation

[TABLE]

with the zeroth-adiabatic order spinors (II.60,II.61), the projector operators are given by equations (II.62). We find

[TABLE]

[TABLE]

where $\mathbb{I}$ is the $2\times 2$ identity matrix. These expressions can be written more compactly introducing the following functions (suppressing the momentum and conformal time arguments),

[TABLE]

as

[TABLE]

Two relevant cases: 1:) Equal time, $\eta=\eta^{\prime}$,

[TABLE]

[TABLE]

2:) Massless fermions,

[TABLE]

Appendix B Minkowski space-time: $m_{\psi}=0$

The Minkowski space-time limit is obtained by replacing $\eta\rightarrow t$ and the frequencies are time independent. The self-energy in this case becomeszenoqft

[TABLE]

Replacing $k^{2}_{0}\rightarrow-d^{2}/dt^{\,{}^{\prime}\,2}$, and introducing a convergence factor $\epsilon\rightarrow 0^{+}$ yields

[TABLE]

and the decay function

[TABLE]

Integrating by parts twice the derivative term in (B.2) we find

[TABLE]

where

[TABLE]

[TABLE]

[TABLE]

where $Si(t)=\int^{t}_{0}\sin(x)/xdx$. Taking $\epsilon\rightarrow 0^{+}$ and keeping $\Lambda$ large but finite yields

[TABLE]

[TABLE]

and

[TABLE]

where for $\Lambda t\rightarrow\infty$ it follows that Si\big{[}(\Lambda-\omega^{\phi}_{k})\,t\big{]}\rightarrow\pi/2, and Ci\big{[}(\Lambda-\omega^{\phi}_{k})\,t\big{]}\rightarrow 0. Taking the limit $\Lambda\rightarrow\infty$ yields

[TABLE]

This is exactly the same result as obtained in ref.zenoqft integrating in $k_{0}$ first.

Appendix C Useful identities:

In this appendix we gather some useful identities valid during the radiation dominated stage (see also appendix (D)).

[TABLE]

The local Lorentz factor in conformal time is given by

[TABLE]

yielding the identity

[TABLE]

The relationship with comoving time $t$ is obtained via eqn. (II.39), namely

[TABLE]

The conformal and comoving time scales $\eta_{nr}\,,t_{nr}$ respectively, determine the scale at which the decaying particle transitions from relativistic with $\gamma(\eta)\gg 1$ for $\eta\ll\eta_{nr}$ or $t\ll t_{nr}$ to non-relativistic with $\gamma(\eta)\simeq 1$ for $\eta\geq\eta_{nr}$ or $t\geq t_{nr}$. In terms of $\eta,\eta_{i}$ we find,

[TABLE]

Appendix D Analysis of $I_{2}$: eqn. (IV.7)

Consider the first term in $I_{2}$ eqn. (IV.7),

[TABLE]

For $\Lambda\gg k,m_{\phi}$ the argument of the cosine becomes simply $\Lambda\,(\eta^{\prime}-\eta_{i})$. Define $x=\Lambda\,(\eta^{\prime}-\eta_{i})$, and change integration variable to $x$, with $x_{f}=\Lambda\,(\eta-\eta_{i})$, yielding

[TABLE]

In the limit $\Lambda\rightarrow\infty$ we find

[TABLE]

where $\gamma_{E}=0.577\cdots$ is Euler’s constant and for $x_{f}\gg 1$ we find $Ci(x_{f})=\sin(x_{f})/x_{f}+\cdots$. We confirmed the result (D.3) numerically. Therefore for $\Lambda\gg k,m_{\phi},1/(\eta-\eta_{i})$ we find

[TABLE]

Let us now consider

[TABLE]

Using the identities obtained in appendix (C) for a radiation dominated cosmology, we write

[TABLE]

where we introduced the definitions

[TABLE]

[TABLE]

and $\gamma_{i}$ is the local Lorentz factor at time $\eta_{i}$.

In terms of these variables we find

[TABLE]

We note that the fulfillment of the adiabatic condition at all times implies that

[TABLE]

For $\xi^{\prime}\ll 1$ it is straightforward to find that $J[\xi^{\prime}]$ features the expansion

[TABLE]

In terms of these variables we find that the subtracted integral $I_{2b}(k,\eta,\eta_{b})$ defined by eqn. (V.4) is given by

[TABLE]

Consider the two contributions to this function,

[TABLE]

During the time scale when $J(\xi^{\prime})\ll 1$ the term $\cos[J(\xi^{\prime})]\simeq 1$ therefore $F_{2}(\xi)\simeq F_{1}(\xi)$ and $I_{2b}\simeq 0$. Figures (6, 7) display $F_{1,2}(\xi)$ for $\omega_{i}\eta_{i}=100$ and $\gamma_{i}=2,10$ respectively for $\xi_{b}=1/\omega_{i}\eta_{i}$. $F_{2}(\xi)$ grows up to a maximum at $\xi_{m}$ at which $J(\xi_{m})=\pi/2$ and begins damped oscillations reaching a plateau. During the rise-time up to the maximum $F_{2}(\xi)\simeq F_{1}(\xi)$ thereby yielding $I_{2b}(k,\eta,\eta_{b})\simeq 0$ during the interval $\xi_{b}\leq\xi\lesssim\xi_{m}$.

Although in general the value of $\xi_{m}$ must be obtained numerically, for $\omega_{i}\eta_{i}\gg 1$ there are two limits that afford an analytic estimate: a:) for $\omega_{i}\eta_{i}\gg 1$ and $\gamma_{i}\simeq 1$, we assume, self consistently that $\xi_{m}\ll 1$, therefore from eqn. (D.11) we obtain

[TABLE]

this expression confirms the assumption that $\xi_{m}\ll 1$ for $\gamma_{i}\simeq 1$. b:) for $\gamma_{i}\gg 1$, it is convenient to carry out the integral (D.9) by expanding $\omega^{\phi}_{k}(\eta_{1})\simeq k+m^{2}_{\phi}\,C^{2}(\eta_{1})/k+\cdots$ and keeping the leading order term, we find

[TABLE]

This latter expression is fairly accurate even for $\gamma_{i}\simeq 2,3$. We have confirmed numerically the validity of these approximate values of the maxima of $F_{2}(\xi)$ (see figures (6,7)). In both cases we find that for $\omega_{i}\eta_{i}\gg 1$ the value at the maxima fulfill $\xi_{m}/\gamma_{i}\ll 1$. In summary, we find that during the time interval $\xi_{b}<\xi<\xi_{m}$ $F_{1}(\xi)\simeq F_{2}(\xi)\simeq 2\ln[\xi/\xi_{b}]$ and $I_{2b}[k,\eta,\eta_{b}]\simeq 0$. For $\xi>\xi_{m}$ the contribution $F_{2}(\xi)\simeq F_{2}(\xi_{m})\simeq 2\ln[\xi_{m}/\xi_{b}]$ remaining constant while $F_{1}(\xi)$ increases monotonically. The above analysis shows that for $\omega_{i}\eta_{i}\gg 1$ it follows that $\xi_{m}\ll\gamma_{i}$ in the whole range of $\gamma_{i}$, therefore during the interval $\xi_{m}<\xi<\gamma_{i}$ and $W[\xi^{\prime}]\simeq 1$, hence

[TABLE]

with $F_{1}(\xi_{m})\simeq F_{2}(\xi_{m})\simeq 2\ln[\xi_{m}/\xi_{b}]$. For $\xi\gg\gamma_{i}$ the function $\sqrt{W[\xi]}+\frac{1}{\sqrt{W[\xi]}}\geq 2$ as shown in fig. (8) hence $F_{1}[\xi]>2\,\ln[\xi]$, with asymptotic behavior

[TABLE]

The behavior of $F_{1,2}$ in the ultrarelativistic case $\gamma_{i}\gg 1$ is summarized as follows:

[TABLE]

Appendix E Analysis of $I_{3b}(k,\eta,\eta_{b})$: eqn. (IV.8)

Let us now consider the following integral in $I_{3b}$, namely the second contribution to eqn. (IV.8):

[TABLE]

this integral is similar to the case of decay into bosonic particles studied in ref.herring .

Following the treatment in this reference we introduce the following definitions

[TABLE]

along with

[TABLE]

In terms of these variables, it follows that

[TABLE]

with

[TABLE]

where there is an implicit $\eta^{\prime}$ dependence in $z$ and $\gamma$.

The argument of the sine function in (E.1) becomes

[TABLE]

with

[TABLE]

where $z\equiv z(\eta^{\prime})\leavevmode\nobreak\ ;\leavevmode\nobreak\ \gamma\equiv\gamma_{k}(\eta^{\prime})$ . Writing

[TABLE]

and using (E.4) it is straightforward to find

[TABLE]

We finally obtain

[TABLE]

where $\widetilde{z}=\omega^{\phi}_{k}(\eta^{\prime})\,(\eta^{\prime}-\eta_{i})$. Gathering together this result with eqn. (IV.13) for $I_{3a}$ we find

[TABLE]

where

[TABLE]

For $\eta^{\prime}\gg\eta_{i}$ and $z(\eta^{\prime})\gg 1$ the integral in (E.10) has an adiabatic expansion, for $\tau\ll z$ we find

[TABLE]

therefore $\delta_{k}$ is of adiabatic order one and higher, furthermore,

[TABLE]

and to leading (zeroth) adiabatic order we can replace $\mathcal{P}[\tau;\eta^{\prime}]=1$. The $\tau$ integral in (E.10) is dominated by the region $\tau\simeq 0$, and the region for which $\tau\simeq z$ yields a contribution $\propto 1/z$, hence of first adiabatic order or smaller. Therefore to leading (zeroth) adiabatic order we neglect $\delta_{k}$ in (E.6) and replace ${P}\rightarrow 1$ in (E.10).

Although the variables (E.2,E.3) allow an explicit identification of the nature of the adiabatic expansion, the most suitable variables to merge the results for $I_{3b}$ with those of the contributions from $I_{2b}$ are those introduced in appendices (C) and (D). We now recast the results for $I_{3b}$ in terms of these variables. Introducing

[TABLE]

in terms of which we find using (C.1)

[TABLE]

Similarly, using (C.3) we obtain

[TABLE]

and the variables $z,\tau$ introduced in eqns. (E.2)) and (E.3) respectively are given by

[TABLE]

which fulfill the identity

[TABLE]

Using these results we find

[TABLE]

and the ratio (E.9) becomes

[TABLE]

and $\delta(\tau,\eta^{\prime})$ in eqn. (E.7) becomes

[TABLE]

Finally the function $A(\tau,\eta^{\prime})$ given by eqn. (E.6), becomes

[TABLE]

with

[TABLE]

and the integral (E.12) becomes

[TABLE]

We have argued above in this appendix that for $\omega^{\phi}_{k}(\eta^{\prime})\,\eta^{\prime}\gg 1$ the term $\delta\equiv\Delta$ is higher order in the adiabatic approximation and can be neglected, and that to leading order in this approximation we can set $P\equiv\mathcal{P}\rightarrow 1$. We now test this assertion numerically in terms of the new variables. Since in the new variables the product $\omega^{\phi}_{k}(\eta^{\prime})\,\eta^{\prime}=\frac{\omega_{i}\,\eta_{i}}{\gamma_{i}}\,y\,f(y)$ it follows that $\omega^{\phi}_{k}(\eta^{\prime})\,\eta^{\prime}\gg 1$ at all times implies that $\omega_{i}\,\eta_{i}\gg 1$ which is precisely the statement of the validity of the adiabatic approximation at the initial time. Fig. (9) compares $\frac{\sin[\mathcal{A}_{0}[y,t]]}{t}$ and $\mathcal{P}[y,t]\,\frac{\sin[\mathcal{A}[y,t]]}{t}$ for $\gamma_{i}=5,\omega_{i}\,\eta_{i}=100,y=10$, confirming the validity of the adiabatic approximation. We have explored a wide range of parameters with similar results.

Therefore to leading order in the adiabatic approximation we can replace the argument of the integral in (E.25) by $\frac{\sin[\mathcal{A}_{0}[y,t]]}{t}$, yielding to lowest adiabatic order

[TABLE]

where $Si[x]$ is the sine-integral function, with $Si[x]\simeq x$ as $x\rightarrow 0$, reaches a maximum at $x=\pi$ and $Si[x]\rightarrow\pi/2$ for $x\gtrsim\pi$. The maximum of $\mathcal{S}(\xi)$ occurs when

[TABLE]

beyond which $\mathcal{S}(\xi)\simeq 1$.

In particular for $\gamma_{i}=1$ (the particle decaying at rest) and $\omega_{i}\eta_{i}\gg 1$, $\mathcal{S}(\xi^{\prime})$ reaches a maximum at $\xi^{\prime}=\xi_{s}\simeq\pi/\omega_{i}\eta_{i}+\mathcal{O}(1/(\omega_{i}\eta_{i})^{2})$ with $\mathcal{S}(\xi_{s})\simeq 1$. In the opposite limit, for an ultrarelativistic particle with $\gamma_{i}\gg 1$ and $\omega_{i}\eta_{i}\gg 1$ we find self-consistently that $S(\xi^{\prime})$ reaches a maximum at $\xi_{s}$ with $\alpha(\xi_{s})=\pi$, where $\xi_{s}$ is a solution of

[TABLE]

For $2\pi\,\gamma^{2}_{i}/\omega_{i}\eta_{i}\ll 1$ we find

[TABLE]

and for $2\pi\,\gamma^{2}_{i}/\omega_{i}\eta_{i}\gg 1$

[TABLE]

In both cases we find that $\frac{\xi_{s}}{\gamma_{i}}\ll 1$ whenever $\gamma_{i}\gg 1$. Fig. (10) displays the behavior of $\mathcal{S}(\xi)$ for $\omega_{i}=100;\gamma_{i}=2,10$.

Using the relations derived in appendix (C) along with the identities $C(\eta^{\prime})=C(\eta_{i})\,(\eta^{\prime}/\eta_{i})=C(\eta_{i})\,(1+\xi^{\prime})$ and $m_{\phi}C(\eta_{i})=\omega_{i}/\gamma_{i}$, it follows that $I_{3}(k,\eta)$ given by (E.11) can be written in terms of the same variables as $I_{2}$, namely $\xi=\eta/\eta_{i}-1$ and $\eta_{b}=\eta_{b}/\eta_{i}-1$. We find

[TABLE]

The contribution from the term with $\mathcal{S}$ in the integrand must be done numerically, however, the first term can be done analytically, yielding

[TABLE]

where we have set $\xi_{b}=0$ to leading adiabatic order.

The function $\mathcal{S}[\alpha(\xi)]$ has the following behavior for $\xi\ll\xi_{s}$ and $\xi\gg\xi_{s}$, corresponding to $\alpha(\xi)\ll\pi$ and $\alpha(\xi)\gg\pi$ respectively:

[TABLE]

[TABLE]